SUBROUTINE CUCNVC C ****************************************************************** C$$$ SUBPROGRAM DOCUMENTATION BLOCK C . . . C SUBPROGRAM: CUCNVC CONVECTIVE PRECIPITATION PARAMETERIZATION C PRGRMMR: JANJIC ORG: W/NP2 DATE: 93-11-02 C C ABSTRACT: C CUCNVC CALCULATES THE SUB-GRID SCALE CONVECTION INCLUDING C DEEP AND SHALLOW CONVECTIVE CLOUDS FOLLOWING THE SCHEME DESCRIBED C BY JANJIC (1994) BUT WITH SIGNIFICANT MODIFICATIONS. IN ADDITION, C THE LATENT HEAT RELEASE AND MOISTURE CHANGE DUE TO PRECIPITATING C AND NON-PRECIPITATING CLOUDS ARE COMPUTED. C C C PROGRAM HISTORY LOG: C 87-09-?? JANJIC - ORIGINATOR C 90-11-21 JANJIC - TWO SETS OF DSP PROFILES (FAST AND SLOW) C REPLACE THE ORIGINAL ONE SET C 95-03-25 BLACK - CONVERSION FROM 1-D TO 2-D IN HORIZONTAL C 96-03-28 BLACK - ADDED EXTERNAL EDGE C 98-11-02 BLACK - MODIFIED FOR DISTRIBUTED MEMORY C C USAGE: CALL CUCNVC FROM MAIN PROGRAM EBU C C INPUT ARGUMENT LIST: C NONE C C OUTPUT ARGUMENT LIST: C NONE C C OUTPUT FILES: C NONE C C SUBPROGRAMS CALLED: C C UNIQUE: C TTBLEX C C LIBRARY: C NONE C C COMMON BLOCKS: CTLBLK C LOOPS C MASKS C PHYS C VRBLS C CNVCLD C PVRBLS C ACMCLH C INDX C C ATTRIBUTES: C LANGUAGE: FORTRAN 90 C MACHINE : IBM SP C$$$ C ****************************************************************** C * * C * REFERENCES: * C * * C * JANJIC, Z.I., 1994: THE STEP-MOUNTAIN ETA COORDINATE MODEL: * C * FURTHER DEVELOPMENTS OF THE CONVECTION, VISCOUS SUBLAYER * C * AND TURBULENCE CLOSURE SCHEMES. MONTHLY WEATHER REVIEW, * C * VOL. 122, 927-945. * C * * C ****************************************************************** C *** WARNING: THIS SUBROUTINE WILL NOT WORK IF LM.LT.12; C MUST BE CALLED IN THE SAME STEP WITH PROFQ2 BECAUSE PROFQ C DEFINES APE; C---------------------------------------------------------------------- INCLUDE "cuparm" C---------------------------------------------------------------------- INCLUDE "parmeta" INCLUDE "parm.tbl" INCLUDE "mpp.h" C---------------------------------------------------------------------- P A R A M E T E R CVVVVVVVVVV INSTABILITY FOR TOO LARGE LSH VVVVVVVVVVVVVVVVVVVVVVVVVVVVV & (KSMUD=0,NROW= 0) C & (KSMUD=0,NROW= 5) CAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA P A R A M E T E R & (IMJM=IM*JM-JM/2,JAM=6+2*(JM-10) &, IMJM_LOC=IDIM2*JDIM2 CVVVVVVVVVV INSTABILITY FOR TOO LARGE LSH VVVVVVVVVVVVVVVVVVVVVVVVVVVVV C &, LP1=LM+1,LM1=LM-1,LNO=1,LSH=LM/3-1,LSHU=LM/2-1,LQM=LM/5,KBM=3 C &, LP1=LM+1,LM1=LM-1,LNO=1,LSH=LM/3 ,LSHU=LM/2-1,LQM=LM/5,KBM=3 C &, LP1=LM+1,LM1=LM-1,LNO=3,LSH=LM/3 ,LSHU=LM/2-1,LQM=LM/5,KBM=3 C &, LP1=LM+1,LM1=LM-1,LNO=2,LSH=LM/3-1,LSHU=LM/2-1,LQM=LM/5,KBM=3 C &, LP1=LM+1,LM1=LM-1,LNO=2,LSH=LM/3-2,LSHU=LM/2-1,LQM=LM/5,KBM=3 &, LP1=LM+1,LM1=LM-1) CAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA C---------------------------------------------------------------------- L O G I C A L & RUN,FIRST,RESTRT,SIGMA C---------------------------------------------------------------------- INCLUDE "CTLBLK.comm" C----------------------------------------------------------------------- INCLUDE "LOOPS.comm" C----------------------------------------------------------------------- INCLUDE "MASKS.comm" C----------------------------------------------------------------------- INCLUDE "PHYS.comm" C----------------------------------------------------------------------- INCLUDE "VRBLS.comm" C----------------------------------------------------------------------- INCLUDE "CNVCLD.comm" C----------------------------------------------------------------------- INCLUDE "PVRBLS.comm" C----------------------------------------------------------------------- INCLUDE "ACMCLH.comm" C----------------------------------------------------------------------- INCLUDE "INDX.comm" CYL--------------------------------------------------------------------- INCLUDE "PPTASM.comm" CYL--------------------------------------------------------------------- D I M E N S I O N & TREFK (LM),QREFK (LM),PK (LM),APEK (LM),TK (LM) &,THSK (LM),PSK (LM),APESK (LM),QK (LM),THERK (LM) &,THVREF(LM),THEVRF(LM),THVMOD(LM),DIFT (LM),DIFQ (LM) &,QSATK (LM),FPK (LM) &,NTOPD (LM),NBOTD (LM),NTOPS (LM),NBOTS (LM) &,NDPTHD(LM),NDPTHS(LM) C D I M E N S I O N & LTOP (IDIM1:IDIM2,JDIM1:JDIM2),LBOT (IDIM1:IDIM2,JDIM1:JDIM2) &,PTOP (IDIM1:IDIM2,JDIM1:JDIM2),PBOT (IDIM1:IDIM2,JDIM1:JDIM2) &,IPTB (IDIM1:IDIM2,JDIM1:JDIM2),ITHTB (IDIM1:IDIM2,JDIM1:JDIM2) &,PDSL (IDIM1:IDIM2,JDIM1:JDIM2),APEBT (IDIM1:IDIM2,JDIM1:JDIM2) &,TBT (IDIM1:IDIM2,JDIM1:JDIM2),Q2BT (IDIM1:IDIM2,JDIM1:JDIM2) &,QQ (IDIM1:IDIM2,JDIM1:JDIM2),PP (IDIM1:IDIM2,JDIM1:JDIM2) &,PSP (IDIM1:IDIM2,JDIM1:JDIM2),THBT (IDIM1:IDIM2,JDIM1:JDIM2) &,THESP (IDIM1:IDIM2,JDIM1:JDIM2),P (IDIM1:IDIM2,JDIM1:JDIM2) &,BTH (IDIM1:IDIM2,JDIM1:JDIM2),STH (IDIM1:IDIM2,JDIM1:JDIM2) &,T00 (IDIM1:IDIM2,JDIM1:JDIM2),T10 (IDIM1:IDIM2,JDIM1:JDIM2) &,T01 (IDIM1:IDIM2,JDIM1:JDIM2),T11 (IDIM1:IDIM2,JDIM1:JDIM2) &,WF1 (IDIM1:IDIM2,JDIM1:JDIM2),WF2 (IDIM1:IDIM2,JDIM1:JDIM2) &,WF3 (IDIM1:IDIM2,JDIM1:JDIM2),WF4 (IDIM1:IDIM2,JDIM1:JDIM2) &,PRECOL(IDIM1:IDIM2,JDIM1:JDIM2) C &,IBUOY (IMJM_LOC),JBUOY (IMJM_LOC) &,IDEEP (IMJM_LOC),JDEEP (IMJM_LOC) &,ISHAL (IMJM_LOC),JSHAL (IMJM_LOC) &,ILRES (IMJM_LOC),JLRES (IMJM_LOC) &,IHRES (IMJM_LOC),JHRES (IMJM_LOC) C &,APE (IDIM1:IDIM2,JDIM1:JDIM2,LM) &,TREF (IDIM1:IDIM2,JDIM1:JDIM2,LM) &,TMOD (IDIM1:IDIM2,JDIM1:JDIM2,LM) &,QMOD (IDIM1:IDIM2,JDIM1:JDIM2,LM) C D I M E N S I O N & DSPB (IDIM1:IDIM2,JDIM1:JDIM2),DSP0 (IDIM1:IDIM2,JDIM1:JDIM2) &,DSPT (IDIM1:IDIM2,JDIM1:JDIM2) &,TL (IDIM1:IDIM2,JDIM1:JDIM2),QL (IDIM1:IDIM2,JDIM1:JDIM2) &,TNE (IDIM1:IDIM2,JDIM1:JDIM2),TSE (IDIM1:IDIM2,JDIM1:JDIM2) &,QNE (IDIM1:IDIM2,JDIM1:JDIM2),QSE (IDIM1:IDIM2,JDIM1:JDIM2) C----------------------------------------------------------------------- CALL ZERO2(DSP0) CALL ZERO2(DSPB) CALL ZERO2(DSPT) CALL ZERO2(PSP) C---------------------------------------------------------------------- AVCNVC=AVCNVC+1. ACUTIM=ACUTIM+1. DTCNVC=NCNVC*DT RDTCNVC=1./DTCNVC TAUK=DTCNVC/TREL c possible future change for shallow convection c TAUKSC=DTCNVC/(5.*TREL) C----------------------------------------------------------------------- C************************** DIAGNOSTICS ******************************** C----------------------------------------------------------------------- DO L=1,LM NTOPD (L)=0 NBOTD (L)=0 NTOPS (L)=0 NBOTS (L)=0 NDPTHS(L)=0 NDPTHD(L)=0 ENDDO C*********************************************************************** C---------------------------PREPARATIONS-------------------------------- C----------------------------------------------------------------------- !$omp parallel do DO 120 J=MYJS,MYJE DO 120 I=MYIS,MYIE LBOT (I,J)=LMH(I,J) THESP(I,J)=0. PDSL (I,J)=RES(I,J)*PD(I,J) PRECOL(I,J)=0. TREF(I,J,1)=T(I,J,1) 120 CONTINUE C----------------------------------------------------------------------- C--------------PADDING SPECIFIC HUMIDITY IF TOO SMALL------------------- C RESTORE APE TO SCRATCH ARRAY C----------------------------------------------------------------------- !$omp parallel do private(apests) DO 130 L=1,LM DO J=MYJS,MYJE DO I=MYIS,MYIE APESTS=PDSL(I,J)*AETA(L)+PT APE(I,J,L)=(1.E5/APESTS)**CAPA IF(Q(I,J,L).LT.EPSQ)Q(I,J,L)=HTM(I,J,L)*EPSQ ENDDO ENDDO 130 CONTINUE C----------------------------------------------------------------------- C--------------SEARCH FOR MAXIMUM BUOYANCY LEVEL------------------------ C----------------------------------------------------------------------- DO 170 KB=1,LM C----------------------------------------------------------------------- C--------------TRIAL MAXIMUM BUOYANCY LEVEL VARIABLES------------------- C----------------------------------------------------------------------- !$omp parallel do !$omp& private(apesp,bq,bqs00k,bqs10k,iq,it,ittb,ittbk,iqtb, !$omp& lmhk,p00k,p01k,p10k,p11k,pkl,pp1,psfck,qbt,qq1, !$omp& sq,sqs00k,sqs10k,tpsp,tq,tth,tthbt,tthes) DO 155 J=MYJS2,MYJE2 DO 150 I=MYIS1,MYIE1 C PKL=AETA(KB)*PDSL(I,J)+PT LMHK=LMH(I,J) PSFCK=AETA(LMHK)*PDSL(I,J)+PT c now searching over a scaled depth in finding the parcel c with the max theta-e instead of the old 130 mb c IF(PKL.LT.PSFCK-PBM)GO TO 150 c IF(KB.GT.LMHK)GO TO 150 IF(KB.LE.LMHK .AND. PKL.GE.0.80*PSFCK) THEN QBT=Q(I,J,KB) TTHBT=T(I,J,KB)*APE(I,J,KB) TTH=(TTHBT-THL)*RDTH QQ1=TTH-AINT(TTH) ITTB=INT(TTH)+1 C--------------KEEPING INDICES WITHIN THE TABLE------------------------- IF(ITTB.LT.1)THEN ITTB=1 QQ1=0. ENDIF IF(ITTB.GE.JTB)THEN ITTB=JTB-1 QQ1=0. ENDIF CONTINUE C--------------BASE AND SCALING FACTOR FOR SPEC. HUMIDITY--------------- ITTBK=ITTB BQS00K=QS0(ITTBK) SQS00K=SQS(ITTBK) BQS10K=QS0(ITTBK+1) SQS10K=SQS(ITTBK+1) C--------------SCALING SPEC. HUMIDITY & TABLE INDEX--------------------- BQ = (BQS10K-BQS00K)*QQ1+BQS00K SQ = (SQS10K-SQS00K)*QQ1+SQS00K TQ = (QBT-BQ)/SQ*RDQ PP1 =TQ - AINT(TQ) IQTB=INT(TQ)+1 C--------------KEEPING INDICES WITHIN THE TABLE------------------------- IF(IQTB.LT.1)THEN IQTB=1 PP1=0. ENDIF IF(IQTB.GE.ITB)THEN IQTB=ITB-1 PP1=0. ENDIF C--------------SATURATION PRESSURE AT FOUR SURROUNDING TABLE PTS.------- IQ=IQTB IT=ITTB P00K=PTBL(IQ ,IT ) P10K=PTBL(IQ+1,IT ) P01K=PTBL(IQ ,IT+1) P11K=PTBL(IQ+1,IT+1) C--------------SATURATION POINT VARIABLES AT THE BOTTOM----------------- TPSP=P00K+(P10K-P00K)*PP1+(P01K-P00K)*QQ1 1 +(P00K-P10K-P01K+P11K)*PP1*QQ1 APESP=(1.E5/TPSP)**CAPA TTHES=TTHBT*EXP(ELOCP*QBT*APESP/TTHBT) C--------------CHECK FOR MAXIMUM BUOYANCY------------------------------- IF(TTHES.GT.THESP(I,J))THEN PSP (I,J)=TPSP THBT (I,J)=TTHBT THESP(I,J)=TTHES ENDIF ENDIF 150 CONTINUE 155 CONTINUE C----------------------------------------------------------------------- 170 CONTINUE C----------------------------------------------------------------------- C---------CHOOSE CLOUD BASE AS MODEL LEVEL JUST BELOW PSP-------------- C----------------------------------------------------------------------- DO 240 L=1,LM1 AETAL=AETA(L) C !$omp parallel do DO J=MYJS2,MYJE2 DO I=MYIS,MYIE P(I,J)=PDSL(I,J)*AETAL+PT IF(P(I,J).LT.PSP(I,J).AND.P(I,J).GE.PQM)LBOT(I,J)=L+1 ENDDO ENDDO C 240 CONTINUE C*** C*** WARNING: LBOT MUST NOT BE GT LMH(I,J)-1 IN SHALLOW CONVECTION C*** MAKE SURE CLOUD BASE IS AT LEAST PONE ABOVE THE SURFACE C*** !$omp parallel do private(lmhij,psfck) DO 250 J=MYJS2,MYJE2 DO 250 I=MYIS,MYIE LMHIJ=LMH(I,J) PBOT(I,J)=AETA(LBOT(I,J))*PDSL(I,J)+PT PSFCK=AETA(LMHIJ)*PDSL(I,J)+PT IF(PBOT(I,J).GE.PSFCK-PONE.OR.LBOT(I,J).GE.LMHIJ)THEN C DO L=1,LMHIJ-1 P(I,J)=AETA(L)*PDSL(I,J)+PT IF(P(I,J).LT.PSFCK-PONE)LBOT(I,J)=L ENDDO C PBOT(I,J)=AETA(LBOT(I,J))*PDSL(I,J)+PT ENDIF 250 CONTINUE C--------------CLOUD TOP COMPUTATION------------------------------------ !$omp parallel do DO J=MYJS,MYJE DO I=MYIS,MYIE LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) ENDDO ENDDO C----------------------------------------------------------------------- !$omp parallel do !$omp& private(ihres,ilres,iptb,ithtb,jhres,jlres, !$omp& knumh,knuml,pp,presk,qq) C DO 290 L=LM,1,-1 C C--------------SCALING PRESSURE & TT TABLE INDEX------------------------ KNUML=0 KNUMH=0 C DO 270 J=MYJS2,MYJE2 DO 270 I=MYIS1,MYIE1 PRESK=PDSL(I,J)*AETA(L)+PT IF(PRESK.LT.PLQ)THEN KNUML=KNUML+1 ILRES(KNUML)=I JLRES(KNUML)=J ELSE KNUMH=KNUMH+1 IHRES(KNUMH)=I JHRES(KNUMH)=J ENDIF 270 CONTINUE C*** C*** COMPUTE PARCEL TEMPERATURE ALONG MOIST ADIABAT FOR PRESSUREPL C** IF(KNUMH.GT.0)THEN CALL TTBLEX(TREF(IDIM1,JDIM1,L),TTBLQ,ITBQ,JTBQ,KNUMH 1, IHRES,JHRES,PDSL,AETA(L),HTM(IDIM1,JDIM1,L) 2, PT,PLQ,QQ(IDIM1,JDIM1),PP(IDIM1,JDIM1) 3, RDPQ,THE0Q,STHEQ,RDTHEQ 4, THESP(IDIM1,JDIM1),IPTB(IDIM1,JDIM1) 5, ITHTB(IDIM1,JDIM1)) ENDIF 290 CONTINUE C--------------BUOYANCY CHECK------------------------------------------- DO 295 L=LM,1,-1 !$omp parallel do DO J=MYJS2,MYJE2 DO I=MYIS1,MYIE1 IF(TREF(I,J,L).GT.T(I,J,L)-DTTOP)LTOP(I,J)=L ENDDO ENDDO C----------------------------------------------------------------------- 295 CONTINUE C-----------------CLOUD TOP PRESSURE------------------------------------ !$omp parallel do DO J=MYJS2,MYJE2 DO I=MYIS1,MYIE1 PTOP(I,J)=AETA(LTOP(I,J))*PDSL(I,J)+PT ENDDO ENDDO C----------------------------------------------------------------------- C--------------DEFINE AND SMOOTH DSPS AND CLDEFI------------------------ C ************ UNIFIED OR SEPARATE LAND/SEA CONV *********************** C----------------------------------------------------------------------- C IF(UNIS)THEN !$omp parallel do private(efi) DO J=MYJS,MYJE DO I=MYIS,MYIE EFI=CLDEFI(I,J) DSPB(I,J)=(EFI-EFIMN)*SLOPBS+DSPBSS DSP0(I,J)=(EFI-EFIMN)*SLOP0S+DSP0SS DSPT(I,J)=(EFI-EFIMN)*SLOPTS+DSPTSS ENDDO ENDDO ELSEIF(.NOT.UNIL)THEN !$omp parallel do private(efi) DO J=MYJS,MYJE DO I=MYIS,MYIE EFI=CLDEFI(I,J) DSPB(I,J)=((EFI-EFIMN)*SLOPBS+DSPBSS)*SM(I,J) 1 +((EFI-EFIMN)*SLOPBL+DSPBSL)*(1.-SM(I,J)) DSP0(I,J)=((EFI-EFIMN)*SLOP0S+DSP0SS)*SM(I,J) 1 +((EFI-EFIMN)*SLOP0L+DSP0SL)*(1.-SM(I,J)) DSPT(I,J)=((EFI-EFIMN)*SLOPTS+DSPTSS)*SM(I,J) 1 +((EFI-EFIMN)*SLOPTL+DSPTSL)*(1.-SM(I,J)) ENDDO ENDDO ELSE !$omp parallel do private(efi) DO J=MYJS,MYJE DO I=MYIS,MYIE EFI=CLDEFI(I,J) DSPB(I,J)=((EFI-EFIMN)*SLOPBL+DSPBSL) DSP0(I,J)=((EFI-EFIMN)*SLOP0L+DSP0SL) DSPT(I,J)=((EFI-EFIMN)*SLOPTL+DSPTSL) ENDDO ENDDO ENDIF C C--------------EXTENDING SEA STRUCTURES INLAND ALONG COASTLINE---------- IF(NROW.GT.0.AND..NOT.UNIS.AND..NOT.UNIL)THEN !$omp parallel do DO J=MYJS,MYJE DO I=MYIS,MYIE WF1(I,J)=0. WF2(I,J)=0. WF3(I,J)=0. WF4(I,J)=0. ENDDO ENDDO C KROW=NROW C DO 350 KOUNT=1,KROW !$omp parallel do DO 345 J=MYJS2,MYJE2 DO 345 I=MYIS1,MYIE1 WF1(I,J)=(DSPB(I+IHW(J),J-1)+DSPB(I+IHE(J),J-1) 1 +DSPB(I+IHW(J),J+1)+DSPB(I+IHE(J),J+1)+4.*DSPB(I,J)) 2 *HBM2(I,J)*0.125 WF2(I,J)=(DSP0(I+IHW(J),J-1)+DSP0(I+IHE(J),J-1) 1 +DSP0(I+IHW(J),J+1)+DSP0(I+IHE(J),J+1)+4.*DSP0(I,J)) 2 *HBM2(I,J)*0.125 WF3(I,J)=(DSPT(I+IHW(J),J-1)+DSPT(I+IHE(J),J-1) 1 +DSPT(I+IHW(J),J+1)+DSPT(I+IHE(J),J+1)+4.*DSPT(I,J)) 2 *HBM2(I,J)*0.125 WF4(I,J)=(CLDEFI(I+IHW(J),J-1)+CLDEFI(I+IHE(J),J-1) 1 +CLDEFI(I+IHW(J),J+1)+CLDEFI(I+IHE(J),J+1)+4.*CLDEFI(I,J)) 2 *HBM2(I,J)*0.125 345 CONTINUE C !$omp parallel do private(rsmk,smk) DO J=MYJS,MYJE DO I=MYIS,MYIE SMK =SM(I,J) RSMK=1.-SMK DSPB (I,J)=DSPB (I,J)*SMK+WF1(I,J)*RSMK DSP0 (I,J)=DSP0 (I,J)*SMK+WF2(I,J)*RSMK DSPT (I,J)=DSPT (I,J)*SMK+WF3(I,J)*RSMK CLDEFI(I,J)=CLDEFI(I,J)*SMK+WF4(I,J)*RSMK ENDDO ENDDO C 350 CONTINUE C----------------------------------------------------------------------- ENDIF C--------------INITIALIZE CHANGES OF T AND Q DUE TO CONVECTION---------- !$omp parallel do DO 360 L=1,LM DO J=MYJS,MYJE DO I=MYIS,MYIE TMOD(I,J,L)=0. QMOD(I,J,L)=0. ENDDO ENDDO 360 CONTINUE C--------------CLEAN UP AND GATHER DEEP CONVECTION POINTS--------------- !$omp parallel do DO 380 J=MYJS2,MYJE2 DO 380 I=MYIS1,MYIE1 IF(LTOP(I,J).GE.LBOT(I,J))THEN LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) ENDIF IF(HBM2(I,J).LT.0.90)THEN LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) ENDIF IF(PTOP(I,J).GT.PBOT(I,J)-PNO.OR.LTOP(I,J).GT.LBOT(I,J)-2) 1 CLDEFI(I,J)=AVGEFI*SM(I,J)+STEFI*(1.-SM(I,J)) 380 CONTINUE C KHDEEP=0 PSHNEW=20000. DO J=MYJS2,MYJE2 DO I=MYIS1,MYIE1 PSFCIJ=PD(I,J)+PT c depth of cloud required to make the point a deep convection point c is now a scaled value of the PSFC instead of 290 mb everywhere DEPMIN=PSHNEW*PSFCIJ*1.E-5 DEPTH=PBOT(I,J)-PTOP(I,J) IF(DEPTH .GE. DEPMIN) THEN c IF(PTOP(I,J).LT.PBOT(I,J)-PSH)THEN KHDEEP=KHDEEP+1 IDEEP(KHDEEP)=I JDEEP(KHDEEP)=J ENDIF ENDDO ENDDO C************* HORIZONTAL LOOP FOR DEEP CONVECTION ********************* !$omp parallel do !$omp& private(apek,apekl,apekxx,apekxy,apesk,avrgt,avtgtl,dentpy, !$omp& depmin,depth,depwl,dhdt,difq,difql,dift,diftl,drheat, !$omp& drhdp,dsp,dsp0k,dspbk,dsptk,dthem,efi,fefi,hcorr, !$omp& i,j,l0,l0m1,lb,lbm1,lbtk,lcor,lqm,lshu,ltp1,ltp2, !$omp& ltpk,ltsh,pbtk,pk,pk0,pkb,pkl,pkt,preck,psfcij,psk, !$omp& pthrs,ptpk,qk,qkl,qrefk,qsatk,rdp0t,rhh,rhl,rhmax, !$omp& sumde,sumdp,therk,therkx,therky,thsk,thskl,tk,tkl, !$omp& trefk,trefkx,tskl) C*********************************************************************** DO 600 N=1,KHDEEP C*********************************************************************** I=IDEEP(N) J=JDEEP(N) PSFCIJ=PD(I,J)+PT LTPK=LTOP(I,J) LBTK=LBOT(I,J) CDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCD CDCDCDCDCDCDC DEEP CONVECTION DCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCD CDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCD LB =LBTK EFI =CLDEFI(I,J) DSPBK=DSPB (I,J) DSP0K=DSP0 (I,J) DSPTK=DSPT (I,J) C--------------INITIALIZE VARIABLES IN THE CONVECTIVE COLUMN------------ C*** C*** ONE SHOULD NOTE THAT THE VALUES ASSIGNED TO THE ARRAY TREFK C*** IN THE 410 LOOP ARE REALLY ONLY RELEVANT IN ANCHORING THE C*** REFERENCE TEMPERATURE PROFILE AT LEVEL LB. WHEN BUILDING THE C*** REFERENCE PROFILE FROM CLOUD BASE, THEN ASSIGNING THE C*** AMBIENT TEMPERATURE TO TREFK IS ACCEPTABLE. HOWEVER, WHEN C*** BUILDING THE REFERENCE PROFILE FROM SOME OTHER LEVEL (SUCH AS C*** ONE LEVEL ABOVE THE GROUND), THEN TREFK SHOULD BE FILLED WITH C*** THE TEMPERATURES IN TREF(I,J,L) WHICH ARE THE TEMPERATURES OF C*** THE MOIST ADIABAT THROUGH CLOUD BASE. BY THE TIME THE LINE C*** NUMBERED 450 HAS BEEN REACHED, TREFK ACTUALLY DOES HOLD THE C*** REFERENCE TEMPERATURE PROFILE. C*** DO 410 L=1,LM DIFT (L)=0. DIFQ (L)=0. TKL =T(I,J,L) TK (L)=TKL TREFK (L)=TKL QKL =Q(I,J,L) QK (L)=QKL QREFK (L)=QKL PKL =AETA(L)*PDSL(I,J)+PT PK (L)=PKL PSK (L)=PKL APEKL =APE(I,J,L) APEK (L)=APEKL THERK (L)=TREF(I,J,L)*APEKL C 410 THVMOD(L)=(QKL*0.608+1.)*TKL*APEKL 410 CONTINUE C--------------DEEP CONVECTION REFERENCE TEMPERATURE PROFILE------------ LTP1=LTPK+1 LBM1=LB-1 PKB=PK(LB) PKT=PK(LTPK) C--------------TEMPERATURE REFERENCE PROFILE BELOW FREEZING LEVEL------- L0=LB PK0=PK(LB) TREFKX=TREFK(LB) THERKX=THERK(LB) APEKXX=APEK(LB) THERKY=THERK(LBM1) APEKXY=APEK(LBM1) C DO 420 L=LBM1,LTPK,-1 IF(T(I,J,L+1).LT.TFRZ)GO TO 430 STABDL=STABD TREFKX=((THERKY-THERKX)*STABDL 1 +TREFKX*APEKXX)/APEKXY TREFK(L)=TREFKX APEKXX=APEKXY THERKX=THERKY APEKXY=APEK(L-1) THERKY=THERK(L-1) L0=L PK0=PK(L0) 420 CONTINUE C--------------FREEZING LEVEL AT OR ABOVE THE CLOUD TOP----------------- L0M1=L0-1 GO TO 450 C--------------TEMPERATURE REFERENCE PROFILE ABOVE FREEZING LEVEL------- 430 L0M1=L0-1 RDP0T=1./(PK0-PKT) DTHEM=THERK(L0)-TREFK(L0)*APEK(L0) C DO L=LTPK,L0M1 TREFK(L)=(THERK(L)-(PK(L)-PKT)*DTHEM*RDP0T)/APEK(L) ENDDO C----------------------------------------------------------------------- C--------------DEEP CONVECTION REFERENCE HUMIDITY PROFILE--------------- C----------------------------------------------------------------------- c reference profile had been too dry in the case where the cloud c base was close to the freezing level - this is now changed 450 DEPTH=PFRZ*PSFCIJ*1.E-5 DEPWL=PKB-PK0 DO 460 L=LTPK,LB C----------------------------------------------------------------------- C--------------SATURATION PRESSURE DIFFERENCE--------------------------- C----------------------------------------------------------------------- c IF(PKB-PK0.GE.PFRZ)THEN IF(DEPWL .GE. DEPTH) THEN IF(L.LT.L0)THEN DSP=((PK0-PK(L))*DSPTK+(PK(L)-PKT)*DSP0K)/(PK0-PKT) ELSE DSP=((PKB-PK(L))*DSP0K+(PK(L)-PK0)*DSPBK)/(PKB-PK0) ENDIF ELSE c DSP=DSPC DSP=DSP0K IF(L.LT.L0) 1 DSP=( (PK0-PK(L))*DSPTK+(PK(L)-PKT)*DSP0K)/(PK0-PKT) ENDIF C--------------HUMIDITY PROFILE----------------------------------------- IF(PK(L).GT.PQM)THEN PSK(L)=PK(L)+DSP APESK(L)=(1.E5/PSK(L))**CAPA THSK(L)=TREFK(L)*APEK(L) QREFK(L)=PQ0/PSK(L)*EXP(A2*(THSK(L)-A3*APESK(L)) 1 /(THSK(L)-A4*APESK(L))) ELSE QREFK(L)=Q(I,J,L) ENDIF 460 CONTINUE C--------------ENTHALPY CONSERVATION INTEGRAL-------------------------- DO 520 ITER=1,2 C----------------------------------------------------------------------- SUMDE=0. SUMDP=0. C DO L=LTPK,LB SUMDE=((TK(L)-TREFK(L))*CP+(QK(L)-QREFK(L))*ELWV)*DETA(L) 1 +SUMDE SUMDP=SUMDP+DETA(L) ENDDO C HCORR=SUMDE/(SUMDP-DETA(LTPK)) LCOR=LTPK+1 C-----------------------FIND LQM---------------------------------------- DO L=1,LB IF(PK(L).LE.PQM)LQM=L ENDDO C--------------ABOVE LQM CORRECT TEMPERATURE ONLY----------------------- IF(LCOR.LE.LQM)THEN DO L=LCOR,LQM TREFK(L)=TREFK(L)+HCORR*RCP ENDDO LCOR=LQM+1 ENDIF C--------------BELOW LQM CORRECT BOTH TEMPERATURE AND MOISTURE---------- C DO 510 L=LCOR,LB TSKL=TREFK(L)*APEK(L)/APESK(L) DHDT=QREFK(L)*A23M4L/(TSKL-A4)**2+CP TREFK(L)=HCORR/DHDT+TREFK(L) THSKL=TREFK(L)*APEK(L) QREFK(L)=PQ0/PSK(L)*EXP(A2*(THSKL-A3*APESK(L)) 1 /(THSKL-A4*APESK(L))) 510 CONTINUE C----------------------------------------------------------------------- 520 CONTINUE C--------------HEATING, MOISTENING, PRECIPITATION----------------------- DENTPY=0. AVRGT =0. PRECK =0. C DO 530 L=LTPK,LB TKL =TK(L) DIFTL =(TREFK(L)-TKL )*TAUK DIFQL =(QREFK(L)-QK(L))*TAUK AVRGTL =(TKL+TKL+DIFTL) DENTPY =(DIFTL*CP+DIFQL*ELWV)*DETA(L)/AVRGTL+DENTPY AVRGT =AVRGTL*DETA(L)+AVRGT PRECK =DETA(L)*DIFTL+PRECK DIFT(L)=DIFTL DIFQ(L)=DIFQL 530 CONTINUE C DENTPY=DENTPY+DENTPY AVRGT =AVRGT/(SUMDP+SUMDP) C-------------SWAP IF ENTROPY AND/OR PRECIP .LT. 0 ...------------------ IF(DENTPY.LT.EPSNTP.OR.PRECK.LT.0.)THEN IF(OCT90)THEN CLDEFI(I,J)=EFIMN ELSE CLDEFI(I,J)=EFIMN*SM(I,J)+STEFI*(1.-SM(I,J)) ENDIF C C--------------SEARCH FOR SHALLOW CLOUD TOP----------------------------- LBTK=LBOT(I,J) LTSH=LBTK LBM1=LBTK-1 PBTK=PK(LBTK) c use new threshold for cloud depth c PTPK=PBTK-PSH PSFCIJ=PD(I,J)+PT DEPMIN=PSHNEW*PSFCIJ*1.E-5 PTPK=PBTK-DEPMIN C-------------CLOUD TOP IS THE LEVEL JUST BELOW PBTK-PSH---------------- DO L=1,LM IF(PK(L).LE.PTPK)LTPK=L+1 ENDDO PTPK=PK(LTPK) C--------------HIGHEST LEVEL ALLOWED IS LEVEL JUST BELOW PSHU----------- IF(PTPK.LE.PSHU)THEN DO L=1,LM IF(PK(L).LE.PSHU)LSHU=L+1 ENDDO LTPK=LSHU PTPK=PK(LTPK) ENDIF C LTP1=LTPK+1 LTP2=LTPK+2 C----------------------------------------------------------------------- DO L=LTPK,LBTK QSATK(L)=PQ0/PK(L)*EXP(A2*(TK(L)-A3)/(TK(L)-A4)) ENDDO C----------------------------------------------------------------------- RHH=QK(LTPK)/QSATK(LTPK) RHMAX=0. C DO 570 L=LTP1,LBM1 RHL=QK(L)/QSATK(L) DRHDP=(RHH-RHL)/(PK(L-1)-PK(L)) IF(DRHDP.GT.RHMAX)THEN LTSH=L-1 RHMAX=DRHDP ENDIF RHH=RHL 570 CONTINUE C LTOP(I,J)=LTSH C---------------CLOUD MUST BE AT LEAST TWO LAYERS THICK----------------- IF(LBOT(I,J)-LTOP(I,J).LT.2)LTOP(I,J)=LBOT(I,J)-2 C PTOP(I,J)=PK(LTOP(I,J)) GO TO 600 ENDIF C----------------------------------------------------------------------- C--------------... DEEP CONVECTION OTHERWISE---------------------------- C----------------------------------------------------------------------- DRHEAT=(PRECK*SM(I,J)+AMAX1(EPSP,PRECK)*(1.-SM(I,J))) 1 *CP/AVRGT EFI=EFIFC*DENTPY/DRHEAT C C************** UNIFIED OR SEPARATE LAND/SEA CONV. ************** C IF(OCT90)THEN IF(UNIS)THEN EFI=CLDEFI(I,J)*FCB+EFI*FCC ELSEIF(.NOT.UNIL)THEN EFI=(CLDEFI(I,J)*FCB+EFI*FCC)*SM(I,J)+1.-SM(I,J) ELSE EFI=1. ENDIF ELSE EFI=CLDEFI(I,J)*FCB+EFI*FCC ENDIF C IF(EFI.GT.1. ) EFI=1. IF(EFI.LT.EFIMN) EFI=EFIMN IF(PRECK.EQ.0.) EFI=1. CLDEFI(I,J)=EFI C FEFI=EFMNT+SLOPE*(EFI-EFIMN) C FEFI=AMAX1(EFI,EFMNT) C PRECK=PRECK*FEFI C C--------------UPDATE PRECIPITATION, TEMPERATURE & MOISTURE------------- C PRECOL(I,J)=PDSL(I,J)*PRECK*CPRLG PREC (I,J)=PDSL(I,J)*PRECK*CPRLG+PREC (I,J) CUPREC(I,J)=PDSL(I,J)*PRECK*CPRLG+CUPREC(I,J) ACPREC(I,J)=PDSL(I,J)*PRECK*CPRLG+ACPREC(I,J) CUPPT(I,J) =PDSL(I,J)*PRECK*CPRLG+CUPPT(I,J) C DO L=LTPK,LB TMOD(I,J,L)=DIFT(L)*FEFI QMOD(I,J,L)=DIFQ(L)*FEFI CYL TLAT(I,J,L)=DIFT(L)*FEFI + TLAT(I,J,L) CYL ENDDO C CDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCD CDCDCDCDCDCDC END OF DEEP CONVECTION DCDCDCDCDCDCDCD CDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCDCD C C----------------------------------------------------------------------- 600 CONTINUE C----------------------------------------------------------------------- NDEEP=0 C DO 620 J=MYJS2,MYJE2 DO 620 I=MYIS,MYIE LTPK=LTOP(I,J) LBTK=LBOT(I,J) LB =LMH(I,J)-1 PSFCIJ=PD(I,J)+PT DEPMIN=PSHNEW*PSFCIJ*1.E-5 IF(PTOP(I,J).LT.PBOT(I,J)-DEPMIN)THEN NDEEP=NDEEP+1 NDEPTH=LB-LTPK NTOPD (LTPK )=NTOPD (LTPK )+1 NBOTD (LB )=NBOTD (LB )+1 NDPTHD(NDEPTH)=NDPTHD(NDEPTH)+1 ENDIF 620 CONTINUE NNEG=KHDEEP-NDEEP C C--------------GATHER SHALLOW CONVECTION POINTS------------------------- C KHSHAL=0 NDSTN =0 NDSTP =0 C DO 630 J=MYJS2,MYJE2 DO 630 I=MYIS,MYIE IF(PTOP(I,J).GT.PBOT(I,J)-PNO.OR. 1 LTOP(I,J).GT.LBOT(I,J)-2)GO TO 630 PSFCIJ=PD(I,J)+PT DEPMIN=PSHNEW*PSFCIJ*1.E-5 IF(PTOP(I,J)+1..GE.PBOT(I,J)-DEPMIN)THEN KHSHAL=KHSHAL+1 ISHAL(KHSHAL)=I JSHAL(KHSHAL)=J ENDIF 630 CONTINUE C C************* HORIZONTAL LOOP FOR SHALLOW CONVECTION ****************** !$omp parallel do !$omp& private(apek,apekl,apekxx,apekxy,bqk,bqs00k,bqs10k,den,dentpy, !$omp& dpkl,dpmix,dqref,dst,dstq,dtdeta,fpk,fptk,i,iq,it,j, !$omp& lbm1,lbtk,ltp1,ltpk,otsum,part1,part2,part3,pk,pkl, !$omp& pkxxxx,pkxxxy,potsum,ppk,psum,ptpk,pz0,qk,qkl,qnew, !$omp& qotsum,qqk,qrefk,qrfkl,qrftp,qsatk,qsum,rdpsum,rtbar, !$omp& smix,sqk,sqs00k,sqs10k,sumdp,sumdt,tcorr,thvmkl, !$omp& thvref,tk,tkl,tqk,trefk,trefkx,trfkl,tthk) C DO 800 N=1,KHSHAL C*********************************************************************** I=ISHAL(N) J=JSHAL(N) C----------------------------------------------------------------------- CSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCS CSCSCSCSCSCSC SHALLOW CONVECTION CSCSCSCSCSCSCSCSCSCSCS CSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCS C----------------------------------------------------------------------- PZ0=PD(I,J)+PT LLMH=LMH(I,J) C DO 650 L=1,LLMH TKL =T (I,J,L) TK (L) =TKL TREFK(L) =TKL QKL =Q (I,J,L) QK (L) =QKL QREFK(L) =QKL PKL =AETA(L)*PDSL(I,J)+PT PK (L) =PKL QSATK(L) =PQ0/PK(L)*EXP(A2*(TK(L)-A3)/(TK(L)-A4)) APEKL =APE(I,J,L) CVVVVVVVVVVVV CHOOSE THE PRESSURE FUNCTION VVVVVVVVVVVVVVVVVVVVVVVVVVVVV C FPK (L) =ALOG(PKL) C FPK (L) =PKL C FPK (L) =-1./PKL CAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA APEK (L) =APEKL THVMKL =TKL*APEKL*(QKL*0.608+1.) C THVMOD(L)=THVMKL THVREF(L)=THVMKL 650 CONTINUE C--------------SHALLOW CLOUD TOP----------------------------- LBTK=LBOT(I,J) LBM1=LBTK-1 PTPK=PTOP(I,J) LTP1=LTOP(I,J) LTPK=LTOP(I,J)-1 C----------------------------------------------------------------------- IF(PTOP(I,J).GT.PBOT(I,J)-PNO.OR.LTOP(I,J).GT.LBOT(I,J)-2)THEN LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) GO TO 800 ENDIF C--------------SCALING POTENTIAL TEMPERATURE & TABLE INDEX AT TOP------- C C THTPK=T(I,J,LTP1)*APE(I,J,LTP1) THTPK=T(I,J,LTPK)*APE(I,J,LTPK) C TTHK =(THTPK-THL)*RDTH QQK =TTHK-AINT(TTHK) IT =INT(TTHK)+1 IF(IT.LT.1)THEN IT=1 QQK=0. ENDIF IF(IT.GE.JTB)THEN IT=JTB-1 QQK=0. ENDIF C--------------BASE AND SCALING FACTOR FOR SPEC. HUMIDITY AT TOP-------- BQS00K=QS0(IT) SQS00K=SQS(IT) BQS10K=QS0(IT+1) SQS10K=SQS(IT+1) C--------------SCALING SPEC. HUMIDITY & TABLE INDEX AT TOP-------------- BQK=(BQS10K-BQS00K)*QQK+BQS00K SQK=(SQS10K-SQS00K)*QQK+SQS00K C C TQK=(Q(I,J,LTP1)-BQK)/SQK*RDQ TQK=(Q(I,J,LTPK)-BQK)/SQK*RDQ C PPK=TQK-AINT(TQK) IQ =INT(TQK)+1 IF(IQ.LT.1)THEN IQ=1 PPK=0. ENDIF IF(IQ.GE.ITB)THEN IQ=ITB-1 PPK=0. ENDIF C--------------CLOUD TOP SATURATION POINT PRESSURE---------------------- PART1=(PTBL(IQ+1,IT)-PTBL(IQ,IT))*PPK PART2=(PTBL(IQ,IT+1)-PTBL(IQ,IT))*QQK PART3=(PTBL(IQ ,IT )-PTBL(IQ+1,IT ) 1 -PTBL(IQ ,IT+1)+PTBL(IQ+1,IT+1))*PPK*QQK PTPK=PTBL(IQ,IT)+PART1+PART2+PART3 C----------------------------------------------------------------------- DPMIX=PTPK-PSP(I,J) IF(ABS(DPMIX).LT.3000.)DPMIX=-3000. C--------------TEMPERATURE PROFILE SLOPE-------------------------------- SMIX=(THTPK-THBT(I,J))/DPMIX*STABS C TREFKX=TREFK(LBTK+1) PKXXXX=PK(LBTK+1) PKXXXY=PK(LBTK) APEKXX=APEK(LBTK+1) APEKXY=APEK(LBTK) C DO 670 L=LBTK,LTP1,-1 TREFKX=((PKXXXY-PKXXXX)*SMIX 1 +TREFKX*APEKXX)/APEKXY TREFK(L)=TREFKX APEKXX=APEKXY PKXXXX=PKXXXY APEKXY=APEK(L-1) PKXXXY=PK(L-1) 670 CONTINUE C--------------TEMPERATURE REFERENCE PROFILE CORRECTION----------------- SUMDT=0. SUMDP=0. C DO L=LTP1,LBTK SUMDT=(TK(L)-TREFK(L))*DETA(L)+SUMDT SUMDP=SUMDP+DETA(L) ENDDO C RDPSUM=1./SUMDP FPK(LBTK)=TREFK(LBTK) C TCORR=SUMDT*RDPSUM C DO L=LTP1,LBTK C TCORR=SUMDT/(SUMDP-DETA(LBTK)) C DO L=LTP1,LBM1 C TRFKL =TREFK(L)+TCORR TREFK(L)=TRFKL FPK (L)=TRFKL ENDDO C--------------HUMIDITY PROFILE EQUATIONS------------------------------- PSUM =0. QSUM =0. POTSUM=0. QOTSUM=0. OTSUM =0. DST =0. FPTK =FPK(LTP1) C DO 700 L=LTP1,LBTK DPKL =FPK(L)-FPTK PSUM =DPKL *DETA(L)+PSUM QSUM =QK(L)*DETA(L)+QSUM RTBAR =2./(TREFK(L)+TK(L)) OTSUM =DETA(L)*RTBAR+OTSUM POTSUM=DPKL *RTBAR*DETA(L)+POTSUM QOTSUM=QK(L) *RTBAR*DETA(L)+QOTSUM DST =(TREFK(L)-TK(L))*RTBAR*DETA(L)+DST 700 CONTINUE C PSUM =PSUM*RDPSUM QSUM =QSUM*RDPSUM ROTSUM=1./OTSUM POTSUM=POTSUM*ROTSUM QOTSUM=QOTSUM*ROTSUM DST =DST *ROTSUM*CP/ELWV C--------------ENSURE POSITIVE ENTROPY CHANGE--------------------------- IF(DST.GT.0.)THEN C C DSTQ=DST*EPSUP LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) cccc NDSTP=NDSTP+1 GO TO 800 C ELSE DSTQ=DST*EPSDN ENDIF C--------------CHECK FOR ISOTHERMAL ATMOSPHERE-------------------------- DEN=POTSUM-PSUM C IF(-DEN/PSUM.LT.5.E-5)THEN LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) GO TO 800 C C--------------SLOPE OF THE REFERENCE HUMIDITY PROFILE------------------ ELSE DQREF=(QOTSUM-DSTQ-QSUM)/DEN ENDIF C------------ HUMIDITY DOES NOT INCREASE WITH HEIGHT-------------------- IF(DQREF.LT.0.)THEN LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) GO TO 800 ENDIF C--------------HUMIDITY AT THE CLOUD TOP-------------------------------- QRFTP=QSUM-DQREF*PSUM C--------------HUMIDITY PROFILE----------------------------------------- C DO 720 L=LTP1,LBTK QRFKL=(FPK(L)-FPTK)*DQREF+QRFTP C C--------------SUPERSATURATION OR NEGATIVE Q NOT ALLOWED---------------- C QNEW =(QRFKL-QK(L))*TAUK+QK(L) IF(QNEW.GT.QSATK(L)*STRESH.OR.QNEW.LT.0.)THEN LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) GO TO 800 ENDIF C C----------------------------------------------------------------------- THVREF(L)=TREFK(L)*APEK(L)*(QRFKL*0.608+1.) QREFK(L)=QRFKL 720 CONTINUE C--------------ELIMINATE IMPOSSIBLE SLOPES (BETTS, DTHETA/DQ)----------- DO 730 L=LTP1,LBTK DTDETA=(THVREF(L-1)-THVREF(L))/(AETA(L)-AETA(L-1)) IF(DTDETA.LT.EPSTH)THEN LBOT(I,J)=0 LTOP(I,J)=LBOT(I,J) PTOP(I,J)=PBOT(I,J) GO TO 800 ENDIF 730 CONTINUE C*************************** DIAGNOSTICS ****************************** cccc IF(DST.GT.0.)THEN cccc NDSTP=NDSTP+1 cccc ELSE cccc NDSTN=NDSTN+1 cccc ENDIF C********************************************************************** DENTPY=0. C DO L=LTP1,LBTK DENTPY=((TREFK(L)-TK(L))*CP+(QREFK(L)-QK(L))*ELWV) 1 /(TK(L)+TREFK(L))*DETA(L)+DENTPY ENDDO C C----------------------------------------------------------------------- C--------------RELAXATION TOWARDS REFERENCE PROFILES-------------------- C----------------------------------------------------------------------- C C DO 750 L=LTP1,LBTK TMOD(I,J,L)=(TREFK(L)-TK(L))*TAUK QMOD(I,J,L)=(QREFK(L)-QK(L))*TAUK 750 CONTINUE C----------------------------------------------------------------------- CSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCS CSCSCSCSCSCSC END OF SHALLOW CONVECTION SCSCSCSCSCSCSCSCS CSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCSCS C----------------------------------------------------------------------- 800 CONTINUE C----------------------------------------------------------------------- C C************************** DIAGNOSTICS ******************************** NSHAL=0 C DO 820 J=MYJS2,MYJE2 DO 820 I=MYIS,MYIE LTPK=LTOP(I,J) LBTK=LBOT(I,J) PTPK=PTOP(I,J) PBTK=PBOT(I,J) IF(PTPK.GT.PBTK-PNO.OR.LTPK.GT.LBTK-2)GO TO 820 PSFCIJ=PD(I,J)+PT DEPMIN=PSHNEW*PSFCIJ*1.E-5 C IF(PTPK.GE.PBTK-PSH)THEN IF(PTPK.GE.PBTK-DEPMIN)THEN NSHAL=NSHAL+1 NTOPS(LTPK)=NTOPS(LTPK)+1 NBOTS(LBTK)=NBOTS(LBTK)+1 NDEPTH=LBTK-LTPK NDPTHS(NDEPTH)=NDPTHS(NDEPTH)+1 ENDIF 820 CONTINUE NEGDS=KHSHAL-NSHAL C*********************************************************************** C C--------------SMOOTHING TEMPERATURE & HUMIDITY CORRECTIONS------------- IF(KSMUD.EQ.0)THEN !$omp parallel do C DO 900 L=1,LM C*** C*** UPDATE THE FUNDAMENTAL PROGNOSTIC ARRAYS C*** DO 830 J=MYJS,MYJE DO 830 I=MYIS,MYIE T(I,J,L)=T(I,J,L)+TMOD(I,J,L) Q(I,J,L)=Q(I,J,L)+QMOD(I,J,L) C*** C*** ACCUMULATE LATENT HEATING DUE TO CONVECTION. C*** SCALE BY THE RECIPROCAL OF THE PERIOD AT WHICH THIS ROUTINE C*** IS CALLED. THIS PERIOD IS THE CONVECTION TIMESTEP. C*** TCUCN(I,J,L)=TCUCN(I,J,L)+TMOD(I,J,L)*RDTCNVC 830 CONTINUE 900 CONTINUE ELSE !$omp parallel do private(ql,qne,qse,tl,tne,tse) C DO 910 L=1,LM C CALL ZERO2(QL) CALL ZERO2(QNE) CALL ZERO2(QSE) CALL ZERO2(TL) CALL ZERO2(TNE) CALL ZERO2(TSE) C----------------------------------------------------------------------- DO J=MYJS,MYJE DO I=MYIS,MYIE TL(I,J)=TMOD(I,J,L) QL(I,J)=QMOD(I,J,L) ENDDO ENDDO C----------------------------------------------------------------------- C----------------------------------------------------------------------- NSMUD=KSMUD C DO 870 KS=1,NSMUD C DO J=MYJS,MYJE1 DO I=MYIS,MYIE TNE(I,J)=(TL(I+IHE(J),J+1)-TL(I,J)) 1 *HTM(I,J,L)*HTM(I+IHE(J),J+1,L) QNE(I,J)=(QL(I+IHE(J),J+1)-QL(I,J)) 1 *HTM(I,J,L)*HTM(I+IHE(J),J+1,L) ENDDO ENDDO C DO J=MYJS1,MYJE DO I=MYIS,MYIE TSE(I,J)=(TL(I+IHE(J),J-1)-TL(I,J)) 1 *HTM(I+IHE(J),J-1,L)*HTM(I,J,L) QSE(I,J)=(QL(I+IHE(J),J-1)-QL(I,J)) 1 *HTM(I+IHE(J),J-1,L)*HTM(I,J,L) ENDDO ENDDO C DO J=MYJS2,MYJE2 DO I=MYIS,MYIE TL(I,J)=(TNE(I,J)-TNE(I+IHW(J),J-1)+TSE(I,J)-TSE(I+IHW(J),J+1)) 1 *HBM2(I,J)*0.125+TL(I,J) QL(I,J)=(QNE(I,J)-QNE(I+IHW(J),J-1)+QSE(I,J)-QSE(I+IHW(J),J+1)) 1 *HBM2(I,J)*0.125+QL(I,J) ENDDO ENDDO C 870 CONTINUE C----------------------------------------------------------------------- C*** C*** UPDATE THE FUNDAMENTAL PROGNOSTIC ARRAYS C*** DO J=MYJS,MYJE DO I=MYIS,MYIE T(I,J,L)=T(I,J,L)+TL(I,J) Q(I,J,L)=Q(I,J,L)+QL(I,J) ENDDO ENDDO C*** C*** ACCUMULATE LATENT HEATING DUE TO CONVECTION. C*** SCALE BY THE RECIPROCAL OF THE PERIOD AT WHICH THIS ROUTINE C*** IS CALLED. THIS PERIOD IS THE CONVECTION TIMESTEP. C*** DO J=MYJS,MYJE DO I=MYIS,MYIE TCUCN(I,J,L)=TCUCN(I,J,L)+TL(I,J)*RDTCNVC ENDDO ENDDO C 910 CONTINUE C ENDIF C--------------SAVE CLOUD TOP AND BOTTOM FOR RADIATION------------------ !$omp parallel do DO J=MYJS,MYJE DO I=MYIS,MYIE IF (LTOP(I,J).GT.0 .AND. LTOP(I,J).LT.LBOT(I,J)) THEN CUtop=FLOAT(LTOP(I,J)) CUbot=FLOAT(LBOT(I,J)) HTOP(I,J)=MIN(CUtop,HTOP(I,J)) HBOT(I,J)=MAX(CUbot,HBOT(I,J)) CNVTOP(I,J)=MIN(CUtop,CNVTOP(I,J)) CNVBOT(I,J)=MAX(CUbot,CNVBOT(I,J)) ENDIF ENDDO ENDDO C----------------------------------------------------------------------- C************************* DIAGNOSTICS ********************************* C C WRITE(LIST,950)NTSD,NSHAL,NDEEP,NNEG,NEGDS,NDSTN,NDSTP C DO 940 L=1,LM C WRITE(LIST,952)L C WRITE(LIST,954)NBOTS(L),NTOPS(L),NDPTHS(L) C 1, NBOTD(L),NTOPD(L),NDPTHD(L) C 940 CONTINUE C 950 FORMAT(' NTSD=',I3,I8,' SHALLOW, ',I4,' DEEP, ', C 1 I4,' NEG., ',I4,' NEG. SHALL.,',I4,' DST.LT.0, ',I4,' DST.GT.0') C 952 FORMAT(' LAYER (FROM TOP),',I2) C 954 FORMAT(' NBOTS=',I4,' NTOPS=',I4,' NDPTHS=',I4, C 1 ' NBOTD=',I4,' NTOPD=',I4,' NDPTHD=',I4) C*********************************************************************** RETURN END