1 !-----------------------------------------------------------------------
3 !WRF:MODEL_LAYER:PHYSICS
5 !####################TIEDTKE SCHEME#########################
6 ! Taken from the IPRC iRAM - Yuqing Wang, University of Hawaii
7 ! Added by Chunxi Zhang and Yuqing Wang to WRF3.2, May, 2010
8 ! refenrence: Tiedtke (1989, MWR, 117, 1779-1800)
9 ! Nordeng, T.E., (1995), CAPE closure and organized entrainment/detrainment
10 ! Yuqing Wang et al. (2003,J. Climate, 16, 1721-1738) for improvements
11 ! for cloud top detrainment
12 ! (2004, Mon. Wea. Rev., 132, 274-296), improvements for PBL clouds
13 ! (2007,Mon. Wea. Rev., 135, 567-585), diurnal cycle of precipitation
14 ! This scheme is on testing
15 !###########################################################
16 MODULE module_cu_tiedtke
18 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
19 ! epsl--- allowed minimum value for floating calculation
20 !---------------------------------------------------------------
21 real,parameter :: epsl = 1.0e-20
22 real,parameter :: t000 = 273.15
23 real,parameter :: hgfr = 233.15 ! defined in param.f in explct
24 !-------------------------------------------------------------
25 ! Ends the parameters set
26 !++++++++++++++++++++++++++++
28 REAL :: API,A,EOMEGA,RD,RV,CPD,RCPD,VTMPC1,VTMPC2, &
29 RHOH2O,ALV,ALS,ALF,CLW,TMELT,SOLC,STBO,DAYL,YEARL, &
30 C1ES,C2ES,C3LES,C3IES,C4LES,C4IES,C5LES,C5IES,ZRG
32 REAL :: ENTRPEN,ENTRSCV,ENTRMID,ENTRDD,CMFCTOP,RHM,RHC, &
33 CMFCMAX,CMFCMIN,CMFDEPS,RHCDD,CPRCON,CRIRH,ZBUO0, &
36 INTEGER :: orgen,nturben,cutrigger
38 REAL :: CVDIFTS, CEVAPCU1, CEVAPCU2,ZDNOPRC
41 PARAMETER(A=6371.22E03, &
46 CPV=1869.46, & ! CPV in module is 1846.4
58 C5LES=C3LES*(TMELT-C4LES), &
61 C5IES=C3IES*(TMELT-C4IES), &
62 API=3.141593, & ! API=2.0*ASIN(1.)
66 CEVAPCU1=1.93E-6*261.0*0.5/G, &
67 CEVAPCU2=1.E3/(38.3*0.293) )
70 ! SPECIFY PARAMETERS FOR MASSFLUX-SCHEME
71 ! --------------------------------------
72 ! These are tunable parameters
74 ! ENTRPEN: AVERAGE ENTRAINMENT RATE FOR PENETRATIVE CONVECTION
77 PARAMETER(ENTRPEN=1.0E-4)
79 ! ENTRSCV: AVERAGE ENTRAINMENT RATE FOR SHALLOW CONVECTION
82 PARAMETER(ENTRSCV=1.2E-3)
84 ! ENTRMID: AVERAGE ENTRAINMENT RATE FOR MIDLEVEL CONVECTION
87 PARAMETER(ENTRMID=1.0E-4)
89 ! ENTRDD: AVERAGE ENTRAINMENT RATE FOR DOWNDRAFTS
92 PARAMETER(ENTRDD =2.0E-4)
94 ! CMFCTOP: RELATIVE CLOUD MASSFLUX AT LEVEL ABOVE NONBUOYANCY LEVEL
97 PARAMETER(CMFCTOP=0.30)
99 ! CMFCMAX: MAXIMUM MASSFLUX VALUE ALLOWED FOR UPDRAFTS ETC
102 PARAMETER(CMFCMAX=1.0)
104 ! CMFCMIN: MINIMUM MASSFLUX VALUE (FOR SAFETY)
107 PARAMETER(CMFCMIN=1.E-10)
109 ! CMFDEPS: FRACTIONAL MASSFLUX FOR DOWNDRAFTS AT LFS
112 PARAMETER(CMFDEPS=0.30)
114 ! CPRCON: COEFFICIENTS FOR DETERMINING CONVERSION FROM CLOUD WATER
116 PARAMETER(CPRCON = 1.1E-3/G)
118 ! ZDNOPRC: The pressure depth below which no precipitation
120 PARAMETER(ZDNOPRC =1.5E4)
121 !--------------------
122 PARAMETER(orgen=1) ! Old organized entrainment rate
123 ! PARAMETER(orgen=2) ! New organized entrainment rate
125 PARAMETER(nturben=1) ! old deep turburent entrainment/detrainment rate
126 ! PARAMETER(nturben=2) ! New deep turburent entrainment/detrainment rate
128 PARAMETER(cutrigger=1) ! Old trigger function
129 ! PARAMETER(cutrigger=2) ! New trigger function
131 !--------------------
132 PARAMETER(RHC=0.80,RHM=1.0,ZBUO0=0.50)
133 !--------------------
134 PARAMETER(CRIRH=0.70,fdbk = 1.0,ZTAU = 1800.0)
135 !--------------------
136 LOGICAL :: LMFPEN,LMFMID,LMFSCV,LMFDD,LMFDUDV
137 PARAMETER(LMFPEN=.TRUE.,LMFMID=.TRUE.,LMFSCV=.TRUE.,LMFDD=.TRUE.,LMFDUDV=.TRUE.)
138 !--------------------
139 !#################### END of Variables definition##########################
140 !-----------------------------------------------------------------------
143 !-----------------------------------------------------------------------
144 SUBROUTINE CU_TIEDTKE( &
145 DT,ITIMESTEP,STEPCU &
146 ,RAINCV,PRATEC,QFX,HFX,ZNU &
147 ,U3D,V3D,W,T3D,QV3D,QC3D,QI3D,PI3D,RHO3D &
149 ,DZ8W,PCPS,P8W,XLAND,CU_ACT_FLAG &
150 ,CUDT, CURR_SECS, ADAPT_STEP_FLAG &
152 ,ids,ide, jds,jde, kds,kde &
153 ,ims,ime, jms,jme, kms,kme &
154 ,its,ite, jts,jte, kts,kte &
155 ,RTHCUTEN,RQVCUTEN,RQCCUTEN,RQICUTEN &
157 ,F_QV ,F_QC ,F_QR ,F_QI ,F_QS &
159 !-------------------------------------------------------------------
161 !-------------------------------------------------------------------
162 !-- U3D 3D u-velocity interpolated to theta points (m/s)
163 !-- V3D 3D v-velocity interpolated to theta points (m/s)
164 !-- TH3D 3D potential temperature (K)
165 !-- T3D temperature (K)
166 !-- QV3D 3D water vapor mixing ratio (Kg/Kg)
167 !-- QC3D 3D cloud mixing ratio (Kg/Kg)
168 !-- QI3D 3D ice mixing ratio (Kg/Kg)
169 !-- RHO3D 3D air density (kg/m^3)
170 !-- P8w 3D hydrostatic pressure at full levels (Pa)
171 !-- Pcps 3D hydrostatic pressure at half levels (Pa)
172 !-- PI3D 3D exner function (dimensionless)
173 !-- RTHCUTEN Theta tendency due to
174 ! cumulus scheme precipitation (K/s)
175 !-- RUCUTEN U wind tendency due to
176 ! cumulus scheme precipitation (K/s)
177 !-- RVCUTEN V wind tendency due to
178 ! cumulus scheme precipitation (K/s)
179 !-- RQVCUTEN Qv tendency due to
180 ! cumulus scheme precipitation (kg/kg/s)
181 !-- RQRCUTEN Qr tendency due to
182 ! cumulus scheme precipitation (kg/kg/s)
183 !-- RQCCUTEN Qc tendency due to
184 ! cumulus scheme precipitation (kg/kg/s)
185 !-- RQSCUTEN Qs tendency due to
186 ! cumulus scheme precipitation (kg/kg/s)
187 !-- RQICUTEN Qi tendency due to
188 ! cumulus scheme precipitation (kg/kg/s)
189 !-- RAINC accumulated total cumulus scheme precipitation (mm)
190 !-- RAINCV cumulus scheme precipitation (mm)
191 !-- PRATEC precipitiation rate from cumulus scheme (mm/s)
192 !-- dz8w dz between full levels (m)
193 !-- QFX upward moisture flux at the surface (kg/m^2/s)
195 !-- ids start index for i in domain
196 !-- ide end index for i in domain
197 !-- jds start index for j in domain
198 !-- jde end index for j in domain
199 !-- kds start index for k in domain
200 !-- kde end index for k in domain
201 !-- ims start index for i in memory
202 !-- ime end index for i in memory
203 !-- jms start index for j in memory
204 !-- jme end index for j in memory
205 !-- kms start index for k in memory
206 !-- kme end index for k in memory
207 !-- its start index for i in tile
208 !-- ite end index for i in tile
209 !-- jts start index for j in tile
210 !-- jte end index for j in tile
211 !-- kts start index for k in tile
212 !-- kte end index for k in tile
213 !-------------------------------------------------------------------
214 INTEGER, INTENT(IN) :: ids,ide, jds,jde, kds,kde, &
215 ims,ime, jms,jme, kms,kme, &
216 its,ite, jts,jte, kts,kte, &
220 REAL, INTENT(IN) :: &
224 REAL, DIMENSION(ims:ime, jms:jme), INTENT(IN) :: &
227 REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: &
230 LOGICAL, DIMENSION(IMS:IME,JMS:JME), INTENT(INOUT) :: &
234 REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: &
250 !--------------------------- OPTIONAL VARS ----------------------------
252 REAL, DIMENSION(ims:ime, kms:kme, jms:jme), &
253 OPTIONAL, INTENT(INOUT) :: &
262 ! Flags relating to the optional tendency arrays declared above
263 ! Models that carry the optional tendencies will provdide the
264 ! optional arguments at compile time; these flags all the model
265 ! to determine at run-time whether a particular tracer is in
268 LOGICAL, OPTIONAL :: &
275 ! Adaptive time-step variables
276 REAL, INTENT(IN ) :: CUDT
277 REAL, INTENT(IN ) :: CURR_SECS
278 LOGICAL,INTENT(IN ) , OPTIONAL :: ADAPT_STEP_FLAG
279 REAL, INTENT (INOUT) :: CUDTACTTIME
281 !--------------------------- LOCAL VARS ------------------------------
283 REAL, DIMENSION(ims:ime, jms:jme) :: &
291 REAL , DIMENSION(its:ite) :: &
297 INTEGER , DIMENSION(its:ite) :: SLIMSK
300 REAL , DIMENSION(its:ite, kts:kte+1) :: &
303 REAL , DIMENSION(its:ite, kts:kte) :: &
323 INTEGER, DIMENSION(its:ite) :: &
337 LOGICAL :: run_param, doing_adapt_dt , decided
339 !-------other local variables----
340 INTEGER,DIMENSION( its:ite ) :: KTYPE
341 REAL, DIMENSION( kts:kte ) :: sig1 ! half sigma levels
342 REAL, DIMENSION( kms:kme ) :: ZNU
344 !-----------------------------------------------------------------------
347 !*** CHECK TO SEE IF THIS IS A CONVECTION TIMESTEP
350 ! Initialization for adaptive time step.
352 doing_adapt_dt = .FALSE.
353 IF ( PRESENT(adapt_step_flag) ) THEN
354 IF ( adapt_step_flag ) THEN
355 doing_adapt_dt = .TRUE.
356 IF ( cudtacttime .EQ. 0. ) THEN
357 cudtacttime = curr_secs + cudt*60.
362 ! Do we run through this scheme or not?
364 ! Test 1: If this is the initial model time, then yes.
366 ! Test 2: If the user asked for the cumulus to be run every time step, then yes.
368 ! Test 3: If not adaptive dt, and this is on the requested cumulus frequency, then yes.
369 ! MOD(ITIMESTEP,STEPCU)=0
370 ! Test 4: If using adaptive dt and the current time is past the last requested activate cumulus time, then yes.
371 ! CURR_SECS >= CUDTACTTIME
373 ! If we do run through the scheme, we set the flag run_param to TRUE and we set the decided flag
374 ! to TRUE. The decided flag says that one of these tests was able to say "yes", run the scheme.
375 ! We only proceed to other tests if the previous tests all have left decided as FALSE.
377 ! If we set run_param to TRUE and this is adaptive time stepping, we set the time to the next
382 IF ( ( .NOT. decided ) .AND. &
383 ( itimestep .EQ. 1 ) ) THEN
388 IF ( ( .NOT. decided ) .AND. &
389 ( ( cudt .EQ. 0. ) .OR. ( stepcu .EQ. 1 ) ) ) THEN
394 IF ( ( .NOT. decided ) .AND. &
395 ( .NOT. doing_adapt_dt ) .AND. &
396 ( MOD(itimestep,stepcu) .EQ. 0 ) ) THEN
401 IF ( ( .NOT. decided ) .AND. &
402 ( doing_adapt_dt ) .AND. &
403 ( curr_secs .GE. cudtacttime ) ) THEN
406 cudtacttime = curr_secs + cudt*60
409 !-----------------------------------------------------------------------
414 CU_ACT_FLAG(I,J)=.TRUE.
423 !------------- J LOOP (OUTER) --------------------------------------------------
427 ! --------------- compute zi and zl -----------------------------------------
435 ZI(I,K)=ZI(I,KM)+dz8w(i,km,j)
442 ZL(I,KM)=(ZI(I,K)+ZI(I,KM))*0.5
447 ZL(I,KTE)=2.*ZI(I,KTE)-ZL(I,KTE-1)
450 ! --------------- end compute zi and zl -------------------------------------
452 SLIMSK(i)=int(ABS(XLAND(i,j)-2.))
458 DOT(i,k)=-0.5*g*rho3d(i,k,j)*(w(i,k,j)+w(i,kp,j))
468 Q1(i,zz)= QV3D(i,k,j)
469 if(itimestep == 1) then
473 Q1B(i,zz)=QVFTEN(i,k,j)
474 Q1BL(i,zz)=QVPBLTEN(i,k,j)
480 PRSL(i,zz) = Pcps(i,k,j)
487 PRSI(i,zz) = P8w(i,k,j)
496 !###############before call TIECNV, we need EVAP########################
497 ! EVAP is the vapor flux at the surface
498 !########################################################################
503 rho2d(i) = rho3d(i,1,j)
505 !########################################################################
506 CALL TIECNV(U1,V1,T1,Q1,Q2,Q3,Q1B,Q1BL,GHT,OMG,PRSL,PRSI,EVAP,heatflux,rho2d, &
507 RN,SLIMSK,KTYPE,IM,KX,KX+1,sig1,DELT)
510 RAINCV(I,J)=RN(I)/STEPCU
511 PRATEC(I,J)=RN(I)/(STEPCU * DT)
517 RTHCUTEN(I,K,J)=(T1(I,zz)-T3D(I,K,J))/PI3D(I,K,J)*RDELT
518 RQVCUTEN(I,K,J)=(Q1(I,zz)-QV3D(I,K,J))*RDELT
519 RUCUTEN(I,K,J) =(U1(I,zz)-U3D(I,K,J))*RDELT
520 RVCUTEN(I,K,J) =(V1(I,zz)-V3D(I,K,J))*RDELT
524 IF(PRESENT(RQCCUTEN))THEN
529 RQCCUTEN(I,K,J)=(Q2(I,zz)-QC3D(I,K,J))*RDELT
535 IF(PRESENT(RQICUTEN))THEN
540 RQICUTEN(I,K,J)=(Q3(I,zz)-QI3D(I,K,J))*RDELT
551 END SUBROUTINE CU_TIEDTKE
553 !====================================================================
554 SUBROUTINE tiedtkeinit(RTHCUTEN,RQVCUTEN,RQCCUTEN,RQICUTEN, &
556 RESTART,P_QC,P_QI,P_FIRST_SCALAR, &
558 ids, ide, jds, jde, kds, kde, &
559 ims, ime, jms, jme, kms, kme, &
560 its, ite, jts, jte, kts, kte)
561 !--------------------------------------------------------------------
563 !--------------------------------------------------------------------
564 LOGICAL , INTENT(IN) :: allowed_to_read,restart
565 INTEGER , INTENT(IN) :: ids, ide, jds, jde, kds, kde, &
566 ims, ime, jms, jme, kms, kme, &
567 its, ite, jts, jte, kts, kte
568 INTEGER , INTENT(IN) :: P_FIRST_SCALAR, P_QI, P_QC
570 REAL, DIMENSION( ims:ime , kms:kme , jms:jme ) , INTENT(OUT) :: &
577 INTEGER :: i, j, k, itf, jtf, ktf
595 IF (P_QC .ge. P_FIRST_SCALAR) THEN
605 IF (P_QI .ge. P_FIRST_SCALAR) THEN
616 END SUBROUTINE tiedtkeinit
618 ! ------------------------------------------------------------------------
620 !------------This is the combined version for tiedtke---------------
621 !----------------------------------------------------------------
622 ! In this module only the mass flux convection scheme of the ECMWF is included
623 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
624 !#############################################################
626 ! LEVEL 1 SUBROUTINEs
628 !#############################################################
629 !********************************************************
631 !********************************************************
632 SUBROUTINE TIECNV(pu,pv,pt,pqv,pqc,pqi,pqvf,pqvbl,poz,pomg, &
633 pap,paph,evap,hfx,rho,zprecc,lndj,KTYPE,lq,km,km1,sig1,dt)
634 !-----------------------------------------------------------------
635 ! This is the interface between the meso-scale model and the mass
636 ! flux convection module
637 !-----------------------------------------------------------------
640 real pu(lq,km),pv(lq,km),pt(lq,km),pqv(lq,km),pqvf(lq,km)
641 real poz(lq,km),pomg(lq,km),evap(lq),zprecc(lq),pqvbl(lq,km)
642 real PHHFL(lq),RHO(lq),hfx(lq)
643 REAL PUM1(lq,km), PVM1(lq,km), &
644 PTTE(lq,km), PQTE(lq,km), PVOM(lq,km), PVOL(lq,km), &
645 PVERV(lq,km), PGEO(lq,km), PAP(lq,km), PAPH(lq,km1)
646 REAL PQHFL(lq), ZQQ(lq,km), PAPRC(lq), PAPRS(lq), &
647 PRSFC(lq), PSSFC(lq), PAPRSM(lq), PCTE(lq,km)
648 REAL ZTP1(lq,km), ZQP1(lq,km), ZTU(lq,km), ZQU(lq,km), &
649 ZLU(lq,km), ZLUDE(lq,km), ZMFU(lq,km), ZMFD(lq,km), &
650 ZQSAT(lq,km), pqc(lq,km), pqi(lq,km), ZRAIN(lq)
652 REAL sig(km1),sig1(km)
653 INTEGER ICBOT(lq), ICTOP(lq), KTYPE(lq), lndj(lq)
657 real PSHEAT,PSRAIN,PSEVAP,PSMELT,PSDISS,TT
658 real ZTMST,ZTPP1,fliq,fice,ZTC,ZALF
659 integer i,j,k,lq,lp,km,km1
664 ! Masv flux diagnostics.
683 ! CONVERT MODEL VARIABLES FOR MFLUX SCHEME
692 ZQP1(j,k)=pqv(j,k)/(1.0+pqv(j,k))
698 ZQSAT(j,k)=TLUCUA(TT)/PAP(j,k)
699 ZQSAT(j,k)=MIN(0.5,ZQSAT(j,k))
700 ZQSAT(j,k)=ZQSAT(j,k)/(1.-VTMPC1*ZQSAT(j,k))
701 PQTE(j,k)=pqvf(j,k)+pqvbl(j,k)
705 !-----------------------------------------------------------------------
706 !* 2. CALL 'CUMASTR'(MASTER-ROUTINE FOR CUMULUS PARAMETERIZATION)
709 (lq, km, km1, km-1, ZTP1, &
710 ZQP1, PUM1, PVM1, PVERV, ZQSAT, &
711 PQHFL, ZTMST, PAP, PAPH, PGEO, &
712 PTTE, PQTE, PVOM, PVOL, PRSFC, &
713 PSSFC, PAPRC, PAPRSM, PAPRS, LOCUM, &
714 KTYPE, ICBOT, ICTOP, ZTU, ZQU, &
715 ZLU, ZLUDE, ZMFU, ZMFD, ZRAIN, &
716 PSRAIN, PSEVAP, PSHEAT, PSDISS, PSMELT, &
717 PCTE, PHHFL, RHO, sig1, lndj)
719 ! TO INCLUDE THE CLOUD WATER AND CLOUD ICE DETRAINED FROM CONVECTION
721 IF(fdbk.ge.1.0e-9) THEN
724 If(PCTE(j,k).GT.0.0) then
725 ZTPP1=pt(j,k)+PTTE(j,k)*ZTMST
726 if(ZTPP1.ge.t000) then
729 else if(ZTPP1.le.hgfr) then
734 fliq=0.0059+0.9941*exp(-0.003102*ZTC*ZTC)
738 pqc(j,k)=pqc(j,k)+fliq*PCTE(j,k)*ZTMST
739 pqi(j,k)=pqi(j,k)+fice*PCTE(j,k)*ZTMST
740 PTTE(j,k)=PTTE(j,k)-ZALF*RCPD*fliq*PCTE(j,k)
747 pt(j,k)=ZTP1(j,k)+PTTE(j,k)*ZTMST
748 ZQP1(j,k)=ZQP1(j,k)+(PQTE(j,k)-ZQQ(j,k))*ZTMST
749 pqv(j,k)=ZQP1(j,k)/(1.0-ZQP1(j,k))
752 zprecc(j)=amax1(0.0,(PRSFC(j)+PSSFC(j))*ZTMST)
757 pu(j,k)=pu(j,k)+PVOM(j,k)*ZTMST
758 pv(j,k)=pv(j,k)+PVOL(j,k)*ZTMST
763 END SUBROUTINE TIECNV
765 !#############################################################
767 ! LEVEL 2 SUBROUTINEs
769 !#############################################################
770 !***********************************************************
771 ! SUBROUTINE CUMASTR_NEW
772 !***********************************************************
773 SUBROUTINE CUMASTR_NEW &
774 (KLON, KLEV, KLEVP1, KLEVM1, PTEN, &
775 PQEN, PUEN, PVEN, PVERV, PQSEN, &
776 PQHFL, ZTMST, PAP, PAPH, PGEO, &
777 PTTE, PQTE, PVOM, PVOL, PRSFC, &
778 PSSFC, PAPRC, PAPRSM, PAPRS, LDCUM, &
779 KTYPE, KCBOT, KCTOP, PTU, PQU, &
780 PLU, PLUDE, PMFU, PMFD, PRAIN, &
781 PSRAIN, PSEVAP, PSHEAT, PSDISS, PSMELT,&
782 PCTE, PHHFL, RHO, sig1, lndj)
784 !***CUMASTR* MASTER ROUTINE FOR CUMULUS MASSFLUX-SCHEME
785 ! M.TIEDTKE E.C.M.W.F. 1986/1987/1989
788 ! THIS ROUTINE COMPUTES THE PHYSICAL TENDENCIES OF THE
789 ! PROGNOSTIC VARIABLES T,Q,U AND V DUE TO CONVECTIVE PROCESSES.
790 ! PROCESSES CONSIDERED ARE: CONVECTIVE FLUXES, FORMATION OF
791 ! PRECIPITATION, EVAPORATION OF FALLING RAIN BELOW CLOUD BASE,
792 ! SATURATED CUMULUS DOWNDRAFTS.
795 ! *CUMASTR* IS CALLED FROM *MSSFLX*
796 ! THE ROUTINE TAKES ITS INPUT FROM THE LONG-TERM STORAGE
797 ! T,Q,U,V,PHI AND P AND MOISTURE TENDENCIES.
798 ! IT RETURNS ITS OUTPUT TO THE SAME SPACE
799 ! 1.MODIFIED TENDENCIES OF MODEL VARIABLES
800 ! 2.RATES OF CONVECTIVE PRECIPITATION
801 ! (USED IN SUBROUTINE SURF)
802 ! 3.CLOUD BASE, CLOUD TOP AND PRECIP FOR RADIATION
803 ! (USED IN SUBROUTINE CLOUD)
806 ! PARAMETERIZATION IS DONE USING A MASSFLUX-SCHEME.
807 ! (1) DEFINE CONSTANTS AND PARAMETERS
808 ! (2) SPECIFY VALUES (T,Q,QS...) AT HALF LEVELS AND
809 ! INITIALIZE UPDRAFT- AND DOWNDRAFT-VALUES IN 'CUINI'
810 ! (3) CALCULATE CLOUD BASE IN 'CUBASE'
811 ! AND SPECIFY CLOUD BASE MASSFLUX FROM PBL MOISTURE BUDGET
812 ! (4) DO CLOUD ASCENT IN 'CUASC' IN ABSENCE OF DOWNDRAFTS
813 ! (5) DO DOWNDRAFT CALCULATIONS:
814 ! (A) DETERMINE VALUES AT LFS IN 'CUDLFS'
815 ! (B) DETERMINE MOIST DESCENT IN 'CUDDRAF'
816 ! (C) RECALCULATE CLOUD BASE MASSFLUX CONSIDERING THE
817 ! EFFECT OF CU-DOWNDRAFTS
818 ! (6) DO FINAL CLOUD ASCENT IN 'CUASC'
819 ! (7) DO FINAL ADJUSMENTS TO CONVECTIVE FLUXES IN 'CUFLX',
820 ! DO EVAPORATION IN SUBCLOUD LAYER
821 ! (8) CALCULATE INCREMENTS OF T AND Q IN 'CUDTDQ'
822 ! (9) CALCULATE INCREMENTS OF U AND V IN 'CUDUDV'
825 ! CUINI: INITIALIZES VALUES AT VERTICAL GRID USED IN CU-PARAMETR.
826 ! CUBASE: CLOUD BASE CALCULATION FOR PENETR.AND SHALLOW CONVECTION
827 ! CUASC: CLOUD ASCENT FOR ENTRAINING PLUME
828 ! CUDLFS: DETERMINES VALUES AT LFS FOR DOWNDRAFTS
829 ! CUDDRAF:DOES MOIST DESCENT FOR CUMULUS DOWNDRAFTS
830 ! CUFLX: FINAL ADJUSTMENTS TO CONVECTIVE FLUXES (ALSO IN PBL)
831 ! CUDQDT: UPDATES TENDENCIES FOR T AND Q
832 ! CUDUDV: UPDATES TENDENCIES FOR U AND V
835 ! LMFPEN=.T. PENETRATIVE CONVECTION IS SWITCHED ON
836 ! LMFSCV=.T. SHALLOW CONVECTION IS SWITCHED ON
837 ! LMFMID=.T. MIDLEVEL CONVECTION IS SWITCHED ON
838 ! LMFDD=.T. CUMULUS DOWNDRAFTS SWITCHED ON
839 ! LMFDUDV=.T. CUMULUS FRICTION SWITCHED ON
841 ! MODEL PARAMETERS (DEFINED IN SUBROUTINE CUPARAM)
842 ! ------------------------------------------------
843 ! ENTRPEN ENTRAINMENT RATE FOR PENETRATIVE CONVECTION
844 ! ENTRSCV ENTRAINMENT RATE FOR SHALLOW CONVECTION
845 ! ENTRMID ENTRAINMENT RATE FOR MIDLEVEL CONVECTION
846 ! ENTRDD ENTRAINMENT RATE FOR CUMULUS DOWNDRAFTS
847 ! CMFCTOP RELATIVE CLOUD MASSFLUX AT LEVEL ABOVE NONBUOYANCY
849 ! CMFCMAX MAXIMUM MASSFLUX VALUE ALLOWED FOR
850 ! CMFCMIN MINIMUM MASSFLUX VALUE (FOR SAFETY)
851 ! CMFDEPS FRACTIONAL MASSFLUX FOR DOWNDRAFTS AT LFS
852 ! CPRCON COEFFICIENT FOR CONVERSION FROM CLOUD WATER TO RAIN
855 ! PAPER ON MASSFLUX SCHEME (TIEDTKE,1989)
856 !-----------------------------------------------------------------
857 !-------------------------------------------------------------------
859 !-------------------------------------------------------------------
860 INTEGER KLON, KLEV, KLEVP1
863 REAL PSRAIN, PSEVAP, PSHEAT, PSDISS, PSMELT, ZCONS2
865 REAL ZQUMQE, ZDQMIN, ZMFMAX, ZALVDCP, ZQALV
866 REAL ZHSAT, ZGAM, ZZZ, ZHHAT, ZBI, ZRO, ZDZ, ZDHDZ, ZDEPTH
867 REAL ZFAC, ZRH, ZPBMPT, DEPT, ZHT, ZEPS
869 REAL PTEN(KLON,KLEV), PQEN(KLON,KLEV), &
870 PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
871 PTTE(KLON,KLEV), PQTE(KLON,KLEV), &
872 PVOM(KLON,KLEV), PVOL(KLON,KLEV), &
873 PQSEN(KLON,KLEV), PGEO(KLON,KLEV), &
874 PAP(KLON,KLEV), PAPH(KLON,KLEVP1),&
875 PVERV(KLON,KLEV), PQHFL(KLON), &
876 PHHFL(KLON), RHO(KLON)
877 REAL PTU(KLON,KLEV), PQU(KLON,KLEV), &
878 PLU(KLON,KLEV), PLUDE(KLON,KLEV), &
879 PMFU(KLON,KLEV), PMFD(KLON,KLEV), &
880 PAPRC(KLON), PAPRS(KLON), &
881 PAPRSM(KLON), PRAIN(KLON), &
882 PRSFC(KLON), PSSFC(KLON)
883 REAL ZTENH(KLON,KLEV), ZQENH(KLON,KLEV),&
884 ZGEOH(KLON,KLEV), ZQSENH(KLON,KLEV),&
885 ZTD(KLON,KLEV), ZQD(KLON,KLEV), &
886 ZMFUS(KLON,KLEV), ZMFDS(KLON,KLEV), &
887 ZMFUQ(KLON,KLEV), ZMFDQ(KLON,KLEV), &
888 ZDMFUP(KLON,KLEV), ZDMFDP(KLON,KLEV),&
889 ZMFUL(KLON,KLEV), ZRFL(KLON), &
890 ZUU(KLON,KLEV), ZVU(KLON,KLEV), &
891 ZUD(KLON,KLEV), ZVD(KLON,KLEV)
892 REAL ZENTR(KLON), ZHCBASE(KLON), &
893 ZMFUB(KLON), ZMFUB1(KLON), &
894 ZDQPBL(KLON), ZDQCV(KLON)
895 REAL ZSFL(KLON), ZDPMEL(KLON,KLEV), &
896 PCTE(KLON,KLEV), ZCAPE(KLON), &
897 ZHEAT(KLON), ZHHATT(KLON,KLEV), &
898 ZHMIN(KLON), ZRELH(KLON)
900 INTEGER ILAB(KLON,KLEV), IDTOP(KLON), &
901 ICTOP0(KLON), ILWMIN(KLON)
902 INTEGER KCBOT(KLON), KCTOP(KLON), &
903 KTYPE(KLON), IHMIN(KLON), &
906 LOGICAL LODDRAF(KLON), LLO1
908 !-------------------------------------------
909 ! 1. SPECIFY CONSTANTS AND PARAMETERS
910 !-------------------------------------------
913 !--------------------------------------------------------------
914 !* 2. INITIALIZE VALUES AT VERTICAL GRID POINTS IN 'CUINI'
915 !--------------------------------------------------------------
918 (KLON, KLEV, KLEVP1, KLEVM1, PTEN, &
919 PQEN, PQSEN, PUEN, PVEN, PVERV, &
920 PGEO, PAPH, ZGEOH, ZTENH, ZQENH, &
921 ZQSENH, ILWMIN, PTU, PQU, ZTD, &
922 ZQD, ZUU, ZVU, ZUD, ZVD, &
923 PMFU, PMFD, ZMFUS, ZMFDS, ZMFUQ, &
924 ZMFDQ, ZDMFUP, ZDMFDP, ZDPMEL, PLU, &
926 !----------------------------------
927 !* 3.0 CLOUD BASE CALCULATIONS
928 !----------------------------------
930 !* (A) DETERMINE CLOUD BASE VALUES IN 'CUBASE'
931 ! -------------------------------------------
933 (KLON, KLEV, KLEVP1, KLEVM1, ZTENH, &
934 ZQENH, ZGEOH, PAPH, PTU, PQU, &
935 PLU, PUEN, PVEN, ZUU, ZVU, &
937 !* (B) DETERMINE TOTAL MOISTURE CONVERGENCE AND
938 !* THEN DECIDE ON TYPE OF CUMULUS CONVECTION
939 ! -----------------------------------------
942 ZDQCV(JL) =PQTE(JL,JK)*(PAPH(JL,JK+1)-PAPH(JL,JK))
948 ZDQCV(JL)=ZDQCV(JL)+PQTE(JL,JK)*(PAPH(JL,JK+1)-PAPH(JL,JK))
949 IF(JK.GE.KCBOT(JL)) ZDQPBL(JL)=ZDQPBL(JL)+PQTE(JL,JK) &
950 *(PAPH(JL,JK+1)-PAPH(JL,JK))
954 if(cutrigger .eq. 1) then
957 IF(ZDQCV(JL).GT.MAX(0.,1.1*PQHFL(JL)*G)) THEN
963 else if(cutrigger .eq. 2) then
965 ( KLON, KLEV, KLEVP1, KLEVM1, &
966 ZTENH, ZQENH, ZQSENH, ZGEOH, PAPH, &
967 RHO, PHHFL, PQHFL, KTYPE, lndj )
969 !* (C) DETERMINE MOISTURE SUPPLY FOR BOUNDARY LAYER
970 !* AND DETERMINE CLOUD BASE MASSFLUX IGNORING
971 !* THE EFFECTS OF DOWNDRAFTS AT THIS STAGE
972 ! ------------------------------------------
974 ! if(ktype(jl) .ge. 1 ) then
975 ! write(6,*)"ktype=", KTYPE(jl)
981 ZQUMQE=PQU(JL,IKB)+PLU(JL,IKB)-ZQENH(JL,IKB)
982 ZDQMIN=MAX(0.01*ZQENH(JL,IKB),1.E-10)
983 IF(ZDQPBL(JL).GT.0..AND.ZQUMQE.GT.ZDQMIN.AND.LDCUM(JL)) THEN
984 ZMFUB(JL)=ZDQPBL(JL)/(G*MAX(ZQUMQE,ZDQMIN))
989 ZMFMAX=(PAPH(JL,IKB)-PAPH(JL,IKB-1))*ZCONS2
990 ZMFUB(JL)=MIN(ZMFUB(JL),ZMFMAX)
991 !------------------------------------------------------
992 !* 4.0 DETERMINE CLOUD ASCENT FOR ENTRAINING PLUME
993 !------------------------------------------------------
995 !* (A) ESTIMATE CLOUD HEIGHT FOR ENTRAINMENT/DETRAINMENT
996 !* CALCULATIONS IN CUASC (MAX.POSSIBLE CLOUD HEIGHT
997 !* FOR NON-ENTRAINING PLUME, FOLLOWING A.-S.,1974)
998 ! -------------------------------------------------------------
1000 ZHCBASE(JL)=CPD*PTU(JL,IKB)+ZGEOH(JL,IKB)+ALV*PQU(JL,IKB)
1001 ICTOP0(JL)=KCBOT(JL)-1
1005 DO 420 JK=KLEVM1,3,-1
1007 ZHSAT=CPD*ZTENH(JL,JK)+ZGEOH(JL,JK)+ALV*ZQSENH(JL,JK)
1008 ZGAM=C5LES*ZALVDCP*ZQSENH(JL,JK)/ &
1009 ((1.-VTMPC1*ZQSENH(JL,JK))*(ZTENH(JL,JK)-C4LES)**2)
1010 ZZZ=CPD*ZTENH(JL,JK)*0.608
1011 ZHHAT=ZHSAT-(ZZZ+ZGAM*ZZZ)/(1.+ZGAM*ZZZ*ZQALV)* &
1012 MAX(ZQSENH(JL,JK)-ZQENH(JL,JK),0.)
1014 IF(JK.LT.ICTOP0(JL).AND.ZHCBASE(JL).GT.ZHHAT) ICTOP0(JL)=JK
1018 ZHSAT=CPD*ZTENH(JL,JK)+ZGEOH(JL,JK)+ALV*ZQSENH(JL,JK)
1019 ZGAM=C5LES*ZALVDCP*ZQSENH(JL,JK)/ &
1020 ((1.-VTMPC1*ZQSENH(JL,JK))*(ZTENH(JL,JK)-C4LES)**2)
1021 ZZZ=CPD*ZTENH(JL,JK)*0.608
1022 ZHHAT=ZHSAT-(ZZZ+ZGAM*ZZZ)/(1.+ZGAM*ZZZ*ZQALV)* &
1023 MAX(ZQSENH(JL,JK)-ZQENH(JL,JK),0.)
1027 ! Find lowest possible org. detrainment level
1031 IF( LDCUM(JL).AND.KTYPE(JL).EQ.1 ) THEN
1032 IHMIN(JL) = KCBOT(JL)
1039 DO 450 JK = KLEV, 1, -1
1041 LLO1 = LDCUM(JL).AND.KTYPE(JL).EQ.1.AND.IHMIN(JL).EQ.KCBOT(JL)
1042 IF (LLO1.AND.JK.LT.KCBOT(JL).AND.JK.GE.ICTOP0(JL)) THEN
1044 ZRO = RD*ZTENH(JL,JK)/(G*PAPH(JL,JK))
1045 ZDZ = (PAPH(JL,JK)-PAPH(JL,JK-1))*ZRO
1046 ZDHDZ=(CPD*(PTEN(JL,JK-1)-PTEN(JL,JK))+ALV*(PQEN(JL,JK-1)- &
1047 PQEN(JL,JK))+(PGEO(JL,JK-1)-PGEO(JL,JK)))*G/(PGEO(JL, &
1049 ZDEPTH = ZGEOH(JL,JK) - ZGEOH(JL,IKB)
1050 ZFAC = SQRT(1.+ZDEPTH*ZBI)
1051 ZHMIN(JL) = ZHMIN(JL) + ZDHDZ*ZFAC*ZDZ
1052 ZRH = -ALV*(ZQSENH(JL,JK)-ZQENH(JL,JK))*ZFAC
1053 IF (ZHMIN(JL).GT.ZRH) IHMIN(JL) = JK
1057 IF (LDCUM(JL).AND.KTYPE(JL).EQ.1) THEN
1058 IF (IHMIN(JL).LT.ICTOP0(JL)) IHMIN(JL) = ICTOP0(JL)
1060 IF(KTYPE(JL).EQ.1) THEN
1065 if(lndj(JL).eq.1) ZENTR(JL)=ZENTR(JL)*1.05
1067 !* (B) DO ASCENT IN 'CUASC'IN ABSENCE OF DOWNDRAFTS
1068 !----------------------------------------------------------
1070 (KLON, KLEV, KLEVP1, KLEVM1, ZTENH, &
1071 ZQENH, PUEN, PVEN, PTEN, PQEN, &
1072 PQSEN, PGEO, ZGEOH, PAP, PAPH, &
1073 PQTE, PVERV, ILWMIN, LDCUM, ZHCBASE, &
1074 KTYPE, ILAB, PTU, PQU, PLU, &
1075 ZUU, ZVU, PMFU, ZMFUB, ZENTR, &
1076 ZMFUS, ZMFUQ, ZMFUL, PLUDE, ZDMFUP, &
1077 KCBOT, KCTOP, ICTOP0, ICUM, ZTMST, &
1078 IHMIN, ZHHATT, ZQSENH)
1079 IF(ICUM.EQ.0) GO TO 1000
1080 !* (C) CHECK CLOUD DEPTH AND CHANGE ENTRAINMENT RATE ACCORDINGLY
1081 ! CALCULATE PRECIPITATION RATE (FOR DOWNDRAFT CALCULATION)
1082 !------------------------------------------------------------------
1084 ZPBMPT=PAPH(JL,KCBOT(JL))-PAPH(JL,KCTOP(JL))
1085 IF(LDCUM(JL)) ICTOP0(JL)=KCTOP(JL)
1086 IF(LDCUM(JL).AND.KTYPE(JL).EQ.1.AND.ZPBMPT.LT.ZDNOPRC) KTYPE(JL)=2
1087 IF(KTYPE(JL).EQ.2) then
1089 if(lndj(JL).eq.1) ZENTR(JL)=ZENTR(JL)*1.05
1091 ZRFL(JL)=ZDMFUP(JL,1)
1095 ZRFL(JL)=ZRFL(JL)+ZDMFUP(JL,JK)
1097 !-----------------------------------------
1098 !* 5.0 CUMULUS DOWNDRAFT CALCULATIONS
1099 !-----------------------------------------
1102 !* (A) DETERMINE LFS IN 'CUDLFS'
1103 !--------------------------------------
1105 (KLON, KLEV, KLEVP1, ZTENH, ZQENH, &
1106 PUEN, PVEN, ZGEOH, PAPH, PTU, &
1107 PQU, ZUU, ZVU, LDCUM, KCBOT, &
1108 KCTOP, ZMFUB, ZRFL, ZTD, ZQD, &
1109 ZUD, ZVD, PMFD, ZMFDS, ZMFDQ, &
1110 ZDMFDP, IDTOP, LODDRAF)
1111 !* (B) DETERMINE DOWNDRAFT T,Q AND FLUXES IN 'CUDDRAF'
1112 !------------------------------------------------------------
1114 (KLON, KLEV, KLEVP1, ZTENH, ZQENH, &
1115 PUEN, PVEN, ZGEOH, PAPH, ZRFL, &
1116 LODDRAF, ZTD, ZQD, ZUD, ZVD, &
1117 PMFD, ZMFDS, ZMFDQ, ZDMFDP)
1118 !* (C) RECALCULATE CONVECTIVE FLUXES DUE TO EFFECT OF
1119 ! DOWNDRAFTS ON BOUNDARY LAYER MOISTURE BUDGET
1120 !-----------------------------------------------------------
1123 !-- 5.1 Recalculate cloud base massflux from a cape closure
1124 ! for deep convection (ktype=1) and by PBL equilibrium
1125 ! taking downdrafts into account for shallow convection
1127 ! implemented by Y. WANG based on ECHAM4 in Nov. 2001.
1133 ZMFUB1(JL)=ZMFUB(JL)
1137 IF(LDCUM(JL).AND.KTYPE(JL).EQ.1) THEN
1139 if(abs(paph(jl,jk)*0.01 - 300) .lt. 50.) then
1140 KTOP0=MAX(jk,KCTOP(JL))
1144 ! KTOP0=MAX(12,KCTOP(JL))
1146 IF(JK.LE.KCBOT(JL).AND.JK.GT.KCTOP(JL)) THEN
1147 ZRO=PAPH(JL,JK)/(RD*ZTENH(JL,JK))
1148 ZDZ=(PAPH(JL,JK)-PAPH(JL,JK-1))/(G*ZRO)
1149 ZHEAT(JL)=ZHEAT(JL)+((PTEN(JL,JK-1)-PTEN(JL,JK) &
1150 +G*ZDZ/CPD)/ZTENH(JL,JK)+0.608*(PQEN(JL,JK-1)- &
1151 PQEN(JL,JK)))*(PMFU(JL,JK)+PMFD(JL,JK))*G/ZRO
1152 ZCAPE(JL)=ZCAPE(JL)+G*((PTU(JL,JK)*(1.+.608*PQU(JL,JK) &
1153 -PLU(JL,JK)))/(ZTENH(JL,JK)*(1.+.608*ZQENH(JL,JK))) &
1156 IF(JK.LE.KCBOT(JL).AND.JK.GT.KTOP0) THEN
1157 dept=(PAPH(JL,JK)-PAPH(JL,JK-1))/(PAPH(JL,KCBOT(JL))- &
1159 ZRELH(JL)=ZRELH(JL)+dept*PQEN(JL,JK)/PQSEN(JL,JK)
1164 if(cutrigger .eq. 1 ) then
1165 IF(lndj(JL).EQ.1) then
1174 IF(ZRELH(JL).GE.CRIRH1 .AND. ZCAPE(JL) .GT. 100.) THEN
1176 ZHT=ZCAPE(JL)/(ZTAU*ZHEAT(JL))
1177 ZMFUB1(JL)=MAX(ZMFUB(JL)*ZHT,0.01)
1178 ZMFMAX=(PAPH(JL,IKB)-PAPH(JL,IKB-1))*ZCONS2
1179 ZMFUB1(JL)=MIN(ZMFUB1(JL),ZMFMAX)
1188 !* 5.2 RECALCULATE CONVECTIVE FLUXES DUE TO EFFECT OF
1189 ! DOWNDRAFTS ON BOUNDARY LAYER MOISTURE BUDGET
1190 !--------------------------------------------------------
1192 IF(KTYPE(JL).NE.1) THEN
1194 IF(PMFD(JL,IKB).LT.0.0.AND.LODDRAF(JL)) THEN
1199 ZQUMQE=PQU(JL,IKB)+PLU(JL,IKB)- &
1200 ZEPS*ZQD(JL,IKB)-(1.-ZEPS)*ZQENH(JL,IKB)
1201 ZDQMIN=MAX(0.01*ZQENH(JL,IKB),1.E-10)
1202 ZMFMAX=(PAPH(JL,IKB)-PAPH(JL,IKB-1))*ZCONS2
1203 IF(ZDQPBL(JL).GT.0..AND.ZQUMQE.GT.ZDQMIN.AND.LDCUM(JL) &
1204 .AND.ZMFUB(JL).LT.ZMFMAX) THEN
1205 ZMFUB1(JL)=ZDQPBL(JL)/(G*MAX(ZQUMQE,ZDQMIN))
1207 ZMFUB1(JL)=ZMFUB(JL)
1209 LLO1=(KTYPE(JL).EQ.2).AND.ABS(ZMFUB1(JL) &
1210 -ZMFUB(JL)).LT.0.2*ZMFUB(JL)
1211 IF(.NOT.LLO1) ZMFUB1(JL)=ZMFUB(JL)
1212 ZMFUB1(JL)=MIN(ZMFUB1(JL),ZMFMAX)
1218 ZFAC=ZMFUB1(JL)/MAX(ZMFUB(JL),1.E-10)
1219 PMFD(JL,JK)=PMFD(JL,JK)*ZFAC
1220 ZMFDS(JL,JK)=ZMFDS(JL,JK)*ZFAC
1221 ZMFDQ(JL,JK)=ZMFDQ(JL,JK)*ZFAC
1222 ZDMFDP(JL,JK)=ZDMFDP(JL,JK)*ZFAC
1232 ZMFUB(JL)=ZMFUB1(JL)
1238 !---------------------------------------------------------------
1239 !* 6.0 DETERMINE FINAL CLOUD ASCENT FOR ENTRAINING PLUME
1240 !* FOR PENETRATIVE CONVECTION (TYPE=1),
1241 !* FOR SHALLOW TO MEDIUM CONVECTION (TYPE=2)
1242 !* AND FOR MID-LEVEL CONVECTION (TYPE=3).
1243 !---------------------------------------------------------------
1246 (KLON, KLEV, KLEVP1, KLEVM1, ZTENH, &
1247 ZQENH, PUEN, PVEN, PTEN, PQEN, &
1248 PQSEN, PGEO, ZGEOH, PAP, PAPH, &
1249 PQTE, PVERV, ILWMIN, LDCUM, ZHCBASE,&
1250 KTYPE, ILAB, PTU, PQU, PLU, &
1251 ZUU, ZVU, PMFU, ZMFUB, ZENTR, &
1252 ZMFUS, ZMFUQ, ZMFUL, PLUDE, ZDMFUP, &
1253 KCBOT, KCTOP, ICTOP0, ICUM, ZTMST, &
1254 IHMIN, ZHHATT, ZQSENH)
1255 !----------------------------------------------------------
1256 !* 7.0 DETERMINE FINAL CONVECTIVE FLUXES IN 'CUFLX'
1257 !----------------------------------------------------------
1260 (KLON, KLEV, KLEVP1, PQEN, PQSEN, &
1261 ZTENH, ZQENH, PAPH, ZGEOH, KCBOT, &
1262 KCTOP, IDTOP, KTYPE, LODDRAF, LDCUM, &
1263 PMFU, PMFD, ZMFUS, ZMFDS, ZMFUQ, &
1264 ZMFDQ, ZMFUL, PLUDE, ZDMFUP, ZDMFDP, &
1265 ZRFL, PRAIN, PTEN, ZSFL, ZDPMEL, &
1266 ITOPM2, ZTMST, sig1)
1267 !----------------------------------------------------------------
1268 !* 8.0 UPDATE TENDENCIES FOR T AND Q IN SUBROUTINE CUDTDQ
1269 !----------------------------------------------------------------
1272 (KLON, KLEV, KLEVP1, ITOPM2, PAPH, &
1273 LDCUM, PTEN, PTTE, PQTE, ZMFUS, &
1274 ZMFDS, ZMFUQ, ZMFDQ, ZMFUL, ZDMFUP, &
1275 ZDMFDP, ZTMST, ZDPMEL, PRAIN, ZRFL, &
1276 ZSFL, PSRAIN, PSEVAP, PSHEAT, PSMELT, &
1277 PRSFC, PSSFC, PAPRC, PAPRSM, PAPRS, &
1278 PQEN, PQSEN, PLUDE, PCTE)
1279 !----------------------------------------------------------------
1280 !* 9.0 UPDATE TENDENCIES FOR U AND U IN SUBROUTINE CUDUDV
1281 !----------------------------------------------------------------
1285 (KLON, KLEV, KLEVP1, ITOPM2, KTYPE, &
1286 KCBOT, PAPH, LDCUM, PUEN, PVEN, &
1287 PVOM, PVOL, ZUU, ZUD, ZVU, &
1288 ZVD, PMFU, PMFD, PSDISS)
1292 END SUBROUTINE CUMASTR_NEW
1295 !#############################################################
1297 ! LEVEL 3 SUBROUTINEs
1299 !#############################################################
1300 !**********************************************
1302 !**********************************************
1305 (KLON, KLEV, KLEVP1, KLEVM1, PTEN, &
1306 PQEN, PQSEN, PUEN, PVEN, PVERV, &
1307 PGEO, PAPH, PGEOH, PTENH, PQENH, &
1308 PQSENH, KLWMIN, PTU, PQU, PTD, &
1309 PQD, PUU, PVU, PUD, PVD, &
1310 PMFU, PMFD, PMFUS, PMFDS, PMFUQ, &
1311 PMFDQ, PDMFUP, PDMFDP, PDPMEL, PLU, &
1313 ! M.TIEDTKE E.C.M.W.F. 12/89
1316 ! THIS ROUTINE INTERPOLATES LARGE-SCALE FIELDS OF T,Q ETC.
1317 ! TO HALF LEVELS (I.E. GRID FOR MASSFLUX SCHEME),
1318 ! AND INITIALIZES VALUES FOR UPDRAFTS AND DOWNDRAFTS
1321 ! THIS ROUTINE IS CALLED FROM *CUMASTR*.
1324 ! FOR EXTRAPOLATION TO HALF LEVELS SEE TIEDTKE(1989)
1327 ! *CUADJTQ* TO SPECIFY QS AT HALF LEVELS
1328 ! ----------------------------------------------------------------
1329 !-------------------------------------------------------------------
1331 !-------------------------------------------------------------------
1332 INTEGER KLON, KLEV, KLEVP1
1334 INTEGER JK,JL,IK, ICALL
1336 REAL PTEN(KLON,KLEV), PQEN(KLON,KLEV), &
1337 PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
1338 PQSEN(KLON,KLEV), PVERV(KLON,KLEV), &
1339 PGEO(KLON,KLEV), PGEOH(KLON,KLEV), &
1340 PAPH(KLON,KLEVP1), PTENH(KLON,KLEV), &
1341 PQENH(KLON,KLEV), PQSENH(KLON,KLEV)
1342 REAL PTU(KLON,KLEV), PQU(KLON,KLEV), &
1343 PTD(KLON,KLEV), PQD(KLON,KLEV), &
1344 PUU(KLON,KLEV), PUD(KLON,KLEV), &
1345 PVU(KLON,KLEV), PVD(KLON,KLEV), &
1346 PMFU(KLON,KLEV), PMFD(KLON,KLEV), &
1347 PMFUS(KLON,KLEV), PMFDS(KLON,KLEV), &
1348 PMFUQ(KLON,KLEV), PMFDQ(KLON,KLEV), &
1349 PDMFUP(KLON,KLEV), PDMFDP(KLON,KLEV), &
1350 PLU(KLON,KLEV), PLUDE(KLON,KLEV)
1351 REAL ZWMAX(KLON), ZPH(KLON), &
1353 INTEGER KLAB(KLON,KLEV), KLWMIN(KLON)
1354 LOGICAL LOFLAG(KLON)
1355 !------------------------------------------------------------
1356 !* 1. SPECIFY LARGE SCALE PARAMETERS AT HALF LEVELS
1357 !* ADJUST TEMPERATURE FIELDS IF STATICLY UNSTABLE
1358 !* FIND LEVEL OF MAXIMUM VERTICAL VELOCITY
1359 ! -----------------------------------------------------------
1364 PGEOH(JL,JK)=PGEO(JL,JK)+(PGEO(JL,JK-1)-PGEO(JL,JK))*ZDP
1365 PTENH(JL,JK)=(MAX(CPD*PTEN(JL,JK-1)+PGEO(JL,JK-1), &
1366 CPD*PTEN(JL,JK)+PGEO(JL,JK))-PGEOH(JL,JK))*RCPD
1367 PQSENH(JL,JK)=PQSEN(JL,JK-1)
1373 CALL CUADJTQ(KLON,KLEV,IK,ZPH,PTENH,PQSENH,LOFLAG,ICALL)
1375 PQENH(JL,JK)=MIN(PQEN(JL,JK-1),PQSEN(JL,JK-1)) &
1376 +(PQSENH(JL,JK)-PQSEN(JL,JK-1))
1377 PQENH(JL,JK)=MAX(PQENH(JL,JK),0.)
1381 PTENH(JL,KLEV)=(CPD*PTEN(JL,KLEV)+PGEO(JL,KLEV)- &
1382 PGEOH(JL,KLEV))*RCPD
1383 PQENH(JL,KLEV)=PQEN(JL,KLEV)
1384 PTENH(JL,1)=PTEN(JL,1)
1385 PQENH(JL,1)=PQEN(JL,1)
1386 PGEOH(JL,1)=PGEO(JL,1)
1390 DO 160 JK=KLEVM1,2,-1
1392 ZZS=MAX(CPD*PTENH(JL,JK)+PGEOH(JL,JK), &
1393 CPD*PTENH(JL,JK+1)+PGEOH(JL,JK+1))
1394 PTENH(JL,JK)=(ZZS-PGEOH(JL,JK))*RCPD
1399 IF(PVERV(JL,JK).LT.ZWMAX(JL)) THEN
1400 ZWMAX(JL)=PVERV(JL,JK)
1405 !-----------------------------------------------------------
1406 !* 2.0 INITIALIZE VALUES FOR UPDRAFTS AND DOWNDRAFTS
1407 !-----------------------------------------------------------
1413 PTU(JL,JK)=PTENH(JL,JK)
1414 PTD(JL,JK)=PTENH(JL,JK)
1415 PQU(JL,JK)=PQENH(JL,JK)
1416 PQD(JL,JK)=PQENH(JL,JK)
1418 PUU(JL,JK)=PUEN(JL,IK)
1419 PUD(JL,JK)=PUEN(JL,IK)
1420 PVU(JL,JK)=PVEN(JL,IK)
1421 PVD(JL,JK)=PVEN(JL,IK)
1436 END SUBROUTINE CUINI
1438 !**********************************************
1440 !**********************************************
1442 (KLON, KLEV, KLEVP1, KLEVM1, PTENH, &
1443 PQENH, PGEOH, PAPH, PTU, PQU, &
1444 PLU, PUEN, PVEN, PUU, PVU, &
1446 ! THIS ROUTINE CALCULATES CLOUD BASE VALUES (T AND Q)
1447 ! FOR CUMULUS PARAMETERIZATION
1448 ! M.TIEDTKE E.C.M.W.F. 7/86 MODIF. 12/89
1451 ! TO PRODUCE CLOUD BASE VALUES FOR CU-PARAMETRIZATION
1454 ! THIS ROUTINE IS CALLED FROM *CUMASTR*.
1455 ! INPUT ARE ENVIRONM. VALUES OF T,Q,P,PHI AT HALF LEVELS.
1456 ! IT RETURNS CLOUD BASE VALUES AND FLAGS AS FOLLOWS;
1457 ! KLAB=1 FOR SUBCLOUD LEVELS
1458 ! KLAB=2 FOR CONDENSATION LEVEL
1461 ! LIFT SURFACE AIR DRY-ADIABATICALLY TO CLOUD BASE
1462 ! (NON ENTRAINING PLUME,I.E.CONSTANT MASSFLUX)
1465 ! *CUADJTQ* FOR ADJUSTING T AND Q DUE TO CONDENSATION IN ASCENT
1466 ! ----------------------------------------------------------------
1467 !-------------------------------------------------------------------
1469 !-------------------------------------------------------------------
1470 INTEGER KLON, KLEV, KLEVP1
1472 INTEGER JL,JK,IS,IK,ICALL,IKB
1474 REAL PTENH(KLON,KLEV), PQENH(KLON,KLEV), &
1475 PGEOH(KLON,KLEV), PAPH(KLON,KLEVP1)
1476 REAL PTU(KLON,KLEV), PQU(KLON,KLEV), &
1478 REAL PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
1479 PUU(KLON,KLEV), PVU(KLON,KLEV)
1480 REAL ZQOLD(KLON,KLEV), ZPH(KLON)
1481 INTEGER KLAB(KLON,KLEV), KCBOT(KLON)
1482 LOGICAL LDCUM(KLON), LOFLAG(KLON)
1483 !***INPUT VARIABLES:
1484 ! PTENH [ZTENH] - Environment Temperature on half levels. (CUINI)
1485 ! PQENH [ZQENH] - Env. specific humidity on half levels. (CUINI)
1486 ! PGEOH [ZGEOH] - Geopotential on half levels, (MSSFLX)
1487 ! PAPH - Pressure of half levels. (MSSFLX)
1488 !***VARIABLES MODIFIED BY CUBASE:
1489 ! LDCUM - Logical denoting profiles. (CUBASE)
1490 ! KTYPE - Convection type - 1: Penetrative (CUMASTR)
1491 ! 2: Stratocumulus (CUMASTR)
1492 ! 3: Mid-level (CUASC)
1493 ! PTU - Cloud Temperature.
1494 ! PQU - Cloud specific Humidity.
1495 ! PLU - Cloud Liquid Water (Moisture condensed out)
1496 ! KCBOT - Cloud Base Level. (CUBASE)
1497 ! KLAB [ILAB] - Level Label - 1: Sub-cloud layer (CUBASE)
1498 !------------------------------------------------
1499 ! 1. INITIALIZE VALUES AT LIFTING LEVEL
1500 !------------------------------------------------
1506 PUU(JL,KLEV)=PUEN(JL,KLEV)*(PAPH(JL,KLEVP1)-PAPH(JL,KLEV))
1507 PVU(JL,KLEV)=PVEN(JL,KLEV)*(PAPH(JL,KLEVP1)-PAPH(JL,KLEV))
1509 !-------------------------------------------------------
1510 ! 2.0 DO ASCENT IN SUBCLOUD LAYER,
1511 ! CHECK FOR EXISTENCE OF CONDENSATION LEVEL,
1512 ! ADJUST T,Q AND L ACCORDINGLY IN *CUADJTQ*,
1513 ! CHECK FOR BUOYANCY AND SET FLAGS
1514 !-------------------------------------------------------
1519 DO 290 JK=KLEVM1,2,-1
1522 IF(KLAB(JL,JK+1).EQ.1) THEN
1530 IF(IS.EQ.0) GO TO 290
1533 PQU(JL,JK)=PQU(JL,JK+1)
1534 PTU(JL,JK)=(CPD*PTU(JL,JK+1)+PGEOH(JL,JK+1) &
1536 ZBUO=PTU(JL,JK)*(1.+VTMPC1*PQU(JL,JK))- &
1537 PTENH(JL,JK)*(1.+VTMPC1*PQENH(JL,JK))+ZBUO0
1538 IF(ZBUO.GT.0.) KLAB(JL,JK)=1
1539 ZQOLD(JL,JK)=PQU(JL,JK)
1544 CALL CUADJTQ(KLON,KLEV,IK,ZPH,PTU,PQU,LOFLAG,ICALL)
1546 IF(LOFLAG(JL).AND.PQU(JL,JK).NE.ZQOLD(JL,JK)) THEN
1548 PLU(JL,JK)=PLU(JL,JK)+ZQOLD(JL,JK)-PQU(JL,JK)
1549 ZBUO=PTU(JL,JK)*(1.+VTMPC1*PQU(JL,JK))- &
1550 PTENH(JL,JK)*(1.+VTMPC1*PQENH(JL,JK))+ZBUO0
1557 ! CALCULATE AVERAGES OF U AND V FOR SUBCLOUD ARA,.
1558 ! THE VALUES WILL BE USED TO DEFINE CLOUD BASE VALUES.
1561 IF(JK.GE.KCBOT(JL)) THEN
1562 PUU(JL,KLEV)=PUU(JL,KLEV)+ &
1563 PUEN(JL,JK)*(PAPH(JL,JK+1)-PAPH(JL,JK))
1564 PVU(JL,KLEV)=PVU(JL,KLEV)+ &
1565 PVEN(JL,JK)*(PAPH(JL,JK+1)-PAPH(JL,JK))
1574 ZZ=1./(PAPH(JL,KLEVP1)-PAPH(JL,IKB))
1575 PUU(JL,KLEV)=PUU(JL,KLEV)*ZZ
1576 PVU(JL,KLEV)=PVU(JL,KLEV)*ZZ
1578 PUU(JL,KLEV)=PUEN(JL,KLEVM1)
1579 PVU(JL,KLEV)=PVEN(JL,KLEVM1)
1584 END SUBROUTINE CUBASE
1586 !**********************************************
1588 !**********************************************
1590 ( KLON, KLEV, KLEVP1, KLEVM1,&
1591 PTENH, PQENH, PQSENH, PGEOH, PAPH,&
1592 RHO, HFX, QFX, KTYPE, lndj )
1593 ! THIS ROUTINE CALCULATES CLOUD BASE and TOP
1594 ! AND RETURN CLOUD TYPES
1595 ! ZHANG & WANG IPRC 12/2010
1598 ! TO PRODUCE CLOUD TYPE for CU-PARAMETERIZATIONS
1601 ! THIS ROUTINE IS CALLED FROM *CUMASTR*.
1602 ! INPUT ARE ENVIRONM. VALUES OF T,Q,P,PHI AT HALF LEVELS.
1603 ! IT RETURNS CLOUD TYPES AS FOLLOWS;
1604 ! KTYPE=1 FOR deep cumulus
1605 ! KTYPE=2 FOR shallow cumulus
1608 ! based on a simplified updraught equation
1609 ! partial(Hup)/partial(z)=eta(H - Hup)
1610 ! eta is the entrainment rate for test parcel
1611 ! H stands for dry static energy or the total water specific humidity
1612 ! references: Christian Jakob, 2003: A new subcloud model for mass-flux convection schemes
1613 ! influence on triggering, updraft properties, and model climate, Mon.Wea.Rev.
1616 ! IFS Documentation - Cy33r1
1620 ! *CUADJTQ* FOR ADJUSTING T AND Q DUE TO CONDENSATION IN ASCENT
1621 ! ----------------------------------------------------------------
1622 !-------------------------------------------------------------------
1624 !-------------------------------------------------------------------
1625 INTEGER KLON, KLEV, KLEVP1
1627 INTEGER JL,JK,IS,IK,ICALL,IKB,LEVELS
1628 REAL PTENH(KLON,KLEV), PQENH(KLON,KLEV), &
1630 PGEOH(KLON,KLEV), PAPH(KLON,KLEVP1)
1632 REAL QFX(KLON),RHO(KLON),HFX(KLON)
1633 REAL ZQOLD(KLON,KLEV), ZPH(KLON)
1634 INTEGER KCTOP(KLON),KCBOT(KLON)
1635 INTEGER KTYPE(KLON),LCLFLAG(KLON)
1636 LOGICAL TOPFLAG(KLON),DEEPFLAG(KLON),MYFLAG(KLON)
1638 REAL part1(klon), part2(klon), root(klon)
1639 REAL conw(klon),deltT(klon),deltQ(klon)
1640 REAL eta(klon),dz(klon),coef(klon)
1641 REAL dhen(KLON,KLEV), dh(KLON,KLEV),qh(KLON,KLEV)
1642 REAL Tup(KLON,KLEV),Qup(KLON,KLEV),ql(KLON,KLEV)
1643 REAL ww(KLON,KLEV),Kup(KLON,KLEV)
1644 REAL Vtup(KLON,KLEV),Vten(KLON,KLEV),buoy(KLON,KLEV)
1648 !***INPUT VARIABLES:
1649 ! PTENH [ZTENH] - Environment Temperature on half levels. (CUINI)
1650 ! PQENH [ZQENH] - Env. specific humidity on half levels. (CUINI)
1651 ! PGEOH [ZGEOH] - Geopotential on half levels, (MSSFLX)
1652 ! PAPH - Pressure of half levels. (MSSFLX)
1653 ! RHO - Density of the lowest Model level
1654 ! QFX - net upward moisture flux at the surface (kg/m^2/s)
1655 ! HFX - net upward heat flux at the surface (W/m^2)
1656 !***VARIABLES OUTPUT BY CUTYPE:
1657 ! KTYPE - Convection type - 1: Penetrative (CUMASTR)
1658 ! 2: Stratocumulus (CUMASTR)
1659 ! 3: Mid-level (CUASC)
1660 !--------------------------------------------------------------
1666 !-----------------------------------------------------------
1667 ! let's do test,and check the shallow convection first
1668 ! the first level is JK+1
1669 ! define deltaT and deltaQ
1670 !-----------------------------------------------------------
1674 ql(jl,jk)=0.0 ! parcel liquid water
1675 Tup(jl,jk)=0.0 ! parcel temperature
1676 Qup(jl,jk)=0.0 ! parcel specific humidity
1677 dh(jl,jk)=0.0 ! parcel dry static energy
1678 qh(jl,jk)=0.0 ! parcel total water specific humidity
1679 ww(jl,jk)=0.0 ! parcel vertical speed (m/s)
1680 dhen(jl,jk)=0.0 ! environment dry static energy
1681 Kup(jl,jk)=0.0 ! updraught kinetic energy for parcel
1682 Vtup(jl,jk)=0.0 ! parcel virtual temperature considering water-loading
1683 Vten(jl,jk)=0.0 ! environment virtual temperature
1684 buoy(jl,jk)=0.0 ! parcel buoyancy
1689 lclflag(jl) = 0 ! flag for the condensation level
1690 conw(jl) = 0.0 ! convective-scale velocity,also used for the vertical speed at the first level
1691 myflag(jl) = .true. ! just as input for cuadjqt subroutine
1692 topflag(jl) = .false.! flag for whether the cloud top is found
1695 ! check the levels from lowest level to second top level
1701 ! define the variables at the first level
1702 if(jk .eq. KLEVM1) then
1704 part1(jl) = 1.5*0.4*pgeoh(jl,jk+1)/(rho(jl)*ptenh(jl,jk+1))
1705 part2(jl) = hfx(jl)/cpd+0.61*ptenh(jl,jk+1)*qfx(jl)
1706 root(jl) = 0.001-part1(jl)*part2(jl)
1707 if(root(jl) .gt. 0) then
1708 conw(jl) = 1.2*(root(jl))**(1.0/3.0)
1710 conw(jl) = -1.2*(-root(jl))**(1.0/3.0)
1712 deltT(jl) = -1.5*hfx(jl)/(rho(jl)*cpd*conw(jl))
1713 deltQ(jl) = -1.5*qfx(jl)/(rho(jl)*conw(jl))
1715 Tup(jl,jk+1) = ptenh(jl,jk+1) + deltT(jl)
1716 Qup(jl,jk+1) = pqenh(jl,jk+1) + deltQ(jl)
1718 dh(jl,jk+1) = pgeoh(jl,jk+1) + Tup(jl,jk+1)*cpd
1719 qh(jl,jk+1) = pqenh(jl,jk+1) + deltQ(jl) + ql(jl,jk+1)
1720 ww(jl,jk+1) = conw(jl)
1724 ! the next levels, we use the variables at the first level as initial values
1726 if(.not. topflag(jl)) then
1727 eta(jl) = 0.5*(0.55/(pgeoh(jl,jk)*zrg)+1.0e-3)
1728 dz(jl) = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg
1729 coef(jl)= eta(jl)*dz(jl)
1730 dhen(jl,jk) = pgeoh(jl,jk) + cpd*ptenh(jl,jk)
1731 dh(jl,jk) = (coef(jl)*dhen(jl,jk) + dh(jl,jk+1))/(1+coef(jl))
1732 qh(jl,jk) = (coef(jl)*pqenh(jl,jk)+ qh(jl,jk+1))/(1+coef(jl))
1733 Tup(jl,jk) = (dh(jl,jk)-pgeoh(jl,jk))*RCPD
1734 Qup(jl,jk) = qh(jl,jk) - ql(jl,jk+1)
1735 zqold(jl,jk) = Qup(jl,jk)
1738 ! check if the parcel is saturated
1741 call CUADJTQ(klon,klev,ik,zph,Tup,Qup,myflag,icall)
1743 if( .not. topflag(jl) .and. zqold(jl,jk) .ne. Qup(jl,jk) ) then
1744 lclflag(jl) = lclflag(jl) + 1
1745 ql(jl,jk) = ql(jl,jk+1) + zqold(jl,jk) - Qup(jl,jk)
1746 dh(jl,jk) = pgeoh(jl,jk) + cpd*Tup(jl,jk)
1750 ! compute the updraft speed
1752 if(.not. topflag(jl))then
1753 Kup(jl,jk+1) = 0.5*ww(jl,jk+1)**2
1754 Vtup(jl,jk) = Tup(jl,jk)*(1.+VTMPC1*Qup(jl,jk)-ql(jl,jk))
1755 Vten(jl,jk) = ptenh(jl,jk)*(1.+VTMPC1*pqenh(jl,jk))
1756 buoy(jl,jk) = (Vtup(jl,jk) - Vten(jl,jk))/Vten(jl,jk)*g
1757 Kup(jl,jk) = (Kup(jl,jk+1) + 0.333*dz(jl)*buoy(jl,jk))/ &
1758 (1+2*2*eta(jl)*dz(jl))
1759 if(Kup(jl,jk) .gt. 0 ) then
1760 ww(jl,jk) = sqrt(2*Kup(jl,jk))
1761 if(lclflag(jl) .eq. 1 ) kcbot(jl) = jk
1769 topflag(jl) = .true.
1773 end do ! end all the levels
1776 if(paph(jl,kcbot(jl)) - paph(jl,kctop(jl)) .lt. ZDNOPRC .and. &
1777 paph(jl,kcbot(jl)) - paph(jl,kctop(jl)) .gt. 0 &
1778 .and. lclflag(jl) .gt. 0) then
1783 !-----------------------------------------------------------
1784 ! Next, let's check the deep convection
1785 ! the first level is JK
1786 ! define deltaT and deltaQ
1787 !----------------------------------------------------------
1788 ! we check the parcel starting level by level (from the second lowest level to the next 12th level,
1789 ! usually, the 12th level around 700 hPa for common eta levels)
1790 do levels=KLEVM1-1,KLEVM1-12,-1
1794 ql(jl,jk)=0.0 ! parcel liquid water
1795 Tup(jl,jk)=0.0 ! parcel temperature
1796 Qup(jl,jk)=0.0 ! parcel specific humidity
1797 dh(jl,jk)=0.0 ! parcel dry static energy
1798 qh(jl,jk)=0.0 ! parcel total water specific humidity
1799 ww(jl,jk)=0.0 ! parcel vertical speed (m/s)
1800 dhen(jl,jk)=0.0 ! environment dry static energy
1801 Kup(jl,jk)=0.0 ! updraught kinetic energy for parcel
1802 Vtup(jl,jk)=0.0 ! parcel virtual temperature considering water-loading
1803 Vten(jl,jk)=0.0 ! environment virtual temperature
1804 buoy(jl,jk)=0.0 ! parcel buoyancy
1809 lclflag(jl) = 0 ! flag for the condensation level
1812 myflag(jl) = .true. ! just as input for cuadjqt subroutine
1813 topflag(jl) = .false.! flag for whether the cloud top is found
1816 ! check the levels from lowest level to second top level
1822 ! define the variables at the first level
1823 if(jk .eq. levels) then
1828 if(paph(jl,KLEVM1-1)-paph(jl,jk) .le. 6.e3) then
1830 Tup(jl,jk+1) = 0.25*(ptenh(jl,jk+1)+ptenh(jl,jk)+ &
1831 ptenh(jl,jk-1)+ptenh(jl,jk-2)) + &
1833 dh(jl,jk+1) = 0.25*(pgeoh(jl,jk+1)+pgeoh(jl,jk)+ &
1834 pgeoh(jl,jk-1)+pgeoh(jl,jk-2)) + &
1836 qh(jl,jk+1) = 0.25*(pqenh(jl,jk+1)+pqenh(jl,jk)+ &
1837 pqenh(jl,jk-1)+pqenh(jl,jk-2))+ &
1838 deltQ(jl) + ql(jl,jk+1)
1839 Qup(jl,jk+1) = qh(jl,jk+1) - ql(jl,jk+1)
1842 Tup(jl,jk+1) = ptenh(jl,jk+1) + deltT(jl)
1843 dh(jl,jk+1) = pgeoh(jl,jk+1) + Tup(jl,jk+1)*cpd
1844 qh(jl,jk+1) = pqenh(jl,jk+1) + deltQ(jl)
1845 Qup(jl,jk+1) = qh(jl,jk+1) - ql(jl,jk+1)
1852 ! the next levels, we use the variables at the first level as initial values
1854 if(.not. topflag(jl)) then
1856 dz(jl) = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg
1857 coef(jl)= eta(jl)*dz(jl)
1858 dhen(jl,jk) = pgeoh(jl,jk) + cpd*ptenh(jl,jk)
1859 dh(jl,jk) = (coef(jl)*dhen(jl,jk) + dh(jl,jk+1))/(1+coef(jl))
1860 qh(jl,jk) = (coef(jl)*pqenh(jl,jk)+ qh(jl,jk+1))/(1+coef(jl))
1861 Tup(jl,jk) = (dh(jl,jk)-pgeoh(jl,jk))*RCPD
1862 Qup(jl,jk) = qh(jl,jk) - ql(jl,jk+1)
1863 zqold(jl,jk) = Qup(jl,jk)
1866 ! check if the parcel is saturated
1869 call CUADJTQ(klon,klev,ik,zph,Tup,Qup,myflag,icall)
1871 if( .not. topflag(jl) .and. zqold(jl,jk) .ne. Qup(jl,jk) ) then
1872 lclflag(jl) = lclflag(jl) + 1
1873 ql(jl,jk) = ql(jl,jk+1) + zqold(jl,jk) - Qup(jl,jk)
1874 dh(jl,jk) = pgeoh(jl,jk) + cpd*Tup(jl,jk)
1878 ! compute the updraft speed
1880 if(.not. topflag(jl))then
1881 Kup(jl,jk+1) = 0.5*ww(jl,jk+1)**2
1882 Vtup(jl,jk) = Tup(jl,jk)*(1.+VTMPC1*Qup(jl,jk)-ql(jl,jk))
1883 Vten(jl,jk) = ptenh(jl,jk)*(1.+VTMPC1*pqenh(jl,jk))
1884 buoy(jl,jk) = (Vtup(jl,jk) - Vten(jl,jk))/Vten(jl,jk)*g
1885 Kup(jl,jk) = (Kup(jl,jk+1) + 0.333*dz(jl)*buoy(jl,jk))/ &
1886 (1+2*2*eta(jl)*dz(jl))
1887 if(Kup(jl,jk) .gt. 0 ) then
1888 ww(jl,jk) = sqrt(2*Kup(jl,jk))
1889 if(lclflag(jl) .eq. 1 ) kcbot(jl) = jk
1897 topflag(jl) = .true.
1901 end do ! end all the levels
1904 if(paph(jl,kcbot(jl)) - paph(jl,kctop(jl)) .gt. ZDNOPRC .and. &
1905 lclflag(jl) .gt. 0 ) then
1907 do jk=kcbot(jl),kctop(jl),-1
1908 ZRELH(JL)=ZRELH(JL)+ PQENH(JL,JK)/PQSENH(JL,JK)
1910 ZRELH(JL) = ZRELH(JL)/(kcbot(jl)-kctop(jl)+1)
1912 if(lndj(JL) .eq. 1) then
1917 if(ZRELH(JL) .ge. CRIRH1) ktype(jl) = 1
1921 end do ! end all cycles
1923 END SUBROUTINE CUTYPE
1926 !**********************************************
1927 ! SUBROUTINE CUASC_NEW
1928 !**********************************************
1929 SUBROUTINE CUASC_NEW &
1930 (KLON, KLEV, KLEVP1, KLEVM1, PTENH, &
1931 PQENH, PUEN, PVEN, PTEN, PQEN, &
1932 PQSEN, PGEO, PGEOH, PAP, PAPH, &
1933 PQTE, PVERV, KLWMIN, LDCUM, PHCBASE,&
1934 KTYPE, KLAB, PTU, PQU, PLU, &
1935 PUU, PVU, PMFU, PMFUB, PENTR, &
1936 PMFUS, PMFUQ, PMFUL, PLUDE, PDMFUP, &
1937 KCBOT, KCTOP, KCTOP0, KCUM, ZTMST, &
1938 KHMIN, PHHATT, PQSENH)
1939 ! THIS ROUTINE DOES THE CALCULATIONS FOR CLOUD ASCENTS
1940 ! FOR CUMULUS PARAMETERIZATION
1941 ! M.TIEDTKE E.C.M.W.F. 7/86 MODIF. 12/89
1942 ! Y.WANG IPRC 11/01 MODIF.
1945 ! TO PRODUCE CLOUD ASCENTS FOR CU-PARAMETRIZATION
1946 ! (VERTICAL PROFILES OF T,Q,L,U AND V AND CORRESPONDING
1947 ! FLUXES AS WELL AS PRECIPITATION RATES)
1950 ! THIS ROUTINE IS CALLED FROM *CUMASTR*.
1953 ! LIFT SURFACE AIR DRY-ADIABATICALLY TO CLOUD BASE
1954 ! AND THEN CALCULATE MOIST ASCENT FOR
1955 ! ENTRAINING/DETRAINING PLUME.
1956 ! ENTRAINMENT AND DETRAINMENT RATES DIFFER FOR
1957 ! SHALLOW AND DEEP CUMULUS CONVECTION.
1958 ! IN CASE THERE IS NO PENETRATIVE OR SHALLOW CONVECTION
1959 ! CHECK FOR POSSIBILITY OF MID LEVEL CONVECTION
1960 ! (CLOUD BASE VALUES CALCULATED IN *CUBASMC*)
1963 ! *CUADJTQ* ADJUST T AND Q DUE TO CONDENSATION IN ASCENT
1964 ! *CUENTR_NEW* CALCULATE ENTRAINMENT/DETRAINMENT RATES
1965 ! *CUBASMC* CALCULATE CLOUD BASE VALUES FOR MIDLEVEL CONVECTION
1969 !***INPUT VARIABLES:
1970 ! PTENH [ZTENH] - Environ Temperature on half levels. (CUINI)
1971 ! PQENH [ZQENH] - Env. specific humidity on half levels. (CUINI)
1972 ! PUEN - Environment wind u-component. (MSSFLX)
1973 ! PVEN - Environment wind v-component. (MSSFLX)
1974 ! PTEN - Environment Temperature. (MSSFLX)
1975 ! PQEN - Environment Specific Humidity. (MSSFLX)
1976 ! PQSEN - Environment Saturation Specific Humidity. (MSSFLX)
1977 ! PGEO - Geopotential. (MSSFLX)
1978 ! PGEOH [ZGEOH] - Geopotential on half levels, (MSSFLX)
1979 ! PAP - Pressure in Pa. (MSSFLX)
1980 ! PAPH - Pressure of half levels. (MSSFLX)
1981 ! PQTE - Moisture convergence (Delta q/Delta t). (MSSFLX)
1982 ! PVERV - Large Scale Vertical Velocity (Omega). (MSSFLX)
1983 ! KLWMIN [ILWMIN] - Level of Minimum Omega. (CUINI)
1984 ! KLAB [ILAB] - Level Label - 1: Sub-cloud layer.
1985 ! 2: Condensation Level (Cloud Base)
1986 ! PMFUB [ZMFUB] - Updraft Mass Flux at Cloud Base. (CUMASTR)
1987 !***VARIABLES MODIFIED BY CUASC:
1988 ! LDCUM - Logical denoting profiles. (CUBASE)
1989 ! KTYPE - Convection type - 1: Penetrative (CUMASTR)
1990 ! 2: Stratocumulus (CUMASTR)
1991 ! 3: Mid-level (CUASC)
1992 ! PTU - Cloud Temperature.
1993 ! PQU - Cloud specific Humidity.
1994 ! PLU - Cloud Liquid Water (Moisture condensed out)
1995 ! PUU [ZUU] - Cloud Momentum U-Component.
1996 ! PVU [ZVU] - Cloud Momentum V-Component.
1997 ! PMFU - Updraft Mass Flux.
1998 ! PENTR [ZENTR] - Entrainment Rate. (CUMASTR ) (CUBASMC)
1999 ! PMFUS [ZMFUS] - Updraft Flux of Dry Static Energy. (CUBASMC)
2000 ! PMFUQ [ZMFUQ] - Updraft Flux of Specific Humidity.
2001 ! PMFUL [ZMFUL] - Updraft Flux of Cloud Liquid Water.
2002 ! PLUDE - Liquid Water Returned to Environment by Detrainment.
2003 ! PDMFUP [ZMFUP] - FLUX DIFFERENCE OF PRECIP. IN UPDRAFTS
2004 ! KCBOT - Cloud Base Level. (CUBASE)
2006 ! KCTOP0 [ICTOP0] - Estimate of Cloud Top. (CUMASTR)
2008 !-------------------------------------------------------------------
2010 !-------------------------------------------------------------------
2011 INTEGER KLON, KLEV, KLEVP1
2013 REAL ZTMST,ZCONS2,ZDZ,ZDRODZ
2014 INTEGER JL,JK,IKB,IK,IS,IKT,ICALL
2015 REAL ZMFMAX,ZFAC,ZMFTEST,ZDPRHO,ZMSE,ZNEVN,ZODMAX
2016 REAL ZQEEN,ZSEEN,ZSCDE,ZGA,ZDT,ZSCOD
2017 REAL ZQUDE,ZQCOD, ZMFUSK, ZMFUQK,ZMFULK
2018 REAL ZBUO, ZPRCON, ZLNEW, ZZ, ZDMFEU, ZDMFDU
2020 REAL PTENH(KLON,KLEV), PQENH(KLON,KLEV), &
2021 PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
2022 PTEN(KLON,KLEV), PQEN(KLON,KLEV), &
2023 PGEO(KLON,KLEV), PGEOH(KLON,KLEV), &
2024 PAP(KLON,KLEV), PAPH(KLON,KLEVP1), &
2025 PQSEN(KLON,KLEV), PQTE(KLON,KLEV), &
2026 PVERV(KLON,KLEV), PQSENH(KLON,KLEV)
2027 REAL PTU(KLON,KLEV), PQU(KLON,KLEV), &
2028 PUU(KLON,KLEV), PVU(KLON,KLEV), &
2029 PMFU(KLON,KLEV), ZPH(KLON), &
2030 PMFUB(KLON), PENTR(KLON), &
2031 PMFUS(KLON,KLEV), PMFUQ(KLON,KLEV), &
2032 PLU(KLON,KLEV), PLUDE(KLON,KLEV), &
2033 PMFUL(KLON,KLEV), PDMFUP(KLON,KLEV)
2034 REAL ZDMFEN(KLON), ZDMFDE(KLON), &
2035 ZMFUU(KLON), ZMFUV(KLON), &
2036 ZPBASE(KLON), ZQOLD(KLON), &
2037 PHHATT(KLON,KLEV), ZODETR(KLON,KLEV), &
2038 ZOENTR(KLON,KLEV), ZBUOY(KLON)
2040 INTEGER KLWMIN(KLON), KTYPE(KLON), &
2041 KLAB(KLON,KLEV), KCBOT(KLON), &
2042 KCTOP(KLON), KCTOP0(KLON), &
2044 LOGICAL LDCUM(KLON), LOFLAG(KLON)
2045 integer leveltop,levelbot
2046 real tt(klon),ttb(klon)
2047 real zqsat(klon), zqsatb(klon)
2050 !--------------------------------
2051 !* 1. SPECIFY PARAMETERS
2052 !--------------------------------
2055 !---------------------------------
2056 ! 2. SET DEFAULT VALUES
2057 !---------------------------------
2063 IF(.NOT.LDCUM(JL)) KTYPE(JL)=0
2076 IF(.NOT.LDCUM(JL).OR.KTYPE(JL).EQ.3) KLAB(JL,JK)=0
2077 IF(.NOT.LDCUM(JL).AND.PAPH(JL,JK).LT.4.E4) KCTOP0(JL)=JK
2079 !------------------------------------------------
2080 ! 3.0 INITIALIZE VALUES AT LIFTING LEVEL
2081 !------------------------------------------------
2084 IF(.NOT.LDCUM(JL)) THEN
2089 PMFU(JL,KLEV)=PMFUB(JL)
2090 PMFUS(JL,KLEV)=PMFUB(JL)*(CPD*PTU(JL,KLEV)+PGEOH(JL,KLEV))
2091 PMFUQ(JL,KLEV)=PMFUB(JL)*PQU(JL,KLEV)
2093 ZMFUU(JL)=PMFUB(JL)*PUU(JL,KLEV)
2094 ZMFUV(JL)=PMFUB(JL)*PVU(JL,KLEV)
2098 !-- 3.1 Find organized entrainment at cloud base
2102 IF (KTYPE(JL).EQ.1) THEN
2104 if(orgen .eq. 1 ) then
2106 ZBUOY(JL)=G*((PTU(JL,IKB)-PTENH(JL,IKB))/PTENH(JL,IKB)+ &
2107 0.608*(PQU(JL,IKB)-PQENH(JL,IKB)))
2108 IF (ZBUOY(JL).GT.0.) THEN
2109 ZDZ = (PGEO(JL,IKB-1)-PGEO(JL,IKB))*ZRG
2110 ZDRODZ = -LOG(PTEN(JL,IKB-1)/PTEN(JL,IKB))/ZDZ - &
2111 G/(RD*PTENH(JL,IKB))
2112 ZOENTR(JL,IKB-1)=ZBUOY(JL)*0.5/(1.+ZBUOY(JL)*ZDZ) &
2114 ZOENTR(JL,IKB-1) = MIN(ZOENTR(JL,IKB-1),1.E-3)
2115 ZOENTR(JL,IKB-1) = MAX(ZOENTR(JL,IKB-1),0.)
2118 ! Let's define the fscale
2119 else if(orgen .eq. 2 ) then
2120 tt(jl) = ptenh(jl,ikb)
2121 zqsat(jl) = TLUCUA(tt(jl))/paph(jl,ikb-1)
2122 zqsat(jl) = zqsat(jl)/(1.-VTMPC1*zqsat(jl))
2123 ttb(jl) = ptenh(jl,ikb)
2124 zqsatb(jl) = TLUCUA(ttb(jl))/paph(jl,ikb)
2125 zqsatb(jl) = zqsatb(jl)/(1.-VTMPC1*zqsatb(jl))
2126 fscale(jl) = (zqsat(jl)/zqsatb(jl))**3
2127 ! end of defining the fscale
2128 zoentr(jl,ikb-1) = 1.E-3*(1.3-PQEN(jl,ikb-1)/PQSEN(jl,ikb-1))*fscale(jl)
2129 zoentr(jl,ikb-1) = MIN(zoentr(jl,ikb-1),1.E-3)
2130 zoentr(jl,ikb-1) = MAX(zoentr(jl,ikb-1),0.)
2135 !-----------------------------------------------------------------
2136 ! 4. DO ASCENT: SUBCLOUD LAYER (KLAB=1) ,CLOUDS (KLAB=2)
2137 ! BY DOING FIRST DRY-ADIABATIC ASCENT AND THEN
2138 ! BY ADJUSTING T,Q AND L ACCORDINGLY IN *CUADJTQ*,
2139 ! THEN CHECK FOR BUOYANCY AND SET FLAGS ACCORDINGLY
2140 !-----------------------------------------------------------------
2143 ! let's define the levels in which the middle level convection could be activated
2145 if(abs(paph(1,jk)*0.01 - 250) .lt. 50.) then
2150 leveltop = min(KLEV-15,leveltop)
2151 levelbot = KLEVM1 - 4
2153 DO 480 JK=KLEVM1,2,-1
2154 ! SPECIFY CLOUD BASE VALUES FOR MIDLEVEL CONVECTION
2155 ! IN *CUBASMC* IN CASE THERE IS NOT ALREADY CONVECTION
2156 ! ---------------------------------------------------------------------
2158 IF(LMFMID.AND.IK.LT.levelbot.AND.IK.GT.leveltop) THEN
2160 (KLON, KLEV, KLEVM1, IK, PTEN, &
2161 PQEN, PQSEN, PUEN, PVEN, PVERV, &
2162 PGEO, PGEOH, LDCUM, KTYPE, KLAB, &
2163 PMFU, PMFUB, PENTR, KCBOT, PTU, &
2164 PQU, PLU, PUU, PVU, PMFUS, &
2165 PMFUQ, PMFUL, PDMFUP, ZMFUU, ZMFUV)
2171 IF(KLAB(JL,JK+1).EQ.0) KLAB(JL,JK)=0
2172 LOFLAG(JL)=KLAB(JL,JK+1).GT.0
2174 IF(KTYPE(JL).EQ.3.AND.JK.EQ.KCBOT(JL)) THEN
2175 ZMFMAX=(PAPH(JL,JK)-PAPH(JL,JK-1))*ZCONS2
2176 IF(PMFUB(JL).GT.ZMFMAX) THEN
2177 ZFAC=ZMFMAX/PMFUB(JL)
2178 PMFU(JL,JK+1)=PMFU(JL,JK+1)*ZFAC
2179 PMFUS(JL,JK+1)=PMFUS(JL,JK+1)*ZFAC
2180 PMFUQ(JL,JK+1)=PMFUQ(JL,JK+1)*ZFAC
2181 ZMFUU(JL)=ZMFUU(JL)*ZFAC
2182 ZMFUV(JL)=ZMFUV(JL)*ZFAC
2187 IF(IS.EQ.0) GO TO 480
2189 !* SPECIFY ENTRAINMENT RATES IN *CUENTR_NEW*
2190 ! -------------------------------------
2193 (KLON, KLEV, KLEVP1, IK, PTENH,&
2194 PAPH, PAP, PGEOH, KLWMIN, LDCUM,&
2195 KTYPE, KCBOT, KCTOP0, ZPBASE, PMFU, &
2196 PENTR, ZDMFEN, ZDMFDE, ZODETR, KHMIN)
2198 ! DO ADIABATIC ASCENT FOR ENTRAINING/DETRAINING PLUME
2199 ! -------------------------------------------------------
2200 ! Do adiabatic ascent for entraining/detraining plume
2201 ! the cloud ensemble entrains environmental values
2202 ! in turbulent detrainment cloud ensemble values are detrained
2203 ! in organized detrainment the dry static energy and
2204 ! moisture that are neutral compared to the
2205 ! environmental air are detrained
2209 IF(JK.LT.KCBOT(JL)) THEN
2210 ZMFTEST=PMFU(JL,JK+1)+ZDMFEN(JL)-ZDMFDE(JL)
2211 ZMFMAX=MIN(ZMFTEST,(PAPH(JL,JK)-PAPH(JL,JK-1))*ZCONS2)
2212 ZDMFEN(JL)=MAX(ZDMFEN(JL)-MAX(ZMFTEST-ZMFMAX,0.),0.)
2214 ZDMFDE(JL)=MIN(ZDMFDE(JL),0.75*PMFU(JL,JK+1))
2215 PMFU(JL,JK)=PMFU(JL,JK+1)+ZDMFEN(JL)-ZDMFDE(JL)
2216 IF (JK.LT.kcbot(jl)) THEN
2217 zdprho = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg
2218 zoentr(jl,jk) = zoentr(jl,jk)*zdprho*pmfu(jl,jk+1)
2219 zmftest = pmfu(jl,jk) + zoentr(jl,jk)-zodetr(jl,jk)
2220 zmfmax = MIN(zmftest,(paph(jl,jk)-paph(jl,jk-1))*zcons2)
2221 zoentr(jl,jk) = MAX(zoentr(jl,jk)-MAX(zmftest-zmfmax,0.),0.)
2224 ! limit organized detrainment to not allowing for too deep clouds
2226 IF (ktype(jl).EQ.1.AND.jk.LT.kcbot(jl).AND.jk.LE.khmin(jl)) THEN
2227 zmse = cpd*ptu(jl,jk+1) + alv*pqu(jl,jk+1) + pgeoh(jl,jk+1)
2229 znevn=(pgeoh(jl,ikt)-pgeoh(jl,jk+1))*(zmse-phhatt(jl, &
2231 IF (znevn.LE.0.) znevn = 1.
2232 zdprho = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg
2233 zodmax = ((phcbase(jl)-zmse)/znevn)*zdprho*pmfu(jl,jk+1)
2234 zodmax = MAX(zodmax,0.)
2235 zodetr(jl,jk) = MIN(zodetr(jl,jk),zodmax)
2237 zodetr(jl,jk) = MIN(zodetr(jl,jk),0.75*pmfu(jl,jk))
2238 pmfu(jl,jk) = pmfu(jl,jk) + zoentr(jl,jk) - zodetr(jl,jk)
2239 ZQEEN=PQENH(JL,JK+1)*ZDMFEN(JL)
2240 zqeen=zqeen + pqenh(jl,jk+1)*zoentr(jl,jk)
2241 ZSEEN=(CPD*PTENH(JL,JK+1)+PGEOH(JL,JK+1))*ZDMFEN(JL)
2242 zseen=zseen+(cpd*ptenh(jl,jk+1)+pgeoh(jl,jk+1))* &
2244 ZSCDE=(CPD*PTU(JL,JK+1)+PGEOH(JL,JK+1))*ZDMFDE(JL)
2245 ! find moist static energy that give nonbuoyant air
2246 zga = alv*pqsenh(jl,jk+1)/(rv*(ptenh(jl,jk+1)**2))
2247 zdt = (plu(jl,jk+1)-0.608*(pqsenh(jl,jk+1)-pqenh(jl, &
2248 jk+1)))/(1./ptenh(jl,jk+1)+0.608*zga)
2249 zscod = cpd*ptenh(jl,jk+1) + pgeoh(jl,jk+1) + cpd*zdt
2250 zscde = zscde + zodetr(jl,jk)*zscod
2251 zqude = pqu(jl,jk+1)*zdmfde(jl)
2252 zqcod = pqsenh(jl,jk+1) + zga*zdt
2253 zqude = zqude + zodetr(jl,jk)*zqcod
2254 plude(jl,jk) = plu(jl,jk+1)*zdmfde(jl)
2255 plude(jl,jk) = plude(jl,jk)+plu(jl,jk+1)*zodetr(jl,jk)
2256 zmfusk = pmfus(jl,jk+1) + zseen - zscde
2257 zmfuqk = pmfuq(jl,jk+1) + zqeen - zqude
2258 zmfulk = pmful(jl,jk+1) - plude(jl,jk)
2259 plu(jl,jk) = zmfulk*(1./MAX(cmfcmin,pmfu(jl,jk)))
2260 pqu(jl,jk) = zmfuqk*(1./MAX(cmfcmin,pmfu(jl,jk)))
2261 ptu(jl,jk)=(zmfusk*(1./MAX(cmfcmin,pmfu(jl,jk)))- &
2263 ptu(jl,jk) = MAX(100.,ptu(jl,jk))
2264 ptu(jl,jk) = MIN(400.,ptu(jl,jk))
2265 zqold(jl) = pqu(jl,jk)
2268 !* DO CORRECTIONS FOR MOIST ASCENT
2269 !* BY ADJUSTING T,Q AND L IN *CUADJTQ*
2270 !------------------------------------------------
2274 CALL CUADJTQ(KLON,KLEV,IK,ZPH,PTU,PQU,LOFLAG,ICALL)
2277 IF(LOFLAG(JL).AND.PQU(JL,JK).NE.ZQOLD(JL)) THEN
2279 PLU(JL,JK)=PLU(JL,JK)+ZQOLD(JL)-PQU(JL,JK)
2280 ZBUO=PTU(JL,JK)*(1.+VTMPC1*PQU(JL,JK)-PLU(JL,JK))- &
2281 PTENH(JL,JK)*(1.+VTMPC1*PQENH(JL,JK))
2282 IF(KLAB(JL,JK+1).EQ.1) ZBUO=ZBUO+ZBUO0
2283 IF(ZBUO.GT.0..AND.PMFU(JL,JK).GT.0.01*PMFUB(JL).AND. &
2284 JK.GE.KCTOP0(JL)) THEN
2287 IF(ZPBASE(JL)-PAPH(JL,JK).GE.ZDNOPRC) THEN
2292 ZLNEW=PLU(JL,JK)/(1.+ZPRCON*(PGEOH(JL,JK)-PGEOH(JL,JK+1)))
2293 PDMFUP(JL,JK)=MAX(0.,(PLU(JL,JK)-ZLNEW)*PMFU(JL,JK))
2301 PMFUL(JL,JK)=PLU(JL,JK)*PMFU(JL,JK)
2302 PMFUS(JL,JK)=(CPD*PTU(JL,JK)+PGEOH(JL,JK))*PMFU(JL,JK)
2303 PMFUQ(JL,JK)=PQU(JL,JK)*PMFU(JL,JK)
2310 zdmfen(jl) = zdmfen(jl) + zoentr(jl,jk)
2311 zdmfde(jl) = zdmfde(jl) + zodetr(jl,jk)
2313 IF(KTYPE(JL).EQ.1.OR.KTYPE(JL).EQ.3) THEN
2314 IF(ZDMFEN(JL).LE.1.E-20) THEN
2320 IF(ZDMFEN(JL).LE.1.0E-20) THEN
2326 ZDMFEU=ZDMFEN(JL)+ZZ*ZDMFDE(JL)
2327 ZDMFDU=ZDMFDE(JL)+ZZ*ZDMFDE(JL)
2328 ZDMFDU=MIN(ZDMFDU,0.75*PMFU(JL,JK+1))
2329 ZMFUU(JL)=ZMFUU(JL)+ &
2330 ZDMFEU*PUEN(JL,JK)-ZDMFDU*PUU(JL,JK+1)
2331 ZMFUV(JL)=ZMFUV(JL)+ &
2332 ZDMFEU*PVEN(JL,JK)-ZDMFDU*PVU(JL,JK+1)
2333 IF(PMFU(JL,JK).GT.0.) THEN
2334 PUU(JL,JK)=ZMFUU(JL)*(1./PMFU(JL,JK))
2335 PVU(JL,JK)=ZMFUV(JL)*(1./PMFU(JL,JK))
2342 ! Compute organized entrainment
2343 ! for use at next level
2346 IF (loflag(jl).AND.ktype(jl).EQ.1) THEN
2348 if(orgen .eq. 1 ) then
2349 zbuoyz=g*((ptu(jl,jk)-ptenh(jl,jk))/ptenh(jl,jk)+ &
2350 0.608*(pqu(jl,jk)-pqenh(jl,jk))-plu(jl,jk))
2351 zbuoyz = MAX(zbuoyz,0.0)
2352 zdz = (pgeo(jl,jk-1)-pgeo(jl,jk))*zrg
2353 zdrodz = -LOG(pten(jl,jk-1)/pten(jl,jk))/zdz - &
2355 zbuoy(jl) = zbuoy(jl) + zbuoyz*zdz
2356 zoentr(jl,jk-1) = zbuoyz*0.5/(1.+zbuoy(jl))+zdrodz
2357 zoentr(jl,jk-1) = MIN(zoentr(jl,jk-1),1.E-3)
2358 zoentr(jl,jk-1) = MAX(zoentr(jl,jk-1),0.)
2359 else if(orgen .eq. 2 ) then
2360 ! Let's define the fscale
2361 tt(jl) = ptenh(jl,jk-1)
2362 zqsat(jl) = TLUCUA(tt(jl))/paph(jl,jk-1)
2363 zqsat(jl) = zqsat(jl)/(1.-VTMPC1*zqsat(jl))
2364 ttb(jl) = ptenh(jl,kcbot(jl))
2365 zqsatb(jl) = TLUCUA(ttb(jl))/paph(jl,kcbot(jl))
2366 zqsatb(jl) = zqsatb(jl)/(1.-VTMPC1*zqsatb(jl))
2367 fscale(jl) = (zqsat(jl)/zqsatb(jl))**3
2368 ! end of defining the fscale
2369 zoentr(jl,jk-1) = 1.E-3*(1.3-PQEN(jl,jk-1)/PQSEN(jl,jk-1))*fscale(jl)
2370 zoentr(jl,jk-1) = MIN(zoentr(jl,jk-1),1.E-3)
2371 zoentr(jl,jk-1) = MAX(zoentr(jl,jk-1),0.)
2372 ! write(6,*) "zoentr=",zoentr(jl,jk-1)
2378 ! -----------------------------------------------------------------
2379 ! 5. DETERMINE CONVECTIVE FLUXES ABOVE NON-BUOYANCY LEVEL
2380 ! -----------------------------------------------------------------
2381 ! (NOTE: CLOUD VARIABLES LIKE T,Q AND L ARE NOT
2382 ! AFFECTED BY DETRAINMENT AND ARE ALREADY KNOWN
2383 ! FROM PREVIOUS CALCULATIONS ABOVE)
2386 IF(KCTOP(JL).EQ.KLEVM1) LDCUM(JL)=.FALSE.
2387 KCBOT(JL)=MAX(KCBOT(JL),KCTOP(JL))
2396 IF(IS.EQ.0) GO TO 800
2401 ZDMFDE(JL)=(1.-ZZDMF)*PMFU(JL,JK+1)
2402 PLUDE(JL,JK)=ZDMFDE(JL)*PLU(JL,JK+1)
2403 PMFU(JL,JK)=PMFU(JL,JK+1)-ZDMFDE(JL)
2404 PMFUS(JL,JK)=(CPD*PTU(JL,JK)+PGEOH(JL,JK))*PMFU(JL,JK)
2405 PMFUQ(JL,JK)=PQU(JL,JK)*PMFU(JL,JK)
2406 PMFUL(JL,JK)=PLU(JL,JK)*PMFU(JL,JK)
2407 PLUDE(JL,JK-1)=PMFUL(JL,JK)
2415 PUU(JL,JK)=PUU(JL,JK+1)
2416 PVU(JL,JK)=PVU(JL,JK+1)
2422 END SUBROUTINE CUASC_NEW
2425 !**********************************************
2427 !**********************************************
2429 (KLON, KLEV, KLEVP1, PTENH, PQENH, &
2430 PUEN, PVEN, PGEOH, PAPH, PTU, &
2431 PQU, PUU, PVU, LDCUM, KCBOT, &
2432 KCTOP, PMFUB, PRFL, PTD, PQD, &
2433 PUD, PVD, PMFD, PMFDS, PMFDQ, &
2434 PDMFDP, KDTOP, LDDRAF)
2435 ! THIS ROUTINE CALCULATES LEVEL OF FREE SINKING FOR
2436 ! CUMULUS DOWNDRAFTS AND SPECIFIES T,Q,U AND V VALUES
2437 ! M.TIEDTKE E.C.M.W.F. 12/86 MODIF. 12/89
2440 ! TO PRODUCE LFS-VALUES FOR CUMULUS DOWNDRAFTS
2441 ! FOR MASSFLUX CUMULUS PARAMETERIZATION
2444 ! THIS ROUTINE IS CALLED FROM *CUMASTR*.
2445 ! INPUT ARE ENVIRONMENTAL VALUES OF T,Q,U,V,P,PHI
2446 ! AND UPDRAFT VALUES T,Q,U AND V AND ALSO
2447 ! CLOUD BASE MASSFLUX AND CU-PRECIPITATION RATE.
2448 ! IT RETURNS T,Q,U AND V VALUES AND MASSFLUX AT LFS.
2451 ! CHECK FOR NEGATIVE BUOYANCY OF AIR OF EQUAL PARTS OF
2452 ! MOIST ENVIRONMENTAL AIR AND CLOUD AIR.
2455 ! *CUADJTQ* FOR CALCULATING WET BULB T AND Q AT LFS
2456 ! ----------------------------------------------------------------
2457 !-------------------------------------------------------------------
2459 !-------------------------------------------------------------------
2460 INTEGER KLON, KLEV, KLEVP1
2461 INTEGER JL,KE,JK,IS,IK,ICALL
2462 REAL ZTTEST, ZQTEST, ZBUO, ZMFTOP
2463 REAL PTENH(KLON,KLEV), PQENH(KLON,KLEV), &
2464 PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
2465 PGEOH(KLON,KLEV), PAPH(KLON,KLEVP1), &
2466 PTU(KLON,KLEV), PQU(KLON,KLEV), &
2467 PUU(KLON,KLEV), PVU(KLON,KLEV), &
2468 PMFUB(KLON), PRFL(KLON)
2469 REAL PTD(KLON,KLEV), PQD(KLON,KLEV), &
2470 PUD(KLON,KLEV), PVD(KLON,KLEV), &
2471 PMFD(KLON,KLEV), PMFDS(KLON,KLEV), &
2472 PMFDQ(KLON,KLEV), PDMFDP(KLON,KLEV)
2473 REAL ZTENWB(KLON,KLEV), ZQENWB(KLON,KLEV), &
2474 ZCOND(KLON), ZPH(KLON)
2475 INTEGER KCBOT(KLON), KCTOP(KLON), &
2477 LOGICAL LDCUM(KLON), LLo2(KLON), &
2479 !-----------------------------------------------
2480 ! 1. SET DEFAULT VALUES FOR DOWNDRAFTS
2481 !-----------------------------------------------
2487 IF(.NOT.LMFDD) GO TO 300
2488 !------------------------------------------------------------
2489 ! 2. DETERMINE LEVEL OF FREE SINKING BY
2490 ! DOING A SCAN FROM TOP TO BASE OF CUMULUS CLOUDS
2491 ! FOR EVERY POINT AND PROCEED AS FOLLOWS:
2492 ! (1) DETEMINE WET BULB ENVIRONMENTAL T AND Q
2493 ! (2) DO MIXING WITH CUMULUS CLOUD AIR
2494 ! (3) CHECK FOR NEGATIVE BUOYANCY
2495 ! THE ASSUMPTION IS THAT AIR OF DOWNDRAFTS IS MIXTURE
2496 ! OF 50% CLOUD AIR + 50% ENVIRONMENTAL AIR AT WET BULB
2497 ! TEMPERATURE (I.E. WHICH BECAME SATURATED DUE TO
2498 ! EVAPORATION OF RAIN AND CLOUD WATER)
2499 !------------------------------------------------------------------
2503 ! 2.1 CALCULATE WET-BULB TEMPERATURE AND MOISTURE
2504 ! FOR ENVIRONMENTAL AIR IN *CUADJTQ*
2505 ! -----------------------------------------------------
2509 ZTENWB(JL,JK)=PTENH(JL,JK)
2510 ZQENWB(JL,JK)=PQENH(JL,JK)
2512 LLO2(JL)=LDCUM(JL).AND.PRFL(JL).GT.0..AND..NOT.LDDRAF(JL).AND. &
2513 (JK.LT.KCBOT(JL).AND.JK.GT.KCTOP(JL))
2518 IF(IS.EQ.0) GO TO 290
2521 CALL CUADJTQ(KLON,KLEV,IK,ZPH,ZTENWB,ZQENWB,LLO2,ICALL)
2522 ! 2.2 DO MIXING OF CUMULUS AND ENVIRONMENTAL AIR
2523 ! AND CHECK FOR NEGATIVE BUOYANCY.
2524 ! THEN SET VALUES FOR DOWNDRAFT AT LFS.
2525 ! -----------------------------------------------------
2529 ZTTEST=0.5*(PTU(JL,JK)+ZTENWB(JL,JK))
2530 ZQTEST=0.5*(PQU(JL,JK)+ZQENWB(JL,JK))
2531 ZBUO=ZTTEST*(1.+VTMPC1*ZQTEST)- &
2532 PTENH(JL,JK)*(1.+VTMPC1*PQENH(JL,JK))
2533 ZCOND(JL)=PQENH(JL,JK)-ZQENWB(JL,JK)
2534 ZMFTOP=-CMFDEPS*PMFUB(JL)
2535 IF(ZBUO.LT.0..AND.PRFL(JL).GT.10.*ZMFTOP*ZCOND(JL)) THEN
2541 PMFDS(JL,JK)=PMFD(JL,JK)*(CPD*PTD(JL,JK)+PGEOH(JL,JK))
2542 PMFDQ(JL,JK)=PMFD(JL,JK)*PQD(JL,JK)
2543 PDMFDP(JL,JK-1)=-0.5*PMFD(JL,JK)*ZCOND(JL)
2544 PRFL(JL)=PRFL(JL)+PDMFDP(JL,JK-1)
2550 IF(PMFD(JL,JK).LT.0.) THEN
2551 PUD(JL,JK)=0.5*(PUU(JL,JK)+PUEN(JL,JK-1))
2552 PVD(JL,JK)=0.5*(PVU(JL,JK)+PVEN(JL,JK-1))
2559 END SUBROUTINE CUDLFS
2562 !**********************************************
2563 ! SUBROUTINE CUDDRAF
2564 !**********************************************
2565 SUBROUTINE CUDDRAF &
2566 (KLON, KLEV, KLEVP1, PTENH, PQENH, &
2567 PUEN, PVEN, PGEOH, PAPH, PRFL, &
2568 LDDRAF, PTD, PQD, PUD, PVD, &
2569 PMFD, PMFDS, PMFDQ, PDMFDP)
2570 ! THIS ROUTINE CALCULATES CUMULUS DOWNDRAFT DESCENT
2571 ! M.TIEDTKE E.C.M.W.F. 12/86 MODIF. 12/89
2574 ! TO PRODUCE THE VERTICAL PROFILES FOR CUMULUS DOWNDRAFTS
2575 ! (I.E. T,Q,U AND V AND FLUXES)
2578 ! THIS ROUTINE IS CALLED FROM *CUMASTR*.
2579 ! INPUT IS T,Q,P,PHI,U,V AT HALF LEVELS.
2580 ! IT RETURNS FLUXES OF S,Q AND EVAPORATION RATE
2581 ! AND U,V AT LEVELS WHERE DOWNDRAFT OCCURS
2584 ! CALCULATE MOIST DESCENT FOR ENTRAINING/DETRAINING PLUME BY
2585 ! A) MOVING AIR DRY-ADIABATICALLY TO NEXT LEVEL BELOW AND
2586 ! B) CORRECTING FOR EVAPORATION TO OBTAIN SATURATED STATE.
2589 ! *CUADJTQ* FOR ADJUSTING T AND Q DUE TO EVAPORATION IN
2594 ! ----------------------------------------------------------------
2595 !-------------------------------------------------------------------
2597 !-------------------------------------------------------------------
2598 INTEGER KLON, KLEV, KLEVP1
2599 INTEGER JK,IS,JL,ITOPDE, IK, ICALL
2600 REAL ZENTR,ZSEEN, ZQEEN, ZSDDE, ZQDDE,ZMFDSK, ZMFDQK
2601 REAL ZBUO, ZDMFDP, ZMFDUK, ZMFDVK
2602 REAL PTENH(KLON,KLEV), PQENH(KLON,KLEV), &
2603 PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
2604 PGEOH(KLON,KLEV), PAPH(KLON,KLEVP1)
2605 REAL PTD(KLON,KLEV), PQD(KLON,KLEV), &
2606 PUD(KLON,KLEV), PVD(KLON,KLEV), &
2607 PMFD(KLON,KLEV), PMFDS(KLON,KLEV), &
2608 PMFDQ(KLON,KLEV), PDMFDP(KLON,KLEV), &
2610 REAL ZDMFEN(KLON), ZDMFDE(KLON), &
2611 ZCOND(KLON), ZPH(KLON)
2612 LOGICAL LDDRAF(KLON), LLO2(KLON)
2613 !--------------------------------------------------------------
2614 ! 1. CALCULATE MOIST DESCENT FOR CUMULUS DOWNDRAFT BY
2615 ! (A) CALCULATING ENTRAINMENT RATES, ASSUMING
2616 ! LINEAR DECREASE OF MASSFLUX IN PBL
2617 ! (B) DOING MOIST DESCENT - EVAPORATIVE COOLING
2618 ! AND MOISTENING IS CALCULATED IN *CUADJTQ*
2619 ! (C) CHECKING FOR NEGATIVE BUOYANCY AND
2620 ! SPECIFYING FINAL T,Q,U,V AND DOWNWARD FLUXES
2621 ! ----------------------------------------------------------------
2627 LLO2(JL)=LDDRAF(JL).AND.PMFD(JL,JK-1).LT.0.
2632 IF(IS.EQ.0) GO TO 180
2635 ZENTR=ENTRDD*PMFD(JL,JK-1)*RD*PTENH(JL,JK-1)/ &
2636 (G*PAPH(JL,JK-1))*(PAPH(JL,JK)-PAPH(JL,JK-1))
2642 IF(JK.GT.ITOPDE) THEN
2646 ZDMFDE(JL)=PMFD(JL,ITOPDE)* &
2647 (PAPH(JL,JK)-PAPH(JL,JK-1))/ &
2648 (PAPH(JL,KLEVP1)-PAPH(JL,ITOPDE))
2654 PMFD(JL,JK)=PMFD(JL,JK-1)+ZDMFEN(JL)-ZDMFDE(JL)
2655 ZSEEN=(CPD*PTENH(JL,JK-1)+PGEOH(JL,JK-1))*ZDMFEN(JL)
2656 ZQEEN=PQENH(JL,JK-1)*ZDMFEN(JL)
2657 ZSDDE=(CPD*PTD(JL,JK-1)+PGEOH(JL,JK-1))*ZDMFDE(JL)
2658 ZQDDE=PQD(JL,JK-1)*ZDMFDE(JL)
2659 ZMFDSK=PMFDS(JL,JK-1)+ZSEEN-ZSDDE
2660 ZMFDQK=PMFDQ(JL,JK-1)+ZQEEN-ZQDDE
2661 PQD(JL,JK)=ZMFDQK*(1./MIN(-CMFCMIN,PMFD(JL,JK)))
2662 PTD(JL,JK)=(ZMFDSK*(1./MIN(-CMFCMIN,PMFD(JL,JK)))- &
2664 PTD(JL,JK)=MIN(400.,PTD(JL,JK))
2665 PTD(JL,JK)=MAX(100.,PTD(JL,JK))
2666 ZCOND(JL)=PQD(JL,JK)
2671 CALL CUADJTQ(KLON,KLEV,IK,ZPH,PTD,PQD,LLO2,ICALL)
2674 ZCOND(JL)=ZCOND(JL)-PQD(JL,JK)
2675 ZBUO=PTD(JL,JK)*(1.+VTMPC1*PQD(JL,JK))- &
2676 PTENH(JL,JK)*(1.+VTMPC1*PQENH(JL,JK))
2677 IF(ZBUO.GE.0..OR.PRFL(JL).LE.(PMFD(JL,JK)*ZCOND(JL))) THEN
2680 PMFDS(JL,JK)=(CPD*PTD(JL,JK)+PGEOH(JL,JK))*PMFD(JL,JK)
2681 PMFDQ(JL,JK)=PQD(JL,JK)*PMFD(JL,JK)
2682 ZDMFDP=-PMFD(JL,JK)*ZCOND(JL)
2683 PDMFDP(JL,JK-1)=ZDMFDP
2684 PRFL(JL)=PRFL(JL)+ZDMFDP
2689 IF(LLO2(JL).AND.PMFD(JL,JK).LT.0.) THEN
2690 ZMFDUK=PMFD(JL,JK-1)*PUD(JL,JK-1)+ &
2691 ZDMFEN(JL)*PUEN(JL,JK-1)-ZDMFDE(JL)*PUD(JL,JK-1)
2692 ZMFDVK=PMFD(JL,JK-1)*PVD(JL,JK-1)+ &
2693 ZDMFEN(JL)*PVEN(JL,JK-1)-ZDMFDE(JL)*PVD(JL,JK-1)
2694 PUD(JL,JK)=ZMFDUK*(1./MIN(-CMFCMIN,PMFD(JL,JK)))
2695 PVD(JL,JK)=ZMFDVK*(1./MIN(-CMFCMIN,PMFD(JL,JK)))
2701 END SUBROUTINE CUDDRAF
2704 !**********************************************
2706 !**********************************************
2708 (KLON, KLEV, KLEVP1, PQEN, PQSEN, &
2709 PTENH, PQENH, PAPH, PGEOH, KCBOT, &
2710 KCTOP, KDTOP, KTYPE, LDDRAF, LDCUM, &
2711 PMFU, PMFD, PMFUS, PMFDS, PMFUQ, &
2712 PMFDQ, PMFUL, PLUDE, PDMFUP, PDMFDP, &
2713 PRFL, PRAIN, PTEN, PSFL, PDPMEL, &
2714 KTOPM2, ZTMST, sig1)
2715 ! M.TIEDTKE E.C.M.W.F. 7/86 MODIF. 12/89
2718 ! THIS ROUTINE DOES THE FINAL CALCULATION OF CONVECTIVE
2719 ! FLUXES IN THE CLOUD LAYER AND IN THE SUBCLOUD LAYER
2722 ! THIS ROUTINE IS CALLED FROM *CUMASTR*.
2726 ! ----------------------------------------------------------------
2727 !-------------------------------------------------------------------
2729 !-------------------------------------------------------------------
2730 INTEGER KLON, KLEV, KLEVP1
2731 INTEGER KTOPM2, ITOP, JL, JK, IKB
2732 REAL ZTMST, ZCONS1, ZCONS2, ZCUCOV, ZTMELP2
2733 REAL ZZP, ZFAC, ZSNMLT, ZRFL, CEVAPCU, ZRNEW
2734 REAL ZRMIN, ZRFLN, ZDRFL, ZDPEVAP
2735 REAL PQEN(KLON,KLEV), PQSEN(KLON,KLEV), &
2736 PTENH(KLON,KLEV), PQENH(KLON,KLEV), &
2737 PAPH(KLON,KLEVP1), PGEOH(KLON,KLEV)
2738 REAL PMFU(KLON,KLEV), PMFD(KLON,KLEV), &
2739 PMFUS(KLON,KLEV), PMFDS(KLON,KLEV), &
2740 PMFUQ(KLON,KLEV), PMFDQ(KLON,KLEV), &
2741 PDMFUP(KLON,KLEV), PDMFDP(KLON,KLEV), &
2742 PMFUL(KLON,KLEV), PLUDE(KLON,KLEV), &
2743 PRFL(KLON), PRAIN(KLON)
2744 REAL PTEN(KLON,KLEV), PDPMEL(KLON,KLEV), &
2745 PSFL(KLON), ZPSUBCL(KLON)
2747 INTEGER KCBOT(KLON), KCTOP(KLON), &
2748 KDTOP(KLON), KTYPE(KLON)
2749 LOGICAL LDDRAF(KLON), LDCUM(KLON)
2750 !* SPECIFY CONSTANTS
2751 ZCONS1=CPD/(ALF*G*ZTMST)
2755 !* 1.0 DETERMINE FINAL CONVECTIVE FLUXES
2756 !---------------------------------------------
2763 ! SWITCH OFF SHALLOW CONVECTION
2764 IF(.NOT.LMFSCV.AND.KTYPE(JL).EQ.2)THEN
2768 ITOP=MIN(ITOP,KCTOP(JL))
2769 IF(.NOT.LDCUM(JL).OR.KDTOP(JL).LT.KCTOP(JL)) LDDRAF(JL)=.FALSE.
2770 IF(.NOT.LDCUM(JL)) KTYPE(JL)=0
2773 DO 120 JK=KTOPM2,KLEV
2775 IF(LDCUM(JL).AND.JK.GE.KCTOP(JL)-1) THEN
2776 PMFUS(JL,JK)=PMFUS(JL,JK)-PMFU(JL,JK)* &
2777 (CPD*PTENH(JL,JK)+PGEOH(JL,JK))
2778 PMFUQ(JL,JK)=PMFUQ(JL,JK)-PMFU(JL,JK)*PQENH(JL,JK)
2779 IF(LDDRAF(JL).AND.JK.GE.KDTOP(JL)) THEN
2780 PMFDS(JL,JK)=PMFDS(JL,JK)-PMFD(JL,JK)* &
2781 (CPD*PTENH(JL,JK)+PGEOH(JL,JK))
2782 PMFDQ(JL,JK)=PMFDQ(JL,JK)-PMFD(JL,JK)*PQENH(JL,JK)
2803 DO 130 JK=KTOPM2,KLEV
2805 IF(LDCUM(JL).AND.JK.GT.KCBOT(JL)) THEN
2807 ZZP=((PAPH(JL,KLEVP1)-PAPH(JL,JK))/ &
2808 (PAPH(JL,KLEVP1)-PAPH(JL,IKB)))
2809 IF(KTYPE(JL).EQ.3) THEN
2812 PMFU(JL,JK)=PMFU(JL,IKB)*ZZP
2813 PMFUS(JL,JK)=PMFUS(JL,IKB)*ZZP
2814 PMFUQ(JL,JK)=PMFUQ(JL,IKB)*ZZP
2815 PMFUL(JL,JK)=PMFUL(JL,IKB)*ZZP
2817 !* 2. CALCULATE RAIN/SNOW FALL RATES
2818 !* CALCULATE MELTING OF SNOW
2819 !* CALCULATE EVAPORATION OF PRECIP
2820 !----------------------------------------------
2822 PRAIN(JL)=PRAIN(JL)+PDMFUP(JL,JK)
2823 IF(PTEN(JL,JK).GT.TMELT) THEN
2824 PRFL(JL)=PRFL(JL)+PDMFUP(JL,JK)+PDMFDP(JL,JK)
2825 IF(PSFL(JL).GT.0..AND.PTEN(JL,JK).GT.ZTMELP2) THEN
2826 ZFAC=ZCONS1*(PAPH(JL,JK+1)-PAPH(JL,JK))
2827 ZSNMLT=MIN(PSFL(JL),ZFAC*(PTEN(JL,JK)-ZTMELP2))
2828 PDPMEL(JL,JK)=ZSNMLT
2829 PSFL(JL)=PSFL(JL)-ZSNMLT
2830 PRFL(JL)=PRFL(JL)+ZSNMLT
2833 PSFL(JL)=PSFL(JL)+PDMFUP(JL,JK)+PDMFDP(JL,JK)
2839 PRFL(JL)=MAX(PRFL(JL),0.)
2840 PSFL(JL)=MAX(PSFL(JL),0.)
2841 ZPSUBCL(JL)=PRFL(JL)+PSFL(JL)
2843 DO 240 JK=KTOPM2,KLEV
2845 IF(LDCUM(JL).AND.JK.GE.KCBOT(JL).AND. &
2846 ZPSUBCL(JL).GT.1.E-20) THEN
2848 CEVAPCU=CEVAPCU1*SQRT(CEVAPCU2*SQRT(sig1(JK)))
2849 ZRNEW=(MAX(0.,SQRT(ZRFL/ZCUCOV)- &
2850 CEVAPCU*(PAPH(JL,JK+1)-PAPH(JL,JK))* &
2851 MAX(0.,PQSEN(JL,JK)-PQEN(JL,JK))))**2*ZCUCOV
2852 ZRMIN=ZRFL-ZCUCOV*MAX(0.,0.8*PQSEN(JL,JK)-PQEN(JL,JK)) &
2853 *ZCONS2*(PAPH(JL,JK+1)-PAPH(JL,JK))
2854 ZRNEW=MAX(ZRNEW,ZRMIN)
2856 ZDRFL=MIN(0.,ZRFLN-ZRFL)
2857 PDMFUP(JL,JK)=PDMFUP(JL,JK)+ZDRFL
2863 ZDPEVAP=ZPSUBCL(JL)-(PRFL(JL)+PSFL(JL))
2864 PRFL(JL)=PRFL(JL)+ZDPEVAP*PRFL(JL)* &
2865 (1./MAX(1.E-20,PRFL(JL)+PSFL(JL)))
2866 PSFL(JL)=PSFL(JL)+ZDPEVAP*PSFL(JL)* &
2867 (1./MAX(1.E-20,PRFL(JL)+PSFL(JL)))
2870 END SUBROUTINE CUFLX
2873 !**********************************************
2875 !**********************************************
2877 (KLON, KLEV, KLEVP1, KTOPM2, PAPH, &
2878 LDCUM, PTEN, PTTE, PQTE, PMFUS, &
2879 PMFDS, PMFUQ, PMFDQ, PMFUL, PDMFUP, &
2880 PDMFDP, ZTMST, PDPMEL, PRAIN, PRFL, &
2881 PSFL, PSRAIN, PSEVAP, PSHEAT, PSMELT, &
2882 PRSFC, PSSFC, PAPRC, PAPRSM, PAPRS, &
2883 PQEN, PQSEN, PLUDE, PCTE)
2884 !**** *CUDTDQ* - UPDATES T AND Q TENDENCIES, PRECIPITATION RATES
2885 ! DOES GLOBAL DIAGNOSTICS
2886 ! M.TIEDTKE E.C.M.W.F. 7/86 MODIF. 12/89
2889 ! *CUDTDQ* IS CALLED FROM *CUMASTR*
2890 ! ----------------------------------------------------------------
2891 !-------------------------------------------------------------------
2893 !-------------------------------------------------------------------
2894 INTEGER KLON, KLEV, KLEVP1
2895 INTEGER KTOPM2,JL, JK
2896 REAL ZTMST, PSRAIN, PSEVAP, PSHEAT, PSMELT, ZDIAGT, ZDIAGW
2897 REAL ZALV, RHK, RHCOE, PLDFD, ZDTDT, ZDQDT
2898 REAL PTTE(KLON,KLEV), PQTE(KLON,KLEV), &
2899 PTEN(KLON,KLEV), PLUDE(KLON,KLEV), &
2900 PGEO(KLON,KLEV), PAPH(KLON,KLEVP1), &
2901 PAPRC(KLON), PAPRS(KLON), &
2902 PAPRSM(KLON), PCTE(KLON,KLEV), &
2903 PRSFC(KLON), PSSFC(KLON)
2904 REAL PMFUS(KLON,KLEV), PMFDS(KLON,KLEV), &
2905 PMFUQ(KLON,KLEV), PMFDQ(KLON,KLEV), &
2906 PMFUL(KLON,KLEV), PQSEN(KLON,KLEV), &
2907 PDMFUP(KLON,KLEV), PDMFDP(KLON,KLEV),&
2908 PRFL(KLON), PRAIN(KLON), &
2910 REAL PDPMEL(KLON,KLEV), PSFL(KLON)
2911 REAL ZSHEAT(KLON), ZMELT(KLON)
2913 !--------------------------------
2914 !* 1.0 SPECIFY PARAMETERS
2915 !--------------------------------
2918 ZDIAGW=ZDIAGT/RHOH2O
2919 !--------------------------------------------------
2920 !* 2.0 INCREMENTATION OF T AND Q TENDENCIES
2921 !--------------------------------------------------
2927 DO 250 JK=KTOPM2,KLEV
2931 IF(PTEN(JL,JK).GT.TMELT) THEN
2936 RHK=MIN(1.0,PQEN(JL,JK)/PQSEN(JL,JK))
2937 RHCOE=MAX(0.0,(RHK-RHC)/(RHM-RHC))
2938 pldfd=MAX(0.0,RHCOE*fdbk*PLUDE(JL,JK))
2939 ZDTDT=(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))*RCPD* &
2940 (PMFUS(JL,JK+1)-PMFUS(JL,JK)+ &
2941 PMFDS(JL,JK+1)-PMFDS(JL,JK)-ALF*PDPMEL(JL,JK) &
2942 -ZALV*(PMFUL(JL,JK+1)-PMFUL(JL,JK)-pldfd- &
2943 (PDMFUP(JL,JK)+PDMFDP(JL,JK))))
2944 PTTE(JL,JK)=PTTE(JL,JK)+ZDTDT
2945 ZDQDT=(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))*&
2946 (PMFUQ(JL,JK+1)-PMFUQ(JL,JK)+ &
2947 PMFDQ(JL,JK+1)-PMFDQ(JL,JK)+ &
2948 PMFUL(JL,JK+1)-PMFUL(JL,JK)-pldfd- &
2949 (PDMFUP(JL,JK)+PDMFDP(JL,JK)))
2950 PQTE(JL,JK)=PQTE(JL,JK)+ZDQDT
2951 PCTE(JL,JK)=(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))*pldfd
2952 ZSHEAT(JL)=ZSHEAT(JL)+ZALV*(PDMFUP(JL,JK)+PDMFDP(JL,JK))
2953 ZMELT(JL)=ZMELT(JL)+PDPMEL(JL,JK)
2959 IF(PTEN(JL,JK).GT.TMELT) THEN
2964 RHK=MIN(1.0,PQEN(JL,JK)/PQSEN(JL,JK))
2965 RHCOE=MAX(0.0,(RHK-RHC)/(RHM-RHC))
2966 pldfd=MAX(0.0,RHCOE*fdbk*PLUDE(JL,JK))
2967 ZDTDT=-(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))*RCPD* &
2968 (PMFUS(JL,JK)+PMFDS(JL,JK)+ALF*PDPMEL(JL,JK)-ZALV* &
2969 (PMFUL(JL,JK)+PDMFUP(JL,JK)+PDMFDP(JL,JK)+pldfd))
2970 PTTE(JL,JK)=PTTE(JL,JK)+ZDTDT
2971 ZDQDT=-(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))* &
2972 (PMFUQ(JL,JK)+PMFDQ(JL,JK)+pldfd+ &
2973 (PMFUL(JL,JK)+PDMFUP(JL,JK)+PDMFDP(JL,JK)))
2974 PQTE(JL,JK)=PQTE(JL,JK)+ZDQDT
2975 PCTE(JL,JK)=(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))*pldfd
2976 ZSHEAT(JL)=ZSHEAT(JL)+ZALV*(PDMFUP(JL,JK)+PDMFDP(JL,JK))
2977 ZMELT(JL)=ZMELT(JL)+PDPMEL(JL,JK)
2982 !---------------------------------------------------------
2983 ! 3. UPDATE SURFACE FIELDS AND DO GLOBAL BUDGETS
2984 !---------------------------------------------------------
2989 PAPRC(JL)=PAPRC(JL)+ZDIAGW*(PRFL(JL)+PSFL(JL))
2990 PAPRS(JL)=PAPRSM(JL)+ZDIAGW*PSFL(JL)
2991 PSHEAT=PSHEAT+ZSHEAT(JL)
2992 PSRAIN=PSRAIN+PRAIN(JL)
2993 PSEVAP=PSEVAP-(PRFL(JL)+PSFL(JL))
2994 PSMELT=PSMELT+ZMELT(JL)
2996 PSEVAP=PSEVAP+PSRAIN
2998 END SUBROUTINE CUDTDQ
3001 !**********************************************
3003 !**********************************************
3005 (KLON, KLEV, KLEVP1, KTOPM2, KTYPE, &
3006 KCBOT, PAPH, LDCUM, PUEN, PVEN, &
3007 PVOM, PVOL, PUU, PUD, PVU, &
3008 PVD, PMFU, PMFD, PSDISS)
3009 !**** *CUDUDV* - UPDATES U AND V TENDENCIES,
3010 ! DOES GLOBAL DIAGNOSTIC OF DISSIPATION
3011 ! M.TIEDTKE E.C.M.W.F. 7/86 MODIF. 12/89
3014 ! *CUDUDV* IS CALLED FROM *CUMASTR*
3015 ! ----------------------------------------------------------------
3016 !-------------------------------------------------------------------
3018 !-------------------------------------------------------------------
3019 INTEGER KLON, KLEV, KLEVP1
3020 INTEGER KTOPM2, JK, IK, JL, IKB
3021 REAL PSDISS,ZZP, ZDUDT ,ZDVDT, ZSUM
3022 REAL PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
3023 PVOL(KLON,KLEV), PVOM(KLON,KLEV), &
3025 REAL PUU(KLON,KLEV), PUD(KLON,KLEV), &
3026 PVU(KLON,KLEV), PVD(KLON,KLEV), &
3027 PMFU(KLON,KLEV), PMFD(KLON,KLEV)
3028 REAL ZMFUU(KLON,KLEV), ZMFDU(KLON,KLEV), &
3029 ZMFUV(KLON,KLEV), ZMFDV(KLON,KLEV), &
3031 INTEGER KTYPE(KLON), KCBOT(KLON)
3033 !------------------------------------------------------------
3034 !* 1.0 CALCULATE FLUXES AND UPDATE U AND V TENDENCIES
3035 ! -----------------------------------------------------------
3037 DO 120 JK=KTOPM2,KLEV
3041 ZMFUU(JL,JK)=PMFU(JL,JK)*(PUU(JL,JK)-PUEN(JL,IK))
3042 ZMFUV(JL,JK)=PMFU(JL,JK)*(PVU(JL,JK)-PVEN(JL,IK))
3043 ZMFDU(JL,JK)=PMFD(JL,JK)*(PUD(JL,JK)-PUEN(JL,IK))
3044 ZMFDV(JL,JK)=PMFD(JL,JK)*(PVD(JL,JK)-PVEN(JL,IK))
3048 DO 140 JK=KTOPM2,KLEV
3050 IF(LDCUM(JL).AND.JK.GT.KCBOT(JL)) THEN
3052 ZZP=((PAPH(JL,KLEVP1)-PAPH(JL,JK))/ &
3053 (PAPH(JL,KLEVP1)-PAPH(JL,IKB)))
3054 IF(KTYPE(JL).EQ.3) THEN
3057 ZMFUU(JL,JK)=ZMFUU(JL,IKB)*ZZP
3058 ZMFUV(JL,JK)=ZMFUV(JL,IKB)*ZZP
3059 ZMFDU(JL,JK)=ZMFDU(JL,IKB)*ZZP
3060 ZMFDV(JL,JK)=ZMFDV(JL,IKB)*ZZP
3067 DO 190 JK=KTOPM2,KLEV
3071 ZDUDT=(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))* &
3072 (ZMFUU(JL,JK+1)-ZMFUU(JL,JK)+ &
3073 ZMFDU(JL,JK+1)-ZMFDU(JL,JK))
3074 ZDVDT=(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))* &
3075 (ZMFUV(JL,JK+1)-ZMFUV(JL,JK)+ &
3076 ZMFDV(JL,JK+1)-ZMFDV(JL,JK))
3077 ZDISS(JL)=ZDISS(JL)+ &
3078 PUEN(JL,JK)*(ZMFUU(JL,JK+1)-ZMFUU(JL,JK)+ &
3079 ZMFDU(JL,JK+1)-ZMFDU(JL,JK))+ &
3080 PVEN(JL,JK)*(ZMFUV(JL,JK+1)-ZMFUV(JL,JK)+ &
3081 ZMFDV(JL,JK+1)-ZMFDV(JL,JK))
3082 PVOM(JL,JK)=PVOM(JL,JK)+ZDUDT
3083 PVOL(JL,JK)=PVOL(JL,JK)+ZDVDT
3089 ZDUDT=-(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))* &
3090 (ZMFUU(JL,JK)+ZMFDU(JL,JK))
3091 ZDVDT=-(G/(PAPH(JL,JK+1)-PAPH(JL,JK)))* &
3092 (ZMFUV(JL,JK)+ZMFDV(JL,JK))
3093 ZDISS(JL)=ZDISS(JL)- &
3094 (PUEN(JL,JK)*(ZMFUU(JL,JK)+ZMFDU(JL,JK))+ &
3095 PVEN(JL,JK)*(ZMFUV(JL,JK)+ZMFDV(JL,JK)))
3096 PVOM(JL,JK)=PVOM(JL,JK)+ZDUDT
3097 PVOL(JL,JK)=PVOL(JL,JK)+ZDVDT
3102 ZSUM=SSUM(KLON,ZDISS(1),1)
3105 END SUBROUTINE CUDUDV
3108 !#################################################################
3110 ! LEVEL 4 SUBROUTINES
3112 !#################################################################
3113 !**************************************************************
3114 ! SUBROUTINE CUBASMC
3115 !**************************************************************
3116 SUBROUTINE CUBASMC &
3117 (KLON, KLEV, KLEVM1, KK, PTEN, &
3118 PQEN, PQSEN, PUEN, PVEN, PVERV, &
3119 PGEO, PGEOH, LDCUM, KTYPE, KLAB, &
3120 PMFU, PMFUB, PENTR, KCBOT, PTU, &
3121 PQU, PLU, PUU, PVU, PMFUS, &
3122 PMFUQ, PMFUL, PDMFUP, PMFUU, PMFUV)
3123 ! M.TIEDTKE E.C.M.W.F. 12/89
3126 ! THIS ROUTINE CALCULATES CLOUD BASE VALUES
3127 ! FOR MIDLEVEL CONVECTION
3130 ! THIS ROUTINE IS CALLED FROM *CUASC*.
3131 ! INPUT ARE ENVIRONMENTAL VALUES T,Q ETC
3132 ! IT RETURNS CLOUDBASE VALUES FOR MIDLEVEL CONVECTION
3139 ! ----------------------------------------------------------------
3140 !-------------------------------------------------------------------
3142 !-------------------------------------------------------------------
3143 INTEGER KLON, KLEV, KLEVP1
3144 INTEGER KLEVM1,KK, JL
3146 REAL PTEN(KLON,KLEV), PQEN(KLON,KLEV), &
3147 PUEN(KLON,KLEV), PVEN(KLON,KLEV), &
3148 PQSEN(KLON,KLEV), PVERV(KLON,KLEV), &
3149 PGEO(KLON,KLEV), PGEOH(KLON,KLEV)
3150 REAL PTU(KLON,KLEV), PQU(KLON,KLEV), &
3151 PUU(KLON,KLEV), PVU(KLON,KLEV), &
3152 PLU(KLON,KLEV), PMFU(KLON,KLEV), &
3153 PMFUB(KLON), PENTR(KLON), &
3154 PMFUS(KLON,KLEV), PMFUQ(KLON,KLEV), &
3155 PMFUL(KLON,KLEV), PDMFUP(KLON,KLEV), &
3156 PMFUU(KLON), PMFUV(KLON)
3157 INTEGER KTYPE(KLON), KCBOT(KLON), &
3160 !--------------------------------------------------------
3161 !* 1. CALCULATE ENTRAINMENT AND DETRAINMENT RATES
3162 ! -------------------------------------------------------
3165 IF( .NOT. LDCUM(JL).AND.KLAB(JL,KK+1).EQ.0.0.AND. &
3166 PQEN(JL,KK).GT.0.80*PQSEN(JL,KK)) THEN
3167 PTU(JL,KK+1)=(CPD*PTEN(JL,KK)+PGEO(JL,KK)-PGEOH(JL,KK+1)) &
3169 PQU(JL,KK+1)=PQEN(JL,KK)
3171 ZZZMB=MAX(CMFCMIN,-PVERV(JL,KK)/G)
3172 ZZZMB=MIN(ZZZMB,CMFCMAX)
3174 PMFU(JL,KK+1)=PMFUB(JL)
3175 PMFUS(JL,KK+1)=PMFUB(JL)*(CPD*PTU(JL,KK+1)+PGEOH(JL,KK+1))
3176 PMFUQ(JL,KK+1)=PMFUB(JL)*PQU(JL,KK+1)
3184 PUU(JL,KK+1)=PUEN(JL,KK)
3185 PVU(JL,KK+1)=PVEN(JL,KK)
3186 PMFUU(JL)=PMFUB(JL)*PUU(JL,KK+1)
3187 PMFUV(JL)=PMFUB(JL)*PVU(JL,KK+1)
3192 END SUBROUTINE CUBASMC
3195 !**************************************************************
3196 ! SUBROUTINE CUADJTQ
3197 !**************************************************************
3198 SUBROUTINE CUADJTQ(KLON,KLEV,KK,PP,PT,PQ,LDFLAG,KCALL)
3199 ! M.TIEDTKE E.C.M.W.F. 12/89
3200 ! D.SALMOND CRAY(UK)) 12/8/91
3203 ! TO PRODUCE T,Q AND L VALUES FOR CLOUD ASCENT
3206 ! THIS ROUTINE IS CALLED FROM SUBROUTINES:
3207 ! *CUBASE* (T AND Q AT CONDENSTION LEVEL)
3208 ! *CUASC* (T AND Q AT CLOUD LEVELS)
3209 ! *CUINI* (ENVIRONMENTAL T AND QS VALUES AT HALF LEVELS)
3210 ! INPUT ARE UNADJUSTED T AND Q VALUES,
3211 ! IT RETURNS ADJUSTED VALUES OF T AND Q
3212 ! NOTE: INPUT PARAMETER KCALL DEFINES CALCULATION AS
3213 ! KCALL=0 ENV. T AND QS IN*CUINI*
3214 ! KCALL=1 CONDENSATION IN UPDRAFTS (E.G. CUBASE, CUASC)
3215 ! KCALL=2 EVAPORATION IN DOWNDRAFTS (E.G. CUDLFS,CUDDRAF
3218 ! 3 LOOKUP TABLES ( TLUCUA, TLUCUB, TLUCUC )
3219 ! FOR CONDENSATION CALCULATIONS.
3220 ! THE TABLES ARE INITIALISED IN *SETPHYS*.
3221 ! ----------------------------------------------------------------
3222 !-------------------------------------------------------------------
3224 !-------------------------------------------------------------------
3226 INTEGER KK, KCALL, ISUM, JL
3227 REAL ZQSAT, ZCOR, ZCOND1, TT
3228 REAL PT(KLON,KLEV), PQ(KLON,KLEV), &
3229 ZCOND(KLON), ZQP(KLON), &
3231 LOGICAL LDFLAG(KLON)
3232 !------------------------------------------------------------------
3233 ! 2. CALCULATE CONDENSATION AND ADJUST T AND Q ACCORDINGLY
3234 !------------------------------------------------------------------
3236 IF (KCALL.EQ.1 ) THEN
3243 ZQSAT=TLUCUA(TT)*ZQP(JL)
3244 ZQSAT=MIN(0.5,ZQSAT)
3245 ZCOR=1./(1.-VTMPC1*ZQSAT)
3247 ZCOND(JL)=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3248 ZCOND(JL)=MAX(ZCOND(JL),0.)
3249 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND(JL)
3250 PQ(JL,KK)=PQ(JL,KK)-ZCOND(JL)
3251 IF(ZCOND(JL).NE.0.0) ISUM=ISUM+1
3254 IF(ISUM.EQ.0) GO TO 230
3256 IF(LDFLAG(JL).AND.ZCOND(JL).NE.0.) THEN
3258 ZQSAT=TLUCUA(TT)*ZQP(JL)
3259 ZQSAT=MIN(0.5,ZQSAT)
3260 ZCOR=1./(1.-VTMPC1*ZQSAT)
3262 ZCOND1=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3263 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND1
3264 PQ(JL,KK)=PQ(JL,KK)-ZCOND1
3276 ZQSAT=TLUCUA(TT)*ZQP(JL)
3277 ZQSAT=MIN(0.5,ZQSAT)
3278 ZCOR=1./(1.-VTMPC1*ZQSAT)
3280 ZCOND(JL)=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3281 ZCOND(JL)=MIN(ZCOND(JL),0.)
3282 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND(JL)
3283 PQ(JL,KK)=PQ(JL,KK)-ZCOND(JL)
3284 IF(ZCOND(JL).NE.0.0) ISUM=ISUM+1
3287 IF(ISUM.EQ.0) GO TO 330
3289 IF(LDFLAG(JL).AND.ZCOND(JL).NE.0.) THEN
3291 ZQSAT=TLUCUA(TT)*ZQP(JL)
3292 ZQSAT=MIN(0.5,ZQSAT)
3293 ZCOR=1./(1.-VTMPC1*ZQSAT)
3295 ZCOND1=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3296 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND1
3297 PQ(JL,KK)=PQ(JL,KK)-ZCOND1
3307 ZQSAT=TLUCUA(TT)*ZQP(JL)
3308 ZQSAT=MIN(0.5,ZQSAT)
3309 ZCOR=1./(1.-VTMPC1*ZQSAT)
3311 ZCOND(JL)=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3312 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND(JL)
3313 PQ(JL,KK)=PQ(JL,KK)-ZCOND(JL)
3314 IF(ZCOND(JL).NE.0.0) ISUM=ISUM+1
3316 IF(ISUM.EQ.0) GO TO 430
3319 ZQSAT=TLUCUA(TT)*ZQP(JL)
3320 ZQSAT=MIN(0.5,ZQSAT)
3321 ZCOR=1./(1.-VTMPC1*ZQSAT)
3323 ZCOND1=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3324 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND1
3325 PQ(JL,KK)=PQ(JL,KK)-ZCOND1
3333 ZQSAT=TLUCUA(TT)*ZQP(JL)
3334 ZQSAT=MIN(0.5,ZQSAT)
3335 ZCOR=1./(1.-VTMPC1*ZQSAT)
3337 ZCOND(JL)=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3338 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND(JL)
3339 PQ(JL,KK)=PQ(JL,KK)-ZCOND(JL)
3343 ZQSAT=TLUCUA(TT)*ZQP(JL)
3344 ZQSAT=MIN(0.5,ZQSAT)
3345 ZCOR=1./(1.-VTMPC1*ZQSAT)
3347 ZCOND1=(PQ(JL,KK)-ZQSAT)/(1.+ZQSAT*ZCOR*TLUCUB(TT))
3348 PT(JL,KK)=PT(JL,KK)+TLUCUC(TT)*ZCOND1
3349 PQ(JL,KK)=PQ(JL,KK)-ZCOND1
3353 END SUBROUTINE CUADJTQ
3356 !**********************************************************
3357 ! SUBROUTINE CUENTR_NEW
3358 !**********************************************************
3359 SUBROUTINE CUENTR_NEW &
3360 (KLON, KLEV, KLEVP1, KK, PTENH, &
3361 PAPH, PAP, PGEOH, KLWMIN, LDCUM, &
3362 KTYPE, KCBOT, KCTOP0, ZPBASE, PMFU, &
3363 PENTR, ZDMFEN, ZDMFDE, ZODETR, KHMIN)
3364 ! M.TIEDTKE E.C.M.W.F. 12/89
3368 ! THIS ROUTINE CALCULATES ENTRAINMENT/DETRAINMENT RATES
3369 ! FOR UPDRAFTS IN CUMULUS PARAMETERIZATION
3372 ! THIS ROUTINE IS CALLED FROM *CUASC*.
3373 ! INPUT ARE ENVIRONMENTAL VALUES T,Q ETC
3374 ! AND UPDRAFT VALUES T,Q ETC
3375 ! IT RETURNS ENTRAINMENT/DETRAINMENT RATES
3378 ! S. TIEDTKE (1989), NORDENG(1996)
3382 ! ----------------------------------------------------------------
3383 !-------------------------------------------------------------------
3385 !-------------------------------------------------------------------
3386 INTEGER KLON, KLEV, KLEVP1
3387 INTEGER KK, JL, IKLWMIN,IKB, IKT, IKH
3388 REAL ZRRHO, ZDPRHO, ZPMID, ZENTR, ZZMZK, ZTMZK, ARG, ZORGDE
3389 REAL PTENH(KLON,KLEV), &
3390 PAP(KLON,KLEV), PAPH(KLON,KLEVP1), &
3391 PMFU(KLON,KLEV), PGEOH(KLON,KLEV), &
3392 PENTR(KLON), ZPBASE(KLON), &
3393 ZDMFEN(KLON), ZDMFDE(KLON), &
3395 INTEGER KLWMIN(KLON), KTYPE(KLON), &
3396 KCBOT(KLON), KCTOP0(KLON), &
3398 LOGICAL LDCUM(KLON),LLO1,LLO2
3400 real tt(klon),ttb(klon)
3401 real zqsat(klon), zqsatb(klon)
3403 !---------------------------------------------------------
3404 !* 1. CALCULATE ENTRAINMENT AND DETRAINMENT RATES
3405 !---------------------------------------------------------
3406 !* 1.1 SPECIFY ENTRAINMENT RATES FOR SHALLOW CLOUDS
3407 !----------------------------------------------------------
3408 !* 1.2 SPECIFY ENTRAINMENT RATES FOR DEEP CLOUDS
3409 !-------------------------------------------------------
3411 zpbase(jl) = paph(jl,kcbot(jl))
3412 zrrho = (rd*ptenh(jl,kk+1))/paph(jl,kk+1)
3413 zdprho = (paph(jl,kk+1)-paph(jl,kk))*zrg
3415 zpmid = 0.5*(zpbase(jl)+paph(jl,kctop0(jl)))
3416 zentr = pentr(jl)*pmfu(jl,kk+1)*zdprho*zrrho
3417 llo1 = kk.LT.kcbot(jl).AND.ldcum(jl)
3420 if(nturben.eq.1) zdmfde(jl) = zentr
3421 if(nturben.eq.2) zdmfde(jl) = zentr*1.2
3426 if(nturben .eq. 1) then
3428 elseif (nturben .eq. 2) then
3429 ! defining the facale
3430 tt(jl) = ptenh(jl,kk+1)
3431 zqsat(jl) = TLUCUA(tt(jl))/paph(jl,kk+1)
3432 zqsat(jl) = zqsat(jl)/(1.-VTMPC1*zqsat(jl))
3433 ttb(jl) = ptenh(jl,kcbot(jl))
3434 zqsatb(jl) = TLUCUA(ttb(jl))/zpbase(jl)
3435 zqsatb(jl) = zqsatb(jl)/(1.-VTMPC1*zqsatb(jl))
3436 fscale(jl) = 4.0*(zqsat(jl)/zqsatb(jl))**2
3438 ! end of defining the fscale
3439 llo2 = llo1.AND.ktype(jl).EQ.2.AND.((zpbase(jl)-paph(jl,kk)) &
3440 .LT.ZDNOPRC.OR.paph(jl,kk).GT.zpmid)
3442 zdmfen(jl) = zentr*fscale(jl)
3446 iklwmin = MAX(klwmin(jl),kctop0(jl)+2)
3447 llo2 = llo1.AND.ktype(jl).EQ.3.AND.(kk.GE.iklwmin.OR.pap(jl,kk) &
3449 IF (llo2) zdmfen(jl) = zentr*fscale(jl)
3450 llo2 = llo1.AND.ktype(jl).EQ.1
3451 ! Turbulent entrainment
3452 IF (llo2) zdmfen(jl) = zentr*fscale(jl)
3453 ! Organized detrainment, detrainment starts at khmin
3456 IF (llo2.AND.kk.LE.khmin(jl).AND.kk.GE.kctop0(jl)) THEN
3459 IF (ikh.GT.ikt) THEN
3460 zzmzk = -(pgeoh(jl,ikh)-pgeoh(jl,kk))*zrg
3461 ztmzk = -(pgeoh(jl,ikh)-pgeoh(jl,ikt))*zrg
3462 arg = 3.1415*(zzmzk/ztmzk)*0.5
3463 zorgde = TAN(arg)*3.1415*0.5/ztmzk
3464 zdprho = (paph(jl,kk+1)-paph(jl,kk))*(zrg*zrrho)
3465 zodetr(jl,kk) = MIN(zorgde,1.E-3)*pmfu(jl,kk+1)*zdprho
3471 END SUBROUTINE CUENTR_NEW
3473 !**********************************************************
3474 ! FUNCTION SSUM, TLUCUA, TLUCUB, TLUCUC
3475 !**********************************************************
3476 REAL FUNCTION SSUM ( N, X, IX )
3478 ! COMPUTES SSUM = SUM OF [X(I)]
3479 ! FOR N ELEMENTS OF X WITH SKIP INCREMENT IX FOR VECTOR X
3484 INTEGER N, IX, JX, JL
3498 REAL FUNCTION TLUCUA(TT)
3500 ! Set up lookup tables for cloud ascent calculations.
3505 IF(TT-TMELT.GT.0.) THEN
3512 TLUCUA=C2ES*EXP(ZCVM3*(TT-TMELT)*(1./(TT-ZCVM4)))
3517 REAL FUNCTION TLUCUB(TT)
3519 ! Set up lookup tables for cloud ascent calculations.
3522 REAL Z5ALVCP,Z5ALSCP,ZCVM4,ZCVM5,TT
3524 Z5ALVCP=C5LES*ALV/CPD
3525 Z5ALSCP=C5IES*ALS/CPD
3526 IF(TT-TMELT.GT.0.) THEN
3533 TLUCUB=ZCVM5*(1./(TT-ZCVM4))**2
3538 REAL FUNCTION TLUCUC(TT)
3540 ! Set up lookup tables for cloud ascent calculations.
3543 REAL ZALVDCP,ZALSDCP,TT,ZLDCP
3547 IF(TT-TMELT.GT.0.) THEN
3558 END MODULE module_cu_tiedtke