Index: LMDZ6/trunk/libf/phylmd/lmdz_lscp_main.f90 =================================================================== --- LMDZ6/trunk/libf/phylmd/lmdz_lscp_main.f90 (revision 5630) +++ LMDZ6/trunk/libf/phylmd/lmdz_lscp_main.f90 (revision 5633) @@ -289,5 +289,5 @@ 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,5 +536,5 @@ ENDIF ! (ok_poprecip) - ! Calculation of qsat, L/Cp*dqsat/dT and ncoreczq counter + ! Calculation of qsat,L/cp*dqsat/dT and ncoreczq counter !------------------------------------------------------- @@ -653,4 +653,9 @@ 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,7 +822,7 @@ 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,17 +851,23 @@ 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 Index: LMDZ6/trunk/libf/phylmd/lmdz_lscp_tools.f90 =================================================================== --- LMDZ6/trunk/libf/phylmd/lmdz_lscp_tools.f90 (revision 5630) +++ LMDZ6/trunk/libf/phylmd/lmdz_lscp_tools.f90 (revision 5633) @@ -241,5 +241,5 @@ -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,5 +256,5 @@ 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,4 +275,5 @@ 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,7 +372,13 @@ 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,5 +452,4 @@ 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,5 +496,5 @@ 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) ) ) )