Ignore:
Timestamp:
May 7, 2019, 10:45:07 AM (5 years ago)
Author:
idelkadi
Message:

Amelioration de la representation des nuages bas (Jean Jouhaud) :
Calcul des nouveaux ecarts types sigmas de la bi-Gaussienne en prenant en compte l'epaisseur des mailles, CF_vol avec ces sigmas, et CF_surf avec la methode de Brooks

Location:
LMDZ6/trunk/libf/phylmd
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/cloudth_mod.F90

    r2960 r3493  
    77       SUBROUTINE cloudth(ngrid,klev,ind2,  &
    88     &           ztv,po,zqta,fraca, &
    9      &           qcloud,ctot,zpspsk,paprs,ztla,zthl, &
     9     &           qcloud,ctot,zpspsk,paprs,pplay,ztla,zthl, &
    1010     &           ratqs,zqs,t)
    1111
    … …  
    3838      REAL zpspsk(ngrid,klev)
    3939      REAL paprs(ngrid,klev+1)
     40      REAL pplay(ngrid,klev)
    4041      REAL ztla(ngrid,klev)
    4142      REAL zthl(ngrid,klev)
    … …  
    7879      CALL cloudth_vert(ngrid,klev,ind2,  &
    7980     &           ztv,po,zqta,fraca, &
    80      &           qcloud,ctot,zpspsk,paprs,ztla,zthl, &
     81     &           qcloud,ctot,zpspsk,paprs,pplay,ztla,zthl, &
    8182     &           ratqs,zqs,t)
    8283      RETURN
    … …  
    251252     SUBROUTINE cloudth_vert(ngrid,klev,ind2,  &
    252253     &           ztv,po,zqta,fraca, &
    253      &           qcloud,ctot,zpspsk,paprs,ztla,zthl, &
     254     &           qcloud,ctot,zpspsk,paprs,pplay,ztla,zthl, &
    254255     &           ratqs,zqs,t)
    255256
    … …  
    282283      REAL zpspsk(ngrid,klev)
    283284      REAL paprs(ngrid,klev+1)
     285      REAL pplay(ngrid,klev)
    284286      REAL ztla(ngrid,klev)
    285287      REAL zthl(ngrid,klev)
    … …  
    585587END SUBROUTINE cloudth_vert
    586588
     589
     590
     591
    587592       SUBROUTINE cloudth_v3(ngrid,klev,ind2,  &
    588593     &           ztv,po,zqta,fraca, &
    589      &           qcloud,ctot,ctot_vol,zpspsk,paprs,ztla,zthl, &
     594     &           qcloud,ctot,ctot_vol,zpspsk,paprs,pplay,ztla,zthl, &
    590595     &           ratqs,zqs,t)
    591596
    … …  
    618623      REAL zpspsk(ngrid,klev)
    619624      REAL paprs(ngrid,klev+1)
     625      REAL pplay(ngrid,klev)
    620626      REAL ztla(ngrid,klev)
    621627      REAL zthl(ngrid,klev)
    … …  
    641647      REAL alth,alenv,ath,aenv
    642648      REAL sth,senv,sigma1s,sigma2s,xth,xenv, exp_xenv1, exp_xenv2,exp_xth1,exp_xth2
     649      REAL inverse_rho,beta,a_Brooks,b_Brooks,A_Maj_Brooks,Dx_Brooks,f_Brooks
    643650      REAL Tbef,zdelta,qsatbef,zcor
    644651      REAL qlbef 
    … …  
    654661      CALL cloudth_vert_v3(ngrid,klev,ind2,  &
    655662     &           ztv,po,zqta,fraca, &
    656      &           qcloud,ctot,ctot_vol,zpspsk,paprs,ztla,zthl, &
     663     &           qcloud,ctot,ctot_vol,zpspsk,paprs,pplay,ztla,zthl, &
    657664     &           ratqs,zqs,t)
    658665      RETURN
    … …  
    808815     SUBROUTINE cloudth_vert_v3(ngrid,klev,ind2,  &
    809816     &           ztv,po,zqta,fraca, &
    810      &           qcloud,ctot,ctot_vol,zpspsk,paprs,ztla,zthl, &
     817     &           qcloud,ctot,ctot_vol,zpspsk,paprs,pplay,ztla,zthl, &
    811818     &           ratqs,zqs,t)
    812819
    … …  
    841848      REAL zpspsk(ngrid,klev)
    842849      REAL paprs(ngrid,klev+1)
     850      REAL pplay(ngrid,klev)
    843851      REAL ztla(ngrid,klev)
    844852      REAL zthl(ngrid,klev)
    … …  
    864872      REAL alth,alenv,ath,aenv
    865873      REAL sth,senv,sigma1s,sigma2s,sigma1s_fraca,sigma1s_ratqs
     874      REAL inverse_rho,beta,a_Brooks,b_Brooks,A_Maj_Brooks,Dx_Brooks,f_Brooks
    866875      REAL xth,xenv,exp_xenv1,exp_xenv2,exp_xth1,exp_xth2
    867876      REAL xth1,xth2,xenv1,xenv2,deltasth, deltasenv
    … …  
    888897      REAL zqs(ngrid), qcloud(ngrid)
    889898      REAL erf
     899
     900      REAL rhodz(ngrid,klev)
     901      REAL zrho(ngrid,klev)
     902      REAL dz(ngrid,klev)
     903
     904      DO ind1 = 1, ngrid
     905        !Layer calculation
     906        rhodz(ind1,ind2) = (paprs(ind1,ind2)-paprs(ind1,ind2+1))/rg !kg/m2
     907        zrho(ind1,ind2) = pplay(ind1,ind2)/t(ind1,ind2)/rd !kg/m3
     908        dz(ind1,ind2) = rhodz(ind1,ind2)/zrho(ind1,ind2) !m : epaisseur de la couche en metre
     909      END DO
     910
    890911
    891912!------------------------------------------------------------------------------
    … …  
    10501071
    10511072      ELSE IF (iflag_cloudth_vert >= 3) THEN
    1052 
     1073      IF (iflag_cloudth_vert < 5) THEN
    10531074!-------------------------------------------------------------------------------
    10541075!  Version 3: Changes by J. Jouhaud; condensation for q > -delta s
    … …  
    11191140      endif
    11201141
    1121 
    11221142      qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2)
    11231143      ctot_vol(ind1,ind2)=fraca(ind1,ind2)*cth_vol(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*cenv_vol(ind1,ind2)
    11241144
     1145
     1146
     1147
     1148
     1149
     1150      ELSE IF (iflag_cloudth_vert == 5) THEN
     1151      sigma1s=(0.71794+0.000498239*dz(ind1,ind2))*(fraca(ind1,ind2)**0.5)/(1-fraca(ind1,ind2))*(((sth-senv)**2)**0.5)+ratqs(ind1,ind2)*po(ind1) !Environment
     1152      sigma2s=(0.03218+0.000092655*dz(ind1,ind2))/((fraca(ind1,ind2)+0.02)**0.5)*(((sth-senv)**2)**0.5)+0.002*zqta(ind1,ind2)                   !Thermals
     1153      !sigma1s=(1.1**0.5)*(fraca(ind1,ind2)**0.6)/(1-fraca(ind1,ind2))*((sth-senv)**2)**0.5+0.002*po(ind1)
     1154      !sigma2s=0.11*((sth-senv)**2)**0.5/(fraca(ind1,ind2)+0.01)**0.4+0.002*zqta(ind1,ind2)
     1155      xth=sth/(sqrt(2.)*sigma2s)
     1156      xenv=senv/(sqrt(2.)*sigma1s)
     1157
     1158      !Volumique
     1159      cth_vol(ind1,ind2)=0.5*(1.+1.*erf(xth))
     1160      cenv_vol(ind1,ind2)=0.5*(1.+1.*erf(xenv))
     1161      ctot_vol(ind1,ind2)=fraca(ind1,ind2)*cth_vol(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*cenv_vol(ind1,ind2)
     1162      !print *,'jeanjean_CV=',ctot_vol(ind1,ind2)
     1163
     1164      qlth(ind1,ind2)=sigma2s*((exp(-1.*xth**2)/sqrt2pi)+xth*sqrt(2.)*cth_vol(ind1,ind2))
     1165      qlenv(ind1,ind2)=sigma1s*((exp(-1.*xenv**2)/sqrt2pi)+xenv*sqrt(2.)*cenv_vol(ind1,ind2)) 
     1166      qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2)
     1167
     1168      !Surfacique
     1169      !Neggers
     1170      !beta=0.0044
     1171      !inverse_rho=1.+beta*dz(ind1,ind2)
     1172      !print *,'jeanjean : beta=',beta
     1173      !cth(ind1,ind2)=cth_vol(ind1,ind2)*inverse_rho
     1174      !cenv(ind1,ind2)=cenv_vol(ind1,ind2)*inverse_rho
     1175      !ctot(ind1,ind2)=fraca(ind1,ind2)*cth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*cenv(ind1,ind2)
     1176
     1177      !Brooks
     1178      a_Brooks=0.6694
     1179      b_Brooks=0.1882
     1180      A_Maj_Brooks=0.1635 !-- sans shear
     1181      !A_Maj_Brooks=0.17   !-- ARM LES
     1182      !A_Maj_Brooks=0.18   !-- RICO LES
     1183      !A_Maj_Brooks=0.19   !-- BOMEX LES
     1184      Dx_Brooks=200000.
     1185      f_Brooks=A_Maj_Brooks*(dz(ind1,ind2)**(a_Brooks))*(Dx_Brooks**(-b_Brooks))
     1186      !print *,'jeanjean_f=',f_Brooks
     1187
     1188      cth(ind1,ind2)=1./(1.+exp(-1.*f_Brooks)*((1./max(1.e-15,min(cth_vol(ind1,ind2),1.)))- 1.))
     1189      cenv(ind1,ind2)=1./(1.+exp(-1.*f_Brooks)*((1./max(1.e-15,min(cenv_vol(ind1,ind2),1.)))- 1.))
     1190      ctot(ind1,ind2)=1./(1.+exp(-1.*f_Brooks)*((1./max(1.e-15,min(ctot_vol(ind1,ind2),1.)))- 1.))
     1191      !print *,'JJ_ctot_1',ctot(ind1,ind2)
     1192
     1193
     1194
     1195
     1196
     1197      ENDIF ! of if (iflag_cloudth_vert<5)
    11251198      ENDIF ! of if (iflag_cloudth_vert==1 or 3 or 4)
    11261199
    1127 
    1128       if (cenv(ind1,ind2).lt.1.e-10.or.cth(ind1,ind2).lt.1.e-10) then
     1200      if (ctot(ind1,ind2).lt.1.e-10) then
    11291201      ctot(ind1,ind2)=0.
    11301202      ctot_vol(ind1,ind2)=0.
    1131       qcloud(ind1)=zqsatenv(ind1,ind2) 
    1132 
    1133       else 
     1203      qcloud(ind1)=zqsatenv(ind1,ind2)
     1204
     1205      else
    11341206               
    11351207      qcloud(ind1)=qltot(ind1,ind2)/ctot(ind1,ind2)+zqs(ind1)
    … …  
    11391211      endif 
    11401212
    1141       else  ! Environment only
     1213      else  ! gaussienne environnement seule
    11421214     
    11431215      zqenv(ind1)=po(ind1)
    … …  
    11511223     
    11521224
    1153 !     qlbef=Max(po(ind1)-zqsatenv(ind1,ind2),0.)
     1225!      qlbef=Max(po(ind1)-zqsatenv(ind1,ind2),0.)
    11541226      zthl(ind1,ind2)=t(ind1,ind2)*(101325/paprs(ind1,ind2))**(rdd/cppd)
    11551227      alenv=(0.622*Lv*zqsatenv(ind1,ind2))/(rdd*zthl(ind1,ind2)**2) 
    11561228      aenv=1./(1.+(alenv*Lv/cppd))
    1157       senv=aenv*(po(ind1)-zqsatenv(ind1,ind2))
    1158       sth=0.
     1229      senv=aenv*(po(ind1)-zqsatenv(ind1,ind2)) 
     1230     
    11591231
    11601232      sigma1s=ratqs(ind1,ind2)*zqenv(ind1)
    1161       sigma2s=0.
    1162 
    1163       xenv=senv/(sqrt2*sigma1s)
     1233
     1234      sqrt2pi=sqrt(2.*pi)
     1235      xenv=senv/(sqrt(2.)*sigma1s)
    11641236      ctot(ind1,ind2)=0.5*(1.+1.*erf(xenv))
    1165       ctot_vol(ind1,ind2)=ctot(ind1,ind2)
    1166       qltot(ind1,ind2)=sigma1s*((exp(-1.*xenv**2)/sqrt2pi)+xenv*sqrt2*cenv(ind1,ind2))
     1237      qltot(ind1,ind2)=sigma1s*((exp(-1.*xenv**2)/sqrt2pi)+xenv*sqrt(2.)*cenv(ind1,ind2))
    11671238     
    11681239      if (ctot(ind1,ind2).lt.1.e-3) then
    … …  
    11721243      else   
    11731244               
     1245      ctot(ind1,ind2)=ctot(ind1,ind2)
    11741246      qcloud(ind1)=qltot(ind1,ind2)/ctot(ind1,ind2)+zqsatenv(ind1,ind2)
    11751247
    1176       endif   
    1177  
     1248      endif 
     1249 
     1250
     1251
     1252
    11781253      endif       ! From the separation (thermal/envrionnement) et (environnement) only, l.335 et l.492
    11791254      ! Outputs used to check the PDFs
    … …  
    11841259
    11851260      enddo       ! from the loop on ngrid l.333
    1186      
    11871261      return
    11881262!     end
    11891263END SUBROUTINE cloudth_vert_v3
    11901264!
     1265
     1266
     1267
     1268
     1269
     1270
     1271
     1272
     1273
     1274
     1275
     1276       SUBROUTINE cloudth_v6(ngrid,klev,ind2,  &
     1277     &           ztv,po,zqta,fraca, &
     1278     &           qcloud,ctot_surf,ctot_vol,zpspsk,paprs,pplay,ztla,zthl, &
     1279     &           ratqs,zqs,T)
     1280
     1281
     1282      USE ioipsl_getin_p_mod, ONLY : getin_p
     1283      USE phys_output_var_mod, ONLY : cloudth_sth,cloudth_senv, &
     1284     &                                cloudth_sigmath,cloudth_sigmaenv
     1285
     1286      IMPLICIT NONE
     1287
     1288#include "YOMCST.h"
     1289#include "YOETHF.h"
     1290#include "FCTTRE.h"
     1291#include "thermcell.h"
     1292#include "nuage.h"
     1293
     1294
     1295        !Domain variables
     1296      INTEGER ngrid !indice Max lat-lon
     1297      INTEGER klev  !indice Max alt
     1298      INTEGER ind1  !indice in [1:ngrid]
     1299      INTEGER ind2  !indice in [1:klev]
     1300        !thermal plume fraction
     1301      REAL fraca(ngrid,klev+1)   !thermal plumes fraction in the gridbox
     1302        !temperatures
     1303      REAL T(ngrid,klev)       !temperature
     1304      REAL zpspsk(ngrid,klev)  !factor (p/p0)**kappa (used for potential variables)
     1305      REAL ztv(ngrid,klev)     !potential temperature (voir thermcell_env.F90)
     1306      REAL ztla(ngrid,klev)    !liquid temperature in the thermals (Tl_th)
     1307      REAL zthl(ngrid,klev)    !liquid temperature in the environment (Tl_env)
     1308        !pressure
     1309      REAL paprs(ngrid,klev+1)   !pressure at the interface of levels
     1310      REAL pplay(ngrid,klev)     !pressure at the middle of the level
     1311        !humidity
     1312      REAL ratqs(ngrid,klev)   !width of the total water subgrid-scale distribution
     1313      REAL po(ngrid)           !total water (qt)
     1314      REAL zqenv(ngrid)        !total water in the environment (qt_env)
     1315      REAL zqta(ngrid,klev)    !total water in the thermals (qt_th)
     1316      REAL zqsatth(ngrid,klev)   !water saturation level in the thermals (q_sat_th)
     1317      REAL zqsatenv(ngrid,klev)  !water saturation level in the environment (q_sat_env)
     1318      REAL qlth(ngrid,klev)    !condensed water in the thermals
     1319      REAL qlenv(ngrid,klev)   !condensed water in the environment
     1320      REAL qltot(ngrid,klev)   !condensed water in the gridbox
     1321        !cloud fractions
     1322      REAL cth_vol(ngrid,klev)   !cloud fraction by volume in the thermals
     1323      REAL cenv_vol(ngrid,klev)  !cloud fraction by volume in the environment
     1324      REAL ctot_vol(ngrid,klev)  !cloud fraction by volume in the gridbox
     1325      REAL cth_surf(ngrid,klev)  !cloud fraction by surface in the thermals
     1326      REAL cenv_surf(ngrid,klev) !cloud fraction by surface in the environment 
     1327      REAL ctot_surf(ngrid,klev) !cloud fraction by surface in the gridbox
     1328        !PDF of saturation deficit variables
     1329      REAL rdd,cppd,Lv
     1330      REAL Tbef,zdelta,qsatbef,zcor
     1331      REAL alth,alenv,ath,aenv
     1332      REAL sth,senv              !saturation deficits in the thermals and environment
     1333      REAL sigma_env,sigma_th    !standard deviations of the biGaussian PDF
     1334        !cloud fraction variables
     1335      REAL xth,xenv
     1336      REAL inverse_rho,beta                                  !Neggers et al. (2011) method
     1337      REAL a_Brooks,b_Brooks,A_Maj_Brooks,Dx_Brooks,f_Brooks !Brooks et al. (2005) method
     1338        !Incloud total water variables
     1339      REAL zqs(ngrid)    !q_sat
     1340      REAL qcloud(ngrid) !eau totale dans le nuage
     1341        !Some arithmetic variables
     1342      REAL erf,pi,sqrt2,sqrt2pi
     1343        !Depth of the layer
     1344      REAL dz(ngrid,klev)    !epaisseur de la couche en metre
     1345      REAL rhodz(ngrid,klev)
     1346      REAL zrho(ngrid,klev)
     1347      DO ind1 = 1, ngrid
     1348        rhodz(ind1,ind2) = (paprs(ind1,ind2)-paprs(ind1,ind2+1))/rg ![kg/m2]
     1349        zrho(ind1,ind2) = pplay(ind1,ind2)/T(ind1,ind2)/rd          ![kg/m3]
     1350        dz(ind1,ind2) = rhodz(ind1,ind2)/zrho(ind1,ind2)            ![m]
     1351      END DO
     1352
     1353!------------------------------------------------------------------------------
     1354! Initialization
     1355!------------------------------------------------------------------------------
     1356      qlth(:,:)=0.
     1357      qlenv(:,:)=0. 
     1358      qltot(:,:)=0.
     1359      cth_vol(:,:)=0.
     1360      cenv_vol(:,:)=0.
     1361      ctot_vol(:,:)=0.
     1362      cth_surf(:,:)=0.
     1363      cenv_surf(:,:)=0.
     1364      ctot_surf(:,:)=0.
     1365      qcloud(:)=0.
     1366      rdd=287.04
     1367      cppd=1005.7
     1368      pi=3.14159
     1369      Lv=2.5e6
     1370      sqrt2=sqrt(2.)
     1371      sqrt2pi=sqrt(2.*pi)
     1372
     1373
     1374      DO ind1=1,ngrid
     1375!-------------------------------------------------------------------------------
     1376!Both thermal and environment in the gridbox
     1377!-------------------------------------------------------------------------------
     1378      IF ((ztv(ind1,1).gt.ztv(ind1,2)).and.(fraca(ind1,ind2).gt.1.e-10)) THEN
     1379        !--------------------------------------------
     1380        !calcul de qsat_env
     1381        !--------------------------------------------
     1382      Tbef=zthl(ind1,ind2)*zpspsk(ind1,ind2)
     1383      zdelta=MAX(0.,SIGN(1.,RTT-Tbef))
     1384      qsatbef= R2ES*FOEEW(Tbef,zdelta)/paprs(ind1,ind2)
     1385      qsatbef=MIN(0.5,qsatbef)
     1386      zcor=1./(1.-retv*qsatbef)
     1387      qsatbef=qsatbef*zcor
     1388      zqsatenv(ind1,ind2)=qsatbef
     1389        !--------------------------------------------
     1390        !calcul de s_env
     1391        !--------------------------------------------
     1392      alenv=(0.622*Lv*zqsatenv(ind1,ind2))/(rdd*zthl(ind1,ind2)**2)     !qsl, p84 these Arnaud Jam
     1393      aenv=1./(1.+(alenv*Lv/cppd))                                      !al, p84 these Arnaud Jam
     1394      senv=aenv*(po(ind1)-zqsatenv(ind1,ind2))                          !s, p84 these Arnaud Jam
     1395        !--------------------------------------------
     1396        !calcul de qsat_th
     1397        !--------------------------------------------
     1398      Tbef=ztla(ind1,ind2)*zpspsk(ind1,ind2)
     1399      zdelta=MAX(0.,SIGN(1.,RTT-Tbef))
     1400      qsatbef= R2ES * FOEEW(Tbef,zdelta)/paprs(ind1,ind2)
     1401      qsatbef=MIN(0.5,qsatbef)
     1402      zcor=1./(1.-retv*qsatbef)
     1403      qsatbef=qsatbef*zcor
     1404      zqsatth(ind1,ind2)=qsatbef
     1405        !--------------------------------------------
     1406        !calcul de s_th 
     1407        !--------------------------------------------
     1408      alth=(0.622*Lv*zqsatth(ind1,ind2))/(rdd*ztla(ind1,ind2)**2)       !qsl, p84 these Arnaud Jam
     1409      ath=1./(1.+(alth*Lv/cppd))                                        !al, p84 these Arnaud Jam
     1410      sth=ath*(zqta(ind1,ind2)-zqsatth(ind1,ind2))                      !s, p84 these Arnaud Jam
     1411        !--------------------------------------------
     1412        !calcul standard deviations bi-Gaussian PDF
     1413        !--------------------------------------------
     1414      sigma_th=(0.03218+0.000092655*dz(ind1,ind2))/((fraca(ind1,ind2)+0.01)**0.5)*(((sth-senv)**2)**0.5)+0.002*zqta(ind1,ind2)
     1415      sigma_env=(0.71794+0.000498239*dz(ind1,ind2))*(fraca(ind1,ind2)**0.5)/(1-fraca(ind1,ind2))*(((sth-senv)**2)**0.5)+ratqs(ind1,ind2)*po(ind1)
     1416      xth=sth/(sqrt2*sigma_th)
     1417      xenv=senv/(sqrt2*sigma_env)
     1418        !--------------------------------------------
     1419        !Cloud fraction by volume CF_vol
     1420        !--------------------------------------------
     1421      cth_vol(ind1,ind2)=0.5*(1.+1.*erf(xth))
     1422      cenv_vol(ind1,ind2)=0.5*(1.+1.*erf(xenv))
     1423      ctot_vol(ind1,ind2)=fraca(ind1,ind2)*cth_vol(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*cenv_vol(ind1,ind2)
     1424        !--------------------------------------------
     1425        !Condensed water qc
     1426        !--------------------------------------------
     1427      qlth(ind1,ind2)=sigma_th*((exp(-1.*xth**2)/sqrt2pi)+xth*sqrt2*cth_vol(ind1,ind2))
     1428      qlenv(ind1,ind2)=sigma_env*((exp(-1.*xenv**2)/sqrt2pi)+xenv*sqrt2*cenv_vol(ind1,ind2)) 
     1429      qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2)
     1430        !--------------------------------------------
     1431        !Cloud fraction by surface CF_surf
     1432        !--------------------------------------------
     1433        !Method Neggers et al. (2011) : ok for cumulus clouds only
     1434      !beta=0.0044 (Jouhaud et al.2018)
     1435      !inverse_rho=1.+beta*dz(ind1,ind2)
     1436      !ctot_surf(ind1,ind2)=ctot_vol(ind1,ind2)*inverse_rho
     1437        !Method Brooks et al. (2005) : ok for all types of clouds
     1438      a_Brooks=0.6694
     1439      b_Brooks=0.1882
     1440      A_Maj_Brooks=0.1635 !-- sans dependence au cisaillement de vent
     1441      Dx_Brooks=200000.   !-- si l'on considere des mailles de 200km de cote
     1442      f_Brooks=A_Maj_Brooks*(dz(ind1,ind2)**(a_Brooks))*(Dx_Brooks**(-b_Brooks))
     1443      ctot_surf(ind1,ind2)=1./(1.+exp(-1.*f_Brooks)*((1./max(1.e-15,min(ctot_vol(ind1,ind2),1.)))- 1.))
     1444        !--------------------------------------------
     1445        !Incloud Condensed water qcloud
     1446        !--------------------------------------------
     1447      if (ctot_surf(ind1,ind2) .lt. 1.e-10) then
     1448      ctot_vol(ind1,ind2)=0.
     1449      ctot_surf(ind1,ind2)=0.
     1450      qcloud(ind1)=zqsatenv(ind1,ind2)
     1451      else
     1452      qcloud(ind1)=qltot(ind1,ind2)/ctot_vol(ind1,ind2)+zqs(ind1)
     1453      endif
     1454
     1455
     1456
     1457!-------------------------------------------------------------------------------
     1458!Environment only in the gridbox
     1459!-------------------------------------------------------------------------------
     1460      ELSE
     1461        !--------------------------------------------
     1462        !calcul de qsat_env
     1463        !--------------------------------------------
     1464      Tbef=zthl(ind1,ind2)*zpspsk(ind1,ind2)
     1465      zdelta=MAX(0.,SIGN(1.,RTT-Tbef))
     1466      qsatbef= R2ES * FOEEW(Tbef,zdelta)/paprs(ind1,ind2)
     1467      qsatbef=MIN(0.5,qsatbef)
     1468      zcor=1./(1.-retv*qsatbef)
     1469      qsatbef=qsatbef*zcor
     1470      zqsatenv(ind1,ind2)=qsatbef
     1471        !--------------------------------------------
     1472        !calcul de s_env
     1473        !--------------------------------------------
     1474      alenv=(0.622*Lv*zqsatenv(ind1,ind2))/(rdd*zthl(ind1,ind2)**2)     !qsl, p84 these Arnaud Jam
     1475      aenv=1./(1.+(alenv*Lv/cppd))                                      !al, p84 these Arnaud Jam
     1476      senv=aenv*(po(ind1)-zqsatenv(ind1,ind2))                          !s, p84 these Arnaud Jam
     1477        !--------------------------------------------
     1478        !calcul standard deviations Gaussian PDF
     1479        !--------------------------------------------
     1480      zqenv(ind1)=po(ind1)
     1481      sigma_env=ratqs(ind1,ind2)*zqenv(ind1)
     1482      xenv=senv/(sqrt2*sigma_env)
     1483        !--------------------------------------------
     1484        !Cloud fraction by volume CF_vol
     1485        !--------------------------------------------
     1486      ctot_vol(ind1,ind2)=0.5*(1.+1.*erf(xenv))
     1487        !--------------------------------------------
     1488        !Condensed water qc
     1489        !--------------------------------------------
     1490      qltot(ind1,ind2)=sigma_env*((exp(-1.*xenv**2)/sqrt2pi)+xenv*sqrt2*ctot_vol(ind1,ind2))
     1491        !--------------------------------------------
     1492        !Cloud fraction by surface CF_surf
     1493        !--------------------------------------------
     1494        !Method Neggers et al. (2011) : ok for cumulus clouds only
     1495      !beta=0.0044 (Jouhaud et al.2018)
     1496      !inverse_rho=1.+beta*dz(ind1,ind2)
     1497      !ctot_surf(ind1,ind2)=ctot_vol(ind1,ind2)*inverse_rho
     1498        !Method Brooks et al. (2005) : ok for all types of clouds
     1499      a_Brooks=0.6694
     1500      b_Brooks=0.1882
     1501      A_Maj_Brooks=0.1635 !-- sans dependence au shear
     1502      Dx_Brooks=200000.
     1503      f_Brooks=A_Maj_Brooks*(dz(ind1,ind2)**(a_Brooks))*(Dx_Brooks**(-b_Brooks))
     1504      ctot_surf(ind1,ind2)=1./(1.+exp(-1.*f_Brooks)*((1./max(1.e-15,min(ctot_vol(ind1,ind2),1.)))- 1.))
     1505        !--------------------------------------------
     1506        !Incloud Condensed water qcloud
     1507        !--------------------------------------------
     1508      if (ctot_surf(ind1,ind2) .lt. 1.e-8) then
     1509      ctot_vol(ind1,ind2)=0.
     1510      ctot_surf(ind1,ind2)=0.
     1511      qcloud(ind1)=zqsatenv(ind1,ind2)
     1512      else
     1513      qcloud(ind1)=qltot(ind1,ind2)/ctot_vol(ind1,ind2)+zqsatenv(ind1,ind2)
     1514      endif
     1515
     1516
     1517      END IF  ! From the separation (thermal/envrionnement) et (environnement only)
     1518
     1519      ! Outputs used to check the PDFs
     1520      cloudth_senv(ind1,ind2) = senv
     1521      cloudth_sth(ind1,ind2) = sth
     1522      cloudth_sigmaenv(ind1,ind2) = sigma_env
     1523      cloudth_sigmath(ind1,ind2) = sigma_th
     1524
     1525      END DO  ! From the loop on ngrid
     1526      return
     1527
     1528END SUBROUTINE cloudth_v6
    11911529END MODULE cloudth_mod
     1530
     1531
     1532
     1533

  • LMDZ6/trunk/libf/phylmd/fisrtilp.F90

    r2969 r3493  
    740740               call cloudth(klon,klev,k,ztv, &
    741741                   zq,zqta,fraca, &
    742                    qcloud,ctot,zpspsk,paprs,ztla,zthl, &
     742                   qcloud,ctot,zpspsk,paprs,pplay,ztla,zthl, &
    743743                   ratqs,zqs,t)
    744               elseif (iflag_cloudth_vert>=3) then
     744              elseif (iflag_cloudth_vert>=3 .and. iflag_cloudth_vert<=5) then
    745745               call cloudth_v3(klon,klev,k,ztv, &
    746746                   zq,zqta,fraca, &
    747                    qcloud,ctot,ctot_vol,zpspsk,paprs,ztla,zthl, &
     747                   qcloud,ctot,ctot_vol,zpspsk,paprs,pplay,ztla,zthl, &
    748748                   ratqs,zqs,t)
     749              !----------------------------------
     750              !Version these Jean Jouhaud, Decembre 2018
     751              !----------------------------------             
     752              elseif (iflag_cloudth_vert==6) then
     753               call cloudth_v6(klon,klev,k,ztv, &
     754                   zq,zqta,fraca, &
     755                   qcloud,ctot,ctot_vol,zpspsk,paprs,pplay,ztla,zthl, &
     756                   ratqs,zqs,t)
     757
    749758              endif
    750759              do i=1,klon
Note: See TracChangeset for help on using the changeset viewer.