Changeset 4690


Ignore:
Timestamp:
Sep 17, 2023, 1:16:26 PM (8 months ago)
Author:
fhourdin
Message:

Preparation/test de convergence num en temps

Location:
LMDZ6/trunk/libf
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • LMDZ6/trunk/libf/phylmd/calltherm.F90

    r4678 r4690  
    44      subroutine calltherm(dtime  &
    55     &      ,pplay,paprs,pphi,weak_inversion  &
    6      &      ,u_seri_,v_seri_,t_seri_,q_seri_,zqsat,debut  &
     6     &      ,u_seri_,v_seri_,t_seri_,q_seri_,t_env,q_env,zqsat,debut  &
    77     &      ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs  &
    88     &      ,fm_therm,entr_therm,detr_therm,zqasc,clwcon0,lmax,ratqscth,  &
     
    5656
    5757      REAL, DIMENSION(klon,klev), INTENT(IN) :: u_seri_,v_seri_
    58       REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri_,q_seri_
     58      REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri_,q_seri_,t_env,q_env
    5959      REAL, DIMENSION(klon,klev) :: u_seri,v_seri
    6060      REAL, DIMENSION(klon,klev) :: t_seri,q_seri
     
    295295            CALL thermcell_main(itap,klon,klev,zdt  &
    296296     &      ,pplay,paprs,pphi,debut  &
    297      &      ,u_seri,v_seri,t_seri,q_seri  &
     297     &      ,u_seri,v_seri,t_seri,q_seri,t_env,q_env  &
    298298     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
    299299     &      ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax  &
  • LMDZ6/trunk/libf/phylmd/lmdz_thermcell_ini.F90

    r4664 r4690  
    55save
    66
    7    integer :: dvdq=1,dqimpl=-1,prt_level=0,lunout
    8    real    :: RG,RD,RCPD,RKAPPA,RLVTT,RLvCp,RETV
    9    real    :: r_aspect_thermals,tau_thermals,fact_thermals_ed_dz
    10    integer :: iflag_thermals_ed,iflag_thermals_optflux,iflag_thermals_closure
    11    integer :: iflag_thermals_down
    12    real    :: fact_thermals_down
     7
     8integer, protected :: dvdq=1,dqimpl=-1,prt_level=0,lunout
     9real   , protected :: RG,RD,RCPD,RKAPPA,RLVTT,RLvCp,RETV
    1310
    1411!$OMP THREADPRIVATE(dvdq,dqimpl,prt_level,lunout)
    1512!$OMP THREADPRIVATE(RG,RD,RCPD,RKAPPA,RLVTT,RLvCp)
     13
     14
     15! Parameters that can be modified directly by a getin call
     16real,    protected :: r_aspect_thermals=2.       ! Aspect ratio for thermal celles
     17real,    protected :: tau_thermals = 0.          ! relaxation time constant
     18real,    protected :: fact_thermals_ed_dz = 0.1  ! bouyancy computed with a delta
     19real,    protected :: betalpha=0.9               !
     20real,    protected :: afact=2./3.                !
     21real,    protected :: fact_shell=1.              !
     22real,    protected :: detr_min=1.e-5             !
     23real,    protected :: entr_min=1.e-5             !
     24real,    protected :: detr_q_coef=0.012          !
     25real,    protected :: detr_q_power=0.5           !
     26real,    protected :: mix0=0.                    !
     27integer, protected :: iflag_thermals_ed = 0      !
     28integer, protected :: iflag_thermals_optflux = 0 !
     29integer, protected :: iflag_thermals_closure = 1 !
     30integer, protected :: iflag_thermals_down = 0    !
     31integer, protected :: fact_thermals_down = 0.5   !
     32integer, protected :: thermals_flag_alim=0       !
     33integer, protected :: iflag_thermals_tenv=0      !
     34
     35   ! WARNING !!! fact_epsilon is not protected. It can be modified in thermcell_plume*
     36   ! depending on other flags.
     37
     38   real               :: fact_epsilon=0.002
     39
    1640!$OMP THREADPRIVATE(r_aspect_thermals,tau_thermals,fact_thermals_ed_dz)
    1741!$OMP THREADPRIVATE(iflag_thermals_ed,iflag_thermals_optflux,iflag_thermals_closure)
    1842!$OMP THREADPRIVATE(iflag_thermals_down)
    1943!$OMP THREADPRIVATE(fact_thermals_down)
    20 
    21 
    22    REAL, SAVE :: fact_epsilon=0.002
    23    REAL, SAVE :: betalpha=0.9
    24    REAL, SAVE :: afact=2./3.
    25    REAL, SAVE :: fact_shell=1.
    26    REAL,SAVE :: detr_min=1.e-5
    27    REAL,SAVE :: entr_min=1.e-5
    28    REAL,SAVE :: detr_q_coef=0.012
    29    REAL,SAVE :: detr_q_power=0.5
    30    REAL,SAVE :: mix0=0.
    31    INTEGER,SAVE :: thermals_flag_alim=0
    32 
    3344!$OMP THREADPRIVATE(fact_epsilon, betalpha, afact, fact_shell)
    3445!$OMP THREADPRIVATE(detr_min, entr_min, detr_q_coef, detr_q_power)
    3546!$OMP THREADPRIVATE( mix0, thermals_flag_alim)
     47!$OMP THREADPRIVATE(iflag_thermals_tenv)
    3648
    3749
     
    7082!=====================================================================
    7183
    72    r_aspect_thermals=2.
    7384   CALL getin_p('r_aspect_thermals',r_aspect_thermals)
    74    
    75    tau_thermals = 0.
    7685   CALL getin_p('tau_thermals',tau_thermals)
    77    
    78    fact_thermals_ed_dz = 0.1
    7986   CALL getin_p('fact_thermals_ed_dz',fact_thermals_ed_dz)
    80    
    81    fact_thermals_ed_dz = 0.1
    82    CALL getin_p('fact_thermals_ed_dz',fact_thermals_ed_dz)
    83    
    84    iflag_thermals_ed = 0
    8587   CALL getin_p('iflag_thermals_ed',iflag_thermals_ed)
    86    
    87    iflag_thermals_optflux = 0
    8888   CALL getin_p('iflag_thermals_optflux',iflag_thermals_optflux)
    89    
    90    iflag_thermals_closure = 1
    9189   CALL getin_p('iflag_thermals_closure',iflag_thermals_closure)
    92 
    93    iflag_thermals_down = 0
    9490   CALL getin_p('iflag_thermals_down',iflag_thermals_down)
    95 
    96    fact_thermals_down = 0.5
    9791   CALL getin_p('fact_thermals_down',fact_thermals_down)
    98 
    99      CALL getin_p('thermals_fact_epsilon',fact_epsilon)
    100      CALL getin_p('thermals_betalpha',betalpha)
    101      CALL getin_p('thermals_afact',afact)
    102      CALL getin_p('thermals_fact_shell',fact_shell)
    103      CALL getin_p('thermals_detr_min',detr_min)
    104      CALL getin_p('thermals_entr_min',entr_min)
    105      CALL getin_p('thermals_detr_q_coef',detr_q_coef)
    106      CALL getin_p('thermals_detr_q_power',detr_q_power)
    107      CALL getin_p('thermals_mix0',mix0)
    108      CALL getin_p('thermals_flag_alim',thermals_flag_alim)
     92   CALL getin_p('thermals_fact_epsilon',fact_epsilon)
     93   CALL getin_p('thermals_betalpha',betalpha)
     94   CALL getin_p('thermals_afact',afact)
     95   CALL getin_p('thermals_fact_shell',fact_shell)
     96   CALL getin_p('thermals_detr_min',detr_min)
     97   CALL getin_p('thermals_entr_min',entr_min)
     98   CALL getin_p('thermals_detr_q_coef',detr_q_coef)
     99   CALL getin_p('thermals_detr_q_power',detr_q_power)
     100   CALL getin_p('thermals_mix0',mix0)
     101   CALL getin_p('thermals_flag_alim',thermals_flag_alim)
     102   CALL getin_p('iflag_thermals_tenv',iflag_thermals_tenv)
    109103
    110104
     
    139133write(lunout,*) 'thermcell_ini ,mix0                     =',  mix0                   
    140134write(lunout,*) 'thermcell_ini ,thermals_flag_alim       =',  thermals_flag_alim
     135write(lunout,*) 'thermcell_ini ,iflag_thermals_tenv      =',  iflag_thermals_tenv
    141136
    142137 RETURN
  • LMDZ6/trunk/libf/phylmd/lmdz_thermcell_main.F90

    r4684 r4690  
    44! A REGARDER !!!!!!!!!!!!!!!!!
    55! ATTENTION : zpspsk est inout et out mais c'est pas forcement pour de bonnes raisons (FH, 2023)
     6! 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
    67CONTAINS
    78
    89      subroutine thermcell_main(itap,ngrid,nlay,ptimestep  &
    910     &                  ,pplay,pplev,pphi,debut  &
    10      &                  ,puwind,pvwind,ptemp,p_o  &
     11     &                  ,puwind,pvwind,ptemp,p_o,ptemp_env, po_env  &
    1112     &                  ,pduadj,pdvadj,pdtadj,pdoadj  &
    1213     &                  ,fm0,entr0,detr0,zqta,zqla,lmax  &
     
    2122
    2223
     24
     25! USE necessaires pour les lignes importees de thermcell_env
    2326      USE lmdz_thermcell_ini, ONLY: thermcell_ini,dqimpl,dvdq,prt_level,lunout,prt_level
    2427      USE lmdz_thermcell_ini, ONLY: iflag_thermals_closure,iflag_thermals_ed,tau_thermals,r_aspect_thermals
    2528      USE lmdz_thermcell_ini, ONLY: iflag_thermals_down,fact_thermals_down
     29      USE lmdz_thermcell_ini, ONLY: iflag_thermals_tenv
    2630      USE lmdz_thermcell_ini, ONLY: RD,RG
    2731
     
    3640      USE lmdz_thermcell_plume, ONLY: thermcell_plume
    3741      USE lmdz_thermcell_plume_6A, ONLY: thermcell_plume_6A,thermcell_plume_5B
     42
     43! USE necessaires pour les lignes importees de thermcell_env
     44   USE lmdz_thermcell_ini, ONLY : RLvCp,RKAPPA,RETV
     45   USE lmdz_thermcell_qsat, ONLY : thermcell_qsat
     46
    3847
    3948#ifdef ISO
     
    105114      integer, intent(in) :: itap,ngrid,nlay
    106115      real, intent(in) ::  ptimestep
    107       real, intent(in), dimension(ngrid,nlay)    :: ptemp,puwind,pvwind,pplay,pphi
     116      real, intent(in), dimension(ngrid,nlay)    :: ptemp,puwind,pvwind,pplay,pphi,ptemp_env,po_env
    108117! ATTENTION : zpspsk est inout et out mais c'est pas forcement pour de bonnes raisons (FH, 2023)
    109118      real, intent(in), dimension(ngrid,nlay)    :: p_o
     
    144153      integer,dimension(ngrid) :: lmin,lmix,lmix_bis,nivcon
    145154      real, dimension(ngrid,nlay) :: ztva_est
    146       real, dimension(ngrid,nlay) :: deltaz,zlay,zh,zdthladj,zu,zv,z_o,zl,zva,zua,z_oa
     155      real, dimension(ngrid,nlay) :: deltaz,zlay,zdthladj,zu,zv,z_o,zl,zva,zua,z_oa
     156      real, dimension(ngrid,nlay) :: ztemp_env ! temperarure liquide de l'environnement
    147157      real, dimension(ngrid,nlay) :: zta,zha,q2,wq,wthl,wthv,thetath2,wth2
    148158      real, dimension(ngrid,nlay) :: rho,masse
     
    154164      real, dimension(ngrid) :: zcon,zcon2,alim_star_tot,f
    155165      real, dimension(ngrid,nlay) :: entrdn,detrdn
     166      logical, dimension(ngrid,nlay) :: mask
    156167
    157168      character (len=20) :: modname='thermcell_main'
     
    199210       enddo
    200211      endif
     212
    201213!-----------------------------------------------------------------------
    202214! Calcul de T,q,ql a partir de Tl et qT dans l environnement
    203215!   --------------------------------------------------------------------
    204216!
    205       CALL thermcell_env(ngrid,nlay,p_o,ptemp,puwind,pvwind,pplay,  &
    206      &           pplev,z_o,zh,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out)
     217          ! On condense l'eau liquide si besoin.
     218          ! En fait on arrive ici d'habitude (jusque 6A) après réévaporation
     219          ! Dans une nouvelle mouture, on passe les profiles
     220          ! avant la couche limite : iflag_thermals_tenv=1
     221          !     dés le début de la physique : iflag_thermals_tenv=2
     222          ! Mais même pour 2) on ne veut sans doute pas réévaporer.
     223          ! On veut comparer thetav dans le thermique, après condensation,
     224          ! avec le theta_v effectif de l'environnement.
     225
     226      if (iflag_thermals_tenv - 10 * ( iflag_thermals_tenv / 10 ) == 0) then
     227
     228          CALL thermcell_env(ngrid,nlay,p_o,ptemp_env,puwind,pvwind,pplay,  &
     229         &           pplev,z_o,ztemp_env,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out)
     230
     231      else
     232
     233        ! Chantier en cours : ne pas effacer (Fredho). 15 septembre 2023
     234        ! Dans la version originale de thermcell_env, on condense l'eau de l'environnement
     235        ! pour calculer une temperature potentielle liquide.
     236        ! On en déduit un Theta v.
     237
     238 ! ...
     239        ! contenu de thermcell_env
     240        !    SUBROUTINE thermcell_env(ngrid,nlay,po,pt,pu,pv,pplay,  &
     241        ! &           pplev,zo,zh,zl,ztv,zthl,zu,zv,zpspsk,pqsat,lev_out)
     242        ! contenu thermcell_env : call thermcell_qsat(ngrid*nlay,mask,pplev,pt,po,pqsat)
     243        ! contenu thermcell_env : do ll=1,nlay
     244        ! contenu thermcell_env :    do ig=1,ngrid
     245        ! contenu thermcell_env :       zl(ig,ll) = max(0.,po(ig,ll)-pqsat(ig,ll))
     246        ! contenu thermcell_env :       zh(ig,ll) = pt(ig,ll)+RLvCp*zl(ig,ll)         !   T = Tl + Lv/Cp ql
     247        ! contenu thermcell_env :       zo(ig,ll) = po(ig,ll)-zl(ig,ll)
     248        ! contenu thermcell_env :    enddo
     249        ! contenu thermcell_env : enddo
     250        ! contenu thermcell_env : do ll=1,nlay
     251        ! contenu thermcell_env :    do ig=1,ngrid
     252        ! contenu thermcell_env :        zpspsk(ig,ll)=(pplay(ig,ll)/100000.)**RKAPPA
     253        ! contenu thermcell_env :        zu(ig,ll)=pu(ig,ll)
     254        ! contenu thermcell_env :        zv(ig,ll)=pv(ig,ll)
     255        ! contenu thermcell_env :          ztv(ig,ll)=zh(ig,ll)/zpspsk(ig,ll)
     256        ! contenu thermcell_env :          ztv(ig,ll)=ztv(ig,ll)*(1.+RETV*(zo(ig,ll))-zl(ig,ll))
     257        ! contenu thermcell_env :          zthl(ig,ll)=pt(ig,ll)/zpspsk(ig,ll)
     258        ! contenu thermcell_env :    enddo
     259        ! contenu thermcell_env : enddo
     260
     261        do l=1,nlay
     262            do ig=1,ngrid
     263                 zl(ig,l)=0.
     264                 zu(ig,l)=puwind(ig,l)
     265                 zv(ig,l)=pvwind(ig,l)
     266                 ztemp_env(ig,l)=ptemp_env(ig,l)
     267                 zpspsk(ig,l)=(pplay(ig,l)/100000.)**RKAPPA
     268                 ztv(ig,l)=ztemp_env(ig,l)/zpspsk(ig,l)
     269                 ztv(ig,l)=ztv(ig,l)*(1.+RETV*po_env(ig,l))
     270                 zthl(ig,l)=ptemp(ig,l)/zpspsk(ig,l)
     271                 mask(ig,l)=.true.
     272            enddo
     273        enddo
     274        call thermcell_qsat(ngrid*nlay,mask,pplev,ptemp_env,p_o,zqsat)
     275         
     276      endif
    207277       
    208278      if (prt_level.ge.1) print*,'thermcell_main apres thermcell_env'
     
    532602      !--------------------------------------------------------------
    533603
     604        ! Temperature potentielle liquide effectivement mélangée par les thermiques
     605        do ll=1,nlay
     606           do ig=1,ngrid
     607              zthl(ig,ll)=ptemp(ig,ll)/zpspsk(ig,ll)
     608           enddo
     609        enddo
    534610        call thermcell_dq(ngrid,nlay,dqimpl,ptimestep,fm0,entr0,masse,  &
    535611        &                    zthl,zdthladj,zta,lev_out)
     
    632708!nouveau calcul
    633709      do ig=1,ngrid
    634       CHI=zh(ig,1)/(1669.0-122.0*z_o(ig,1)/zqsat(ig,1)-zh(ig,1))
     710      ! WARNING !!! verifier que c'est bien ztemp_env qu'on veut là
     711      CHI=ztemp_env(ig,1)/(1669.0-122.0*z_o(ig,1)/zqsat(ig,1)-ztemp_env(ig,1))
    635712      pcon(ig)=pplay(ig,1)*(z_o(ig,1)/zqsat(ig,1))**CHI
    636713      enddo
  • LMDZ6/trunk/libf/phylmd/physiq_mod.F90

    r4687 r4690  
    8484    USE yamada_ini_mod, ONLY : yamada_ini
    8585    USE lmdz_atke_turbulence_ini, ONLY : atke_ini
    86     USE lmdz_thermcell_ini, ONLY : thermcell_ini
     86    USE lmdz_thermcell_ini, ONLY : thermcell_ini, iflag_thermals_tenv
    8787    USE lmdz_thermcell_dtke, ONLY : thermcell_dtke
    8888    USE blowing_snow_ini_mod, ONLY : blowing_snow_ini , qbst_bs
     
    12521252    LOGICAL, PARAMETER :: mass_fixer=.FALSE.
    12531253    REAL qql1(klon),qql2(klon),corrqql
     1254
     1255    REAL, dimension(klon,klev) :: t_env,q_env
    12541256
    12551257    REAL pi
     
    27372739    ! Calcul de l'humidite de saturation au niveau du sol
    27382740
     2741! Tests Fredho, instensibilite au pas de temps -------------------------------
     2742! A detruire en 2024 une fois les tests documentes et les choix faits        !
     2743! Conservation des variables avant l'appel à l a diffusion pour les tehrmic  !
     2744    if (iflag_thermals_tenv / 10 == 1 ) then                                 !
     2745        do k=1,klev                                                          !
     2746           do i=1,klon                                                       !
     2747              t_env(i,k)=t_seri(i,k)                                         !
     2748              q_env(i,k)=q_seri(i,k)                                         !
     2749           enddo                                                             !
     2750        enddo                                                                !
     2751    else if (iflag_thermals_tenv / 10 == 2 ) then                            !
     2752        do k=1,klev                                                          !
     2753           do i=1,klon                                                       !
     2754              t_env(i,k)=t_seri(i,k)                                         !
     2755           enddo                                                             !
     2756        enddo                                                                !
     2757    endif                                                                    !
     2758! Tests Fredho, instensibilite au pas de temps -------------------------------
    27392759
    27402760
     
    35703590
    35713591       IF (iflag_thermals>=1) THEN
     3592
     3593! Tests Fredho, instensibilite au pas de temps -------------------------------
     3594! A detruire en 2024 une fois les tests documentes et les choix faits        !
     3595          if (iflag_thermals_tenv /10 == 0 ) then                            !
     3596            do k=1,klev                                                      !
     3597               do i=1,klon                                                   !
     3598                  t_env(i,k)=t_seri(i,k)                                     !
     3599                  q_env(i,k)=q_seri(i,k)                                     !
     3600               enddo                                                         !
     3601            enddo                                                            !
     3602          else if (iflag_thermals_tenv / 10 == 2 ) then                      !
     3603            do k=1,klev                                                      !
     3604               do i=1,klon                                                   !
     3605                  q_env(i,k)=q_seri(i,k)                                     !
     3606               enddo                                                         !
     3607            enddo                                                            !
     3608          else if (iflag_thermals_tenv / 10 == 3 ) then                      !
     3609            do k=1,klev                                                      !
     3610               do i=1,klon                                                   !
     3611                  t_env(i,k)=t(i,k)                                          !
     3612                  q_env(i,k)=qx(i,k,1)                                       !
     3613               enddo                                                         !
     3614            enddo                                                            !
     3615          endif                                                              !
     3616! Tests Fredho, instensibilite au pas de temps ------------------------------
     3617
    35723618          !jyg<
    35733619!!       IF (mod(iflag_pbl_split/2,2) .EQ. 1) THEN
    3574        IF (mod(iflag_pbl_split/10,10) .GE. 1) THEN
     3620          IF (mod(iflag_pbl_split/10,10) .GE. 1) THEN
    35753621             !  Appel des thermiques avec les profils exterieurs aux poches
    35763622             DO k=1,klev
     
    35783624                   t_therm(i,k) = t_seri(i,k) - wake_s(i)*wake_deltat(i,k)
    35793625                   q_therm(i,k) = q_seri(i,k) - wake_s(i)*wake_deltaq(i,k)
     3626                   t_env(i,k)   = t_env(i,k) - wake_s(i)*wake_deltat(i,k)
     3627                   q_env(i,k)   = q_env(i,k) - wake_s(i)*wake_deltaq(i,k)
    35803628                   u_therm(i,k) = u_seri(i,k)
    35813629                   v_therm(i,k) = v_seri(i,k)
     
    35973645               ,pplay,paprs,pphi,weak_inversion &
    35983646                        ! ,u_seri,v_seri,t_seri,q_seri,zqsat,debut & !jyg
    3599                ,u_therm,v_therm,t_therm,q_therm,zqsat,debut &  !jyg
     3647               ,u_therm,v_therm,t_therm,q_therm,t_env,q_env,zqsat,debut &  !jyg
    36003648               ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs &
    36013649               ,fm_therm,entr_therm,detr_therm &
  • LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90

    r4687 r4690  
    8383    USE yamada_ini_mod, ONLY : yamada_ini
    8484    USE lmdz_atke_turbulence_ini, ONLY : atke_ini
    85     USE lmdz_thermcell_ini, ONLY : thermcell_ini
     85    USE lmdz_thermcell_ini, ONLY : thermcell_ini, iflag_thermals_tenv
    8686    USE lmdz_thermcell_dtke, ONLY : thermcell_dtke
    8787    USE blowing_snow_ini_mod, ONLY : blowing_snow_ini , qbst_bs
     
    44654465#endif
    44664466          !>jyg
     4467          ! WARNING !!! introduction d'une subtilite avec la distinction
     4468          ! entre deux types de profiles pour les thermiques.
     4469          ! Ne change pas les resultats avec les isotopes tant qu'on
     4470          ! tourne avec le defaut : iflag_thermals_tenv=0
     4471          if ( iflag_thermals_tenv == 0 ) then
     4472              abort_message ='iflag_thermals_env/=0 : non prevu avec les isotopes'
     4473              CALL abort_physic (modname,abort_message,1)
     4474          endif
     4475
    44674476          CALL calltherm(pdtphys &
    44684477               ,pplay,paprs,pphi,weak_inversion &
    44694478                        ! ,u_seri,v_seri,t_seri,q_seri,zqsat,debut & !jyg
    4470                ,u_therm,v_therm,t_therm,q_therm,zqsat,debut &  !jyg
     4479               ,u_therm,v_therm,t_therm,q_therm,t_therm,q_therm,zqsat,debut &  !jyg
    44714480               ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs &
    44724481               ,fm_therm,entr_therm,detr_therm &
Note: See TracChangeset for help on using the changeset viewer.