source: LMDZ6/trunk/libf/phylmd/convect1.F90 @ 3639

Last change on this file since 3639 was 2346, checked in by Ehouarn Millour, 9 years ago

Physics/dynamics separation:

  • remove all references to dimensions.h from physics. nbp_lon (==iim) , nbp_lat (==jjm+1) and nbp_lev (==llm) from mod_grid_phy_lmdz should be used instead.
  • added module regular_lonlat_mod in phy_common to store information about the global (lon-lat) grid cell boundaries and centers.

EM

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 19.8 KB
Line 
1
2! $Header$
3
4SUBROUTINE convect1(len, nd, ndp1, noff, minorig, t, q, qs, u, v, p, ph, &
5    iflag, ft, fq, fu, fv, precip, cbmf, delt, ma)
6  ! .............................START PROLOGUE............................
7
8  ! SCCS IDENTIFICATION:  @(#)convect1.f        1.1 04/21/00
9  ! 19:40:52 /h/cm/library/nogaps4/src/sub/fcst/convect1.f_v
10
11  ! CONFIGURATION IDENTIFICATION:  None
12
13  ! MODULE NAME:  convect1
14
15  ! DESCRIPTION:
16
17  ! convect1     The Emanuel Cumulus Convection Scheme
18
19  ! CONTRACT NUMBER AND TITLE:  None
20
21  ! REFERENCES: Programmers  K. Emanuel (MIT), Timothy F. Hogan, M. Peng
22  ! (NRL)
23
24  ! CLASSIFICATION:  Unclassified
25
26  ! RESTRICTIONS: None
27
28  ! COMPILER DEPENDENCIES: FORTRAN 77, FORTRAN 90
29
30  ! COMPILE OPTIONS: Fortran 77: -Zu -Wf"-ei -o aggress"
31  ! Fortran 90: -O vector3,scalar3,task1,aggress,overindex  -ei -r 2
32
33  ! LIBRARIES OF RESIDENCE: /a/ops/lib/libfcst159.a
34
35  ! USAGE: call convect1(len,nd,noff,minorig,
36  ! &                   t,q,qs,u,v,
37  ! &                   p,ph,iflag,ft,
38  ! &                   fq,fu,fv,precip,cbmf,delt)
39
40  ! PARAMETERS:
41  ! Name            Type         Usage            Description
42  ! ----------      ----------     -------  ----------------------------
43
44  ! len           Integer        Input        first (i) dimension
45  ! nd            Integer        Input        vertical (k) dimension
46  ! ndp1          Integer        Input        nd + 1
47  ! noff          Integer        Input        integer limit for convection
48  ! (nd-noff)
49  ! minorig       Integer        Input        First level of convection
50  ! t             Real           Input        temperature
51  ! q             Real           Input        specific hum
52  ! qs            Real           Input        sat specific hum
53  ! u             Real           Input        u-wind
54  ! v             Real           Input        v-wind
55  ! p             Real           Input        full level pressure
56  ! ph            Real           Input        half level pressure
57  ! iflag         Integer        Output       iflag on latitude strip
58  ! ft            Real           Output       temp tend
59  ! fq            Real           Output       spec hum tend
60  ! fu            Real           Output       u-wind tend
61  ! fv            Real           Output       v-wind tend
62  ! cbmf          Real           In/Out       cumulus mass flux
63  ! delt          Real           Input        time step
64  ! iflag         Integer        Output       integer flag for Emanuel
65  ! conditions
66
67  ! COMMON BLOCKS:
68  ! Block      Name     Type    Usage              Notes
69  ! --------  --------   ----    ------   ------------------------
70
71  ! FILES: None
72
73  ! DATA BASES: None
74
75  ! NON-FILE INPUT/OUTPUT: None
76
77  ! ERROR CONDITIONS: None
78
79  ! ADDITIONAL COMMENTS: None
80
81  ! .................MAINTENANCE SECTION................................
82
83  ! MODULES CALLED:
84  ! Name           Description
85  ! convect2        Emanuel cumulus convection tendency calculations
86  ! -------     ----------------------
87  ! LOCAL VARIABLES AND
88  ! STRUCTURES:
89  ! Name     Type    Description
90  ! -------  ------  -----------
91  ! See Comments Below
92
93  ! i        Integer loop index
94  ! k        Integer loop index
95
96  ! METHOD:
97
98  ! See Emanuel, K. and M. Zivkovic-Rothman, 2000: Development and evaluation
99  ! of a
100  ! convective scheme for use in climate models.
101
102  ! FILES: None
103
104  ! INCLUDE FILES: None
105
106  ! MAKEFILE: /a/ops/met/nogaps/src/sub/fcst/fcst159lib.mak
107
108  ! ..............................END PROLOGUE.............................
109
110
111  USE dimphy
112  IMPLICIT NONE
113
114  INTEGER len
115  INTEGER nd
116  INTEGER ndp1
117  INTEGER noff
118  REAL t(len, nd)
119  REAL q(len, nd)
120  REAL qs(len, nd)
121  REAL u(len, nd)
122  REAL v(len, nd)
123  REAL p(len, nd)
124  REAL ph(len, ndp1)
125  INTEGER iflag(len)
126  REAL ft(len, nd)
127  REAL fq(len, nd)
128  REAL fu(len, nd)
129  REAL fv(len, nd)
130  REAL precip(len)
131  REAL cbmf(len)
132  REAL ma(len, nd)
133  INTEGER minorig
134  REAL delt, cpd, cpv, cl, rv, rd, lv0, g
135  REAL sigs, sigd, elcrit, tlcrit, omtsnow, dtmax, damp
136  REAL alpha, entp, coeffs, coeffr, omtrain, cu
137
138  ! -------------------------------------------------------------------
139  ! --- ARGUMENTS
140  ! -------------------------------------------------------------------
141  ! --- On input:
142
143  ! t:   Array of absolute temperature (K) of dimension ND, with first
144  ! index corresponding to lowest model level. Note that this array
145  ! will be altered by the subroutine if dry convective adjustment
146  ! occurs and if IPBL is not equal to 0.
147
148  ! q:   Array of specific humidity (gm/gm) of dimension ND, with first
149  ! index corresponding to lowest model level. Must be defined
150  ! at same grid levels as T. Note that this array will be altered
151  ! if dry convective adjustment occurs and if IPBL is not equal to 0.
152
153  ! qs:  Array of saturation specific humidity of dimension ND, with first
154  ! index corresponding to lowest model level. Must be defined
155  ! at same grid levels as T. Note that this array will be altered
156  ! if dry convective adjustment occurs and if IPBL is not equal to 0.
157
158  ! u:   Array of zonal wind velocity (m/s) of dimension ND, witth first
159  ! index corresponding with the lowest model level. Defined at
160  ! same levels as T. Note that this array will be altered if
161  ! dry convective adjustment occurs and if IPBL is not equal to 0.
162
163  ! v:   Same as u but for meridional velocity.
164
165  ! tra: Array of passive tracer mixing ratio, of dimensions (ND,NTRA),
166  ! where NTRA is the number of different tracers. If no
167  ! convective tracer transport is needed, define a dummy
168  ! input array of dimension (ND,1). Tracers are defined at
169  ! same vertical levels as T. Note that this array will be altered
170  ! if dry convective adjustment occurs and if IPBL is not equal to 0.
171
172  ! p:   Array of pressure (mb) of dimension ND, with first
173  ! index corresponding to lowest model level. Must be defined
174  ! at same grid levels as T.
175
176  ! ph:  Array of pressure (mb) of dimension ND+1, with first index
177  ! corresponding to lowest level. These pressures are defined at
178  ! levels intermediate between those of P, T, Q and QS. The first
179  ! value of PH should be greater than (i.e. at a lower level than)
180  ! the first value of the array P.
181
182  ! nl:  The maximum number of levels to which convection can penetrate, plus
183  ! 1.
184  ! NL MUST be less than or equal to ND-1.
185
186  ! delt: The model time step (sec) between calls to CONVECT
187
188  ! ----------------------------------------------------------------------------
189  ! ---   On Output:
190
191  ! iflag: An output integer whose value denotes the following:
192  ! VALUE   INTERPRETATION
193  ! -----   --------------
194  ! 0     Moist convection occurs.
195  ! 1     Moist convection occurs, but a CFL condition
196  ! on the subsidence warming is violated. This
197  ! does not cause the scheme to terminate.
198  ! 2     Moist convection, but no precip because ep(inb) lt 0.0001
199  ! 3     No moist convection because new cbmf is 0 and old cbmf is 0.
200  ! 4     No moist convection; atmosphere is not
201  ! unstable
202  ! 6     No moist convection because ihmin le minorig.
203  ! 7     No moist convection because unreasonable
204  ! parcel level temperature or specific humidity.
205  ! 8     No moist convection: lifted condensation
206  ! level is above the 200 mb level.
207  ! 9     No moist convection: cloud base is higher
208  ! then the level NL-1.
209
210  ! ft:   Array of temperature tendency (K/s) of dimension ND, defined at
211  ! same
212  ! grid levels as T, Q, QS and P.
213
214  ! fq:   Array of specific humidity tendencies ((gm/gm)/s) of dimension ND,
215  ! defined at same grid levels as T, Q, QS and P.
216
217  ! fu:   Array of forcing of zonal velocity (m/s^2) of dimension ND,
218  ! defined at same grid levels as T.
219
220  ! fv:   Same as FU, but for forcing of meridional velocity.
221
222  ! ftra: Array of forcing of tracer content, in tracer mixing ratio per
223  ! second, defined at same levels as T. Dimensioned (ND,NTRA).
224
225  ! precip: Scalar convective precipitation rate (mm/day).
226
227  ! wd:   A convective downdraft velocity scale. For use in surface
228  ! flux parameterizations. See convect.ps file for details.
229
230  ! tprime: A convective downdraft temperature perturbation scale (K).
231  ! For use in surface flux parameterizations. See convect.ps
232  ! file for details.
233
234  ! qprime: A convective downdraft specific humidity
235  ! perturbation scale (gm/gm).
236  ! For use in surface flux parameterizations. See convect.ps
237  ! file for details.
238
239  ! cbmf: The cloud base mass flux ((kg/m**2)/s). THIS SCALAR VALUE MUST
240  ! BE STORED BY THE CALLING PROGRAM AND RETURNED TO CONVECT AT
241  ! ITS NEXT CALL. That is, the value of CBMF must be "remembered"
242  ! by the calling program between calls to CONVECT.
243
244  ! det:   Array of detrainment mass flux of dimension ND.
245
246  ! -------------------------------------------------------------------
247
248  ! Local arrays
249
250  INTEGER nl
251  INTEGER nlp
252  INTEGER nlm
253  INTEGER i, k, n
254  REAL delti
255  REAL rowl
256  REAL clmcpv
257  REAL clmcpd
258  REAL cpdmcp
259  REAL cpvmcpd
260  REAL eps
261  REAL epsi
262  REAL epsim1
263  REAL ginv
264  REAL hrd
265  REAL prccon1
266  INTEGER icbmax
267  REAL lv(klon, klev)
268  REAL cpn(klon, klev)
269  REAL cpx(klon, klev)
270  REAL tv(klon, klev)
271  REAL gz(klon, klev)
272  REAL hm(klon, klev)
273  REAL h(klon, klev)
274  REAL work(klon)
275  INTEGER ihmin(klon)
276  INTEGER nk(klon)
277  REAL rh(klon)
278  REAL chi(klon)
279  REAL plcl(klon)
280  INTEGER icb(klon)
281  REAL tnk(klon)
282  REAL qnk(klon)
283  REAL gznk(klon)
284  REAL pnk(klon)
285  REAL qsnk(klon)
286  REAL ticb(klon)
287  REAL gzicb(klon)
288  REAL tp(klon, klev)
289  REAL tvp(klon, klev)
290  REAL clw(klon, klev)
291
292  REAL ah0(klon), cpp(klon)
293  REAL tg, qg, s, alv, tc, ahg, denom, es, rg
294
295  INTEGER ncum
296  INTEGER idcum(klon)
297
298  cpd = 1005.7
299  cpv = 1870.0
300  cl = 4190.0
301  rv = 461.5
302  rd = 287.04
303  lv0 = 2.501E6
304  g = 9.8
305
306  ! *** ELCRIT IS THE AUTOCONVERSION THERSHOLD WATER CONTENT (gm/gm) ***
307  ! ***  TLCRIT IS CRITICAL TEMPERATURE BELOW WHICH THE AUTO-        ***
308  ! ***       CONVERSION THRESHOLD IS ASSUMED TO BE ZERO             ***
309  ! ***     (THE AUTOCONVERSION THRESHOLD VARIES LINEARLY            ***
310  ! ***               BETWEEN 0 C AND TLCRIT)                        ***
311  ! ***   ENTP IS THE COEFFICIENT OF MIXING IN THE ENTRAINMENT       ***
312  ! ***                       FORMULATION                            ***
313  ! ***  SIGD IS THE FRACTIONAL AREA COVERED BY UNSATURATED DNDRAFT  ***
314  ! ***  SIGS IS THE FRACTION OF PRECIPITATION FALLING OUTSIDE       ***
315  ! ***                        OF CLOUD                              ***
316  ! ***        OMTRAIN IS THE ASSUMED FALL SPEED (P/s) OF RAIN       ***
317  ! ***     OMTSNOW IS THE ASSUMED FALL SPEED (P/s) OF SNOW          ***
318  ! ***  COEFFR IS A COEFFICIENT GOVERNING THE RATE OF EVAPORATION   ***
319  ! ***                          OF RAIN                             ***
320  ! ***  COEFFS IS A COEFFICIENT GOVERNING THE RATE OF EVAPORATION   ***
321  ! ***                          OF SNOW                             ***
322  ! ***     CU IS THE COEFFICIENT GOVERNING CONVECTIVE MOMENTUM      ***
323  ! ***                         TRANSPORT                            ***
324  ! ***    DTMAX IS THE MAXIMUM NEGATIVE TEMPERATURE PERTURBATION    ***
325  ! ***        A LIFTED PARCEL IS ALLOWED TO HAVE BELOW ITS LFC      ***
326  ! ***    ALPHA AND DAMP ARE PARAMETERS THAT CONTROL THE RATE OF    ***
327  ! ***                 APPROACH TO QUASI-EQUILIBRIUM                ***
328  ! ***   (THEIR STANDARD VALUES ARE  0.20 AND 0.1, RESPECTIVELY)    ***
329  ! ***                   (DAMP MUST BE LESS THAN 1)                 ***
330
331  sigs = 0.12
332  sigd = 0.05
333  elcrit = 0.0011
334  tlcrit = -55.0
335  omtsnow = 5.5
336  dtmax = 0.9
337  damp = 0.1
338  alpha = 0.2
339  entp = 1.5
340  coeffs = 0.8
341  coeffr = 1.0
342  omtrain = 50.0
343
344  cu = 0.70
345  damp = 0.1
346
347
348  ! Define nl, nlp, nlm, and delti
349
350  nl = nd - noff
351  nlp = nl + 1
352  nlm = nl - 1
353  delti = 1.0/delt
354
355  ! -------------------------------------------------------------------
356  ! --- SET CONSTANTS
357  ! -------------------------------------------------------------------
358
359  rowl = 1000.0
360  clmcpv = cl - cpv
361  clmcpd = cl - cpd
362  cpdmcp = cpd - cpv
363  cpvmcpd = cpv - cpd
364  eps = rd/rv
365  epsi = 1.0/eps
366  epsim1 = epsi - 1.0
367  ginv = 1.0/g
368  hrd = 0.5*rd
369  prccon1 = 86400.0*1000.0/(rowl*g)
370
371  ! dtmax is the maximum negative temperature perturbation.
372
373  ! =====================================================================
374  ! --- INITIALIZE OUTPUT ARRAYS AND PARAMETERS
375  ! =====================================================================
376
377  DO k = 1, nd
378    DO i = 1, len
379      ft(i, k) = 0.0
380      fq(i, k) = 0.0
381      fu(i, k) = 0.0
382      fv(i, k) = 0.0
383      tvp(i, k) = 0.0
384      tp(i, k) = 0.0
385      clw(i, k) = 0.0
386      gz(i, k) = 0.
387    END DO
388  END DO
389  DO i = 1, len
390    precip(i) = 0.0
391    iflag(i) = 0
392  END DO
393
394  ! =====================================================================
395  ! --- CALCULATE ARRAYS OF GEOPOTENTIAL, HEAT CAPACITY & STATIC ENERGY
396  ! =====================================================================
397  DO k = 1, nl + 1
398    DO i = 1, len
399      lv(i, k) = lv0 - clmcpv*(t(i,k)-273.15)
400      cpn(i, k) = cpd*(1.0-q(i,k)) + cpv*q(i, k)
401      cpx(i, k) = cpd*(1.0-q(i,k)) + cl*q(i, k)
402      tv(i, k) = t(i, k)*(1.0+q(i,k)*epsim1)
403    END DO
404  END DO
405
406  ! gz = phi at the full levels (same as p).
407
408  DO i = 1, len
409    gz(i, 1) = 0.0
410  END DO
411  DO k = 2, nlp
412    DO i = 1, len
413      gz(i, k) = gz(i, k-1) + hrd*(tv(i,k-1)+tv(i,k))*(p(i,k-1)-p(i,k))/ph(i, &
414        k)
415    END DO
416  END DO
417
418  ! h  = phi + cpT (dry static energy).
419  ! hm = phi + cp(T-Tbase)+Lq
420
421  DO k = 1, nlp
422    DO i = 1, len
423      h(i, k) = gz(i, k) + cpn(i, k)*t(i, k)
424      hm(i, k) = gz(i, k) + cpx(i, k)*(t(i,k)-t(i,1)) + lv(i, k)*q(i, k)
425    END DO
426  END DO
427
428  ! -------------------------------------------------------------------
429  ! --- Find level of minimum moist static energy
430  ! --- If level of minimum moist static energy coincides with
431  ! --- or is lower than minimum allowable parcel origin level,
432  ! --- set iflag to 6.
433  ! -------------------------------------------------------------------
434  DO i = 1, len
435    work(i) = 1.0E12
436    ihmin(i) = nl
437  END DO
438  DO k = 2, nlp
439    DO i = 1, len
440      IF ((hm(i,k)<work(i)) .AND. (hm(i,k)<hm(i,k-1))) THEN
441        work(i) = hm(i, k)
442        ihmin(i) = k
443      END IF
444    END DO
445  END DO
446  DO i = 1, len
447    ihmin(i) = min(ihmin(i), nlm)
448    IF (ihmin(i)<=minorig) THEN
449      iflag(i) = 6
450    END IF
451  END DO
452
453  ! -------------------------------------------------------------------
454  ! --- Find that model level below the level of minimum moist static
455  ! --- energy that has the maximum value of moist static energy
456  ! -------------------------------------------------------------------
457
458  DO i = 1, len
459    work(i) = hm(i, minorig)
460    nk(i) = minorig
461  END DO
462  DO k = minorig + 1, nl
463    DO i = 1, len
464      IF ((hm(i,k)>work(i)) .AND. (k<=ihmin(i))) THEN
465        work(i) = hm(i, k)
466        nk(i) = k
467      END IF
468    END DO
469  END DO
470  ! -------------------------------------------------------------------
471  ! --- Check whether parcel level temperature and specific humidity
472  ! --- are reasonable
473  ! -------------------------------------------------------------------
474  DO i = 1, len
475    IF (((t(i,nk(i))<250.0) .OR. (q(i,nk(i))<=0.0) .OR. (p(i,ihmin(i))< &
476      400.0)) .AND. (iflag(i)==0)) iflag(i) = 7
477  END DO
478  ! -------------------------------------------------------------------
479  ! --- Calculate lifted condensation level of air at parcel origin level
480  ! --- (Within 0.2% of formula of Bolton, MON. WEA. REV.,1980)
481  ! -------------------------------------------------------------------
482  DO i = 1, len
483    tnk(i) = t(i, nk(i))
484    qnk(i) = q(i, nk(i))
485    gznk(i) = gz(i, nk(i))
486    pnk(i) = p(i, nk(i))
487    qsnk(i) = qs(i, nk(i))
488
489    rh(i) = qnk(i)/qsnk(i)
490    rh(i) = min(1.0, rh(i))
491    chi(i) = tnk(i)/(1669.0-122.0*rh(i)-tnk(i))
492    plcl(i) = pnk(i)*(rh(i)**chi(i))
493    IF (((plcl(i)<200.0) .OR. (plcl(i)>=2000.0)) .AND. (iflag(i)==0)) iflag(i &
494      ) = 8
495  END DO
496  ! -------------------------------------------------------------------
497  ! --- Calculate first level above lcl (=icb)
498  ! -------------------------------------------------------------------
499  DO i = 1, len
500    icb(i) = nlm
501  END DO
502
503  DO k = minorig, nl
504    DO i = 1, len
505      IF ((k>=(nk(i)+1)) .AND. (p(i,k)<plcl(i))) icb(i) = min(icb(i), k)
506    END DO
507  END DO
508
509  DO i = 1, len
510    IF ((icb(i)>=nlm) .AND. (iflag(i)==0)) iflag(i) = 9
511  END DO
512
513  ! Compute icbmax.
514
515  icbmax = 2
516  DO i = 1, len
517    icbmax = max(icbmax, icb(i))
518  END DO
519
520  ! -------------------------------------------------------------------
521  ! --- Calculates the lifted parcel virtual temperature at nk,
522  ! --- the actual temperature, and the adiabatic
523  ! --- liquid water content. The procedure is to solve the equation.
524  ! cp*tp+L*qp+phi=cp*tnk+L*qnk+gznk.
525  ! -------------------------------------------------------------------
526
527  DO i = 1, len
528    tnk(i) = t(i, nk(i))
529    qnk(i) = q(i, nk(i))
530    gznk(i) = gz(i, nk(i))
531    ticb(i) = t(i, icb(i))
532    gzicb(i) = gz(i, icb(i))
533  END DO
534
535  ! ***  Calculate certain parcel quantities, including static energy   ***
536
537  DO i = 1, len
538    ah0(i) = (cpd*(1.-qnk(i))+cl*qnk(i))*tnk(i) + qnk(i)*(lv0-clmcpv*(tnk(i)- &
539      273.15)) + gznk(i)
540    cpp(i) = cpd*(1.-qnk(i)) + qnk(i)*cpv
541  END DO
542
543  ! ***   Calculate lifted parcel quantities below cloud base   ***
544
545  DO k = minorig, icbmax - 1
546    DO i = 1, len
547      tp(i, k) = tnk(i) - (gz(i,k)-gznk(i))/cpp(i)
548      tvp(i, k) = tp(i, k)*(1.+qnk(i)*epsi)
549    END DO
550  END DO
551
552  ! ***  Find lifted parcel quantities above cloud base    ***
553
554  DO i = 1, len
555    tg = ticb(i)
556    qg = qs(i, icb(i))
557    alv = lv0 - clmcpv*(ticb(i)-273.15)
558
559    ! First iteration.
560
561    s = cpd + alv*alv*qg/(rv*ticb(i)*ticb(i))
562    s = 1./s
563    ahg = cpd*tg + (cl-cpd)*qnk(i)*ticb(i) + alv*qg + gzicb(i)
564    tg = tg + s*(ah0(i)-ahg)
565    tg = max(tg, 35.0)
566    tc = tg - 273.15
567    denom = 243.5 + tc
568    IF (tc>=0.0) THEN
569      es = 6.112*exp(17.67*tc/denom)
570    ELSE
571      es = exp(23.33086-6111.72784/tg+0.15215*log(tg))
572    END IF
573    qg = eps*es/(p(i,icb(i))-es*(1.-eps))
574
575    ! Second iteration.
576
577    s = cpd + alv*alv*qg/(rv*ticb(i)*ticb(i))
578    s = 1./s
579    ahg = cpd*tg + (cl-cpd)*qnk(i)*ticb(i) + alv*qg + gzicb(i)
580    tg = tg + s*(ah0(i)-ahg)
581    tg = max(tg, 35.0)
582    tc = tg - 273.15
583    denom = 243.5 + tc
584    IF (tc>=0.0) THEN
585      es = 6.112*exp(17.67*tc/denom)
586    ELSE
587      es = exp(23.33086-6111.72784/tg+0.15215*log(tg))
588    END IF
589    qg = eps*es/(p(i,icb(i))-es*(1.-eps))
590
591    alv = lv0 - clmcpv*(ticb(i)-273.15)
592    tp(i, icb(i)) = (ah0(i)-(cl-cpd)*qnk(i)*ticb(i)-gz(i,icb(i))-alv*qg)/cpd
593    clw(i, icb(i)) = qnk(i) - qg
594    clw(i, icb(i)) = max(0.0, clw(i,icb(i)))
595    rg = qg/(1.-qnk(i))
596    tvp(i, icb(i)) = tp(i, icb(i))*(1.+rg*epsi)
597  END DO
598
599  DO k = minorig, icbmax
600    DO i = 1, len
601      tvp(i, k) = tvp(i, k) - tp(i, k)*qnk(i)
602    END DO
603  END DO
604
605  ! -------------------------------------------------------------------
606  ! --- Test for instability.
607  ! --- If there was no convection at last time step and parcel
608  ! --- is stable at icb, then set iflag to 4.
609  ! -------------------------------------------------------------------
610
611  DO i = 1, len
612    IF ((cbmf(i)==0.0) .AND. (iflag(i)==0) .AND. (tvp(i, &
613      icb(i))<=(tv(i,icb(i))-dtmax))) iflag(i) = 4
614  END DO
615
616  ! =====================================================================
617  ! --- IF THIS POINT IS REACHED, MOIST CONVECTIVE ADJUSTMENT IS NECESSARY
618  ! =====================================================================
619
620  ncum = 0
621  DO i = 1, len
622    IF (iflag(i)==0) THEN
623      ncum = ncum + 1
624      idcum(ncum) = i
625    END IF
626  END DO
627
628  ! Call convect2, which compresses the points and computes the heating,
629  ! moistening, velocity mixing, and precipiation.
630
631  ! print*,'cpd avant convect2 ',cpd
632  IF (ncum>0) THEN
633    CALL convect2(ncum, idcum, len, nd, ndp1, nl, minorig, nk, icb, t, q, qs, &
634      u, v, gz, tv, tp, tvp, clw, h, lv, cpn, p, ph, ft, fq, fu, fv, tnk, &
635      qnk, gznk, plcl, precip, cbmf, iflag, delt, cpd, cpv, cl, rv, rd, lv0, &
636      g, sigs, sigd, elcrit, tlcrit, omtsnow, dtmax, damp, alpha, entp, &
637      coeffs, coeffr, omtrain, cu, ma)
638  END IF
639
640  RETURN
641END SUBROUTINE convect1
Note: See TracBrowser for help on using the repository browser.