/[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

revision 62 by guez, Thu Jul 26 14:37:37 2012 UTC revision 71 by guez, Mon Jul 8 18:12:18 2013 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 pap(klon, klev), paph(klon, klev+1)
28        REAL pqte(klon, klev)
29        REAL pvervel(klon, klev) ! vitesse verticale en Pa/s
30        LOGICAL 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    
63        !----------------------------------------------------------------------
64    
65        ztglace = RTT - 13.
66    
67        ! Chercher le niveau où la vitesse verticale est maximale :
68    
69        DO i = 1, klon
70           klwmin(i) = klev
71           zwmax(i) = 0.0
72        ENDDO
73    
74        DO k = klev, 3, -1
75           DO i = 1, klon
76              IF (pvervel(i, k) < zwmax(i)) THEN
77                 zwmax(i) = pvervel(i, k)
78                 klwmin(i) = k
79              ENDIF
80           ENDDO
81        ENDDO
82    
83        ! Set default values:
84    
85        DO i = 1, klon
86           IF (.NOT. ldcum(i)) ktype(i)=0
87        ENDDO
88    
89        DO k=1, klev
90           DO i = 1, klon
91              plu(i, k)=0.
92              pmfu(i, k)=0.
93              pmfus(i, k)=0.
94              pmfuq(i, k)=0.
95              pmful(i, k)=0.
96              plude(i, k)=0.
97              pdmfup(i, k)=0.
98              IF(.NOT. ldcum(i).OR.ktype(i) == 3) klab(i, k)=0
99              IF(.NOT. ldcum(i).AND.paph(i, k) < 4.E4) kctop0(i)=k
100           ENDDO
101        ENDDO
102    
103        DO i = 1, klon
104           IF (ldland(i)) THEN
105              zdland(i)=3.0E4
106              zdphi=pgeoh(i, kctop0(i))-pgeoh(i, kcbot(i))
107              IF (ptu(i, kctop0(i)) >= ztglace) zdland(i)=zdphi
108              zdland(i)=MAX(3.0E4, zdland(i))
109              zdland(i)=MIN(5.0E4, zdland(i))
110           ENDIF
111        ENDDO
112    
113        ! Initialiser les valeurs au niveau d'ascendance
114    
115        DO i = 1, klon
116           kctop(i) = klev-1
117           IF (.NOT. ldcum(i)) THEN
118              kcbot(i) = klev-1
119              pmfub(i) = 0.
120              pqu(i, klev) = 0.
121           ENDIF
122           pmfu(i, klev) = pmfub(i)
123           pmfus(i, klev) = pmfub(i) * (RCPD * ptu(i, klev)+pgeoh(i, klev))
124           pmfuq(i, klev) = pmfub(i) * pqu(i, klev)
125        ENDDO
126    
127        DO i = 1, klon
128           ldcum(i) = .FALSE.
129        ENDDO
130    
131        ! Do ascent: subcloud layer (klab=1), clouds (klab=2) by doing
132        ! first dry-adiabatic ascent and then by adjusting t, q and l
133        ! accordingly in flxadjtq, then check for buoyancy and set flags
134        ! accordingly.
135    
136        DO k = klev - 1, 3, -1
137           IF (LMFMID .AND. k < klev - 1 .AND. k > klev / 2) THEN
138              DO i = 1, klon
139                 IF (.NOT. ldcum(i) .AND. klab(i, k + 1) == 0 .AND. &
140                      pqen(i, k) > 0.9 * pqsen(i, k)) THEN
141                    ptu(i, k+1) = pten(i, k) +(pgeo(i, k)-pgeoh(i, k+1))/RCPD
142                    pqu(i, k+1) = pqen(i, k)
143                    plu(i, k+1) = 0.0
144                    zzzmb = MAX(CMFCMIN, -pvervel(i, k)/RG)
145                    zmfmax = (paph(i, k)-paph(i, k-1))/(RG * pdtime)
146                    pmfub(i) = MIN(zzzmb, zmfmax)
147                    pmfu(i, k+1) = pmfub(i)
148                    pmfus(i, k+1) = pmfub(i) * (RCPD * ptu(i, k+1)+pgeoh(i, k+1))
149                    pmfuq(i, k+1) = pmfub(i) * pqu(i, k+1)
150                    pmful(i, k+1) = 0.0
151                    pdmfup(i, k+1) = 0.0
152                    kcbot(i) = k
153                    klab(i, k+1) = 1
154                    ktype(i) = 3
155                    pentr(i) = ENTRMID
156                 ENDIF
157              ENDDO
158           ENDIF
159    
160           is = 0
161           DO i = 1, klon
162              is = is + klab(i, k+1)
163              IF (klab(i, k+1) == 0) klab(i, k) = 0
164              llflag(i) = .FALSE.
165              IF (klab(i, k+1) > 0) llflag(i) = .TRUE.
166           ENDDO
167           IF (is == 0) cycle
168    
169           ! Calculer le taux d'entraînement et de détraînement :
170    
171           DO i = 1, klon
172              pen_u(i, k) = 0.0
173              pde_u(i, k) = 0.0
174              zrho(i)=paph(i, k+1)/(RD * ptenh(i, k+1))
175              zpbot(i)=paph(i, kcbot(i))
176              zptop(i)=paph(i, kctop0(i))
177           ENDDO
178    
179           DO i = 1, klon
180              IF(ldcum(i)) THEN
181                 zdprho=(paph(i, k+1)-paph(i, k))/(RG * zrho(i))
182                 zentr=pentr(i) * pmfu(i, k+1) * zdprho
183                 llo1=k < kcbot(i)
184                 IF(llo1) pde_u(i, k)=zentr
185                 zpmid=0.5 * (zpbot(i)+zptop(i))
186                 llo2=llo1.AND.ktype(i) == 2.AND. &
187                      (zpbot(i)-paph(i, k) < 0.2E5.OR. &
188                      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                    zentr=zentr * (1.+3. * (1.-MIN(1., (zpbot(i)-pap(i, k))/1.5E4)))
196                    pen_u(i, k)=pen_u(i, k) * (1.+3. * (1.-MIN(1., &
197                         (zpbot(i)-pap(i, k))/1.5E4)))
198                    pde_u(i, k)=pde_u(i, k) * (1.+3. * (1.-MIN(1., &
199                         (zpbot(i)-pap(i, k))/1.5E4)))
200                 ENDIF
201                 IF(llo2.AND.pqenh(i, k+1) > 1.E-5) &
202                      pen_u(i, k)=zentr+MAX(pqte(i, k), 0.)/pqenh(i, k+1) * &
203                      zrho(i) * zdprho
204              ENDIF
205           end DO
206    
207           ! Do adiabatic ascent for entraining/detraining plume
208    
209           DO i = 1, klon
210              IF (llflag(i)) THEN
211                 IF (k < kcbot(i)) THEN
212                    zmftest = pmfu(i, k+1)+pen_u(i, k)-pde_u(i, k)
213                    zmfmax = MIN(zmftest, (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.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    
277        end DO
278    
279        ! Determine convective fluxes above non-buoyancy level (note:
280        ! cloud variables like t, q and l are not affected by detrainment
281        ! and are already known from previous calculations above).
282    
283        DO i = 1, klon
284           IF (kctop(i) == klev-1) ldcum(i) = .FALSE.
285           kcbot(i) = MAX(kcbot(i), kctop(i))
286        ENDDO
287    
288        ldcum(1)=ldcum(1)
289    
290        is = 0
291        DO i = 1, klon
292           if (ldcum(i)) is = is + 1
293        ENDDO
294        kcum = is
295        IF (is /= 0) then
296           DO i = 1, klon
297              IF (ldcum(i)) THEN
298                 k=kctop(i)-1
299                 pde_u(i, k)=(1.-CMFCTOP) * pmfu(i, k+1)
300                 plude(i, k)=pde_u(i, k) * plu(i, k+1)
301                 pmfu(i, k)=pmfu(i, k+1)-pde_u(i, k)
302                 zlnew=plu(i, k)
303                 pdmfup(i, k)=MAX(0., (plu(i, k)-zlnew) * pmfu(i, k))
304                 plu(i, k)=zlnew
305                 pmfus(i, k)=(RCPD * ptu(i, k)+pgeoh(i, k)) * pmfu(i, k)
306                 pmfuq(i, k)=pqu(i, k) * pmfu(i, k)
307                 pmful(i, k)=plu(i, k) * pmfu(i, k)
308                 plude(i, k-1)=pmful(i, k)
309              ENDIF
310           end DO
311        end IF
312    
313      END SUBROUTINE flxasc
314    
315    end module flxasc_m
Legend:
Removed from v.62  
changed lines
  Added in v.71
  ViewVC Help
Powered by ViewVC 1.1.21