Changeset 5633 for LMDZ6/trunk

Timestamp:

Apr 28, 2025, 2:46:35 PM (6 weeks ago)

Author:

evignon

Message:

small changes related to mixed phase clouds treatment
Etienne and Lea

Location:

LMDZ6/trunk/libf/phylmd

Files:

• lmdz_lscp_main.f90 (modified) (5 diffs)
• lmdz_lscp_tools.f90 (modified) (6 diffs)

LMDZ6/trunk/libf/phylmd/lmdz_lscp_main.f90

@@ -289 +289 @@
   REAL, DIMENSION(klon) :: zdqsdT_raw
   REAL, DIMENSION(klon) :: gammasat,dgammasatdt                   ! coefficient to make cold condensation at the correct RH and derivative wrt T
-  REAL, DIMENSION(klon) :: Tbef,Tbefth,qlibef,DT                  ! temperature, humidity and temp. variation during condensation iteration
+  REAL, DIMENSION(klon) :: Tbef,Tbefth,Tbefthm1,qlibef,DT                  ! temperature, humidity and temp. variation during condensation iteration
   REAL :: num,denom
   REAL :: cste
@@ -536 +536 @@
     ENDIF ! (ok_poprecip)
     
-    ! Calculation of qsat, L/Cp*dqsat/dT and ncoreczq counter
+    ! Calculation of qsat,L/cp*dqsat/dT and ncoreczq counter
     !-------------------------------------------------------
 
@@ -653 +653 @@
         qlibef(:)=0.
         Tbefth(:)=tla(:,k)*pspsk(:,k)
+        IF (k .GT. 1) THEN
+         Tbefthm1(:)=tla(:,k-1)*pspsk(:,k-1)
+        ELSE
+         Tbefthm1(:)=Tbefth(:)       
+        ENDIF
         zqth=qta(:,k)
 
@@ -817 +822 @@
                         invtau_e(i) = 0.
                      ENDDO
-                     CALL icefrac_lscp_turb(klon, dtime, pticefracturb, Tbef, zp, paprs(:,k), paprs(:,k+1), wls(:,k), zqi_ini, ziflcld, qincloud, &
-                     zcf, tke(:,k), tke_dissip(:,k), sursat_e, invtau_e, zqliq, zqvapcl, zqice, zficeenv, dzficeenv,                   &
-                     cldfraliq(:,k),sigma2_icefracturb(:,k),mean_icefracturb(:,k))                     
+                     CALL icefrac_lscp_turb(klon, dtime, pticefracturb, Tbef, zp, paprs(:,k), paprs(:,k+1), wls(:,k), zqi_ini,    &
+                     ziflcld, znebprecipcld, qincloud, zcf, tke(:,k), tke_dissip(:,k), sursat_e, invtau_e, zqliq, zqvapcl, zqice, & 
+                     zficeenv, dzficeenv, cldfraliq(:,k),sigma2_icefracturb(:,k),mean_icefracturb(:,k))                     
                      DO i=1,klon
                         IF (pticefracturb(i)) THEN
@@ -846 +851 @@
                         invtau_e(i)=gamma_mixth*MAX(entr_therm(i,k)-detr_therm(i,k),0.)*RD*Tbef(i)/zp(i)
                      ENDDO
-                     CALL icefrac_lscp_turb(klon, dtime, pticefracturb, Tbefth, zp, paprs(:,k), paprs(:,k+1), wth(:,k), zqi_ini, zeroklon, qincloud, &
-                     zcf, zeroklon, zeroklon, sursat_e, invtau_e, zqliqth, zqvapclth, zqiceth, zficeth, dzficeth,                      &
-                     cldfraliqth(:,k), sigma2_icefracturbth(:,k), mean_icefracturbth(:,k))
+                     CALL icefrac_lscp_turb(klon, dtime, pticefracturb, Tbefth, zp, paprs(:,k), paprs(:,k+1), wth(:,k), zqi_ini,  &
+                     zeroklon, znebprecipcld, qincloud, zcf, zeroklon, zeroklon, sursat_e, invtau_e, zqliqth, zqvapclth, zqiceth, &
+                     zficeth, dzficeth,cldfraliqth(:,k), sigma2_icefracturbth(:,k), mean_icefracturbth(:,k))
                      !Environment
                      DO i=1,klon
                         qincloud(i)=zqn(i)/(1.-fraca(i,k))
                         zcf(i)=MIN(MAX(rneb(i,k)-rnebth(i,k), 0.),1.)*(1.-fraca(i,k))
-                        sursat_e(i)=cloudth_sth(i,k)/zqsith(i)
-                        invtau_e(i)=gamma_mixth*MAX(detr_therm(i,k)-entr_therm(i,k),0.)*RD*Tbef(i)/zp(i)
+                        IF (k .GT. 1) THEN
+                           ! evaluate the mixing sursaturation using saturation deficit at level below
+                           sursat_e(i)=cloudth_sth(i,k-1)/(zqsith(i)+zdqsith(i)*RCPD/RLSTT*(Tbefthm1(i)-Tbefth(i)))
+                           invtau_e(i)=gamma_mixth*MAX(detr_therm(i,k)-entr_therm(i,k),0.)*RD*Tbef(i)/zp(i)
+                        ELSE
+                           sursat_e(i)=0.
+                           invtau_e(i)=0.
+                        ENDIF
                      ENDDO
-                     CALL icefrac_lscp_turb(klon, dtime, pticefracturb, Tbef, zp, paprs(:,k), paprs(:,k+1), wls(:,k), zqi_ini, ziflcld, qincloud, &
-                     zcf, tke(:,k), tke_dissip(:,k), sursat_e, invtau_e, zqliq, zqvapcl, zqice, zfice, dzfice,                      &
-                     cldfraliq(:,k),sigma2_icefracturb(:,k), mean_icefracturb(:,k)) 
+                     CALL icefrac_lscp_turb(klon, dtime, pticefracturb, Tbef, zp, paprs(:,k), paprs(:,k+1), wls(:,k), zqi_ini,    &
+                     ziflcld, znebprecipcld, qincloud, zcf, tke(:,k), tke_dissip(:,k), sursat_e, invtau_e, zqliq, zqvapcl, zqice, &
+                     zfice, dzfice, cldfraliq(:,k),sigma2_icefracturb(:,k), mean_icefracturb(:,k)) 
   
                     ! adjust zfice to account for condensates in thermals'fraction

LMDZ6/trunk/libf/phylmd/lmdz_lscp_tools.f90

@@ -241 +241 @@
 
 
-SUBROUTINE ICEFRAC_LSCP_TURB(klon, dtime, pticefracturb, temp, pplay, paprsdn, paprsup, wvel, qice_ini, snowcld, qtot_incl, cldfra, tke,   &
+SUBROUTINE ICEFRAC_LSCP_TURB(klon, dtime, pticefracturb, temp, pplay, paprsdn, paprsup, wvel, qice_ini, snowcld, snowfracld, qtot_incl, cldfra, tke,   &
                              tke_dissip, sursat_e, invtau_e, qliq, qvap_cld, qice, icefrac, dicefracdT, cldfraliq, sigma2_icefracturb, mean_icefracturb)
 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -256 +256 @@
    USE lmdz_lscp_ini, ONLY : seuil_neb, temp_nowater
    USE lmdz_lscp_ini, ONLY : naero5, gamma_snwretro, gamma_taud, capa_crystal, rho_ice
-   USE lmdz_lscp_ini, ONLY : eps
+   USE lmdz_lscp_ini, ONLY : eps, snow_fallspeed
 
    IMPLICIT NONE
@@ -275 +275 @@
    REAL,      INTENT(IN),       DIMENSION(klon)    :: qice_ini            !--initial specific ice content gridbox-mean     [kg/kg]
    REAL,      INTENT(IN),       DIMENSION(klon)    :: snowcld             !--in-cloud snowfall flux                        [kg/m2/s]
+   REAL,      INTENT(IN),       DIMENSION(klon)    :: snowfracld          !--cloudy precip fraction                        [-]
    REAL,      INTENT(IN),       DIMENSION(klon)    :: sursat_e            !--environment supersaturation                   [-]
    REAL,      INTENT(IN),       DIMENSION(klon)    :: invtau_e            !--inverse time-scale of mixing with environment [s-1]
@@ -371 +372 @@
         ELSE
 
-           qiceini_incl  = qice_ini(i) / cldfra1D + gamma_snwretro * snowcld(i) * RG * dtime / ( paprsdn(i) - paprsup(i) ) / cldfra1D
-           ! assuming snowfall velocity = 1 m/s:
-           !qiceini_incl  = qice_ini(i) / cldfra1D + gamma_snwretro * snowcld(i) / pplay(i) * RD * temp(i) / 1. / cldfra1D
+           ! gamma_snwretro controls the contribution of snowflakes to the negative feedback
+           ! note that for reasons related to inetarctions with the condensation iteration in lscp_main
+           ! we consider here the mean snowflake concentration in the mesh (not the in-cloud concentration)
+           ! when poprecip is active, it will be worth testing considering the incloud fraction, dividing 
+           ! by snowfracld     
+           ! qiceini_incl  = qice_ini(i) / cldfra1D + &
+           !              gamma_snwretro * snowcld(i) * RG * dtime / ( paprsdn(i) - paprsup(i) )
+           ! assuming constant snowfall velocity 
+           qiceini_incl  = qice_ini(i) / cldfra1D + gamma_snwretro * snowcld(i) / pplay(i) * RD * temp(i) / snow_fallspeed 
 
            !--1. No preexisting ice and no mixing with environment: if vertical motion, fully liquid
@@ -445 +452 @@
               ai = RG / RD / temp(i) * ( RD * RLSTT / RCPD / RV / temp(i) - 1. )
               sursat_equ    = (ai * wvel(i) + sursat_e(i)*invtau_e(i)) / (invtau_phaserelax + invtau_e(i))
-              sursat_equ    = MAX(-1., MIN(sursat_equ , 1.)) 
               ! as sursaturation is by definition lower than -1 and 
               ! because local supersaturation > 1 are never found in the atmosphere
@@ -490 +496 @@
 
                  liqfra_max = MAX(0., (MIN (1.,( qtot_incl(i) - (qice_ini(i) / cldfra1D) - qsati(i) * (1 + sursat_iceext ) ) / ( qsatl(i) - qsati(i) ) ) ) )
-                 sigma2_pdf = MIN(1., 1./2. * ( ai**2 ) *  2./3. * tke(i) * tau_dissipturb / (invtau_phaserelax + invtau_e(i))) 
+                 sigma2_pdf = 1./2. * ( ai**2 ) *  2./3. * tke(i) * tau_dissipturb / (invtau_phaserelax + invtau_e(i))
                  ! sursat ranges between -1 and 1, so we prevent sigma2 so exceed 1
                  cldfraliq(i) = 0.5 * (1. - erf( ( sursat_iceliq - sursat_equ) / (SQRT(2.* sigma2_pdf) ) ) )