SUBROUTINE PRECPD
C ******************************************************************
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
C . . .
C SUBPROGRAM: PRECPD LARGE SCALE PRECIPITATION
C PRGRMMR: ZHAO ORG: W/NP22 DATE: ??-??-??
C
C ABSTRACT:
C PRECPD COMPUTES THE GRID SCALE PRECIPITATION.
C
C PROGRAM HISTORY LOG:
C 94-??-?? ZHAO - ORIGINATOR
C 95-03-25 BLACK - CONVERSION FROM 1-D TO 2-D IN HORIZONTAL
C 95-11-20 ABELES - PARALLEL OPTIMIZATION
C 96-03-29 BLACK - REMOVED SCRCH COMMON
C 96-07-18 ZHAO - NEW WMIN CALCULATION
C 96-09-25 BALDWIN - NEW SR CALCULATION
C 98-11-02 BLACK - MODIFICATION FOR DISTRIBUTED MEMORY
C
C USAGE: CALL PRECPD 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: NONE
C
C LIBRARY: NONE
C
C COMMON BLOCKS: CTLBLK
C LOOPS
C MASKS
C PHYS
C VRBLS
C CLDWTR
C PVRBLS
C ACMCLH
C
C ATTRIBUTES:
C LANGUAGE: FORTRAN 90
C MACHINE : IBM SP
C$$$
C ******************************************************************
P A R A M E T E R
& (A1=610.78,A2=17.2693882,A3=273.16,A4=35.86
&, PQ0=379.90516,TRESH=.95,R=287.04,C0=0.15
&, CP=1004.6,ELWV=2.50E6,ELIV=2.834E6,ROW=1.E3,G=9.8
&, EPSQ=2.E-12,DLDT=2274.0,ELIW=ELIV-ELWV)
C
P A R A M E T E R
& (ARCP=A2*(A3-A4)/CP,RCP=1./CP,PQ0C=PQ0*TRESH,RROG=1./(ROW*G)
&, RROW=1./ROW)
C----------------------------------------------------------------------
INCLUDE "parmeta"
INCLUDE "parm.tbl"
INCLUDE "mpp.h"
#include "sp.h"
C----------------------------------------------------------------------
P A R A M E T E R
& (IMJM=IM*JM-JM/2,JAM=6+2*(JM-10)
&, LP1=LM+1,LTOP=1,LBOT=LM)
C
P A R A M E T E R
& (IMJM_LOC=IDIM2*JDIM2,LDA=(IDIM2-IDIM1+1)*(JDIM2-JDIM1+1))
C-----------------------------------------------------------------------
L O G I C A L
& RUN,FIRST,RESTRT,SIGMA,NOZ
C----------------------------------------------------------------------
INCLUDE "CTLBLK.comm"
C-----------------------------------------------------------------------
INCLUDE "LOOPS.comm"
C-----------------------------------------------------------------------
INCLUDE "MASKS.comm"
C-----------------------------------------------------------------------
INCLUDE "PHYS.comm"
C-----------------------------------------------------------------------
INCLUDE "VRBLS.comm"
C-----------------------------------------------------------------------
INCLUDE "CLDWTR.comm"
C-----------------------------------------------------------------------
INCLUDE "PVRBLS.comm"
C-----------------------------------------------------------------------
INCLUDE "ACMCLH.comm"
C-----------------------------------------------------------------------
INCLUDE "PPTASM.comm"
C-----------------------------------------------------------------------
D I M E N S I O N
& IPREC(IMJM_LOC),JPREC(IMJM_LOC)
&,T_T(LM,IDIM1:IDIM2,JDIM1:JDIM2)
&,Q_T(LM,IDIM1:IDIM2,JDIM1:JDIM2)
&,TRAIN_T(LM,IDIM1:IDIM2,JDIM1:JDIM2)
&,HTM_T(LM,IDIM1:IDIM2,JDIM1:JDIM2)
&,CWM_T(LM,IDIM1:IDIM2,JDIM1:JDIM2)
&,TLAT_T(LM,IDIM1:IDIM2,JDIM1:JDIM2)
C-----------------------------------------------------------------------
D I M E N S I O N
& PRECRL1(IMJM_LOC),PRECSL1(IMJM_LOC),IWL1(IMJM_LOC)
C-----------------------------------------------------------------------
R E A L
& KE,INIT,MI0
C-----------------------------------------------------------------------
C***********************************************************************
C--------------PREPARATORY CALCULATIONS---------------------------------
DTPH = NPHS * DT
RDTPH = 1. / DTPH
TWODT= DTPH
RTWODT=1./TWODT
KE=2.0E-5
US=1.
EPS=0.622E0
CCLIMIT=1.0E-3
CLIMIT=1.0E-20
CWS=0.025
CSM1=5.0000E-8
CRS1=5.00000E-6
CRS2=6.66600E-10
CR=5.0E-4
MI0=5.0E-4
AA2=1.25E-3
AVRAIN=AVRAIN+1.
ARATIM=ARATIM+1
C-------------------PADDING CLOUD MIXING RATIO IF TOO SMALL-------------
!$omp parallel do
DO 20 L=1,LM
DO J=MYJS,MYJE
DO I=MYIS,MYIE
CWM(I,J,L)=CWM(I,J,L)*HTM(I,J,L)*HBM2(I,J)
IF(CWM(I,J,L).LT.0.)CWM(I,J,L)=0.
C------------------PADDING SPECIFIC HUMIDITY IF TOO SMALL---------------
IF(Q(I,J,L).LT.EPSQ)Q(I,J,L)=EPSQ*HTM(I,J,L)
ENDDO
ENDDO
20 CONTINUE
C
UTIM=1.
C-----------------------------------------------------------------------
!$omp parallel do
DO N=1,IMJM_LOC
IWL1(N)=0
PRECRL1(N)=0.
PRECSL1(N)=0.
ENDDO
C------------CHOOSE THE COLUMNS WHERE PREC CAN BE PRODUCED--------------
NPRE=0
DO 35 J=MYJS2,MYJE2
DO 35 I=MYIS,MYIE
C
DO L=2,LM
TTEMP=0.025*(T(I,J,L)-273.16)
WFIX=0.9814*EXP(0.01873*L)
WMIN=0.1E-3*EXP(TTEMP)*WFIX
IF(CWM(I,J,L).GT.WMIN)GO TO 33
ENDDO
C
GO TO 35
33 NPRE=NPRE+1
IPREC(NPRE)=I
JPREC(NPRE)=J
35 CONTINUE
C------------------------------------------------------------------------
C***
C*** TRANSPOSE ARRAYS
C***
!$omp parallel sections
!$omp section
CALL SGETMO(T,LDA,LDA,LM,T_T,LM)
!$omp section
CALL SGETMO(Q,LDA,LDA,LM,Q_T,LM)
!$omp section
CALL SGETMO(CWM,LDA,LDA,LM,CWM_T,LM)
!$omp section
CALL SGETMO(HTM,LDA,LDA,LM,HTM_T,LM)
!$omp section
CALL SGETMO(TRAIN,LDA,LDA,LM,TRAIN_T,LM)
!$omp section
CALL SGETMO(TLAT,LDA,LDA,LM,TLAT_T,LM)
!$omp end parallel sections
C-----------------------------------------------------------------------
C-----------------BEGINING OF PRECIPITATION CALCULATION-----------------
C-----------------------------------------------------------------------
C***
C*** LOOP OVER ALL POSSIBLE PRECIPITATION POINTS
C***
!$omp parallel do
!$omp& private(aai,aetal,ai,amaxcm,amaxps,amaxrq,bi,ccr,conde,
!$omp& const,cpdr,cs,cwmk,detal,erk,err,errt,ers,erst,expf,fi,
!$omp& fiw,hbm2k,hh,htmk,i,iwl,j,lml,mi0,pdsl,pid,ppr,pps,
!$omp& pracw,praut,precrk,precrl,precsk,precsl,precss,psaci,
!$omp& psaut,psm,psm1,psm2,qc,qct,qi,qint,qintt,qit,qk,qq,
!$omp& qtemp,qw,qwt,rconde,rprs,rq,rqkll,rqt,rqtt,tk,tmt0,
!$omp& tmt0k,tmt0t,tmt15,tmt15t,totppt,tt,ttemp,ttlat,u00ij,
!$omp& u00kl,u00klt,ull,wfix,wmink,ww)
C
DO 300 N=1,NPRE
C
I=IPREC(N)
J=JPREC(N)
HBM2K=HBM2(I,J)
PDSL=RES(I,J)*PD(I,J)
CONST=PDSL/G*TWODT
LML=LM-LMH(I,J)
U00IJ=U00(I,J)
C
DO 180 L=2,LM
C
DETAL=DETA(L)
ULL=UL(L)
AETAL=AETA(L)
WFIX=0.9814*EXP(0.01873*L)
C
PRECRL=0.
PRECSL=0.
PRAUT=0.
PSAUT=0.
PRACW=0.
PSACI=0.
ERR =0.
ERS =0.
PSM =0.
PSM1 =0.
PSM2 =0.
PPR =0.
PPS =0.
CPDR =0.
HH =0.
PID =0.
IWL =0.
CONDE=0.
RCONDE=0.
C-----------------------------------------------------------------------
TT=T_T(L,I,J)
QQ=Q_T(L,I,J)
WW=CWM_T(L,I,J)
HTMK=HTM_T(L,I,J)
TTLAT=TLAT_T(L,I,J)
C-----------------------------------------------------------------------
U00KL=U00IJ+UL(L+LML)*(0.95-U00IJ)*UTIM
TTEMP=0.025*(TT-273.16)
WMINK=0.1E-3*EXP(TTEMP)*WFIX
C-----------------------------------------------------------------------
C----------CHOOSE THE POINTS WHERE PRECIPITATION CAN BE PRODUCED--------
C-----------------------------------------------------------------------
PRECRK=AMAX1(0.,PRECRL1(N))
PRECSK=AMAX1(0.,PRECSL1(N))
HH=HTMK*HBM2K
IF(WW.LT.WMINK.AND.(PRECRK+PRECSK).EQ.0.)THEN
PID=0.
ELSE
PID=HH
ENDIF
C-----------------------------------------------------------------------
C-------------------QW, QI AND QINT-------------------------------------
C-----------------------------------------------------------------------
IF(PID.EQ.1.)THEN
CONDE=CONST*DETAL
RCONDE=1./CONDE
TK=TT
QK=QQ
TMT0=(TK-273.16)*HH
TMT15=AMIN1(TMT0,-15.)*HH
AI=0.008855
BI=1.
C
IF(TMT0.LT.-20.)THEN
AI=0.007225
BI=0.9674
ENDIF
C
QW=HH*PQ0/(PDSL*AETAL+PT)
1 *EXP(HH*A2*(TK-A3)/(TK-A4))
QI=QW*(BI+AI*AMIN1(TMT0,0.))
QINT=QW*(1.-0.00032*TMT15*(TMT15+15.))
IF(TMT0.LE.-40.)QINT=QI
C-------------------ICE-WATER ID NUMBER IW------------------------------
IF(TMT0.LT.-15.)THEN
FI=QK-U00KL*QI
IF(FI.GT.0..OR.WW.GT.CLIMIT) THEN
IWL=1
ELSE
IWL=0
ENDIF
ENDIF
C
IF(TMT0.LT.0.0.AND.TMT0.GE.-15.0)THEN
IWL=0
IF(IWL1(N).EQ.1.AND.WW.GT.CLIMIT)IWL=1
ENDIF
C
IF(TMT0.GE.0.)THEN
IWL=0
ENDIF
C----------------THE SATUATION SPECIFIC HUMIDITY------------------------
FIW=FLOAT(IWL)
QC=(1.-FIW)*QINT+FIW*QI
C----------------THE RELATIVE HUMIDITY----------------------------------
IF(QC.LE.0.)THEN
RQ=1.E-10
ELSE
RQ=QK/QC
ENDIF
C----------------CLOUD COVER RATIO CCR----------------------------------
IF(RQ.LE.U00KL)THEN
CCR=0.
ELSE
RQKLL=AMIN1(US,RQ)
CCR=1.-SQRT((US-RQKLL)/(US-U00KL))
ENDIF
C-----------CORRECT CCR IF IT IS TOO SMALL IN LARGE CWM REGIONS--------
IF(CCR.GE.0.01.AND.CCR.LE.0.2.AND
1 .WW.GE.0.2E-3)THEN
CCR=AMIN1(1.,WW*1.0E3)
ENDIF
ENDIF
60 CONTINUE
C-----------------------------------------------------------------------
C------------------PRECIPITATION PRODUCTION RATES-----------------------
C------------------AUTO-CONVERT RATES-----------------------------------
C-----------------------------------------------------------------------
IF(PID.EQ.1.)THEN
IWK=IWL
CWMK=AMAX1(0.,WW-CLIMIT)
MI0=WMINK
C
IF(IWK.EQ.1)THEN
EXPF=EXP(0.025*TMT0)
AA1=1.E-3*EXPF
PSAUT=AA1*AMAX1(0.,CWMK-MI0)
CPDR=-PSAUT*TWODT
IF(-CPDR.GE.CWMK)THEN
CPDR=-CWMK
PSAUT=-CPDR*RTWODT
ENDIF
ELSE
AMAXCM=AMAX1(0.,CWMK-MI0)
PRAUT=C0*AMAXCM*AMAXCM
CPDR=-PRAUT*TWODT
IF(-CPDR.GE.CWMK)THEN
CPDR=-CWMK
PRAUT=-CPDR*RTWODT
ENDIF
ENDIF
PPR=PRAUT*CONDE
PPS=PSAUT*CONDE
ENDIF
C
IF(PID.EQ.1.)THEN
WW=CPDR*HH+WW
PRECRL=PRECRL1(N)+PPR*HH
PRECSL=PRECSL1(N)+PPS*HH
ENDIF
C-----------------------------------------------------------------------
C-----------------------ACCRETIONS--------------------------------------
C-----------------------------------------------------------------------
IF(PID.EQ.1.)THEN
IWK=IWL
CWMK=WW
PRECRK=AMAX1(0.,PRECRL1(N))
PRECSK=AMAX1(0.,PRECSL1(N))
IF(IWK.EQ.1)THEN
EXPF=EXP(0.025*TMT0)
CS=AA2*EXPF
PSACI=CS*AMAX1(0.,CWMK)*PRECSK
CPDR=-PSACI*TWODT
IF(-CPDR.GE.CWMK)THEN
CPDR=-CWMK
PSACI=-CPDR*RTWODT
ENDIF
ELSE
PRACW=CR*AMAX1(0.,CWMK)*(PRECRK+PRECSK)
CPDR=-PRACW*TWODT
IF(-CPDR.GE.CWMK)THEN
CPDR=-CWMK
PRACW=-CPDR*RTWODT
ENDIF
ENDIF
PPR=PRACW*CONDE
PPS=PSACI*CONDE
ENDIF
C
IF(PID.EQ.1.)THEN
WW=CPDR*HH+WW
PRECRL=PRECRL+PPR*HH
PRECSL=PRECSL+PPS*HH
ENDIF
C-----------------------------------------------------------------------
C-----EVAPORATION/CONDENSATION OF PRECIPITATION-------------------------
C***** ERR & ERS POSITIVE--EVAPORATION
C***** ERR & ERS NEGTIVE---CONDENSATION
C-----------------------------------------------------------------------
IF(PID.EQ.1.0)THEN
QK=QQ
TMT0K=TMT0
IF(TMT0K.LT.-30.)TMT0K=-30.
PRECRK=AMAX1(0.,PRECRL)
PRECSK=AMAX1(0.,PRECSL)
C----------------------------------------------------------------------
C INCREASE THE EVAPORATION/CONDENSATION FOR STRONG/LIGHT PREC
C----------------------------------------------------------------------
U00KLT=U00KL
AMAXRQ=AMAX1(0.,U00KL-RQ)
ERR=KE*AMAXRQ*PRECRK**0.5
C
IF(TMT0.GE.0.)THEN
ERS=0.
ELSE
ERS=(CRS1+CRS2*TMT0K)*AMAXRQ*PRECSK/U00KLT
ENDIF
C
IF(ERR+ERS.LE.1.E-20) GO TO 125
C---------------CORRECT IF OVER-EVAPO./COND. OCCURS--------------------
HHT=HH*TWODT
TTEMP=TT-RCP*(ELWV*ERR+ELIV*ERS)*HHT
QTEMP=QQ+HHT*(ERR+ERS)
TMT0T=(TTEMP-273.16)*HH
IF(TMT0T.LT.-30.)TMT0T=-30.
TMT15T=AMIN1(TMT0T,-15.)*HH
AI=0.008855
BI=1.
C
IF(TMT0T.LT.-20.)THEN
AI=0.007225
BI=0.9674
ENDIF
C
QWT=HH*PQ0/(PDSL*AETAL+PT)
1 *EXP(HH*A2*(TTEMP-A3)/(TTEMP-A4))
QIT=QWT*(BI+AI*AMIN1(TMT0T,0.))
QINTT=QWT*(1.-0.00032*TMT15T*(TMT15T+15.))
IF(TMT0T.LE.-40.)QINTT=QIT
FIW=FLOAT(IWL)
QCT=(1.-FIW)*QINTT+FIW*QIT
C
IF(QCT.LE.1.E-10) THEN
RQT=1.E-10
RQTT=1.E-10
ELSE
RQT=QTEMP/QCT
RQTT=QQ/QCT
ENDIF
C
IF(RQT.LE.U00KL) GO TO 125
C
ERK=(U00KL-RQTT)*QCT*RTWODT
RPRS=ERK/(PRECRK+PRECSK)
ERRT=PRECRK*RPRS
ERST=PRECSK*RPRS
ERR=AMAX1(0.,0.5*(ERR+ERRT))
ERS=AMAX1(0.,0.5*(ERS+ERST))
C
125 CONTINUE
C
PPR=-ERR*CONDE
PPS=-ERS*CONDE
C
IF(-PPR.GE.PRECRK)THEN
PPR=-PRECRK
ERR=-PPR*RCONDE
ENDIF
C
IF(-PPS.GE.PRECSK)THEN
PPS=-PRECSK
ERS=-PPS*RCONDE
ENDIF
C
ENDIF
C
IF(PID.EQ.1.)THEN
PRECRL=PRECRL+PPR*HH
PRECSL=PRECSL+PPS*HH
ENDIF
C-----------------------------------------------------------------------
C--------------------MELTING OF THE SNOW--------------------------------
C-----------------------------------------------------------------------
IF(PID.EQ.1.)THEN
CWMK=WW
AMAXPS=AMAX1(0.,PRECSL)
C
IF(TMT0.GT.0.)THEN
PSM1=CSM1*TMT0*TMT0*AMAXPS
PSM2=CWS*CR*CWMK*AMAXPS
PSM=PSM1+PSM2
ELSE
PSM1=0.
PSM2=0.
PSM=0.
ENDIF
C
PPR=PSM*CONDE
PPS=-PSM*CONDE
C
IF(-PPS.GE.AMAXPS)THEN
PPS=-AMAXPS
PPR=AMAXPS
PSM1=-PPS*RCONDE
PSM2=0.
PSM=PSM1
ENDIF
C
ENDIF
C
IF(PID.EQ.1.)THEN
PRECRL=PRECRL+PPR*HH
PRECSL=PRECSL+PPS*HH
ENDIF
C-----------------------------------------------------------------------
C---------------UPDATE T AND Q------------------------------------------
C-----------------------------------------------------------------------
IF(PID.EQ.1.)THEN
HHT=HH*TWODT
TT=-RCP*(ELWV*ERR+ELIV*ERS+ELIW*PSM1)
1 *HHT+TT
QQ=(ERR+ERS)*HHT+QQ
TTLAT=-RCP*(ELWV*ERR+ELIV*ERS+ELIW*PSM1
& -ELWV*(PRAUT+PRACW)
& -ELIV*(PSAUT+PSACI))*HHT+TTLAT
ENDIF
C
IF(HH.EQ.1.)THEN
IWL1(N)=IWL
PRECRL1(N)=PRECRL
PRECSL1(N)=PRECSL
ENDIF
C
C ACCUMULATE LATENT HEATING DUE TO GRID-SCALE PRECIP/EVAP.
C SCALE BY THE RECIPROCAL OF THE PERIOD AT WHICH THIS ROUTINE
C IS CALLED. THIS PERIOD IS THE PHYSICS TIMESTEP.
C
TRAIN_T(L,I,J)=TRAIN_T(L,I,J)+(TT-T_T(L,I,J))*RDTPH
T_T(L,I,J)=TT
Q_T(L,I,J)=QQ
CWM_T(L,I,J)=WW
TLAT_T(L,I,J)=TTLAT
180 CONTINUE
C-----------------------------------------------------------------------
C-------------------THE PRECIPITATION ON SFC----------------------------
C-----------------------------------------------------------------------
PRECRS=PRECRL1(N)*RROW
PRECSS=PRECSL1(N)*RROW
C
APREC(I,J)=PRECRS+PRECSS
PREC(I,J)=PREC(I,J)+PRECRS+PRECSS
ACPREC(I,J)=ACPREC(I,J)+APREC(I,J)
C-----------------------------------------------------------------------
C---------------THE SNOW AND RAIN RATIO OF SFC PREC---------------------
C----SR IS THE RATIO OF SNOW TO THE TOTAL PRECIP------------------------
C----IF TOTAL PRECIP IS ZERO, SR IS ZERO--------------------------------
C-----------------------------------------------------------------------
TOTPPT=PRECRS+PRECSS
IF (TOTPPT.GT.1.E-8) THEN
SR(I,J)=PRECSS/TOTPPT
ELSE
SR(I,J)=0.
ENDIF
C-----------------------------------------------------------------------
300 CONTINUE
C-----------------------------------------------------------------------
C***
C*** TRANSPOSE BACK
C***
!$omp parallel sections
!$omp section
CALL SGETMO(T_T,LM,LM,LDA,T,LDA)
!$omp section
CALL SGETMO(Q_T,LM,LM,LDA,Q,LDA)
!$omp section
CALL SGETMO(CWM_T,LM,LM,LDA,CWM,LDA)
!$omp section
CALL SGETMO(TRAIN_T,LM,LM,LDA,TRAIN,LDA)
!$omp section
CALL SGETMO(TLAT_T,LM,LM,LDA,TLAT,LDA)
!$omp end parallel sections
C-----------------------------------------------------------------------
RETURN
END