/[lmdze]/trunk/Sources/phylmd/CV30_routines/cv30_yield.f
ViewVC logotype

Diff of /trunk/Sources/phylmd/CV30_routines/cv30_yield.f

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

revision 190 by guez, Thu Apr 14 15:15:56 2016 UTC revision 200 by guez, Thu Jun 2 15:40:30 2016 UTC
# Line 9  contains Line 9  contains
9         qent, uent, vent, nent, elij, sig, tv, tvp, iflag, precip, VPrecip, &         qent, uent, vent, nent, elij, sig, tv, tvp, iflag, precip, VPrecip, &
10         ft, fr, fu, fv, upwd, dnwd, dnwd0, ma, mike, tls, tps, qcondc)         ft, fr, fu, fv, upwd, dnwd, dnwd0, ma, mike, tls, tps, qcondc)
11    
12        ! Tendencies, precipitation, variables of interface with other
13        ! processes, etc.
14    
15      use conema3_m, only: iflag_clw      use conema3_m, only: iflag_clw
16      use cv30_param_m, only: delta, minorig, nl, sigd      use cv30_param_m, only: minorig, nl, sigd
17      use cv_thermo_m, only: cl, cpd, cpv, grav, rowl, rrd, rrv      use cv_thermo_m, only: cl, cpd, cpv, rowl, rrd, rrv
18      USE dimphy, ONLY: klev, klon      USE dimphy, ONLY: klev, klon
19        use SUPHEC_M, only: rg
20    
21      ! inputs:      ! inputs:
22      integer, intent(in):: icb(:), inb(:) ! (ncum)      integer, intent(in):: icb(:), inb(:) ! (ncum)
23      real, intent(in):: delt      real, intent(in):: delt
24      real t(klon, klev), rr(klon, klev), u(klon, klev), v(klon, klev)      real, intent(in):: t(klon, klev), rr(klon, klev)
25        real, intent(in):: u(klon, klev), v(klon, klev)
26      real gz(klon, klev)      real gz(klon, klev)
27      real p(klon, klev)      real p(klon, klev)
28      real ph(klon, klev + 1), h(klon, klev), hp(klon, klev)      real ph(klon, klev + 1), h(klon, klev), hp(klon, klev)
# Line 27  contains Line 32  contains
32      real m(klon, klev)      real m(klon, klev)
33      real tp(klon, klev)      real tp(klon, klev)
34      real mp(klon, klev), rp(klon, klev), up(klon, klev)      real mp(klon, klev), rp(klon, klev), up(klon, klev)
35      real vp(klon, klev), wt(klon, klev)      real, intent(in):: vp(:, 2:) ! (ncum, 2:nl)
36        real, intent(in):: wt(:, :) ! (ncum, nl - 1)
37      real, intent(in):: water(:, :), evap(:, :) ! (ncum, nl)      real, intent(in):: water(:, :), evap(:, :) ! (ncum, nl)
38      real, intent(in):: b(:, :) ! (ncum, nl - 1)      real, intent(in):: b(:, :) ! (ncum, nl - 1)
39      real ment(klon, klev, klev), qent(klon, klev, klev), uent(klon, klev, klev)      real ment(klon, klev, klev), qent(klon, klev, klev), uent(klon, klev, klev)
# Line 37  contains Line 43  contains
43      real sig(klon, klev)      real sig(klon, klev)
44      real tv(klon, klev), tvp(klon, klev)      real tv(klon, klev), tvp(klon, klev)
45    
     ! input / output:  
     integer iflag(klon)  
   
46      ! outputs:      ! outputs:
47        integer, intent(out):: iflag(:) ! (ncum)
48      real precip(klon)      real precip(klon)
49      real VPrecip(klon, klev + 1)      real VPrecip(klon, klev + 1)
50      real ft(klon, klev), fr(klon, klev), fu(klon, klev), fv(klon, klev)      real ft(klon, klev), fr(klon, klev), fu(klon, klev), fv(klon, klev)
# Line 52  contains Line 56  contains
56      real qcondc(klon, klev)      real qcondc(klon, klev)
57    
58      ! Local:      ! Local:
59        real, parameter:: delta = 0.01 ! interface cloud parameterization
60      integer ncum      integer ncum
61      integer i, k, il, n, j, num1      integer i, k, il, n, j
62      real rat, awat, delti      real awat, delti
63      real ax, bx, cx, dx      real ax, bx, cx, dx
64      real cpinv, rdcp, dpinv      real cpinv, rdcp, dpinv
65      real lvcp(klon, klev)      real lvcp(klon, klev)
# Line 67  contains Line 72  contains
72      !-------------------------------------------------------------      !-------------------------------------------------------------
73    
74      ncum = size(icb)      ncum = size(icb)
75        iflag = 0
76    
77      ! initialization:      ! initialization:
78    
# Line 96  contains Line 102  contains
102         enddo         enddo
103      enddo      enddo
104    
105      ! calculate surface precipitation in mm / day      ! calculate surface precipitation in mm / day
106    
107      do il = 1, ncum      do il = 1, ncum
108         if (ep(il, inb(il)) >= 1e-4) precip(il) = wt(il, 1) * sigd &         if (ep(il, inb(il)) >= 1e-4) precip(il) = wt(il, 1) * sigd &
109              * water(il, 1) * 86400. * 1000. / (rowl * grav)              * water(il, 1) * 86400. * 1000. / (rowl * rg)
110      enddo      enddo
111    
112      ! CALCULATE VERTICAL PROFILE OF PRECIPITATIONs IN kg / m2 / s ===      ! CALCULATE VERTICAL PROFILE OF PRECIPITATIONs IN kg / m2 / s ===
# Line 109  contains Line 115  contains
115      do k = 1, nl - 1      do k = 1, nl - 1
116         do il = 1, ncum         do il = 1, ncum
117            if (k <= inb(il)) VPrecip(il, k) = wt(il, k) * sigd * water(il, k) &            if (k <= inb(il)) VPrecip(il, k) = wt(il, k) * sigd * water(il, k) &
118                 / grav                 / rg
119         end do         end do
120      end do      end do
121    
122      ! calculate tendencies of lowest level potential temperature      ! calculate tendencies of lowest level potential temperature
123      ! and mixing ratio      ! and mixing ratio
124    
125      do il = 1, ncum      do il = 1, ncum
126         work(il) = 1.0 / (ph(il, 1) - ph(il, 2))         work(il) = 1.0 / (ph(il, 1) - ph(il, 2))
# Line 128  contains Line 134  contains
134      enddo      enddo
135    
136      do il = 1, ncum      do il = 1, ncum
137         ! Consist vect:         if (0.01 * rg * work(il) * am(il) >= delti) iflag(il) = 1
        if (0.01 * grav * work(il) * am(il) >= delti) iflag(il) = 1  
   
        ft(il, 1) = 0.01 * grav * work(il) * am(il) * (t(il, 2) - t(il, 1) &  
             + (gz(il, 2) - gz(il, 1)) / cpn(il, 1))  
138    
139         ft(il, 1) = ft(il, 1) - 0.5 * lvcp(il, 1) * sigd * (evap(il, 1) &         ft(il, 1) = 0.01 * rg * work(il) * am(il) * (t(il, 2) - t(il, 1) &
140              + evap(il, 2))              + (gz(il, 2) - gz(il, 1)) / cpn(il, 1)) - 0.5 * lvcp(il, 1) &
141                * sigd * (evap(il, 1) + evap(il, 2)) - 0.009 * rg * sigd &
142         ft(il, 1) = ft(il, 1) - 0.009 * grav * sigd * mp(il, 2) &              * mp(il, 2) * t(il, 1) * b(il, 1) * work(il) + 0.01 * sigd &
143              * t(il, 1) * b(il, 1) * work(il)              * wt(il, 1) * (cl - cpd) * water(il, 2) * (t(il, 2) - t(il, 1)) &
144                * work(il) / cpn(il, 1)
        ft(il, 1) = ft(il, 1) + 0.01 * sigd * wt(il, 1) * (cl - cpd) &  
             * water(il, 2) * (t(il, 2) - t(il, 1)) * work(il) / cpn(il, 1)  
145    
146         !jyg1 Correction pour mieux conserver l'eau (conformite avec CONVECT4.3)         !jyg1 Correction pour mieux conserver l'eau (conformite avec CONVECT4.3)
147         ! (sb: pour l'instant, on ne fait que le chgt concernant grav, pas evap)         ! (sb: pour l'instant, on ne fait que le chgt concernant rg, pas evap)
148         fr(il, 1) = 0.01 * grav * mp(il, 2) * (rp(il, 2) - rr(il, 1)) &         fr(il, 1) = 0.01 * rg * mp(il, 2) * (rp(il, 2) - rr(il, 1)) &
149              * work(il) + sigd * 0.5 * (evap(il, 1) + evap(il, 2))              * work(il) + sigd * 0.5 * (evap(il, 1) + evap(il, 2))
150         ! + tard : + sigd * evap(il, 1)         ! + tard : + sigd * evap(il, 1)
151    
152         fr(il, 1) = fr(il, 1) + 0.01 * grav * am(il) * (rr(il, 2) - rr(il, 1)) &         fr(il, 1) = fr(il, 1) + 0.01 * rg * am(il) * (rr(il, 2) - rr(il, 1)) &
153              * work(il)              * work(il)
154    
155         fu(il, 1) = fu(il, 1) + 0.01 * grav * work(il) * (mp(il, 2) &         fu(il, 1) = fu(il, 1) + 0.01 * rg * work(il) * (mp(il, 2) &
156              * (up(il, 2) - u(il, 1)) + am(il) * (u(il, 2) - u(il, 1)))              * (up(il, 2) - u(il, 1)) + am(il) * (u(il, 2) - u(il, 1)))
157         fv(il, 1) = fv(il, 1) + 0.01 * grav * work(il) * (mp(il, 2) &         fv(il, 1) = fv(il, 1) + 0.01 * rg * work(il) * (mp(il, 2) &
158              * (vp(il, 2) - v(il, 1)) + am(il) * (v(il, 2) - v(il, 1)))              * (vp(il, 2) - v(il, 1)) + am(il) * (v(il, 2) - v(il, 1)))
159      enddo ! il      enddo
160    
161      do j = 2, nl      do j = 2, nl
162         do il = 1, ncum         do il = 1, ncum
163            if (j <= inb(il)) then            if (j <= inb(il)) then
164               fr(il, 1) = fr(il, 1) + 0.01 * grav * work(il) * ment(il, j, 1) &               fr(il, 1) = fr(il, 1) + 0.01 * rg * work(il) * ment(il, j, 1) &
165                    * (qent(il, j, 1) - rr(il, 1))                    * (qent(il, j, 1) - rr(il, 1))
166               fu(il, 1) = fu(il, 1) + 0.01 * grav * work(il) * ment(il, j, 1) &               fu(il, 1) = fu(il, 1) + 0.01 * rg * work(il) * ment(il, j, 1) &
167                    * (uent(il, j, 1) - u(il, 1))                    * (uent(il, j, 1) - u(il, 1))
168               fv(il, 1) = fv(il, 1) + 0.01 * grav * work(il) * ment(il, j, 1) &               fv(il, 1) = fv(il, 1) + 0.01 * rg * work(il) * ment(il, j, 1) &
169                    * (vent(il, j, 1) - v(il, 1))                    * (vent(il, j, 1) - v(il, 1))
170            endif            endif
171         enddo         enddo
172      enddo      enddo
173    
174      ! calculate tendencies of potential temperature and mixing ratio      ! calculate tendencies of potential temperature and mixing ratio
175      ! at levels above the lowest level      ! at levels above the lowest level
176    
177      ! first find the net saturated updraft and downdraft mass fluxes      ! first find the net saturated updraft and downdraft mass fluxes
178      ! through each level      ! through each level
179    
180      loop_i: do i = 2, nl - 1      loop_i: do i = 2, nl - 1
181         num1 = 0         if (any(inb >= i)) then
   
        do il = 1, ncum  
           if (i <= inb(il)) num1 = num1 + 1  
        enddo  
   
        if (num1 > 0) then  
182            amp1(:ncum) = 0.            amp1(:ncum) = 0.
183            ad(:ncum) = 0.            ad(:ncum) = 0.
184    
# Line 221  contains Line 215  contains
215                  dpinv = 1.0 / (ph(il, i) - ph(il, i + 1))                  dpinv = 1.0 / (ph(il, i) - ph(il, i + 1))
216                  cpinv = 1.0 / cpn(il, i)                  cpinv = 1.0 / cpn(il, i)
217    
218                  ! Vecto:                  if (0.01 * rg * dpinv * amp1(il) >= delti) iflag(il) = 1
                 if (0.01 * grav * dpinv * amp1(il) >= delti) iflag(il) = 1  
219    
220                  ft(il, i) = 0.01 * grav * dpinv * (amp1(il) * (t(il, i + 1) &                  ft(il, i) = 0.01 * rg * dpinv * (amp1(il) * (t(il, i + 1) &
221                       - t(il, i) + (gz(il, i + 1) - gz(il, i)) * cpinv) &                       - t(il, i) + (gz(il, i + 1) - gz(il, i)) * cpinv) &
222                       - ad(il) * (t(il, i) - t(il, i - 1) + (gz(il, i) &                       - ad(il) * (t(il, i) - t(il, i - 1) + (gz(il, i) &
223                       - gz(il, i - 1)) * cpinv)) - 0.5 * sigd * lvcp(il, i) &                       - gz(il, i - 1)) * cpinv)) - 0.5 * sigd * lvcp(il, i) &
224                       * (evap(il, i) + evap(il, i + 1))                       * (evap(il, i) + evap(il, i + 1)) - 0.009 * rg * sigd &
225                  rat = cpn(il, i - 1) * cpinv                       * (mp(il, i + 1) * t(il, i) * b(il, i) - mp(il, i) &
226                  ft(il, i) = ft(il, i) - 0.009 * grav * sigd * (mp(il, i + 1) &                       * t(il, i - 1) * cpn(il, i - 1) * cpinv * b(il, i - 1)) &
227                       * t(il, i) * b(il, i) - mp(il, i) * t(il, i - 1) * rat &                       * dpinv + 0.01 * rg * dpinv * ment(il, i, i) &
                      * b(il, i - 1)) * dpinv  
                 ft(il, i) = ft(il, i) + 0.01 * grav * dpinv * ment(il, i, i) &  
228                       * (hp(il, i) - h(il, i) + t(il, i) * (cpv - cpd) &                       * (hp(il, i) - h(il, i) + t(il, i) * (cpv - cpd) &
229                       * (rr(il, i) - qent(il, i, i))) * cpinv                       * (rr(il, i) - qent(il, i, i))) * cpinv + 0.01 * sigd &
230                         * wt(il, i) * (cl - cpd) * water(il, i + 1) &
231                  ft(il, i) = ft(il, i) + 0.01 * sigd * wt(il, i) * (cl - cpd) &                       * (t(il, i + 1) - t(il, i)) * dpinv * cpinv
232                       * water(il, i + 1) * (t(il, i + 1) - t(il, i)) * dpinv &                  fr(il, i) = 0.01 * rg * dpinv * (amp1(il) * (rr(il, i + 1) &
                      * cpinv  
   
                 fr(il, i) = 0.01 * grav * dpinv * (amp1(il) * (rr(il, i + 1) &  
233                       - rr(il, i)) - ad(il) * (rr(il, i) - rr(il, i - 1)))                       - rr(il, i)) - ad(il) * (rr(il, i) - rr(il, i - 1)))
234                  fu(il, i) = fu(il, i) + 0.01 * grav * dpinv * (amp1(il) &                  fu(il, i) = fu(il, i) + 0.01 * rg * dpinv * (amp1(il) &
235                       * (u(il, i + 1) - u(il, i)) - ad(il) * (u(il, i) &                       * (u(il, i + 1) - u(il, i)) - ad(il) * (u(il, i) &
236                       - u(il, i - 1)))                       - u(il, i - 1)))
237                  fv(il, i) = fv(il, i) + 0.01 * grav * dpinv * (amp1(il) &                  fv(il, i) = fv(il, i) + 0.01 * rg * dpinv * (amp1(il) &
238                       * (v(il, i + 1) - v(il, i)) - ad(il) * (v(il, i) &                       * (v(il, i + 1) - v(il, i)) - ad(il) * (v(il, i) &
239                       - v(il, i - 1)))                       - v(il, i - 1)))
240               endif               endif
# Line 261  contains Line 249  contains
249                     awat = elij(il, k, i) - (1. - ep(il, i)) * clw(il, i)                     awat = elij(il, k, i) - (1. - ep(il, i)) * clw(il, i)
250                     awat = amax1(awat, 0.0)                     awat = amax1(awat, 0.0)
251    
252                     fr(il, i) = fr(il, i) + 0.01 * grav * dpinv &                     fr(il, i) = fr(il, i) + 0.01 * rg * dpinv &
253                          * ment(il, k, i) * (qent(il, k, i) - awat - rr(il, i))                          * ment(il, k, i) * (qent(il, k, i) - awat - rr(il, i))
254                     fu(il, i) = fu(il, i) + 0.01 * grav * dpinv &                     fu(il, i) = fu(il, i) + 0.01 * rg * dpinv &
255                          * ment(il, k, i) * (uent(il, k, i) - u(il, i))                          * ment(il, k, i) * (uent(il, k, i) - u(il, i))
256                     fv(il, i) = fv(il, i) + 0.01 * grav * dpinv &                     fv(il, i) = fv(il, i) + 0.01 * rg * dpinv &
257                          * ment(il, k, i) * (vent(il, k, i) - v(il, i))                          * ment(il, k, i) * (vent(il, k, i) - v(il, i))
258    
259                     ! (saturated updrafts resulting from mixing)                     ! (saturated updrafts resulting from mixing)
# Line 281  contains Line 269  contains
269                     dpinv = 1.0 / (ph(il, i) - ph(il, i + 1))                     dpinv = 1.0 / (ph(il, i) - ph(il, i + 1))
270                     cpinv = 1.0 / cpn(il, i)                     cpinv = 1.0 / cpn(il, i)
271    
272                     fr(il, i) = fr(il, i) + 0.01 * grav * dpinv &                     fr(il, i) = fr(il, i) + 0.01 * rg * dpinv &
273                          * ment(il, k, i) * (qent(il, k, i) - rr(il, i))                          * ment(il, k, i) * (qent(il, k, i) - rr(il, i))
274                     fu(il, i) = fu(il, i) + 0.01 * grav * dpinv &                     fu(il, i) = fu(il, i) + 0.01 * rg * dpinv &
275                          * ment(il, k, i) * (uent(il, k, i) - u(il, i))                          * ment(il, k, i) * (uent(il, k, i) - u(il, i))
276                     fv(il, i) = fv(il, i) + 0.01 * grav * dpinv &                     fv(il, i) = fv(il, i) + 0.01 * rg * dpinv &
277                          * ment(il, k, i) * (vent(il, k, i) - v(il, i))                          * ment(il, k, i) * (vent(il, k, i) - v(il, i))
278                  endif                  endif
279               end do               end do
# Line 299  contains Line 287  contains
287                  ! sb: on ne fait pas encore la correction permettant de mieux                  ! sb: on ne fait pas encore la correction permettant de mieux
288                  ! conserver l'eau:                  ! conserver l'eau:
289                  fr(il, i) = fr(il, i) + 0.5 * sigd * (evap(il, i) &                  fr(il, i) = fr(il, i) + 0.5 * sigd * (evap(il, i) &
290                       + evap(il, i + 1)) + 0.01 * grav * (mp(il, i + 1) &                       + evap(il, i + 1)) + 0.01 * rg * (mp(il, i + 1) &
291                       * (rp(il, i + 1) - rr(il, i)) - mp(il, i) * (rp(il, i) &                       * (rp(il, i + 1) - rr(il, i)) - mp(il, i) * (rp(il, i) &
292                       - rr(il, i - 1))) * dpinv                       - rr(il, i - 1))) * dpinv
293    
294                  fu(il, i) = fu(il, i) + 0.01 * grav * (mp(il, i + 1) &                  fu(il, i) = fu(il, i) + 0.01 * rg * (mp(il, i + 1) &
295                       * (up(il, i + 1) - u(il, i)) - mp(il, i) * (up(il, i) &                       * (up(il, i + 1) - u(il, i)) - mp(il, i) * (up(il, i) &
296                       - u(il, i - 1))) * dpinv                       - u(il, i - 1))) * dpinv
297                  fv(il, i) = fv(il, i) + 0.01 * grav * (mp(il, i + 1) &                  fv(il, i) = fv(il, i) + 0.01 * rg * (mp(il, i + 1) &
298                       * (vp(il, i + 1) - v(il, i)) - mp(il, i) * (vp(il, i) &                       * (vp(il, i + 1) - v(il, i)) - mp(il, i) * (vp(il, i) &
299                       - v(il, i - 1))) * dpinv                       - v(il, i - 1))) * dpinv
300               endif               endif
# Line 340  contains Line 328  contains
328         end if         end if
329      end do loop_i      end do loop_i
330    
331      ! move the detrainment at level inb down to level inb - 1      ! move the detrainment at level inb down to level inb - 1
332      ! in such a way as to preserve the vertically      ! in such a way as to preserve the vertically
333      ! integrated enthalpy and water tendencies      ! integrated enthalpy and water tendencies
334    
335      do il = 1, ncum      do il = 1, ncum
336         ax = 0.1 * ment(il, inb(il), inb(il)) * (hp(il, inb(il)) &         ax = 0.1 * ment(il, inb(il), inb(il)) * (hp(il, inb(il)) &
# Line 377  contains Line 365  contains
365    
366      end do      end do
367    
368      ! homoginize tendencies below cloud base      ! homoginize tendencies below cloud base
369    
370      do il = 1, ncum      do il = 1, ncum
371         asum(il) = 0.0         asum(il) = 0.0
# Line 409  contains Line 397  contains
397         enddo         enddo
398      enddo      enddo
399    
400      ! reset counter and return      ! reset counter and return
401    
402      do il = 1, ncum      do il = 1, ncum
403         sig(il, klev) = 2.0         sig(il, klev) = 2.0
# Line 475  contains Line 463  contains
463         enddo         enddo
464      enddo      enddo
465    
466      !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc      ! D\'etermination de la variation de flux ascendant entre
467      ! determination de la variation de flux ascendant entre      ! deux niveaux non dilu\'es mike
     ! deux niveau non dilue mike  
     !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc  
468    
469      do i = 1, nl      do i = 1, nl
470         do il = 1, ncum         do il = 1, ncum
# Line 520  contains Line 506  contains
506         enddo         enddo
507      enddo      enddo
508    
509      !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc      ! icb repr\'esente le niveau o\`u se trouve la base du nuage, et
510      ! icb represente de niveau ou se trouve la      ! inb le sommet du nuage
     ! base du nuage, et inb le top du nuage  
     !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc  
511    
512      do i = 1, klev      do i = 1, klev
513         DO il = 1, ncum         DO il = 1, ncum
# Line 534  contains Line 518  contains
518         end DO         end DO
519      enddo      enddo
520    
521      ! diagnose the in-cloud mixing ratio      ! Diagnose the in-cloud mixing ratio of condensed water
     ! of condensed water  
     !  
522    
523      do i = 1, klev      do i = 1, klev
524         do il = 1, ncum         do il = 1, ncum
Legend:
Removed from v.190  
changed lines
  Added in v.200
  ViewVC Help
Powered by ViewVC 1.1.21