Index: LMDZ6/trunk/libf/phylmd/calltherm.F90 =================================================================== --- LMDZ6/trunk/libf/phylmd/calltherm.F90 (revision 4689) +++ LMDZ6/trunk/libf/phylmd/calltherm.F90 (revision 4690) @@ -4,5 +4,5 @@ subroutine calltherm(dtime & & ,pplay,paprs,pphi,weak_inversion & - & ,u_seri_,v_seri_,t_seri_,q_seri_,zqsat,debut & + & ,u_seri_,v_seri_,t_seri_,q_seri_,t_env,q_env,zqsat,debut & & ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs & & ,fm_therm,entr_therm,detr_therm,zqasc,clwcon0,lmax,ratqscth, & @@ -56,5 +56,5 @@ REAL, DIMENSION(klon,klev), INTENT(IN) :: u_seri_,v_seri_ - REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri_,q_seri_ + REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri_,q_seri_,t_env,q_env REAL, DIMENSION(klon,klev) :: u_seri,v_seri REAL, DIMENSION(klon,klev) :: t_seri,q_seri @@ -295,5 +295,5 @@ CALL thermcell_main(itap,klon,klev,zdt & & ,pplay,paprs,pphi,debut & - & ,u_seri,v_seri,t_seri,q_seri & + & ,u_seri,v_seri,t_seri,q_seri,t_env,q_env & & ,d_u_the,d_v_the,d_t_the,d_q_the & & ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax & Index: LMDZ6/trunk/libf/phylmd/lmdz_thermcell_ini.F90 =================================================================== --- LMDZ6/trunk/libf/phylmd/lmdz_thermcell_ini.F90 (revision 4689) +++ LMDZ6/trunk/libf/phylmd/lmdz_thermcell_ini.F90 (revision 4690) @@ -5,33 +5,45 @@ save - integer :: dvdq=1,dqimpl=-1,prt_level=0,lunout - real :: RG,RD,RCPD,RKAPPA,RLVTT,RLvCp,RETV - real :: r_aspect_thermals,tau_thermals,fact_thermals_ed_dz - integer :: iflag_thermals_ed,iflag_thermals_optflux,iflag_thermals_closure - integer :: iflag_thermals_down - real :: fact_thermals_down + +integer, protected :: dvdq=1,dqimpl=-1,prt_level=0,lunout +real , protected :: RG,RD,RCPD,RKAPPA,RLVTT,RLvCp,RETV !$OMP THREADPRIVATE(dvdq,dqimpl,prt_level,lunout) !$OMP THREADPRIVATE(RG,RD,RCPD,RKAPPA,RLVTT,RLvCp) + + +! Parameters that can be modified directly by a getin call +real, protected :: r_aspect_thermals=2. ! Aspect ratio for thermal celles +real, protected :: tau_thermals = 0. ! relaxation time constant +real, protected :: fact_thermals_ed_dz = 0.1 ! bouyancy computed with a delta +real, protected :: betalpha=0.9 ! +real, protected :: afact=2./3. ! +real, protected :: fact_shell=1. ! +real, protected :: detr_min=1.e-5 ! +real, protected :: entr_min=1.e-5 ! +real, protected :: detr_q_coef=0.012 ! +real, protected :: detr_q_power=0.5 ! +real, protected :: mix0=0. ! +integer, protected :: iflag_thermals_ed = 0 ! +integer, protected :: iflag_thermals_optflux = 0 ! +integer, protected :: iflag_thermals_closure = 1 ! +integer, protected :: iflag_thermals_down = 0 ! +integer, protected :: fact_thermals_down = 0.5 ! +integer, protected :: thermals_flag_alim=0 ! +integer, protected :: iflag_thermals_tenv=0 ! + + ! WARNING !!! fact_epsilon is not protected. It can be modified in thermcell_plume* + ! depending on other flags. + + real :: fact_epsilon=0.002 + !$OMP THREADPRIVATE(r_aspect_thermals,tau_thermals,fact_thermals_ed_dz) !$OMP THREADPRIVATE(iflag_thermals_ed,iflag_thermals_optflux,iflag_thermals_closure) !$OMP THREADPRIVATE(iflag_thermals_down) !$OMP THREADPRIVATE(fact_thermals_down) - - - REAL, SAVE :: fact_epsilon=0.002 - REAL, SAVE :: betalpha=0.9 - REAL, SAVE :: afact=2./3. - REAL, SAVE :: fact_shell=1. - REAL,SAVE :: detr_min=1.e-5 - REAL,SAVE :: entr_min=1.e-5 - REAL,SAVE :: detr_q_coef=0.012 - REAL,SAVE :: detr_q_power=0.5 - REAL,SAVE :: mix0=0. - INTEGER,SAVE :: thermals_flag_alim=0 - !$OMP THREADPRIVATE(fact_epsilon, betalpha, afact, fact_shell) !$OMP THREADPRIVATE(detr_min, entr_min, detr_q_coef, detr_q_power) !$OMP THREADPRIVATE( mix0, thermals_flag_alim) +!$OMP THREADPRIVATE(iflag_thermals_tenv) @@ -70,41 +82,23 @@ !===================================================================== - r_aspect_thermals=2. CALL getin_p('r_aspect_thermals',r_aspect_thermals) - - tau_thermals = 0. CALL getin_p('tau_thermals',tau_thermals) - - fact_thermals_ed_dz = 0.1 CALL getin_p('fact_thermals_ed_dz',fact_thermals_ed_dz) - - fact_thermals_ed_dz = 0.1 - CALL getin_p('fact_thermals_ed_dz',fact_thermals_ed_dz) - - iflag_thermals_ed = 0 CALL getin_p('iflag_thermals_ed',iflag_thermals_ed) - - iflag_thermals_optflux = 0 CALL getin_p('iflag_thermals_optflux',iflag_thermals_optflux) - - iflag_thermals_closure = 1 CALL getin_p('iflag_thermals_closure',iflag_thermals_closure) - - iflag_thermals_down = 0 CALL getin_p('iflag_thermals_down',iflag_thermals_down) - - fact_thermals_down = 0.5 CALL getin_p('fact_thermals_down',fact_thermals_down) - - CALL getin_p('thermals_fact_epsilon',fact_epsilon) - CALL getin_p('thermals_betalpha',betalpha) - CALL getin_p('thermals_afact',afact) - CALL getin_p('thermals_fact_shell',fact_shell) - CALL getin_p('thermals_detr_min',detr_min) - CALL getin_p('thermals_entr_min',entr_min) - CALL getin_p('thermals_detr_q_coef',detr_q_coef) - CALL getin_p('thermals_detr_q_power',detr_q_power) - CALL getin_p('thermals_mix0',mix0) - CALL getin_p('thermals_flag_alim',thermals_flag_alim) + CALL getin_p('thermals_fact_epsilon',fact_epsilon) + CALL getin_p('thermals_betalpha',betalpha) + CALL getin_p('thermals_afact',afact) + CALL getin_p('thermals_fact_shell',fact_shell) + CALL getin_p('thermals_detr_min',detr_min) + CALL getin_p('thermals_entr_min',entr_min) + CALL getin_p('thermals_detr_q_coef',detr_q_coef) + CALL getin_p('thermals_detr_q_power',detr_q_power) + CALL getin_p('thermals_mix0',mix0) + CALL getin_p('thermals_flag_alim',thermals_flag_alim) + CALL getin_p('iflag_thermals_tenv',iflag_thermals_tenv) @@ -139,4 +133,5 @@ write(lunout,*) 'thermcell_ini ,mix0 =', mix0 write(lunout,*) 'thermcell_ini ,thermals_flag_alim =', thermals_flag_alim +write(lunout,*) 'thermcell_ini ,iflag_thermals_tenv =', iflag_thermals_tenv RETURN Index: LMDZ6/trunk/libf/phylmd/lmdz_thermcell_main.F90 =================================================================== --- LMDZ6/trunk/libf/phylmd/lmdz_thermcell_main.F90 (revision 4689) +++ LMDZ6/trunk/libf/phylmd/lmdz_thermcell_main.F90 (revision 4690) @@ -4,9 +4,10 @@ ! A REGARDER !!!!!!!!!!!!!!!!! ! ATTENTION : zpspsk est inout et out mais c'est pas forcement pour de bonnes raisons (FH, 2023) +! ATTENTION : dans thermcell_env, on condense potentiellement de l'eau. Mais comme on ne mélange pas l'eau liquide supposant qu'il n'y en n'a pas, c'est potentiellement un souci CONTAINS subroutine thermcell_main(itap,ngrid,nlay,ptimestep & & ,pplay,pplev,pphi,debut & - & ,puwind,pvwind,ptemp,p_o & + & ,puwind,pvwind,ptemp,p_o,ptemp_env, po_env & & ,pduadj,pdvadj,pdtadj,pdoadj & & ,fm0,entr0,detr0,zqta,zqla,lmax & @@ -21,7 +22,10 @@ + +! USE necessaires pour les lignes importees de thermcell_env USE lmdz_thermcell_ini, ONLY: thermcell_ini,dqimpl,dvdq,prt_level,lunout,prt_level USE lmdz_thermcell_ini, ONLY: iflag_thermals_closure,iflag_thermals_ed,tau_thermals,r_aspect_thermals USE lmdz_thermcell_ini, ONLY: iflag_thermals_down,fact_thermals_down + USE lmdz_thermcell_ini, ONLY: iflag_thermals_tenv USE lmdz_thermcell_ini, ONLY: RD,RG @@ -36,4 +40,9 @@ USE lmdz_thermcell_plume, ONLY: thermcell_plume USE lmdz_thermcell_plume_6A, ONLY: thermcell_plume_6A,thermcell_plume_5B + +! USE necessaires pour les lignes importees de thermcell_env + USE lmdz_thermcell_ini, ONLY : RLvCp,RKAPPA,RETV + USE lmdz_thermcell_qsat, ONLY : thermcell_qsat + #ifdef ISO @@ -105,5 +114,5 @@ integer, intent(in) :: itap,ngrid,nlay real, intent(in) :: ptimestep - real, intent(in), dimension(ngrid,nlay) :: ptemp,puwind,pvwind,pplay,pphi + real, intent(in), dimension(ngrid,nlay) :: ptemp,puwind,pvwind,pplay,pphi,ptemp_env,po_env ! ATTENTION : zpspsk est inout et out mais c'est pas forcement pour de bonnes raisons (FH, 2023) real, intent(in), dimension(ngrid,nlay) :: p_o @@ -144,5 +153,6 @@ integer,dimension(ngrid) :: lmin,lmix,lmix_bis,nivcon real, dimension(ngrid,nlay) :: ztva_est - real, dimension(ngrid,nlay) :: deltaz,zlay,zh,zdthladj,zu,zv,z_o,zl,zva,zua,z_oa + real, dimension(ngrid,nlay) :: deltaz,zlay,zdthladj,zu,zv,z_o,zl,zva,zua,z_oa + real, dimension(ngrid,nlay) :: ztemp_env ! temperarure liquide de l'environnement real, dimension(ngrid,nlay) :: zta,zha,q2,wq,wthl,wthv,thetath2,wth2 real, dimension(ngrid,nlay) :: rho,masse @@ -154,4 +164,5 @@ real, dimension(ngrid) :: zcon,zcon2,alim_star_tot,f real, dimension(ngrid,nlay) :: entrdn,detrdn + logical, dimension(ngrid,nlay) :: mask character (len=20) :: modname='thermcell_main' @@ -199,10 +210,69 @@ enddo endif + !----------------------------------------------------------------------- ! Calcul de T,q,ql a partir de Tl et qT dans l environnement ! -------------------------------------------------------------------- ! - CALL thermcell_env(ngrid,nlay,p_o,ptemp,puwind,pvwind,pplay, & - & pplev,z_o,zh,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out) + ! On condense l'eau liquide si besoin. + ! En fait on arrive ici d'habitude (jusque 6A) après réévaporation + ! Dans une nouvelle mouture, on passe les profiles + ! avant la couche limite : iflag_thermals_tenv=1 + ! dés le début de la physique : iflag_thermals_tenv=2 + ! Mais même pour 2) on ne veut sans doute pas réévaporer. + ! On veut comparer thetav dans le thermique, après condensation, + ! avec le theta_v effectif de l'environnement. + + if (iflag_thermals_tenv - 10 * ( iflag_thermals_tenv / 10 ) == 0) then + + CALL thermcell_env(ngrid,nlay,p_o,ptemp_env,puwind,pvwind,pplay, & + & pplev,z_o,ztemp_env,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out) + + else + + ! Chantier en cours : ne pas effacer (Fredho). 15 septembre 2023 + ! Dans la version originale de thermcell_env, on condense l'eau de l'environnement + ! pour calculer une temperature potentielle liquide. + ! On en déduit un Theta v. + + ! ... + ! contenu de thermcell_env + ! SUBROUTINE thermcell_env(ngrid,nlay,po,pt,pu,pv,pplay, & + ! & pplev,zo,zh,zl,ztv,zthl,zu,zv,zpspsk,pqsat,lev_out) + ! contenu thermcell_env : call thermcell_qsat(ngrid*nlay,mask,pplev,pt,po,pqsat) + ! contenu thermcell_env : do ll=1,nlay + ! contenu thermcell_env : do ig=1,ngrid + ! contenu thermcell_env : zl(ig,ll) = max(0.,po(ig,ll)-pqsat(ig,ll)) + ! contenu thermcell_env : zh(ig,ll) = pt(ig,ll)+RLvCp*zl(ig,ll) ! T = Tl + Lv/Cp ql + ! contenu thermcell_env : zo(ig,ll) = po(ig,ll)-zl(ig,ll) + ! contenu thermcell_env : enddo + ! contenu thermcell_env : enddo + ! contenu thermcell_env : do ll=1,nlay + ! contenu thermcell_env : do ig=1,ngrid + ! contenu thermcell_env : zpspsk(ig,ll)=(pplay(ig,ll)/100000.)**RKAPPA + ! contenu thermcell_env : zu(ig,ll)=pu(ig,ll) + ! contenu thermcell_env : zv(ig,ll)=pv(ig,ll) + ! contenu thermcell_env : ztv(ig,ll)=zh(ig,ll)/zpspsk(ig,ll) + ! contenu thermcell_env : ztv(ig,ll)=ztv(ig,ll)*(1.+RETV*(zo(ig,ll))-zl(ig,ll)) + ! contenu thermcell_env : zthl(ig,ll)=pt(ig,ll)/zpspsk(ig,ll) + ! contenu thermcell_env : enddo + ! contenu thermcell_env : enddo + + do l=1,nlay + do ig=1,ngrid + zl(ig,l)=0. + zu(ig,l)=puwind(ig,l) + zv(ig,l)=pvwind(ig,l) + ztemp_env(ig,l)=ptemp_env(ig,l) + zpspsk(ig,l)=(pplay(ig,l)/100000.)**RKAPPA + ztv(ig,l)=ztemp_env(ig,l)/zpspsk(ig,l) + ztv(ig,l)=ztv(ig,l)*(1.+RETV*po_env(ig,l)) + zthl(ig,l)=ptemp(ig,l)/zpspsk(ig,l) + mask(ig,l)=.true. + enddo + enddo + call thermcell_qsat(ngrid*nlay,mask,pplev,ptemp_env,p_o,zqsat) + + endif if (prt_level.ge.1) print*,'thermcell_main apres thermcell_env' @@ -532,4 +602,10 @@ !-------------------------------------------------------------- + ! Temperature potentielle liquide effectivement mélangée par les thermiques + do ll=1,nlay + do ig=1,ngrid + zthl(ig,ll)=ptemp(ig,ll)/zpspsk(ig,ll) + enddo + enddo call thermcell_dq(ngrid,nlay,dqimpl,ptimestep,fm0,entr0,masse, & & zthl,zdthladj,zta,lev_out) @@ -632,5 +708,6 @@ !nouveau calcul do ig=1,ngrid - CHI=zh(ig,1)/(1669.0-122.0*z_o(ig,1)/zqsat(ig,1)-zh(ig,1)) + ! WARNING !!! verifier que c'est bien ztemp_env qu'on veut là + CHI=ztemp_env(ig,1)/(1669.0-122.0*z_o(ig,1)/zqsat(ig,1)-ztemp_env(ig,1)) pcon(ig)=pplay(ig,1)*(z_o(ig,1)/zqsat(ig,1))**CHI enddo Index: LMDZ6/trunk/libf/phylmd/physiq_mod.F90 =================================================================== --- LMDZ6/trunk/libf/phylmd/physiq_mod.F90 (revision 4689) +++ LMDZ6/trunk/libf/phylmd/physiq_mod.F90 (revision 4690) @@ -84,5 +84,5 @@ USE yamada_ini_mod, ONLY : yamada_ini USE lmdz_atke_turbulence_ini, ONLY : atke_ini - USE lmdz_thermcell_ini, ONLY : thermcell_ini + USE lmdz_thermcell_ini, ONLY : thermcell_ini, iflag_thermals_tenv USE lmdz_thermcell_dtke, ONLY : thermcell_dtke USE blowing_snow_ini_mod, ONLY : blowing_snow_ini , qbst_bs @@ -1252,4 +1252,6 @@ LOGICAL, PARAMETER :: mass_fixer=.FALSE. REAL qql1(klon),qql2(klon),corrqql + + REAL, dimension(klon,klev) :: t_env,q_env REAL pi @@ -2737,4 +2739,22 @@ ! Calcul de l'humidite de saturation au niveau du sol +! Tests Fredho, instensibilite au pas de temps ------------------------------- +! A detruire en 2024 une fois les tests documentes et les choix faits ! +! Conservation des variables avant l'appel à l a diffusion pour les tehrmic ! + if (iflag_thermals_tenv / 10 == 1 ) then ! + do k=1,klev ! + do i=1,klon ! + t_env(i,k)=t_seri(i,k) ! + q_env(i,k)=q_seri(i,k) ! + enddo ! + enddo ! + else if (iflag_thermals_tenv / 10 == 2 ) then ! + do k=1,klev ! + do i=1,klon ! + t_env(i,k)=t_seri(i,k) ! + enddo ! + enddo ! + endif ! +! Tests Fredho, instensibilite au pas de temps ------------------------------- @@ -3570,7 +3590,33 @@ IF (iflag_thermals>=1) THEN + +! Tests Fredho, instensibilite au pas de temps ------------------------------- +! A detruire en 2024 une fois les tests documentes et les choix faits ! + if (iflag_thermals_tenv /10 == 0 ) then ! + do k=1,klev ! + do i=1,klon ! + t_env(i,k)=t_seri(i,k) ! + q_env(i,k)=q_seri(i,k) ! + enddo ! + enddo ! + else if (iflag_thermals_tenv / 10 == 2 ) then ! + do k=1,klev ! + do i=1,klon ! + q_env(i,k)=q_seri(i,k) ! + enddo ! + enddo ! + else if (iflag_thermals_tenv / 10 == 3 ) then ! + do k=1,klev ! + do i=1,klon ! + t_env(i,k)=t(i,k) ! + q_env(i,k)=qx(i,k,1) ! + enddo ! + enddo ! + endif ! +! Tests Fredho, instensibilite au pas de temps ------------------------------ + !jyg< !! IF (mod(iflag_pbl_split/2,2) .EQ. 1) THEN - IF (mod(iflag_pbl_split/10,10) .GE. 1) THEN + IF (mod(iflag_pbl_split/10,10) .GE. 1) THEN ! Appel des thermiques avec les profils exterieurs aux poches DO k=1,klev @@ -3578,4 +3624,6 @@ t_therm(i,k) = t_seri(i,k) - wake_s(i)*wake_deltat(i,k) q_therm(i,k) = q_seri(i,k) - wake_s(i)*wake_deltaq(i,k) + t_env(i,k) = t_env(i,k) - wake_s(i)*wake_deltat(i,k) + q_env(i,k) = q_env(i,k) - wake_s(i)*wake_deltaq(i,k) u_therm(i,k) = u_seri(i,k) v_therm(i,k) = v_seri(i,k) @@ -3597,5 +3645,5 @@ ,pplay,paprs,pphi,weak_inversion & ! ,u_seri,v_seri,t_seri,q_seri,zqsat,debut & !jyg - ,u_therm,v_therm,t_therm,q_therm,zqsat,debut & !jyg + ,u_therm,v_therm,t_therm,q_therm,t_env,q_env,zqsat,debut & !jyg ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs & ,fm_therm,entr_therm,detr_therm & Index: LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90 =================================================================== --- LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90 (revision 4689) +++ LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90 (revision 4690) @@ -83,5 +83,5 @@ USE yamada_ini_mod, ONLY : yamada_ini USE lmdz_atke_turbulence_ini, ONLY : atke_ini - USE lmdz_thermcell_ini, ONLY : thermcell_ini + USE lmdz_thermcell_ini, ONLY : thermcell_ini, iflag_thermals_tenv USE lmdz_thermcell_dtke, ONLY : thermcell_dtke USE blowing_snow_ini_mod, ONLY : blowing_snow_ini , qbst_bs @@ -4465,8 +4465,17 @@ #endif !>jyg + ! WARNING !!! introduction d'une subtilite avec la distinction + ! entre deux types de profiles pour les thermiques. + ! Ne change pas les resultats avec les isotopes tant qu'on + ! tourne avec le defaut : iflag_thermals_tenv=0 + if ( iflag_thermals_tenv == 0 ) then + abort_message ='iflag_thermals_env/=0 : non prevu avec les isotopes' + CALL abort_physic (modname,abort_message,1) + endif + CALL calltherm(pdtphys & ,pplay,paprs,pphi,weak_inversion & ! ,u_seri,v_seri,t_seri,q_seri,zqsat,debut & !jyg - ,u_therm,v_therm,t_therm,q_therm,zqsat,debut & !jyg + ,u_therm,v_therm,t_therm,q_therm,t_therm,q_therm,zqsat,debut & !jyg ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs & ,fm_therm,entr_therm,detr_therm &