/[MITgcm]/MITgcm/pkg/fizhi/fizhi_moist.F
ViewVC logotype

Diff of /MITgcm/pkg/fizhi/fizhi_moist.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.3 by molod, Thu Jun 24 19:57:02 2004 UTC revision 1.6 by molod, Thu Jul 8 14:53:45 2004 UTC
# Line 2  C $Header$ Line 2  C $Header$
2  C $Name$  C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5        subroutine moistio (ndmoist,istrip,npcs,pz,tz,qz,bi,bj,        subroutine moistio (ndmoist,istrip,npcs,
6       .   ntracer,ptracer,       .   lowlevel,midlevel,nltop,nsubmin,nsubmax,Lup,
7       .   pkht,qqz,dumoist,dvmoist,dtmoist,dqmoist,       .   pz,plz,plze,dpres,pkht,pkl,tz,qz,bi,bj,ntracer,ptracer,
8       .   im,jm,lm,sige,sig,dsig,ptop,       .   qqz,dumoist,dvmoist,dtmoist,dqmoist,
9         .   im,jm,lm,ptop,
10       .   iras,rainlsp,rainconv,snowfall,       .   iras,rainlsp,rainconv,snowfall,
11       .   nswcld,cldtot_sw,cldras_sw,cldlsp_sw,nswlz,swlz,       .   nswcld,cldtot_sw,cldras_sw,cldlsp_sw,nswlz,swlz,
12       .   nlwcld,cldtot_lw,cldras_lw,cldlsp_lw,nlwlz,lwlz,       .   nlwcld,cldtot_lw,cldras_lw,cldlsp_lw,nlwlz,lwlz,
# Line 17  C $Name$ Line 18  C $Name$
18    
19  c Input Variables  c Input Variables
20  c ---------------  c ---------------
21        integer ndmoist,istrip,npcs,myid,bi,bj        integer ndmoist,istrip,npcs
22          integer lowlevel,midlevel,nltop,nsubmin,nsubmax,Lup
23        integer im,jm,lm                        real pz(im,jm),plz(im,jm,lm),plze(im,jm,lm+1),dpres(im,jm,lm)
24        real  ptop                              real pkht(im,jm,lm+1),pkl(im,jm,lm)
25        real  sige(lm+1)                        real tz(im,jm,lm),qz(im,jm,lm,ntracer)      
26        real   sig(lm)                          integer bi,bj,ntracer,ptracer        
27        real  dsig(lm)                          real qqz(im,jm,lm)
28          real dumoist(im,jm,lm),dvmoist(im,jm,lm)
29        integer ntracer,ptracer                real dtmoist(im,jm,lm),dqmoist(im,jm,lm,ntracer)
30          integer im,jm,lm
31        real pz(im,jm)                          real ptop
32        real tz(im,jm,lm)                      integer iras
33        real qz(im,jm,lm,ntracer)              real rainlsp(im,jm),rainconv(im,jm),snowfall(im,jm)
34          integer nswcld,nswlz
35        real  pkht(im,jm,lm)                    real cldlsp_sw(im,jm,lm),cldras_sw(im,jm,lm)
36          real cldtot_sw(im,jm,lm),swlz(im,jm,lm)
37        real   qqz(im,jm,lm)                    integer nlwcld,nlwlz
38          real  cldlsp_lw(im,jm,lm),cldras_lw(im,jm,lm)
39        real dumoist(im,jm,lm)                  real  cldtot_lw(im,jm,lm),lwlz(im,jm,lm)
40        real dvmoist(im,jm,lm)                  logical lpnt
41        real dtmoist(im,jm,lm)                  integer myid
       real dqmoist(im,jm,lm,ntracer)    
   
       integer iras                      
       real   rainlsp(im,jm)            
       real  rainconv(im,jm)            
       real  snowfall(im,jm)            
   
       integer nswcld,nswlz              
       real  cldlsp_sw(im,jm,lm)        
       real  cldras_sw(im,jm,lm)        
       real  cldtot_sw(im,jm,lm)        
       real       swlz(im,jm,lm)        
   
       integer nlwcld,nlwlz              
       real  cldlsp_lw(im,jm,lm)        
       real  cldras_lw(im,jm,lm)        
       real  cldtot_lw(im,jm,lm)        
       real       lwlz(im,jm,lm)        
   
       logical lpnt                      
42    
43  c Local Variables  c Local Variables
44  c ---------------  c ---------------
45        integer    ncrnd,nsecf,nsubmax        integer    ncrnd,nsecf
46    
47        real       fracqq, rh,temp1,temp2,dum        real       fracqq, rh,temp1,temp2,dum
48        integer    snowcrit, lup        integer    snowcrit
49        parameter (fracqq = 0.1)        parameter (fracqq = 0.1)
50    
51        real   cldsr(im,jm,lm)        real   cldsr(im,jm,lm)
# Line 75  c --------------- Line 56  c ---------------
56        real cldprs(im,jm),cldtmp(im,jm)        real cldprs(im,jm),cldtmp(im,jm)
57        real cldhi (im,jm),cldlow(im,jm)        real cldhi (im,jm),cldlow(im,jm)
58        real cldmid(im,jm),totcld(im,jm)        real cldmid(im,jm),totcld(im,jm)
       integer midlevel,lowlevel  
59    
60        real   CLDLS(im,jm,lm)  , CPEN(im,jm,lm)        real   CLDLS(im,jm,lm)  , CPEN(im,jm,lm)
61        real    tmpimjm(im,jm)        real    tmpimjm(im,jm)
# Line 94  c -------------------------------------- Line 74  c --------------------------------------
74        real     thgather(im*jm,lm)        real     thgather(im*jm,lm)
75        real     shgather(im*jm,lm)        real     shgather(im*jm,lm)
76        real    pkzgather(im*jm,lm)        real    pkzgather(im*jm,lm)
77        real    pkegather(im*jm,lm)        real    pkegather(im*jm,lm+1)
78          real    plzgather(im*jm,lm)
79          real    plegather(im*jm,lm+1)
80          real     dpgather(im*jm,lm)
81        real    tmpgather(im*jm,lm)        real    tmpgather(im*jm,lm)
82        real   deltgather(im*jm,lm)        real   deltgather(im*jm,lm)
83        real   delqgather(im*jm,lm)        real   delqgather(im*jm,lm)
# Line 115  c --------------- Line 98  c ---------------
98        real usubcl (istrip,   ntracer)        real usubcl (istrip,   ntracer)
99    
100        real     ple(istrip,lm+1), gam(istrip,lm)        real     ple(istrip,lm+1), gam(istrip,lm)
101          real      dp(istrip,lm)
102        real      TL(ISTRIP,lm)  , SHL(ISTRIP,lm)        real      TL(ISTRIP,lm)  , SHL(ISTRIP,lm)
103        real      PL(ISTRIP,lm)  , PLK(ISTRIP,lm)        real      PL(ISTRIP,lm)  , PLK(ISTRIP,lm)
104        real    PLKE(ISTRIP,lm+1)        real    PLKE(ISTRIP,lm+1)
# Line 149  c --------------- Line 133  c ---------------
133        real rnd(lm/2)        real rnd(lm/2)
134        DATA      FIRST /.TRUE./        DATA      FIRST /.TRUE./
135    
136        integer imstp,nltop,nsubcl,nlras,nsubmin        integer imstp,nsubcl,nlras
137        integer i,j,iloop,index,l,nn,num,numdeps,nt        integer i,j,iloop,index,l,nn,num,numdeps,nt
138        real tmstp,tminv,sday,grav,alhl,cp,elocp,gamfac        real tmstp,tminv,sday,grav,alhl,cp,elocp,gamfac
139        real rkappa,p0kappa,p0kinv,ptopkap,pcheck        real rkappa,p0kappa,p0kinv,ptopkap,pcheck
# Line 192  C Threshold for Cloud Liquid Water Memor Line 176  C Threshold for Cloud Liquid Water Memor
176        tice     = getcon('FREEZING-POINT')        tice     = getcon('FREEZING-POINT')
177        PI       = 4.*atan(1.)        PI       = 4.*atan(1.)
178    
179  c Determine Upper  Level for Cumulus Convection  c Determine Total number of Random Clouds to Check
 c and Total number of Random Clouds to Check  
180  c ---------------------------------------------  c ---------------------------------------------
   
       NLTOP = 1  
       DO L=1,lm  
        PCHECK = (1000.-ptop)*SIG(L) + PTOP  
        IF (PCHECK.GE.10.) THEN  
         NLTOP = L  
         GO TO 2  
        ENDIF  
       ENDDO  
  2    CONTINUE  
181        ncrnd = (lm-nltop+1)/2        ncrnd = (lm-nltop+1)/2
182    
 c Determine Minimum Number of Levels in Sub-Cloud (50 mb) Layer  
 c -------------------------------------------------------------  
       nsubmin = lm  
       nsubmax = 1  
       DO L=lm-1,1,-1  
           PCHECK = (1000.-ptop)*SIG(L) + PTOP  
       IF( PCHECK.GE.950.0 ) nsubmin = L  
       IF( PCHECK.GE.750.0 ) nsubmax = L  
       ENDDO  
   
183        if(first .and. myid.eq.0) then        if(first .and. myid.eq.0) then
184         print *         print *
185         print *,'Top Level Allowed for Convection : ',nltop,         print *,'Top Level Allowed for Convection : ',nltop
186       .                    ' (',(1000.-ptop)*SIG(nltop) + PTOP,' mb)'         print *,'          Highest Sub-Cloud Level: ',nsubmax
187         print *,'          Highest Sub-Cloud Level: ',nsubmax,         print *,'          Lowest  Sub-Cloud Level: ',nsubmin
      .                    ' (',(1000.-ptop)*SIG(nsubmax) + PTOP,' mb)'  
        print *,'          Lowest  Sub-Cloud Level: ',nsubmin,  
      .                    ' (',(1000.-ptop)*SIG(nsubmin) + PTOP,' mb)'  
188         print *,'    Total Number of Random Clouds: ',ncrnd         print *,'    Total Number of Random Clouds: ',ncrnd
189         print *         print *
190         first = .false.         first = .false.
# Line 279  c -------------------------------------- Line 239  c --------------------------------------
239        do index = 1,im*jm        do index = 1,im*jm
240         levgather(index) = levpbl(pblindex(index),1)         levgather(index) = levpbl(pblindex(index),1)
241          pigather(index) =     pz(pblindex(index),1)          pigather(index) =     pz(pblindex(index),1)
242            pkegather(index,lm+1) = pkht(pblindex(index),1,lm+1)
243            plegather(index,lm+1) = ple(pblindex(index),1,lm+1)
244        enddo        enddo
245    
246        do L = 1,lm        do L = 1,lm
247         do index = 1,im*jm         do index = 1,im*jm
248           thgather(index,L) =   tz(pblindex(index),1,L)           thgather(index,L) = tz(pblindex(index),1,L)
249           shgather(index,L) =   qz(pblindex(index),1,L,1)           shgather(index,L) = qz(pblindex(index),1,L,1)
250          pkegather(index,L) = pkht(pblindex(index),1,L)          pkegather(index,L) = pkht(pblindex(index),1,L)
251            pkzgather(index,L) = pkl(pblindex(index),1,L)
252            plegather(index,L) = plze(pblindex(index),1,L)
253            plzgather(index,L) = plz(pblindex(index),1,L)
254             dpgather(index,L) = dpres(pblindex(index),1,L)
255         enddo         enddo
256        enddo        enddo
257        do nt = 1,ntracer-ptracer        do nt = 1,ntracer-ptracer
# Line 296  c -------------------------------------- Line 262  c --------------------------------------
262        enddo        enddo
263        enddo        enddo
264    
       call pkappa(pigather,pkegather,pkzgather,ptop,sige,dsig,im,jm,lm)  
   
265  c bump the counter for number of calls to convection  c bump the counter for number of calls to convection
266  c --------------------------------------------------  c --------------------------------------------------
267                          iras = iras + 1                          iras = iras + 1
# Line 330  C ************************************** Line 294  C **************************************
294    
295         CALL STRIP (  pigather, SP      ,im*jm,ISTRIP,1 ,NN )         CALL STRIP (  pigather, SP      ,im*jm,ISTRIP,1 ,NN )
296         CALL STRIP ( pkzgather, PLK     ,im*jm,ISTRIP,lm,NN )         CALL STRIP ( pkzgather, PLK     ,im*jm,ISTRIP,lm,NN )
297         CALL STRIP ( pkegather, PLKE    ,im*jm,ISTRIP,lm,NN )         CALL STRIP ( pkegather, PLKE    ,im*jm,ISTRIP,lm+1,NN )
298           CALL STRIP ( plzgather, PL      ,im*jm,ISTRIP,lm,NN )
299           CALL STRIP ( plegather, PLE     ,im*jm,ISTRIP,lm+1,NN )
300           CALL STRIP (  dpgather, dp      ,im*jm,ISTRIP,lm,NN )
301         CALL STRIP (  thgather, TH      ,im*jm,ISTRIP,lm,NN )         CALL STRIP (  thgather, TH      ,im*jm,ISTRIP,lm,NN )
302         CALL STRIP (  shgather, SHL     ,im*jm,ISTRIP,lm,NN )         CALL STRIP (  shgather, SHL     ,im*jm,ISTRIP,lm,NN )
303         CALL STRINT( levgather, pbl     ,im*jm,ISTRIP,1 ,NN )         CALL STRINT( levgather, pbl     ,im*jm,ISTRIP,1 ,NN )
# Line 339  C ************************************** Line 306  C **************************************
306         call strip ( ugather(1,1,nt), ul(1,1,nt),im*jm,istrip,lm,nn )         call strip ( ugather(1,1,nt), ul(1,1,nt),im*jm,istrip,lm,nn )
307         enddo         enddo
308    
       do l = 1,lm  
       do i = 1,istrip  
        PL(I,L) =  SIG(L)*SP(I) + PTOP  
       PLE(I,L) = SIGE(L)*SP(I) + PTOP  
       enddo  
       enddo  
   
       do i = 1,istrip  
       PLE(I,lm+1) = SP(I) + PTOP  
       enddo  
   
309  C **********************************************************************  C **********************************************************************
310  C ****        SETUP FOR RAS CUMULUS PARAMETERIZATION                ****  C ****        SETUP FOR RAS CUMULUS PARAMETERIZATION                ****
311  C **********************************************************************  C **********************************************************************
# Line 445  c ------------------------------- Line 401  c -------------------------------
401         ENDDO         ENDDO
402         DO L=2,lm         DO L=2,lm
403         DO I=num,num+nindeces(nsubcl)-1         DO I=num,num+nindeces(nsubcl)-1
404          TMP5(I,L) = PLKE(I,L-1)*P0KINV          TMP5(I,L) = PLKE(I,L)*P0KINV
405         ENDDO         ENDDO
406         ENDDO         ENDDO
407         DO  I=num,num+nindeces(nsubcl)-1         DO  I=num,num+nindeces(nsubcl)-1
408          TMP4(I,lm+1) = PLE (I,lm+1)          TMP4(I,lm+1) = PLE (I,lm+1)
409          TMP5(I,lm+1) = PLKE(I,lm)*P0KINV          TMP5(I,lm+1) = PLKE(I,lm+1)*P0KINV
410         ENDDO         ENDDO
411         DO 113 I=num,num+nindeces(nsubcl)-1         DO 113 I=num,num+nindeces(nsubcl)-1
412          TMP4(I,NSUBCL+1) = PLE (I,lm+1)          TMP4(I,NSUBCL+1) = PLE (I,lm+1)
413          TMP5(I,NSUBCL+1) = PLKE(I,lm)*P0KINV          TMP5(I,NSUBCL+1) = PLKE(I,lm+1)*P0KINV
414   113   CONTINUE   113   CONTINUE
415    
416        do i=num,num+nindeces(nsubcl)-1        do i=num,num+nindeces(nsubcl)-1
# Line 491  C Top level of atan func above this rh_t Line 447  C Top level of atan func above this rh_t
447           rhcrit(i,L) = 1.           rhcrit(i,L) = 1.
448         enddo         enddo
449         do L = 1, nsubcl-1         do L = 1, nsubcl-1
450          pcheck = (1000.-ptop)*sig(L) + ptop          pcheck = pl(i,L)
451          if (pcheck .le. pup) then          if (pcheck .le. pup) then
452           rhcrit(i,L) = rhmin           rhcrit(i,L) = rhmin
453          else          else
454           ppbl = (1000.-ptop)*sig(nsubcl) + ptop           ppbl = pl(i,nsubcl)
455           rhcrit(i,L) = rhmin + (1.-rhmin)/(19.) *           rhcrit(i,L) = rhmin + (1.-rhmin)/(19.) *
456       .      ((atan( (2.*(pcheck-pup)/(ppbl-pup)-1.) *       .      ((atan( (2.*(pcheck-pup)/(ppbl-pup)-1.) *
457       .       tan(20.*pi/21.-0.5*pi) )       .       tan(20.*pi/21.-0.5*pi) )
# Line 525  c -------------------------------------- Line 481  c --------------------------------------
481  c Compute Diagnostic CLDMAS in RAS Subcloud Layers  c Compute Diagnostic CLDMAS in RAS Subcloud Layers
482  c ------------------------------------------------  c ------------------------------------------------
483         do L=nsubcl,lm         do L=nsubcl,lm
         dum = dsig(L)/(1.0-sige(nsubcl))  
484         do I=num,num+nindeces(nsubcl)-1         do I=num,num+nindeces(nsubcl)-1
485            dum = dp(i,L)/(ple(i,lm+1)-ple(i,nsubcl))
486          cldmas(i,L) = cldmas(i,L-1) - dum*cldmas(i,nsubcl-1)          cldmas(i,L) = cldmas(i,L-1) - dum*cldmas(i,nsubcl-1)
487         enddo         enddo
488         enddo         enddo
# Line 640  C ************************************** Line 596  C **************************************
596        ENDDO        ENDDO
597        ENDDO        ENDDO
598    
599         CALL RNEVP (NN,ISTRIP,lm,TL,SHL,PCPEN,PL,CLFRAC,SP,DSIG,PLKE,         CALL RNEVP (NN,ISTRIP,lm,TL,SHL,PCPEN,PL,CLFRAC,SP,DP,PLKE,
600       .  PLK,TH,TMP1,TMP2,TMP3,ITMP1,ITMP2,PCNET,PRECIP,       .  PLK,TH,TMP1,TMP2,TMP3,ITMP1,ITMP2,PCNET,PRECIP,
601       .  CLSBTH,TMSTP,1.,cp,grav,alhl,gamfac,cldlz,rhcrit,offset,alpha)       .  CLSBTH,TMSTP,1.,cp,grav,alhl,gamfac,cldlz,rhcrit,offset,alpha)
602    
# Line 729  c  snow algorthm: Line 685  c  snow algorthm:
685  c  if temperature profile from the surface level to 700 mb  c  if temperature profile from the surface level to 700 mb
686  c  uniformaly c  below zero, then precipitation (total) is  c  uniformaly c  below zero, then precipitation (total) is
687  c  snowfall.  Else there is no snow.  c  snowfall.  Else there is no snow.
 c  For version of level 70, the sigma level corresponding  
 c  to 700mb (assume the surface pressure is 1000mb) is  
 c  the 13th level from the surface  
 c   Runhua Yang Aug. 24 98  
 c  added pcheck for 700mb - sharon sept 18, 1998  
688  c-------------------------------------------------------  c-------------------------------------------------------
689    
         pup = 700.  
         do L = lm, 1, -1  
           pcheck = (1000.-ptop)*sig(L) + ptop  
           if (pcheck .ge. pup) then  
             lup = L  
           endif  
         enddo  
690          do i = 1,istrip          do i = 1,istrip
691            snowcrit=0            snowcrit=0
692            do l=lup,lm            do l=lup,lm
# Line 826  C ************************************** Line 770  C **************************************
770  C                          BUMP DIAGNOSTICS  C                          BUMP DIAGNOSTICS
771  C **********************************************************************  C **********************************************************************
772    
 c Determine Level Indices for Low-Mid-High Cloud Regions  
 c ------------------------------------------------------  
       lowlevel = lm  
       midlevel = lm  
       do L = lm-1,1,-1  
       pcheck = (1000.-ptop)*sig(l) + ptop  
       if (pcheck.gt.700.0) lowlevel = L  
       if (pcheck.gt.400.0) midlevel = L  
       enddo  
   
   
773  c Clear-Sky (Above 400mb) Temperature  c Clear-Sky (Above 400mb) Temperature
774  c -----------------------------------  c -----------------------------------
775        if( itmpuclr.ne.0 .or. isphuclr.ne.0 ) then        if( itmpuclr.ne.0 .or. isphuclr.ne.0 ) then
# Line 1031  c ------------------------------------ Line 964  c ------------------------------------
964    
965        do L = 1,lm        do L = 1,lm
966        do i = 1,im*jm        do i = 1,im*jm
967         plev = sig(L)*pz(i,1)+ptop         plev = pl(i,L)
968    
969  c Compute Time-averaged Cloud and Water Amounts for Longwave Radiation  c Compute Time-averaged Cloud and Water Amounts for Longwave Radiation
970  c --------------------------------------------------------------------  c --------------------------------------------------------------------
# Line 2263  C*************************************** Line 2196  C***************************************
2196    
2197        RETURN        RETURN
2198        END        END
2199         SUBROUTINE RNEVP(NN,IRUN,NLAY,TL,QL,RAIN,PL,CLFRAC,SP,DSIG,PLKE,         SUBROUTINE RNEVP(NN,IRUN,NLAY,TL,QL,RAIN,PL,CLFRAC,SP,DP,PLKE,
2200       1   PLK,TH,TEMP1,TEMP2,TEMP3,ITMP1,ITMP2,RCON,RLAR,CLSBTH,tmscl,       1   PLK,TH,TEMP1,TEMP2,TEMP3,ITMP1,ITMP2,RCON,RLAR,CLSBTH,tmscl,
2201       2   tmfrc,cp,gravity,alhl,gamfac,cldlz,RHCRIT,offset,alpha)       2   tmfrc,cp,gravity,alhl,gamfac,cldlz,RHCRIT,offset,alpha)
2202    
# Line 2282  C*************************************** Line 2215  C***************************************
2215  C  C
2216         DIMENSION TL(IRUN,NLAY),QL(IRUN,NLAY),RAIN(IRUN,NLAY),         DIMENSION TL(IRUN,NLAY),QL(IRUN,NLAY),RAIN(IRUN,NLAY),
2217       $ PL(IRUN,NLAY),CLFRAC(IRUN,NLAY),SP(IRUN),TEMP1(IRUN,NLAY),       $ PL(IRUN,NLAY),CLFRAC(IRUN,NLAY),SP(IRUN),TEMP1(IRUN,NLAY),
2218       $ TEMP2(IRUN,NLAY),PLKE(IRUN,NLAY),       $ TEMP2(IRUN,NLAY),PLKE(IRUN,NLAY+1),
2219       $ RCON(IRUN),RLAR(IRUN),DSIG(NLAY),PLK(IRUN,NLAY),TH(IRUN,NLAY),       $ RCON(IRUN),RLAR(IRUN),DP(IRUN,NLAY),PLK(IRUN,NLAY),TH(IRUN,NLAY),
2220       $ TEMP3(IRUN,NLAY),ITMP1(IRUN,NLAY),       $ TEMP3(IRUN,NLAY),ITMP1(IRUN,NLAY),
2221       $ ITMP2(IRUN,NLAY),CLSBTH(IRUN,NLAY)       $ ITMP2(IRUN,NLAY),CLSBTH(IRUN,NLAY)
2222  C  C
# Line 2348  C INVERSE OF MASS IN EACH LAYER Line 2281  C INVERSE OF MASS IN EACH LAYER
2281  c -----------------------------  c -----------------------------
2282        DO L = 1,NLAY        DO L = 1,NLAY
2283        DO I = 1,IRUN        DO I = 1,IRUN
2284        TEMP3(I,L) = SP(I) * DSIG(L)        TEMP3(I,L) = GRAVITY*ZP01 / DP(I,L)
       TEMP3(I,L) = GRAVITY*ZP01 / TEMP3(I,L)  
2285        ENDDO        ENDDO
2286        ENDDO        ENDDO
2287    
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.6
  ViewVC Help
Powered by ViewVC 1.1.22