diff --git a/src/common/turb/modd_turbn.F90 b/src/common/turb/modd_turbn.F90 index 867228d2721c5bcb82b64ce6a5440eb332c22cb2..b2a7955e916475ec92b980edb275ceefbccd3481 100644 --- a/src/common/turb/modd_turbn.F90 +++ b/src/common/turb/modd_turbn.F90 @@ -41,6 +41,7 @@ !! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O !! D. Ricard May 2021 add the switches for Leonard terms !! JL Redelsperger 03/2021 Add O-A flux for auto-coupled LES case +!! S. Riette June 2023: add LSMOOTH_PRANDTL, XMINSIGS and XBL89EXP/XUSRBL89 !! !------------------------------------------------------------------------------- ! @@ -105,6 +106,10 @@ TYPE TURB_t LOGICAL :: LROTATE_WIND !< .TRUE. to rotate wind components LOGICAL :: LTKEMINTURB !< set a minimum value for the TKE in the turbulence scheme LOGICAL :: LPROJQITURB !< project the rt tendency on rc/ri + LOGICAL :: LSMOOTH_PRANDTL !< .TRUE. to smooth prandtl functions + REAL :: XMINSIGS !< minimum value for SIGS computed by the turbulence scheme + REAL :: XBL89EXP, XUSRBL89 !< exponent on final BL89 length + LOGICAL :: LNEWBL89EXP !< .TRUE. to use the bug correction aound the expononent in BL89 length computation ! END TYPE TURB_t @@ -152,14 +157,18 @@ REAL, POINTER :: XLINI=>NULL() LOGICAL, POINTER :: LROTATE_WIND=>NULL() LOGICAL, POINTER :: LTKEMINTURB=>NULL() LOGICAL, POINTER :: LPROJQITURB=>NULL() - +LOGICAL, POINTER :: LSMOOTH_PRANDTL=>NULL() +REAL, POINTER :: XMINSIGS=>NULL() +REAL, POINTER :: XBL89EXP=>NULL(), XUSRBL89=>NULL() +LOGICAL, POINTER :: LNEWBL89EXP=>NULL() ! NAMELIST/NAM_TURBn/XIMPL,CTURBLEN,CTURBDIM,LTURB_FLX,LTURB_DIAG, & LSIG_CONV,LRMC01,CTOM,& XTKEMIN,XCED,XCTP,XCADAP,& LLEONARD,XCOEFHGRADTHL, XCOEFHGRADRM, & XALTHGRAD, XCLDTHOLD, XLINI, LHARAT, & - LPROJQITURB + LPROJQITURB, LSMOOTH_PRANDTL, XMINSIGS, & + LNEWBL89EXP ! !------------------------------------------------------------------------------- ! @@ -235,6 +244,11 @@ XLINI=>TURB_MODEL(KTO)%XLINI LROTATE_WIND=>TURB_MODEL(KTO)%LROTATE_WIND LTKEMINTURB=>TURB_MODEL(KTO)%LTKEMINTURB LPROJQITURB=>TURB_MODEL(KTO)%LPROJQITURB +LSMOOTH_PRANDTL=>TURB_MODEL(KTO)%LSMOOTH_PRANDTL +XMINSIGS=>TURB_MODEL(KTO)%XMINSIGS +XBL89EXP=>TURB_MODEL(KTO)%XBL89EXP +XUSRBL89=>TURB_MODEL(KTO)%XUSRBL89 +LNEWBL89EXP=>TURB_MODEL(KTO)%LNEWBL89EXP ! ENDIF ! @@ -335,14 +349,20 @@ IF(LLDEFAULTVAL) THEN LROTATE_WIND=.FALSE. LTKEMINTURB=.TRUE. LPROJQITURB=.TRUE. + LSMOOTH_PRANDTL=.TRUE. + XMINSIGS=0. + LNEWBL89EXP=.TRUE. IF(HPROGRAM=='AROME') THEN XTKEMIN=1.E-6 XLINI=0. LPROJQITURB=.FALSE. + LSMOOTH_PRANDTL=.FALSE. + LNEWBL89EXP=.FALSE. ELSEIF(HPROGRAM=='MESONH') THEN LROTATE_WIND=.TRUE. LTKEMINTURB=.FALSE. + XMINSIGS=1.E-12 ELSEIF(HPROGRAM=='LMDZ') THEN XTKEMIN=1.E-6 XLINI=0. diff --git a/src/common/turb/mode_bl89.F90 b/src/common/turb/mode_bl89.F90 index 80e0dc0db96e00dc1b99e56414c80a10adbf7752..f60bbc580aeae4aa9132cf7238dfb3b347624351 100644 --- a/src/common/turb/mode_bl89.F90 +++ b/src/common/turb/mode_bl89.F90 @@ -112,7 +112,7 @@ INTEGER :: JRR ! moist loop counter REAL :: ZRVORD ! Rv/Rd REAL :: ZPOTE,ZLWORK1,ZLWORK2 REAL :: ZTEST,ZTEST0,ZTESTM ! test for vectorization -REAL :: Z2SQRT2,ZUSRBL89,ZBL89EXP +REAL :: Z2SQRT2 !------------------------------------------------------------------------------- ! REAL(KIND=JPHOOK) :: ZHOOK_HANDLE @@ -157,9 +157,6 @@ END IF ZSQRT_TKE(IIJB:IIJE,1:IKT) = SQRT(PTKEM(IIJB:IIJE,1:IKT)) !$mnh_end_expand_array(JIJ=IIJB:IIJE,JK=1:IKT) ! -!ZBL89EXP is defined here because (and not in ini_cturb) because TURBN%XCED was defined in read_exseg (depending on BL89/RM17) -ZBL89EXP = LOG(16.)/(4.*LOG(CST%XKARMAN)+LOG(TURBN%XCED)-3.*LOG(CSTURB%XCMFS)) -ZUSRBL89 = 1./ZBL89EXP !------------------------------------------------------------------------------- ! !* 2. Virtual potential temperature on the model grid @@ -317,13 +314,13 @@ DO JK=IKTB,IKTE ZLWORK1=MAX(PLMDN(JIJ,JK),1.E-10_MNHREAL) ZLWORK2=MAX(ZLWORK(JIJ),1.E-10_MNHREAL) ZPOTE = ZLWORK1 / ZLWORK2 -#ifdef REPRO48 + IF (.NOT. TURBN%LNEWBL89EXP) THEN ZLWORK2=1.d0 + ZPOTE**(2./3.) PLM(JIJ,JK) = Z2SQRT2*ZLWORK1/(ZLWORK2*SQRT(ZLWORK2)) -#else - ZLWORK2=1.d0 + ZPOTE**ZBL89EXP - PLM(JIJ,JK) = ZLWORK1*(2./ZLWORK2)**ZUSRBL89 -#endif + ELSE + ZLWORK2=1.d0 + ZPOTE**TURBN%XBL89EXP + PLM(JIJ,JK) = ZLWORK1*(2./ZLWORK2)**TURBN%XUSRBL89 + ENDIF PLM(JIJ,JK)=MAX(PLM(JIJ,JK),TURBN%XLINI) END DO diff --git a/src/common/turb/mode_ini_turb.F90 b/src/common/turb/mode_ini_turb.F90 index 321453659fe0918e7e873f1549f251500f444987..fc70cfaf945db2312defa67a60b887f9f3bea4af 100644 --- a/src/common/turb/mode_ini_turb.F90 +++ b/src/common/turb/mode_ini_turb.F90 @@ -52,7 +52,8 @@ CONTAINS ! ------------ ! USE MODD_CST -USE MODD_TURB_n, ONLY : XCTP, XCED, XCSHF, XCHF, XCTV, XCHV, XCHT1, XCHT2, XCPR1, CTURBLEN +USE MODD_TURB_n, ONLY : XCTP, XCED, XCSHF, XCHF, XCTV, XCHV, XCHT1, XCHT2, XCPR1, CTURBLEN, & + & XBL89EXP, XUSRBL89 USE MODD_NEB_n, ONLY: LSTATNW USE MODD_CTURB ! For true constants (not tunable) USE MODD_PARAMETERS, ONLY : XUNDEF @@ -255,6 +256,12 @@ XSBL_O_BL = 0.05 ! SBL height / BL height ratio XFTOP_O_FSURF = 0.05 ! Fraction of surface (heat or momentum) flux used to define top of BL ! ! +! 7. Constants for BL89 computation +! +XBL89EXP=LOG(16.)/(4.*LOG(XKARMAN)+LOG(XCED)-3.*LOG(XCMFS)) +XUSRBL89=1./XBL89EXP +! +! IF (LHOOK) CALL DR_HOOK('INI_TURB',1,ZHOOK_HANDLE) END SUBROUTINE INI_TURB END MODULE MODE_INI_TURB diff --git a/src/common/turb/mode_prandtl.F90 b/src/common/turb/mode_prandtl.F90 index ca06510b284b06d69e8672066c2e240e951c205d..98c2af638ccac5eb401efa917539db095ff9dc62 100644 --- a/src/common/turb/mode_prandtl.F90 +++ b/src/common/turb/mode_prandtl.F90 @@ -11,6 +11,7 @@ USE YOMHOOK , ONLY : LHOOK, DR_HOOK, JPHOOK !* modification 08/2010 V. Masson smoothing of the discontinuity in functions ! used for implicitation of exchange coefficients ! 05/2020 V. Masson and C. Lac : bug in D_PHI3DTDZ2_O_DDTDZ +! 06/2023 S. Riette add the LSMOOTH_PRANDTL key ! USE MODD_CTURB, ONLY : CSTURB_t USE MODD_TURB_n, ONLY : TURB_t @@ -753,12 +754,14 @@ IIJE=D%NIJE IIJB=D%NIJB IKT=D%NKT ! -!$mnh_expand_array(JIJ=IIJB:IIJE,JK=1:IKT) -ZCOEF(IIJB:IIJE,1:IKT) = MAX(MIN(( 10.*(1.-PPHI3(IIJB:IIJE,1:IKT)/CSTURB%XPHI_LIM)) ,1.), 0.) -! -PF(IIJB:IIJE,1:IKT) = ZCOEF(IIJB:IIJE,1:IKT) * PF(IIJB:IIJE,1:IKT) & - + (1.-ZCOEF(IIJB:IIJE,1:IKT)) * PF_LIM(IIJB:IIJE,1:IKT) -!$mnh_end_expand_array(JIJ=IIJB:IIJE,JK=1:IKT) +IF(TURBN%LSMOOTH_PRANDTL) THEN + !$mnh_expand_array(JIJ=IIJB:IIJE,JK=1:IKT) + ZCOEF(IIJB:IIJE,1:IKT) = MAX(MIN(( 10.*(1.-PPHI3(IIJB:IIJE,1:IKT)/CSTURB%XPHI_LIM)) ,1.), 0.) + ! + PF(IIJB:IIJE,1:IKT) = ZCOEF(IIJB:IIJE,1:IKT) * PF(IIJB:IIJE,1:IKT) & + + (1.-ZCOEF(IIJB:IIJE,1:IKT)) * PF_LIM(IIJB:IIJE,1:IKT) + !$mnh_end_expand_array(JIJ=IIJB:IIJE,JK=1:IKT) +ENDIF ! END SUBROUTINE SMOOTH_TURB_FUNCT !---------------------------------------------------------------------------- @@ -974,11 +977,8 @@ DO JK=1,IKT ENDDO END IF ! -#ifdef REPRO48 -#else !* smoothing CALL SMOOTH_TURB_FUNCT(D,CSTURB,TURBN,PPHI3,PPHI3,PD_PHI3DTDZ_O_DDTDZ) -#endif ! PD_PHI3DTDZ_O_DDTDZ(IIJB:IIJE,IKB-1)=PD_PHI3DTDZ_O_DDTDZ(IIJB:IIJE,IKB) PD_PHI3DTDZ_O_DDTDZ(IIJB:IIJE,IKE+1)=PD_PHI3DTDZ_O_DDTDZ(IIJB:IIJE,IKE) @@ -1052,11 +1052,8 @@ ELSE !$mnh_end_expand_where(JIJ=IIJB:IIJE,JK=1:IKT) END IF ! -#ifdef REPRO48 -#else !* smoothing CALL SMOOTH_TURB_FUNCT(D,CSTURB,TURBN,PPHI3,PPHI3,PD_PHI3DRDZ_O_DDRDZ) -#endif ! PD_PHI3DRDZ_O_DDRDZ(IIJB:IIJE,IKB-1)=PD_PHI3DRDZ_O_DDRDZ(IIJB:IIJE,IKB) PD_PHI3DRDZ_O_DDRDZ(IIJB:IIJE,IKE+1)=PD_PHI3DRDZ_O_DDRDZ(IIJB:IIJE,IKE) @@ -1112,11 +1109,8 @@ ELSE !$mnh_end_expand_where(JIJ=IIJB:IIJE,JK=1:IKT) END IF ! -#ifdef REPRO48 -#else !* smoothing CALL SMOOTH_TURB_FUNCT(D,CSTURB,TURBN,PPHI3,PPHI3*2.*PDTDZ,PD_PHI3DTDZ2_O_DDTDZ) -#endif ! ! PD_PHI3DTDZ2_O_DDTDZ(IIJB:IIJE,IKB-1)=PD_PHI3DTDZ2_O_DDTDZ(IIJB:IIJE,IKB) diff --git a/src/common/turb/mode_turb_ver_thermo_corr.F90 b/src/common/turb/mode_turb_ver_thermo_corr.F90 index 6a54ba96effa46a639ab5934f9b87f8dbeba64d6..369137dd7c51d966077d2bdd3f2881d22f638185 100644 --- a/src/common/turb/mode_turb_ver_thermo_corr.F90 +++ b/src/common/turb/mode_turb_ver_thermo_corr.F90 @@ -199,6 +199,7 @@ SUBROUTINE TURB_VER_THERMO_CORR(D,CST,CSTURB,TURBN,NEBN,TLES, & !! Modifications July 2015 (Wim de Rooy) TURBN%LHARAT switch !! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O !! Modifications June 2019 (Wim de Rooy) New set up cloud scheme +!! Modifications: June 2023 (S. Riette) tunable value for SIGS minimum value !!-------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -1253,12 +1254,7 @@ ENDIF PSIGS(IIJB:IIJE,IKU) = PSIGS(IIJB:IIJE,IKE) !$mnh_end_expand_array(JIJ=IIJB:IIJE) !$mnh_expand_array(JIJ=IIJB:IIJE,JK=1:IKT) -#ifdef REPRO48 - PSIGS(IIJB:IIJE,1:IKT) = MAX (PSIGS(IIJB:IIJE,1:IKT) , 0.) - PSIGS(IIJB:IIJE,1:IKT) = SQRT(PSIGS(IIJB:IIJE,1:IKT)) -#else - PSIGS(IIJB:IIJE,1:IKT) = SQRT( MAX (PSIGS(IIJB:IIJE,1:IKT) , 1.E-12) ) -#endif + PSIGS(IIJB:IIJE,1:IKT) = SQRT( MAX (PSIGS(IIJB:IIJE,1:IKT) , TURBN%XMINSIGS) ) !$mnh_end_expand_array(JIJ=IIJB:IIJE,JK=1:IKT) END IF