/[lmdze]/trunk/libf/phylmd/Conflx/flxasc.f90
ViewVC logotype

Diff of /trunk/libf/phylmd/Conflx/flxasc.f90

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

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