/[lmdze]/trunk/Sources/phylmd/Conflx/flxasc.f
ViewVC logotype

Diff of /trunk/Sources/phylmd/Conflx/flxasc.f

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

trunk/libf/phylmd/Conflx/flxasc.f90 revision 62 by guez, Thu Jul 26 14:37:37 2012 UTC trunk/phylmd/Conflx/flxasc.f revision 82 by guez, Wed Mar 5 14:57:53 2014 UTC
# Line 1  Line 1 
1  SUBROUTINE flxasc(pdtime, ptenh, pqenh, pten, pqen, pqsen, &  module flxasc_m
2       pgeo, pgeoh, pap, paph, pqte, pvervel, &  
      ldland, ldcum, ktype, klab, ptu, pqu, plu, &  
      pmfu, pmfub, pentr, pmfus, pmfuq, &  
      pmful, plude, pdmfup, kcbot, kctop, kctop0, kcum, &  
      pen_u, pde_u)  
   use dimens_m  
   use dimphy  
   use SUPHEC_M  
   use yoethf_m  
   use yoecumf  
3    IMPLICIT none    IMPLICIT none
   !----------------------------------------------------------------------  
   ! THIS ROUTINE DOES THE CALCULATIONS FOR CLOUD ASCENTS  
   ! FOR CUMULUS PARAMETERIZATION  
   !----------------------------------------------------------------------  
   !  
   REAL, intent(in):: pdtime  
   REAL pten(klon,klev), ptenh(klon,klev)  
   REAL pqen(klon,klev), pqenh(klon,klev), pqsen(klon,klev)  
   REAL pgeo(klon,klev), pgeoh(klon,klev)  
   REAL pap(klon,klev), paph(klon,klev+1)  
   REAL pqte(klon,klev)  
   REAL pvervel(klon,klev) ! vitesse verticale en Pa/s  
   !  
   REAL pmfub(klon), pentr(klon)  
   REAL ptu(klon,klev), pqu(klon,klev), plu(klon,klev)  
   REAL plude(klon,klev)  
   REAL pmfu(klon,klev), pmfus(klon,klev)  
   REAL pmfuq(klon,klev), pmful(klon,klev)  
   REAL pdmfup(klon,klev)  
   INTEGER ktype(klon), klab(klon,klev), kcbot(klon), kctop(klon)  
   INTEGER kctop0(klon)  
   LOGICAL ldland(klon), ldcum(klon)  
   !  
   REAL pen_u(klon,klev), pde_u(klon,klev)  
   REAL zqold(klon)  
   REAL zdland(klon)  
   LOGICAL llflag(klon)  
   INTEGER k, i, is, icall, kcum  
   REAL ztglace, zdphi, zqeen, zseen, zscde, zqude  
   REAL zmfusk, zmfuqk, zmfulk, zbuo, zdnoprc, zprcon, zlnew  
   !  
   REAL zpbot(klon), zptop(klon), zrho(klon)  
   REAL zdprho, zentr, zpmid, zmftest, zmfmax  
   LOGICAL llo1, llo2  
   !  
   REAL zwmax(klon), zzzmb  
   INTEGER klwmin(klon) ! level of maximum vertical velocity  
   !----------------------------------------------------------------------  
   ztglace = RTT - 13.  
   !  
   ! Chercher le niveau ou la vitesse verticale est maximale:  
   DO i = 1, klon  
      klwmin(i) = klev  
      zwmax(i) = 0.0  
   ENDDO  
   DO k = klev, 3, -1  
      DO i = 1, klon  
         IF (pvervel(i,k).LT.zwmax(i)) THEN  
            zwmax(i) = pvervel(i,k)  
            klwmin(i) = k  
         ENDIF  
      ENDDO  
   ENDDO  
   !----------------------------------------------------------------------  
   ! SET DEFAULT VALUES  
   !----------------------------------------------------------------------  
   DO i = 1, klon  
      IF (.NOT.ldcum(i)) ktype(i)=0  
   ENDDO  
   !  
   DO k=1,klev  
      DO i = 1, klon  
         plu(i,k)=0.  
         pmfu(i,k)=0.  
         pmfus(i,k)=0.  
         pmfuq(i,k)=0.  
         pmful(i,k)=0.  
         plude(i,k)=0.  
         pdmfup(i,k)=0.  
         IF(.NOT.ldcum(i).OR.ktype(i).EQ.3) klab(i,k)=0  
         IF(.NOT.ldcum(i).AND.paph(i,k).LT.4.E4) kctop0(i)=k  
      ENDDO  
   ENDDO  
   !  
   DO i = 1, klon  
      IF (ldland(i)) THEN  
         zdland(i)=3.0E4  
         zdphi=pgeoh(i,kctop0(i))-pgeoh(i,kcbot(i))  
         IF (ptu(i,kctop0(i)).GE.ztglace) zdland(i)=zdphi  
         zdland(i)=MAX(3.0E4,zdland(i))  
         zdland(i)=MIN(5.0E4,zdland(i))  
      ENDIF  
   ENDDO  
   !  
   ! Initialiser les valeurs au niveau d'ascendance  
   !  
   DO i = 1, klon  
      kctop(i) = klev-1  
      IF (.NOT.ldcum(i)) THEN  
         kcbot(i) = klev-1  
         pmfub(i) = 0.  
         pqu(i,klev) = 0.  
      ENDIF  
      pmfu(i,klev) = pmfub(i)  
      pmfus(i,klev) = pmfub(i)*(RCPD*ptu(i,klev)+pgeoh(i,klev))  
      pmfuq(i,klev) = pmfub(i)*pqu(i,klev)  
   ENDDO  
   !  
   DO i = 1, klon  
      ldcum(i) = .FALSE.  
   ENDDO  
   !----------------------------------------------------------------------  
   !  DO ASCENT: SUBCLOUD LAYER (klab=1) ,CLOUDS (klab=2)  
   !  BY DOING FIRST DRY-ADIABATIC ASCENT AND THEN  
   !  BY ADJUSTING T,Q AND L ACCORDINGLY IN *flxadjtq*,  
   !  THEN CHECK FOR BUOYANCY AND SET FLAGS ACCORDINGLY  
   !----------------------------------------------------------------------  
   DO  k = klev-1,3,-1  
      !  
      IF (LMFMID .AND. k.LT.klev-1 .AND. k.GT.klev/2) THEN  
         DO i = 1, klon  
            IF (.NOT.ldcum(i) .AND. klab(i,k+1).EQ.0 .AND. &  
                 pqen(i,k).GT.0.9*pqsen(i,k)) THEN  
               ptu(i,k+1) = pten(i,k) +(pgeo(i,k)-pgeoh(i,k+1))/RCPD  
               pqu(i,k+1) = pqen(i,k)  
               plu(i,k+1) = 0.0  
               zzzmb = MAX(CMFCMIN, -pvervel(i,k)/RG)  
               zmfmax = (paph(i,k)-paph(i,k-1))/(RG*pdtime)  
               pmfub(i) = MIN(zzzmb,zmfmax)  
               pmfu(i,k+1) = pmfub(i)  
               pmfus(i,k+1) = pmfub(i)*(RCPD*ptu(i,k+1)+pgeoh(i,k+1))  
               pmfuq(i,k+1) = pmfub(i)*pqu(i,k+1)  
               pmful(i,k+1) = 0.0  
               pdmfup(i,k+1) = 0.0  
               kcbot(i) = k  
               klab(i,k+1) = 1  
               ktype(i) = 3  
               pentr(i) = ENTRMID  
            ENDIF  
         ENDDO  
      ENDIF  
      !  
      is = 0  
      DO i = 1, klon  
         is = is + klab(i,k+1)  
         IF (klab(i,k+1) .EQ. 0) klab(i,k) = 0  
         llflag(i) = .FALSE.  
         IF (klab(i,k+1) .GT. 0) llflag(i) = .TRUE.  
      ENDDO  
      IF (is .EQ. 0) cycle  
      !  
      ! calculer le taux d'entrainement et de detrainement  
      !  
      DO i = 1, klon  
         pen_u(i,k) = 0.0  
         pde_u(i,k) = 0.0  
         zrho(i)=paph(i,k+1)/(RD*ptenh(i,k+1))  
         zpbot(i)=paph(i,kcbot(i))  
         zptop(i)=paph(i,kctop0(i))  
      ENDDO  
      !  
      DO i = 1, klon  
         IF(ldcum(i)) THEN  
            zdprho=(paph(i,k+1)-paph(i,k))/(RG*zrho(i))  
            zentr=pentr(i)*pmfu(i,k+1)*zdprho  
            llo1=k.LT.kcbot(i)  
            IF(llo1) pde_u(i,k)=zentr  
            zpmid=0.5*(zpbot(i)+zptop(i))  
            llo2=llo1.AND.ktype(i).EQ.2.AND. &  
                 (zpbot(i)-paph(i,k).LT.0.2E5.OR. &  
                 paph(i,k).GT.zpmid)  
            IF(llo2) pen_u(i,k)=zentr  
            llo2=llo1.AND.(ktype(i).EQ.1.OR.ktype(i).EQ.3).AND. &  
                 (k.GE.MAX(klwmin(i),kctop0(i)+2).OR.pap(i,k).GT.zpmid)  
            IF(llo2) pen_u(i,k)=zentr  
            llo1=pen_u(i,k).GT.0..AND.(ktype(i).EQ.1.OR.ktype(i).EQ.2)  
            IF(llo1) THEN  
               zentr=zentr*(1.+3.*(1.-MIN(1.,(zpbot(i)-pap(i,k))/1.5E4)))  
               pen_u(i,k)=pen_u(i,k)*(1.+3.*(1.-MIN(1., &  
                    (zpbot(i)-pap(i,k))/1.5E4)))  
               pde_u(i,k)=pde_u(i,k)*(1.+3.*(1.-MIN(1., &  
                    (zpbot(i)-pap(i,k))/1.5E4)))  
            ENDIF  
            IF(llo2.AND.pqenh(i,k+1).GT.1.E-5) &  
                 pen_u(i,k)=zentr+MAX(pqte(i,k),0.)/pqenh(i,k+1)* &  
                 zrho(i)*zdprho  
         ENDIF  
      end DO  
      !  
      !----------------------------------------------------------------------  
      ! DO ADIABATIC ASCENT FOR ENTRAINING/DETRAINING PLUME  
      !----------------------------------------------------------------------  
      !  
      DO  i = 1, klon  
         IF (llflag(i)) THEN  
            IF (k.LT.kcbot(i)) THEN  
               zmftest = pmfu(i,k+1)+pen_u(i,k)-pde_u(i,k)  
               zmfmax = MIN(zmftest,(paph(i,k)-paph(i,k-1))/(RG*pdtime))  
               pen_u(i,k)=MAX(pen_u(i,k)-MAX(0.0,zmftest-zmfmax),0.0)  
            ENDIF  
            pde_u(i,k)=MIN(pde_u(i,k),0.75*pmfu(i,k+1))  
            ! calculer le flux de masse du niveau k a partir de celui du k+1  
            pmfu(i,k)=pmfu(i,k+1)+pen_u(i,k)-pde_u(i,k)  
            ! calculer les valeurs Su, Qu et l du niveau k dans le panache montant  
            zqeen=pqenh(i,k+1)*pen_u(i,k)  
            zseen=(RCPD*ptenh(i,k+1)+pgeoh(i,k+1))*pen_u(i,k)  
            zscde=(RCPD*ptu(i,k+1)+pgeoh(i,k+1))*pde_u(i,k)  
            zqude=pqu(i,k+1)*pde_u(i,k)  
            plude(i,k)=plu(i,k+1)*pde_u(i,k)  
            zmfusk=pmfus(i,k+1)+zseen-zscde  
            zmfuqk=pmfuq(i,k+1)+zqeen-zqude  
            zmfulk=pmful(i,k+1)    -plude(i,k)  
            plu(i,k)=zmfulk*(1./MAX(CMFCMIN,pmfu(i,k)))  
            pqu(i,k)=zmfuqk*(1./MAX(CMFCMIN,pmfu(i,k)))  
            ptu(i,k)=(zmfusk*(1./MAX(CMFCMIN,pmfu(i,k)))- &  
                 pgeoh(i,k))/RCPD  
            ptu(i,k)=MAX(100.,ptu(i,k))  
            ptu(i,k)=MIN(400.,ptu(i,k))  
            zqold(i)=pqu(i,k)  
         ELSE  
            zqold(i)=0.0  
         ENDIF  
      end DO  
      !  
      !----------------------------------------------------------------------  
      ! DO CORRECTIONS FOR MOIST ASCENT BY ADJUSTING T,Q AND L  
      !----------------------------------------------------------------------  
      !  
      icall = 1  
      CALL flxadjtq(paph(1,k), ptu(1,k), pqu(1,k), llflag, icall)  
      !  
      DO i = 1, klon  
         IF(llflag(i).AND.pqu(i,k).NE.zqold(i)) THEN  
            klab(i,k) = 2  
            plu(i,k) = plu(i,k)+zqold(i)-pqu(i,k)  
            zbuo = ptu(i,k)*(1.+RETV*pqu(i,k))- &  
                 ptenh(i,k)*(1.+RETV*pqenh(i,k))  
            IF (klab(i,k+1).EQ.1) zbuo=zbuo+0.5  
            IF (zbuo.GT.0..AND.pmfu(i,k).GE.0.1*pmfub(i)) THEN  
               kctop(i) = k  
               ldcum(i) = .TRUE.  
               zdnoprc = 1.5E4  
               IF (ldland(i)) zdnoprc = zdland(i)  
               zprcon = CPRCON  
               IF ((zpbot(i)-paph(i,k)).LT.zdnoprc) zprcon = 0.0  
               zlnew=plu(i,k)/(1.+zprcon*(pgeoh(i,k)-pgeoh(i,k+1)))  
               pdmfup(i,k)=MAX(0.,(plu(i,k)-zlnew)*pmfu(i,k))  
               plu(i,k)=zlnew  
            ELSE  
               klab(i,k)=0  
               pmfu(i,k)=0.  
            ENDIF  
         ENDIF  
      end DO  
      DO  i = 1, klon  
         IF (llflag(i)) THEN  
            pmful(i,k)=plu(i,k)*pmfu(i,k)  
            pmfus(i,k)=(RCPD*ptu(i,k)+pgeoh(i,k))*pmfu(i,k)  
            pmfuq(i,k)=pqu(i,k)*pmfu(i,k)  
         ENDIF  
      end DO  
      !  
   end DO  
   !----------------------------------------------------------------------  
   ! DETERMINE CONVECTIVE FLUXES ABOVE NON-BUOYANCY LEVEL  
   !    (NOTE: CLOUD VARIABLES LIKE T,Q AND L ARE NOT  
   !           AFFECTED BY DETRAINMENT AND ARE ALREADY KNOWN  
   !           FROM PREVIOUS CALCULATIONS ABOVE)  
   !----------------------------------------------------------------------  
   DO i = 1, klon  
      IF (kctop(i).EQ.klev-1) ldcum(i) = .FALSE.  
      kcbot(i) = MAX(kcbot(i),kctop(i))  
   ENDDO  
   !  
   ldcum(1)=ldcum(1)  
   !  
   is = 0  
   DO i = 1, klon  
      if (ldcum(i)) is = is + 1  
   ENDDO  
   kcum = is  
   IF (is /= 0) then  
      !  
      DO  i = 1, klon  
         IF (ldcum(i)) THEN  
            k=kctop(i)-1  
            pde_u(i,k)=(1.-CMFCTOP)*pmfu(i,k+1)  
            plude(i,k)=pde_u(i,k)*plu(i,k+1)  
            pmfu(i,k)=pmfu(i,k+1)-pde_u(i,k)  
            zlnew=plu(i,k)  
            pdmfup(i,k)=MAX(0.,(plu(i,k)-zlnew)*pmfu(i,k))  
            plu(i,k)=zlnew  
            pmfus(i,k)=(RCPD*ptu(i,k)+pgeoh(i,k))*pmfu(i,k)  
            pmfuq(i,k)=pqu(i,k)*pmfu(i,k)  
            pmful(i,k)=plu(i,k)*pmfu(i,k)  
            plude(i,k-1)=pmful(i,k)  
         ENDIF  
      end DO  
      !  
   end IF  
4    
5  END SUBROUTINE flxasc  contains
6    
7      SUBROUTINE flxasc(pdtime, ptenh, pqenh, pten, pqen, pqsen, pgeo, pgeoh, &
8           pap, paph, pqte, pvervel, ldland, ldcum, ktype, klab, ptu, pqu, plu, &
9           pmfu, pmfub, pentr, pmfus, pmfuq, pmful, plude, pdmfup, kcbot, kctop, &
10           kctop0, kcum, pen_u, pde_u)
11    
12        ! This routine does the calculations for cloud ascents for cumulus
13        ! parameterization.
14    
15        USE dimphy, ONLY: klev, klon
16        use flxadjtq_m, only: flxadjtq
17        USE suphec_m, ONLY: rcpd, rd, retv, rg, rtt
18        USE yoecumf, ONLY: cmfcmin, cmfctop, cprcon, entrmid, lmfmid
19    
20        REAL, intent(in):: pdtime
21        REAL, intent(in):: ptenh(klon, klev)
22        REAL, intent(in):: pqenh(klon, klev)
23        REAL, intent(in):: pten(klon, klev)
24        REAL, intent(in):: pqen(klon, klev)
25        REAL, intent(in):: pqsen(klon, klev)
26        REAL, intent(in):: pgeo(klon, klev), pgeoh(klon, klev)
27        REAL, intent(in):: pap(klon, klev), paph(klon, klev+1)
28        REAL, intent(in):: pqte(klon, klev)
29        REAL, intent(in):: pvervel(klon, klev) ! vitesse verticale en Pa/s
30        LOGICAL, intent(in):: ldland(klon)
31        LOGICAL, intent(inout):: ldcum(klon)
32        INTEGER, intent(inout):: ktype(klon)
33        integer klab(klon, klev)
34        REAL ptu(klon, klev), pqu(klon, klev), plu(klon, klev)
35        REAL pmfu(klon, klev)
36        REAL, intent(inout):: pmfub(klon)
37        real pentr(klon)
38        real pmfus(klon, klev)
39        REAL pmfuq(klon, klev), pmful(klon, klev)
40        REAL plude(klon, klev)
41        REAL pdmfup(klon, klev)
42        integer kcbot(klon), kctop(klon)
43        INTEGER kctop0(klon)
44        integer, intent(out):: kcum
45        REAL pen_u(klon, klev), pde_u(klon, klev)
46    
47        ! Local:
48    
49        REAL zqold(klon)
50        REAL zdland(klon)
51        LOGICAL llflag(klon)
52        INTEGER k, i, is, icall
53        REAL ztglace, zdphi, zqeen, zseen, zscde, zqude
54        REAL zmfusk, zmfuqk, zmfulk, zbuo, zdnoprc, zprcon, zlnew
55    
56        REAL zpbot(klon), zptop(klon), zrho(klon)
57        REAL zdprho, zentr, zpmid, zmftest, zmfmax
58        LOGICAL llo1, llo2
59    
60        REAL zwmax(klon), zzzmb
61        INTEGER klwmin(klon) ! level of maximum vertical velocity
62        real fact
63    
64        !----------------------------------------------------------------------
65    
66        ztglace = RTT - 13.
67    
68        ! Chercher le niveau où la vitesse verticale est maximale :
69    
70        DO i = 1, klon
71           klwmin(i) = klev
72           zwmax(i) = 0.0
73        ENDDO
74    
75        DO k = klev, 3, -1
76           DO i = 1, klon
77              IF (pvervel(i, k) < zwmax(i)) THEN
78                 zwmax(i) = pvervel(i, k)
79                 klwmin(i) = k
80              ENDIF
81           ENDDO
82        ENDDO
83    
84        ! Set default values:
85    
86        DO i = 1, klon
87           IF (.NOT. ldcum(i)) ktype(i)=0
88        ENDDO
89    
90        DO k=1, klev
91           DO i = 1, klon
92              plu(i, k)=0.
93              pmfu(i, k)=0.
94              pmfus(i, k)=0.
95              pmfuq(i, k)=0.
96              pmful(i, k)=0.
97              plude(i, k)=0.
98              pdmfup(i, k)=0.
99              IF (.NOT. ldcum(i) .OR. ktype(i) == 3) klab(i, k)=0
100              IF (.NOT. ldcum(i) .AND. paph(i, k) < 4e4) kctop0(i) = k
101           ENDDO
102        ENDDO
103    
104        DO i = 1, klon
105           IF (ldland(i)) THEN
106              zdland(i)=3.0E4
107              zdphi=pgeoh(i, kctop0(i))-pgeoh(i, kcbot(i))
108              IF (ptu(i, kctop0(i)) >= ztglace) zdland(i)=zdphi
109              zdland(i)=MAX(3.0E4, zdland(i))
110              zdland(i)=MIN(5.0E4, zdland(i))
111           ENDIF
112        ENDDO
113    
114        ! Initialiser les valeurs au niveau d'ascendance
115    
116        DO i = 1, klon
117           kctop(i) = klev-1
118           IF (.NOT. ldcum(i)) THEN
119              kcbot(i) = klev-1
120              pmfub(i) = 0.
121              pqu(i, klev) = 0.
122           ENDIF
123           pmfu(i, klev) = pmfub(i)
124           pmfus(i, klev) = pmfub(i) * (RCPD * ptu(i, klev)+pgeoh(i, klev))
125           pmfuq(i, klev) = pmfub(i) * pqu(i, klev)
126        ENDDO
127    
128        DO i = 1, klon
129           ldcum(i) = .FALSE.
130        ENDDO
131    
132        ! Do ascent: subcloud layer (klab=1), clouds (klab=2) by doing
133        ! first dry-adiabatic ascent and then by adjusting t, q and l
134        ! accordingly in flxadjtq, then check for buoyancy and set flags
135        ! accordingly.
136    
137        DO k = klev - 1, 3, -1
138           IF (LMFMID .AND. k < klev - 1 .AND. k > klev / 2) THEN
139              DO i = 1, klon
140                 IF (.NOT. ldcum(i) .AND. klab(i, k + 1) == 0 .AND. &
141                      pqen(i, k) > 0.9 * pqsen(i, k)) THEN
142                    ptu(i, k+1) = pten(i, k) +(pgeo(i, k)-pgeoh(i, k+1))/RCPD
143                    pqu(i, k+1) = pqen(i, k)
144                    plu(i, k+1) = 0.0
145                    zzzmb = MAX(CMFCMIN, -pvervel(i, k)/RG)
146                    zmfmax = (paph(i, k) - paph(i, k-1)) / (RG * pdtime)
147                    pmfub(i) = MIN(zzzmb, zmfmax)
148                    pmfu(i, k+1) = pmfub(i)
149                    pmfus(i, k+1) = pmfub(i) * (RCPD * ptu(i, k+1)+pgeoh(i, k+1))
150                    pmfuq(i, k+1) = pmfub(i) * pqu(i, k+1)
151                    pmful(i, k+1) = 0.0
152                    pdmfup(i, k+1) = 0.0
153                    kcbot(i) = k
154                    klab(i, k+1) = 1
155                    ktype(i) = 3
156                    pentr(i) = ENTRMID
157                 ENDIF
158              ENDDO
159           ENDIF
160    
161           is = 0
162           DO i = 1, klon
163              is = is + klab(i, k+1)
164              IF (klab(i, k+1) == 0) klab(i, k) = 0
165              llflag(i) = .FALSE.
166              IF (klab(i, k+1) > 0) llflag(i) = .TRUE.
167           ENDDO
168           IF (is == 0) cycle
169    
170           ! Calculer le taux d'entraînement et de détraînement :
171    
172           DO i = 1, klon
173              pen_u(i, k) = 0.0
174              pde_u(i, k) = 0.0
175              zrho(i) = paph(i, k + 1) / (RD * ptenh(i, k + 1))
176              zpbot(i) = paph(i, kcbot(i))
177              zptop(i) = paph(i, kctop0(i))
178           ENDDO
179    
180           DO i = 1, klon
181              IF (ldcum(i)) THEN
182                 zdprho = (paph(i, k + 1) - paph(i, k)) / (RG * zrho(i))
183                 zentr=pentr(i) * pmfu(i, k+1) * zdprho
184                 llo1=k < kcbot(i)
185                 IF (llo1) pde_u(i, k)=zentr
186                 zpmid=0.5 * (zpbot(i)+zptop(i))
187                 llo2 = llo1 .AND. ktype(i) == 2 &
188                      .AND. (zpbot(i) - paph(i, k) < 0.2E5 .OR. paph(i, k) > zpmid)
189                 IF (llo2) pen_u(i, k)=zentr
190                 llo2 = llo1 .AND. (ktype(i) == 1 .OR. ktype(i) == 3) .AND. &
191                      (k >= MAX(klwmin(i), kctop0(i) + 2) .OR. pap(i, k) > zpmid)
192                 IF (llo2) pen_u(i, k)=zentr
193                 llo1=pen_u(i, k) > 0. .AND. (ktype(i) == 1 .OR. ktype(i) == 2)
194                 IF (llo1) THEN
195                    fact = 1. + 3. * (1. - MIN(1., (zpbot(i) - pap(i, k)) / 1.5E4))
196                    zentr = zentr * fact
197                    pen_u(i, k)=pen_u(i, k) * fact
198                    pde_u(i, k)=pde_u(i, k) * fact
199                 ENDIF
200                 IF (llo2 .AND. pqenh(i, k+1) > 1e-5) &
201                      pen_u(i, k)=zentr+MAX(pqte(i, k), 0.)/pqenh(i, k+1) * &
202                      zrho(i) * zdprho
203              ENDIF
204           end DO
205    
206           ! Do adiabatic ascent for entraining/detraining plume
207    
208           DO i = 1, klon
209              IF (llflag(i)) THEN
210                 IF (k < kcbot(i)) THEN
211                    zmftest = pmfu(i, k+1)+pen_u(i, k)-pde_u(i, k)
212                    zmfmax = MIN(zmftest, &
213                         (paph(i, k) - paph(i, k - 1)) / (RG * pdtime))
214                    pen_u(i, k)=MAX(pen_u(i, k)-MAX(0.0, zmftest-zmfmax), 0.0)
215                 ENDIF
216                 pde_u(i, k)=MIN(pde_u(i, k), 0.75 * pmfu(i, k+1))
217                 ! calculer le flux de masse du niveau k a partir de celui du k+1
218                 pmfu(i, k)=pmfu(i, k+1)+pen_u(i, k)-pde_u(i, k)
219                 ! calculer les valeurs Su, Qu et l du niveau k dans le
220                 ! panache montant
221                 zqeen=pqenh(i, k+1) * pen_u(i, k)
222                 zseen=(RCPD * ptenh(i, k+1)+pgeoh(i, k+1)) * pen_u(i, k)
223                 zscde=(RCPD * ptu(i, k+1)+pgeoh(i, k+1)) * pde_u(i, k)
224                 zqude=pqu(i, k+1) * pde_u(i, k)
225                 plude(i, k)=plu(i, k+1) * pde_u(i, k)
226                 zmfusk=pmfus(i, k+1)+zseen-zscde
227                 zmfuqk=pmfuq(i, k+1)+zqeen-zqude
228                 zmfulk=pmful(i, k+1) -plude(i, k)
229                 plu(i, k)=zmfulk * (1./MAX(CMFCMIN, pmfu(i, k)))
230                 pqu(i, k)=zmfuqk * (1./MAX(CMFCMIN, pmfu(i, k)))
231                 ptu(i, k)=(zmfusk * (1./MAX(CMFCMIN, pmfu(i, k)))- &
232                      pgeoh(i, k))/RCPD
233                 ptu(i, k)=MAX(100., ptu(i, k))
234                 ptu(i, k)=MIN(400., ptu(i, k))
235                 zqold(i)=pqu(i, k)
236              ELSE
237                 zqold(i)=0.0
238              ENDIF
239           end DO
240    
241           ! Do corrections for moist ascent by adjusting t, q and l
242    
243           icall = 1
244           CALL flxadjtq(paph(1, k), ptu(1, k), pqu(1, k), llflag, icall)
245    
246           DO i = 1, klon
247              IF (llflag(i) .AND. pqu(i, k).NE.zqold(i)) THEN
248                 klab(i, k) = 2
249                 plu(i, k) = plu(i, k)+zqold(i)-pqu(i, k)
250                 zbuo = ptu(i, k) * (1.+RETV * pqu(i, k))- &
251                      ptenh(i, k) * (1.+RETV * pqenh(i, k))
252                 IF (klab(i, k+1) == 1) zbuo=zbuo+0.5
253                 IF (zbuo > 0. .AND. pmfu(i, k) >= 0.1 * pmfub(i)) THEN
254                    kctop(i) = k
255                    ldcum(i) = .TRUE.
256                    zdnoprc = 1.5E4
257                    IF (ldland(i)) zdnoprc = zdland(i)
258                    zprcon = CPRCON
259                    IF ((zpbot(i) - paph(i, k)) < zdnoprc) zprcon = 0.
260                    zlnew=plu(i, k)/(1.+zprcon * (pgeoh(i, k)-pgeoh(i, k+1)))
261                    pdmfup(i, k)=MAX(0., (plu(i, k)-zlnew) * pmfu(i, k))
262                    plu(i, k)=zlnew
263                 ELSE
264                    klab(i, k)=0
265                    pmfu(i, k)=0.
266                 ENDIF
267              ENDIF
268           end DO
269           DO i = 1, klon
270              IF (llflag(i)) THEN
271                 pmful(i, k)=plu(i, k) * pmfu(i, k)
272                 pmfus(i, k)=(RCPD * ptu(i, k)+pgeoh(i, k)) * pmfu(i, k)
273                 pmfuq(i, k)=pqu(i, k) * pmfu(i, k)
274              ENDIF
275           end DO
276        end DO
277    
278        ! Determine convective fluxes above non-buoyancy level (note:
279        ! cloud variables like t, q and l are not affected by detrainment
280        ! and are already known from previous calculations above).
281    
282        DO i = 1, klon
283           IF (kctop(i) == klev-1) ldcum(i) = .FALSE.
284           kcbot(i) = MAX(kcbot(i), kctop(i))
285        ENDDO
286    
287        ldcum(1)=ldcum(1)
288    
289        is = 0
290        DO i = 1, klon
291           if (ldcum(i)) is = is + 1
292        ENDDO
293        kcum = is
294        IF (is /= 0) then
295           DO i = 1, klon
296              IF (ldcum(i)) THEN
297                 k=kctop(i)-1
298                 pde_u(i, k)=(1.-CMFCTOP) * pmfu(i, k+1)
299                 plude(i, k)=pde_u(i, k) * plu(i, k+1)
300                 pmfu(i, k)=pmfu(i, k+1)-pde_u(i, k)
301                 zlnew=plu(i, k)
302                 pdmfup(i, k)=MAX(0., (plu(i, k)-zlnew) * pmfu(i, k))
303                 plu(i, k)=zlnew
304                 pmfus(i, k)=(RCPD * ptu(i, k)+pgeoh(i, k)) * pmfu(i, k)
305                 pmfuq(i, k)=pqu(i, k) * pmfu(i, k)
306                 pmful(i, k)=plu(i, k) * pmfu(i, k)
307                 plude(i, k-1)=pmful(i, k)
308              ENDIF
309           end DO
310        end IF
311    
312      END SUBROUTINE flxasc
313    
314    end module flxasc_m
Legend:
Removed from v.62  
changed lines
  Added in v.82
  ViewVC Help
Powered by ViewVC 1.1.21