source: LMDZ5/trunk/libf/phylmd/cva_driver.F90 @ 2079

Last change on this file since 2079 was 2079, checked in by lguez, 10 years ago

Protect against division by 0 in computation of proba_notrig in procedure physiq: initialize s2 at some low value instead of 0.

  • 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: 38.8 KB
Line 
1
2! $Id: cva_driver.F90 2079 2014-07-07 13:38:36Z lguez $
3
4SUBROUTINE cva_driver(len, nd, ndp1, ntra, nloc, &
5                      iflag_con, iflag_mix, iflag_ice_thermo, iflag_clos, ok_conserv_q, &
6                      delt, t1, q1, qs1, t1_wake, q1_wake, qs1_wake, s1_wake, &
7                      u1, v1, tra1, &
8                      p1, ph1, &
9                      Ale1, Alp1, &
10                      sig1feed1, sig2feed1, wght1, &
11                      iflag1, ft1, fq1, fu1, fv1, ftra1, &
12                      precip1, kbas1, ktop1, &
13                      cbmf1, plcl1, plfc1, wbeff1, &
14                      sig1, w01, & !input/output
15                      ptop21, sigd1, &
16                      ma1, mip1, Vprecip1, upwd1, dnwd1, dnwd01, &
17                      qcondc1, wd1, &
18                      cape1, cin1, tvp1, &
19                      ftd1, fqd1, &
20                      Plim11, Plim21, asupmax1, supmax01, asupmaxmin1, &
21                      lalim_conv, &
22!!                      da1,phi1,mp1,phi21,d1a1,dam1,sigij1,clw1, &        ! RomP
23!!                      elij1,evap1,ep1,epmlmMm1,eplaMm1, &                ! RomP
24                      da1, phi1, mp1, phi21, d1a1, dam1, sigij1, wghti1, & ! RomP, RL
25                      clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, &        ! RomP, RL
26                      wdtrainA1, wdtrainM1)                                ! RomP
27! **************************************************************
28! *
29! CV_DRIVER                                                   *
30! *
31! *
32! written by   : Sandrine Bony-Lena , 17/05/2003, 11.19.41    *
33! modified by :                                               *
34! **************************************************************
35! **************************************************************
36
37  USE dimphy
38  IMPLICIT NONE
39
40! .............................START PROLOGUE............................
41
42
43! All argument names (except len,nd,ntra,nloc,delt and the flags) have a "1" appended.
44! The "1" is removed for the corresponding compressed variables.
45! PARAMETERS:
46! Name            Type         Usage            Description
47! ----------      ----------     -------  ----------------------------
48
49! len           Integer        Input        first (i) dimension
50! nd            Integer        Input        vertical (k) dimension
51! ndp1          Integer        Input        nd + 1
52! ntra          Integer        Input        number of tracors
53! iflag_con     Integer        Input        version of convect (3/4)
54! iflag_mix     Integer        Input        version of mixing  (0/1/2)
55! iflag_ice_thermo Integer        Input        accounting for ice thermodynamics (0/1)
56! iflag_clos    Integer        Input        version of closure (0/1)
57! ok_conserv_q  Logical        Input        when true corrections for water conservation are swtiched on
58! delt          Real           Input        time step
59! t1            Real           Input        temperature (sat draught envt)
60! q1            Real           Input        specific hum (sat draught envt)
61! qs1           Real           Input        sat specific hum (sat draught envt)
62! t1_wake       Real           Input        temperature (unsat draught envt)
63! q1_wake       Real           Input        specific hum(unsat draught envt)
64! qs1_wake      Real           Input        sat specific hum(unsat draughts envt)
65! s1_wake       Real           Input        fractionnal area covered by wakes
66! u1            Real           Input        u-wind
67! v1            Real           Input        v-wind
68! tra1          Real           Input        tracors
69! p1            Real           Input        full level pressure
70! ph1           Real           Input        half level pressure
71! ALE1          Real           Input        Available lifting Energy
72! ALP1          Real           Input        Available lifting Power
73! sig1feed1     Real           Input        sigma coord at lower bound of feeding layer
74! sig2feed1     Real           Input        sigma coord at upper bound of feeding layer
75! wght1         Real           Input        weight density determining the feeding mixture
76! iflag1        Integer        Output       flag for Emanuel conditions
77! ft1           Real           Output       temp tend
78! fq1           Real           Output       spec hum tend
79! fu1           Real           Output       u-wind tend
80! fv1           Real           Output       v-wind tend
81! ftra1         Real           Output       tracor tend
82! precip1       Real           Output       precipitation
83! kbas1         Integer        Output       cloud base level
84! ktop1         Integer        Output       cloud top level
85! cbmf1         Real           Output       cloud base mass flux
86! sig1          Real           In/Out       section adiabatic updraft
87! w01           Real           In/Out       vertical velocity within adiab updraft
88! ptop21        Real           In/Out       top of entraining zone
89! Ma1           Real           Output       mass flux adiabatic updraft
90! mip1          Real           Output       mass flux shed by the adiabatic updraft
91! Vprecip1      Real           Output       vertical profile of precipitations
92! upwd1         Real           Output       total upward mass flux (adiab+mixed)
93! dnwd1         Real           Output       saturated downward mass flux (mixed)
94! dnwd01        Real           Output       unsaturated downward mass flux
95! qcondc1       Real           Output       in-cld mixing ratio of condensed water
96! wd1           Real           Output       downdraft velocity scale for sfc fluxes
97! cape1         Real           Output       CAPE
98! cin1          Real           Output       CIN
99! tvp1          Real           Output       adiab lifted parcell virt temp
100! ftd1          Real           Output       precip temp tend
101! fqt1          Real           Output       precip spec hum tend
102! Plim11        Real           Output
103! Plim21        Real           Output
104! asupmax1      Real           Output
105! supmax01      Real           Output
106! asupmaxmin1   Real           Output
107
108! ftd1          Real           Output  Array of temperature tendency due to precipitations (K/s) of dimension ND,
109!                                      defined at same grid levels as T, Q, QS and P.
110
111! fqd1          Real           Output  Array of specific humidity tendencies due to precipitations ((gm/gm)/s)
112!                                      of dimension ND, defined at same grid levels as T, Q, QS and P.
113
114! wdtrainA1     Real           Output   precipitation detrained from adiabatic draught;
115!                                         used in tracer transport (cvltr)
116! wdtrainM1     Real           Output   precipitation detrained from mixed draughts;
117!                                         used in tracer transport (cvltr)
118! da1           Real           Output     used in tracer transport (cvltr)
119! phi1          Real           Output     used in tracer transport (cvltr)
120! mp1           Real           Output     used in tracer transport (cvltr)
121                                         
122! phi21         Real           Output     used in tracer transport (cvltr)
123                                         
124! d1a1          Real           Output     used in tracer transport (cvltr)
125! dam1          Real           Output     used in tracer transport (cvltr)
126                                         
127! epmlmMm1      Real           Output     used in tracer transport (cvltr)
128! eplaMm1       Real           Output     used in tracer transport (cvltr)
129                                         
130! evap1         Real           Output   
131! ep1           Real           Output   
132! sigij1        Real           Output     used in tracer transport (cvltr)
133! elij1         Real           Output
134! wghti1        Real           Output   final weight of the feeding layers,
135!                                         used in tracer transport (cvltr)
136
137
138! S. Bony, Mar 2002:
139! * Several modules corresponding to different physical processes
140! * Several versions of convect may be used:
141!         - iflag_con=3: version lmd  (previously named convect3)
142!         - iflag_con=4: version 4.3b (vect. version, previously convect1/2)
143! + tard: - iflag_con=5: version lmd with ice (previously named convectg)
144! S. Bony, Oct 2002:
145! * Vectorization of convect3 (ie version lmd)
146
147! ..............................END PROLOGUE.............................
148
149
150  include "dimensions.h"
151!!!!!#include "dimphy.h"
152  include 'iniprint.h'
153
154
155! Input
156  INTEGER len
157  INTEGER nd
158  INTEGER ndp1
159  INTEGER ntra
160  INTEGER iflag_con
161  INTEGER iflag_mix
162  INTEGER iflag_ice_thermo
163  INTEGER iflag_clos
164  LOGICAL ok_conserv_q
165  REAL delt
166  REAL t1(len, nd)
167  REAL q1(len, nd)
168  REAL qs1(len, nd)
169  REAL t1_wake(len, nd)
170  REAL q1_wake(len, nd)
171  REAL qs1_wake(len, nd)
172  REAL s1_wake(len)
173  REAL u1(len, nd)
174  REAL v1(len, nd)
175  REAL tra1(len, nd, ntra)
176  REAL p1(len, nd)
177  REAL ph1(len, ndp1)
178  REAL Ale1(len)
179  REAL Alp1(len)
180  REAL sig1feed1 ! pressure at lower bound of feeding layer
181  REAL sig2feed1 ! pressure at upper bound of feeding layer
182  REAL wght1(nd) ! weight density determining the feeding mixture
183
184! Output
185  INTEGER iflag1(len)
186  REAL ft1(len, nd)
187  REAL fq1(len, nd)
188  REAL fu1(len, nd)
189  REAL fv1(len, nd)
190  REAL ftra1(len, nd, ntra)
191  REAL precip1(len)
192  INTEGER kbas1(len)
193  INTEGER ktop1(len)
194  REAL cbmf1(len)
195  REAL plcl1(klon)
196  REAL plfc1(klon)
197  REAL wbeff1(klon)
198  REAL sig1(len, klev) !input/output
199  REAL w01(len, klev) !input/output
200  REAL ptop21(len)
201  REAL sigd1(len)
202  REAL ma1(len, nd)
203  REAL mip1(len, nd)
204! real Vprecip1(len,nd)
205  REAL vprecip1(len, nd+1)
206  REAL upwd1(len, nd)
207  REAL dnwd1(len, nd)
208  REAL dnwd01(len, nd)
209  REAL qcondc1(len, nd) ! cld
210  REAL wd1(len) ! gust
211  REAL cape1(len)
212  REAL cin1(len)
213  REAL tvp1(len, nd)
214
215!AC!
216!!      real da1(len,nd),phi1(len,nd,nd)
217!!      real da(len,nd),phi(len,nd,nd)
218!AC!
219  REAL ftd1(len, nd)
220  REAL fqd1(len, nd)
221  REAL Plim11(len)
222  REAL Plim21(len)
223  REAL asupmax1(len, nd)
224  REAL supmax01(len)
225  REAL asupmaxmin1(len)
226  INTEGER lalim_conv(len)
227! RomP >>>
228  REAL wdtrainA1(len, nd), wdtrainM1(len, nd)
229  REAL da1(len, nd), phi1(len, nd, nd), mp1(len, nd)
230  REAL epmlmMm1(len, nd, nd), eplaMm1(len, nd)
231  REAL evap1(len, nd), ep1(len, nd)
232  REAL sigij1(len, nd, nd), elij1(len, nd, nd)
233!JYG,RL
234  REAL wghti1(len, nd) ! final weight of the feeding layers
235!JYG,RL
236  REAL phi21(len, nd, nd)
237  REAL d1a1(len, nd), dam1(len, nd)
238! RomP <<<
239
240! -------------------------------------------------------------------
241! Prolog by Kerry Emanuel.
242! -------------------------------------------------------------------
243! --- ARGUMENTS
244! -------------------------------------------------------------------
245! --- On input:
246
247! t:   Array of absolute temperature (K) of dimension ND, with first
248! index corresponding to lowest model level. Note that this array
249! will be altered by the subroutine if dry convective adjustment
250! occurs and if IPBL is not equal to 0.
251
252! q:   Array of specific humidity (gm/gm) of dimension ND, with first
253! index corresponding to lowest model level. Must be defined
254! at same grid levels as T. Note that this array will be altered
255! if dry convective adjustment occurs and if IPBL is not equal to 0.
256
257! qs:  Array of saturation specific humidity of dimension ND, with first
258! index corresponding to lowest model level. Must be defined
259! at same grid levels as T. Note that this array will be altered
260! if dry convective adjustment occurs and if IPBL is not equal to 0.
261
262! t_wake: Array of absolute temperature (K), seen by unsaturated draughts,
263! of dimension ND, with first index corresponding to lowest model level.
264
265! q_wake: Array of specific humidity (gm/gm), seen by unsaturated draughts,
266! of dimension ND, with first index corresponding to lowest model level.
267! Must be defined at same grid levels as T.
268
269! qs_wake: Array of saturation specific humidity, seen by unsaturated draughts,
270! of dimension ND, with first index corresponding to lowest model level.
271! Must be defined at same grid levels as T.
272
273! s_wake: Array of fractionnal area occupied by the wakes.
274
275! u:   Array of zonal wind velocity (m/s) of dimension ND, witth first
276! index corresponding with the lowest model level. Defined at
277! same levels as T. Note that this array will be altered if
278! dry convective adjustment occurs and if IPBL is not equal to 0.
279
280! v:   Same as u but for meridional velocity.
281
282! tra: Array of passive tracer mixing ratio, of dimensions (ND,NTRA),
283! where NTRA is the number of different tracers. If no
284! convective tracer transport is needed, define a dummy
285! input array of dimension (ND,1). Tracers are defined at
286! same vertical levels as T. Note that this array will be altered
287! if dry convective adjustment occurs and if IPBL is not equal to 0.
288
289! p:   Array of pressure (mb) of dimension ND, with first
290! index corresponding to lowest model level. Must be defined
291! at same grid levels as T.
292
293! ph:  Array of pressure (mb) of dimension ND+1, with first index
294! corresponding to lowest level. These pressures are defined at
295! levels intermediate between those of P, T, Q and QS. The first
296! value of PH should be greater than (i.e. at a lower level than)
297! the first value of the array P.
298
299! ALE:  Available lifting Energy
300
301! ALP:  Available lifting Power
302
303! nl:  The maximum number of levels to which convection can penetrate, plus 1.
304!       NL MUST be less than or equal to ND-1.
305
306! delt: The model time step (sec) between calls to CONVECT
307
308! ----------------------------------------------------------------------------
309! ---   On Output:
310
311! iflag: An output integer whose value denotes the following:
312!       VALUE   INTERPRETATION
313!       -----   --------------
314!         0     Moist convection occurs.
315!         1     Moist convection occurs, but a CFL condition
316!               on the subsidence warming is violated. This
317!               does not cause the scheme to terminate.
318!         2     Moist convection, but no precip because ep(inb) lt 0.0001
319!         3     No moist convection because new cbmf is 0 and old cbmf is 0.
320!         4     No moist convection; atmosphere is not
321!               unstable
322!         6     No moist convection because ihmin le minorig.
323!         7     No moist convection because unreasonable
324!               parcel level temperature or specific humidity.
325!         8     No moist convection: lifted condensation
326!               level is above the 200 mb level.
327!         9     No moist convection: cloud base is higher
328!               then the level NL-1.
329
330! ft:   Array of temperature tendency (K/s) of dimension ND, defined at same
331!       grid levels as T, Q, QS and P.
332
333! fq:   Array of specific humidity tendencies ((gm/gm)/s) of dimension ND,
334!       defined at same grid levels as T, Q, QS and P.
335
336! fu:   Array of forcing of zonal velocity (m/s^2) of dimension ND,
337!      defined at same grid levels as T.
338
339! fv:   Same as FU, but for forcing of meridional velocity.
340
341! ftra: Array of forcing of tracer content, in tracer mixing ratio per
342!       second, defined at same levels as T. Dimensioned (ND,NTRA).
343
344! precip: Scalar convective precipitation rate (mm/day).
345
346! wd:   A convective downdraft velocity scale. For use in surface
347!       flux parameterizations. See convect.ps file for details.
348
349! tprime: A convective downdraft temperature perturbation scale (K).
350!         For use in surface flux parameterizations. See convect.ps
351!         file for details.
352
353! qprime: A convective downdraft specific humidity
354!         perturbation scale (gm/gm).
355!         For use in surface flux parameterizations. See convect.ps
356!         file for details.
357
358! cbmf: The cloud base mass flux ((kg/m**2)/s). THIS SCALAR VALUE MUST
359!       BE STORED BY THE CALLING PROGRAM AND RETURNED TO CONVECT AT
360!       ITS NEXT CALL. That is, the value of CBMF must be "remembered"
361!       by the calling program between calls to CONVECT.
362
363! det:   Array of detrainment mass flux of dimension ND.
364! -------------------------------------------------------------------
365
366! Local (non compressed) arrays
367
368
369  INTEGER i, k, n, il, j
370  INTEGER nword1, nword2, nword3, nword4
371  INTEGER icbmax
372  INTEGER nk1(klon)
373  INTEGER icb1(klon)
374  INTEGER icbs1(klon)
375
376  LOGICAL ok_inhib ! True => possible inhibition of convection by dryness
377  LOGICAL, SAVE :: debut = .TRUE.
378!$OMP THREADPRIVATE(debut)
379
380  REAL tnk1(klon)
381  REAL thnk1(klon)
382  REAL qnk1(klon)
383  REAL gznk1(klon)
384  REAL pnk1(klon)
385  REAL qsnk1(klon)
386  REAL unk1(klon)
387  REAL vnk1(klon)
388  REAL cpnk1(klon)
389  REAL hnk1(klon)
390  REAL pbase1(klon)
391  REAL buoybase1(klon)
392
393  REAL lf1(klon, klev), lf1_wake(klon, klev)
394  REAL lv1(klon, klev), lv1_wake(klon, klev)
395  REAL cpn1(klon, klev), cpn1_wake(klon, klev)
396  REAL tv1(klon, klev), tv1_wake(klon, klev)
397  REAL gz1(klon, klev), gz1_wake(klon, klev)
398  REAL hm1(klon, klev), hm1_wake(klon, klev)
399  REAL h1(klon, klev), h1_wake(klon, klev)
400  REAL tp1(klon, klev)
401  REAL clw1(klon, klev)
402  REAL th1(klon, klev), th1_wake(klon, klev)
403
404  REAL bid(klon, klev) ! dummy array
405
406  INTEGER ncum
407
408  INTEGER j1feed(klon)
409  INTEGER j2feed(klon)
410  REAL p1feed1(len) ! pressure at lower bound of feeding layer
411  REAL p2feed1(len) ! pressure at upper bound of feeding layer
412!JYG,RL
413!!      real wghti1(len,nd) ! weights of the feeding layers
414!JYG,RL
415
416! (local) compressed fields:
417
418  INTEGER nloc
419! parameter (nloc=klon) ! pour l'instant
420
421  INTEGER idcum(nloc)
422  INTEGER iflag(nloc), nk(nloc), icb(nloc)
423  INTEGER nent(nloc, klev)
424  INTEGER icbs(nloc)
425  INTEGER inb(nloc), inbis(nloc)
426
427  REAL cbmf(nloc), plcl(nloc), plfc(nloc), wbeff(nloc)
428  REAL t(nloc, klev), q(nloc, klev), qs(nloc, klev)
429  REAL t_wake(nloc, klev), q_wake(nloc, klev), qs_wake(nloc, klev)
430  REAL s_wake(nloc)
431  REAL u(nloc, klev), v(nloc, klev)
432  REAL gz(nloc, klev), h(nloc, klev), hm(nloc, klev)
433  REAL h_wake(nloc, klev), hm_wake(nloc, klev)
434  REAL lv(nloc, klev), lf(nloc, klev), cpn(nloc, klev)
435  REAL lv_wake(nloc, klev), lf_wake(nloc, klev), cpn_wake(nloc, klev)
436  REAL p(nloc, klev), ph(nloc, klev+1), tv(nloc, klev), tp(nloc, klev)
437  REAL tv_wake(nloc, klev)
438  REAL clw(nloc, klev)
439  REAL dph(nloc, klev)
440  REAL pbase(nloc), buoybase(nloc), th(nloc, klev)
441  REAL th_wake(nloc, klev)
442  REAL tvp(nloc, klev)
443  REAL sig(nloc, klev), w0(nloc, klev), ptop2(nloc)
444  REAL hp(nloc, klev), ep(nloc, klev), sigp(nloc, klev)
445  REAL buoy(nloc, klev)
446  REAL cape(nloc)
447  REAL cin(nloc)
448  REAL m(nloc, klev)
449  REAL ment(nloc, klev, klev), sigij(nloc, klev, klev)
450  REAL qent(nloc, klev, klev)
451  REAL hent(nloc, klev, klev)
452  REAL uent(nloc, klev, klev), vent(nloc, klev, klev)
453  REAL ments(nloc, klev, klev), qents(nloc, klev, klev)
454  REAL elij(nloc, klev, klev)
455  REAL supmax(nloc, klev)
456  REAL Ale(nloc), Alp(nloc), coef_clos(nloc)
457  REAL sigd(nloc)
458! real mp(nloc,klev), qp(nloc,klev), up(nloc,klev), vp(nloc,klev)
459! real wt(nloc,klev), water(nloc,klev), evap(nloc,klev), ice(nloc,klev)
460! real b(nloc,klev), sigd(nloc)
461! save mp,qp,up,vp,wt,water,evap,b
462  REAL, SAVE, ALLOCATABLE :: mp(:, :), qp(:, :), up(:, :), vp(:, :)
463  REAL, SAVE, ALLOCATABLE :: wt(:, :), water(:, :), evap(:, :)
464  REAL, SAVE, ALLOCATABLE :: ice(:, :), fondue(:, :), b(:, :)
465  REAL, SAVE, ALLOCATABLE :: frac(:, :), faci(:, :)
466!$OMP THREADPRIVATE(mp,qp,up,vp,wt,water,evap,ice,fondue,b,frac,faci)
467  REAL ft(nloc, klev), fq(nloc, klev)
468  REAL ftd(nloc, klev), fqd(nloc, klev)
469  REAL fu(nloc, klev), fv(nloc, klev)
470  REAL upwd(nloc, klev), dnwd(nloc, klev), dnwd0(nloc, klev)
471  REAL ma(nloc, klev), mip(nloc, klev), tls(nloc, klev)
472  REAL tps(nloc, klev), qprime(nloc), tprime(nloc)
473  REAL precip(nloc)
474! real Vprecip(nloc,klev)
475  REAL vprecip(nloc, klev+1)
476  REAL tra(nloc, klev, ntra), trap(nloc, klev, ntra)
477  REAL ftra(nloc, klev, ntra), traent(nloc, klev, klev, ntra)
478  REAL qcondc(nloc, klev)      ! cld
479  REAL wd(nloc)                ! gust
480  REAL Plim1(nloc), plim2(nloc)
481  REAL asupmax(nloc, klev)
482  REAL supmax0(nloc)
483  REAL asupmaxmin(nloc)
484
485  REAL tnk(nloc), qnk(nloc), gznk(nloc)
486  REAL wghti(nloc, nd)
487  REAL hnk(nloc), unk(nloc), vnk(nloc)
488
489! RomP >>>
490  REAL wdtrainA(nloc, klev), wdtrainM(nloc, klev)
491  REAL da(len, nd), phi(len, nd, nd)
492  REAL epmlmMm(nloc, klev, klev), eplaMm(nloc, klev)
493  REAL phi2(len, nd, nd)
494  REAL d1a(len, nd), dam(len, nd)
495! RomP <<<
496
497  LOGICAL, SAVE :: first = .TRUE.
498!$OMP THREADPRIVATE(first)
499  CHARACTER (LEN=20) :: modname = 'cva_driver'
500  CHARACTER (LEN=80) :: abort_message
501
502
503! print *, 't1, t1_wake ',(k,t1(1,k),t1_wake(1,k),k=1,klev)
504! print *, 'q1, q1_wake ',(k,q1(1,k),q1_wake(1,k),k=1,klev)
505
506! -------------------------------------------------------------------
507! --- SET CONSTANTS AND PARAMETERS
508! -------------------------------------------------------------------
509
510  IF (first) THEN
511    ALLOCATE (mp(nloc,klev), qp(nloc,klev), up(nloc,klev))
512    ALLOCATE (vp(nloc,klev), wt(nloc,klev), water(nloc,klev))
513    ALLOCATE (ice(nloc,klev), fondue(nloc,klev))
514    ALLOCATE (evap(nloc,klev), b(nloc,klev))
515    ALLOCATE (frac(nloc,klev), faci(nloc,klev))
516    first = .FALSE.
517  END IF
518! -- set simulation flags:
519! (common cvflag)
520
521  CALL cv_flag(iflag_ice_thermo)
522
523! -- set thermodynamical constants:
524! (common cvthermo)
525
526  CALL cv_thermo(iflag_con)
527
528! -- set convect parameters
529
530! includes microphysical parameters and parameters that
531! control the rate of approach to quasi-equilibrium)
532! (common cvparam)
533
534  IF (iflag_con==3) THEN
535    CALL cv3_param(nd, delt)
536
537  END IF
538
539  IF (iflag_con==4) THEN
540    CALL cv_param(nd)
541  END IF
542
543! ---------------------------------------------------------------------
544! --- INITIALIZE OUTPUT ARRAYS AND PARAMETERS
545! ---------------------------------------------------------------------
546  nword1 = len
547  nword2 = len*nd
548  nword3 = len*nd*ntra
549  nword4 = len*nd*nd
550
551  iflag1(:) = 0
552  ktop1(:) = 0
553  kbas1(:) = 0
554  ft1(:, :) = 0.0
555  fq1(:, :) = 0.0
556  fu1(:, :) = 0.0
557  fv1(:, :) = 0.0
558  ftra1(:, :, :) = 0.
559  precip1(:) = 0.
560  cbmf1(:) = 0.
561  ptop21(:) = 0.
562  sigd1(:) = 0.
563  ma1(:, :) = 0.
564  mip1(:, :) = 0.
565  vprecip1(:, :) = 0.
566  upwd1(:, :) = 0.
567  dnwd1(:, :) = 0.
568  dnwd01(:, :) = 0.
569  qcondc1(:, :) = 0.
570  wd1(:) = 0.
571  cape1(:) = 0.
572  cin1(:) = 0.
573  tvp1(:, :) = 0.
574  ftd1(:, :) = 0.
575  fqd1(:, :) = 0.
576  Plim11(:) = 0.
577  Plim21(:) = 0.
578  asupmax1(:, :) = 0.
579  supmax01(:) = 0.
580  asupmaxmin1(:) = 0.
581
582  DO il = 1, len
583    cin1(il) = -100000.
584    cape1(il) = -1.
585  END DO
586
587  IF (iflag_con==3) THEN
588    DO il = 1, len
589      sig1(il, nd) = sig1(il, nd) + 1.
590      sig1(il, nd) = amin1(sig1(il,nd), 12.1)
591    END DO
592  END IF
593
594! RomP >>>
595  wdtrainA1(:, :) = 0.
596  wdtrainM1(:, :) = 0.
597  da1(:, :) = 0.
598  phi1(:, :, :) = 0.
599  epmlmMm1(:, :, :) = 0.
600  eplaMm1(:, :) = 0.
601  mp1(:, :) = 0.
602  evap1(:, :) = 0.
603  ep1(:, :) = 0.
604  sigij1(:, :, :) = 0.
605  elij1(:, :, :) = 0.
606  phi21(:, :, :) = 0.
607  d1a1(:, :) = 0.
608  dam1(:, :) = 0.
609! RomP <<<
610! ---------------------------------------------------------------------
611! --- INITIALIZE LOCAL ARRAYS AND PARAMETERS
612! ---------------------------------------------------------------------
613
614  DO il = 1, nloc
615    coef_clos(il) = 1.
616  END DO
617
618! --------------------------------------------------------------------
619! --- CALCULATE ARRAYS OF GEOPOTENTIAL, HEAT CAPACITY & STATIC ENERGY
620! --------------------------------------------------------------------
621
622  IF (iflag_con==3) THEN
623
624    IF (debut) THEN
625      PRINT *, 'Emanuel version 3 nouvelle'
626    END IF
627! print*,'t1, q1 ',t1,q1
628    CALL cv3_prelim(len, nd, ndp1, t1, q1, p1, ph1, &           ! nd->na
629                    lv1, lf1, cpn1, tv1, gz1, h1, hm1, th1)
630
631
632    CALL cv3_prelim(len, nd, ndp1, t1_wake, q1_wake, p1, ph1, & ! nd->na
633                    lv1_wake, lf1_wake, cpn1_wake, tv1_wake, gz1_wake, &
634                    h1_wake, bid, th1_wake)
635
636  END IF
637
638  IF (iflag_con==4) THEN
639    PRINT *, 'Emanuel version 4 '
640    CALL cv_prelim(len, nd, ndp1, t1, q1, p1, ph1, &
641                   lv1, cpn1, tv1, gz1, h1, hm1)
642  END IF
643
644! --------------------------------------------------------------------
645! --- CONVECTIVE FEED
646! --------------------------------------------------------------------
647
648! compute feeding layer potential temperature and mixing ratio :
649
650! get bounds of feeding layer
651
652! test niveaux couche alimentation KE
653  IF (sig1feed1==sig2feed1) THEN
654    WRITE (lunout, *) 'impossible de choisir sig1feed=sig2feed'
655    WRITE (lunout, *) 'changer la valeur de sig2feed dans physiq.def'
656    abort_message = ''
657    CALL abort_gcm(modname, abort_message, 1)
658  END IF
659
660  DO i = 1, len
661    p1feed1(i) = sig1feed1*ph1(i, 1)
662    p2feed1(i) = sig2feed1*ph1(i, 1)
663!test maf
664!   p1feed1(i)=ph1(i,1)
665!   p2feed1(i)=ph1(i,2)
666!   p2feed1(i)=ph1(i,3)
667!testCR: on prend la couche alim des thermiques
668!   p2feed1(i)=ph1(i,lalim_conv(i)+1)
669!   print*,'lentr=',lentr(i),ph1(i,lentr(i)+1),ph1(i,2)
670  END DO
671
672  IF (iflag_con==3) THEN
673  END IF
674  DO i = 1, len
675! print*,'avant cv3_feed Plim',p1feed1(i),p2feed1(i)
676  END DO
677  IF (iflag_con==3) THEN
678
679! print*, 'IFLAG1 avant cv3_feed'
680! print*,'len,nd',len,nd
681! write(*,'(64i1)') iflag1(2:klon-1)
682
683    CALL cv3_feed(len, nd, ok_conserv_q, &                 ! nd->na
684                  t1, q1, u1, v1, p1, ph1, hm1, gz1, &
685                  p1feed1, p2feed1, wght1, &
686                  wghti1, tnk1, thnk1, qnk1, qsnk1, unk1, vnk1, &
687                  cpnk1, hnk1, nk1, icb1, icbmax, iflag1, gznk1, plcl1)
688  END IF
689
690! print*, 'IFLAG1 apres cv3_feed'
691! print*,'len,nd',len,nd
692! write(*,'(64i1)') iflag1(2:klon-1)
693
694  IF (iflag_con==4) THEN
695    CALL cv_feed(len, nd, t1, q1, qs1, p1, hm1, gz1, &
696                 nk1, icb1, icbmax, iflag1, tnk1, qnk1, gznk1, plcl1)
697  END IF
698
699! print *, 'cv3_feed-> iflag1, plcl1 ',iflag1(1),plcl1(1)
700
701! --------------------------------------------------------------------
702! --- UNDILUTE (ADIABATIC) UPDRAFT / 1st part
703! (up through ICB for convect4, up through ICB+1 for convect3)
704! Calculates the lifted parcel virtual temperature at nk, the
705! actual temperature, and the adiabatic liquid water content.
706! --------------------------------------------------------------------
707
708  IF (iflag_con==3) THEN
709
710    CALL cv3_undilute1(len, nd, t1, qs1, gz1, plcl1, p1, icb1, tnk1, qnk1, & ! nd->na
711                       gznk1, tp1, tvp1, clw1, icbs1)
712  END IF
713
714
715  IF (iflag_con==4) THEN
716    CALL cv_undilute1(len, nd, t1, q1, qs1, gz1, p1, nk1, icb1, icbmax, &
717                      tp1, tvp1, clw1)
718  END IF
719
720! -------------------------------------------------------------------
721! --- TRIGGERING
722! -------------------------------------------------------------------
723
724! print *,' avant triggering, iflag_con ',iflag_con
725
726  IF (iflag_con==3) THEN
727
728    CALL cv3_trigger(len, nd, icb1, plcl1, p1, th1, tv1, tvp1, thnk1, & ! nd->na
729                      pbase1, buoybase1, iflag1, sig1, w01)
730
731
732! print*, 'IFLAG1 apres cv3_triger'
733! print*,'len,nd',len,nd
734! write(*,'(64i1)') iflag1(2:klon-1)
735
736! call dump2d(iim,jjm-1,sig1(2)
737  END IF
738
739  IF (iflag_con==4) THEN
740    CALL cv_trigger(len, nd, icb1, cbmf1, tv1, tvp1, iflag1)
741  END IF
742
743
744! =====================================================================
745! --- IF THIS POINT IS REACHED, MOIST CONVECTIVE ADJUSTMENT IS NECESSARY
746! =====================================================================
747
748  ncum = 0
749  DO i = 1, len
750    IF (iflag1(i)==0) THEN
751      ncum = ncum + 1
752      idcum(ncum) = i
753    END IF
754  END DO
755
756! print*,'klon, ncum = ',len,ncum
757
758  IF (ncum>0) THEN
759
760! ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
761! --- COMPRESS THE FIELDS
762!       (-> vectorization over convective gridpoints)
763! ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
764
765    IF (iflag_con==3) THEN
766! print*,'ncum tv1 ',ncum,tv1
767! print*,'tvp1 ',tvp1
768      CALL cv3a_compress(len, nloc, ncum, nd, ntra, &
769                         iflag1, nk1, icb1, icbs1, &
770                         plcl1, tnk1, qnk1, gznk1, hnk1, unk1, vnk1, &
771                         wghti1, pbase1, buoybase1, &
772                         t1, q1, qs1, t1_wake, q1_wake, qs1_wake, s1_wake, &
773                         u1, v1, gz1, th1, th1_wake, &
774                         tra1, &
775                         h1, lv1, lf1, cpn1, p1, ph1, tv1, tp1, tvp1, clw1, &
776                         h1_wake, lv1_wake, lf1_wake, cpn1_wake, tv1_wake, &
777                         sig1, w01, ptop21, &
778                         Ale1, Alp1, &
779                         iflag, nk, icb, icbs, &
780                         plcl, tnk, qnk, gznk, hnk, unk, vnk, &
781                         wghti, pbase, buoybase, &
782                         t, q, qs, t_wake, q_wake, qs_wake, s_wake, &
783                         u, v, gz, th, th_wake, &
784                         tra, &
785                         h, lv, lf, cpn, p, ph, tv, tp, tvp, clw, &
786                         h_wake, lv_wake, lf_wake, cpn_wake, tv_wake, &
787                         sig, w0, ptop2, &
788                         Ale, Alp)
789
790! print*,'tv ',tv
791! print*,'tvp ',tvp
792
793    END IF
794
795    IF (iflag_con==4) THEN
796      CALL cv_compress(len, nloc, ncum, nd, &
797                       iflag1, nk1, icb1, &
798                       cbmf1, plcl1, tnk1, qnk1, gznk1, &
799                       t1, q1, qs1, u1, v1, gz1, &
800                       h1, lv1, cpn1, p1, ph1, tv1, tp1, tvp1, clw1, &
801                       iflag, nk, icb, &
802                       cbmf, plcl, tnk, qnk, gznk, &
803                       t, q, qs, u, v, gz, h, lv, cpn, p, ph, tv, tp, tvp, clw, &
804                       dph)
805    END IF
806
807! -------------------------------------------------------------------
808! --- UNDILUTE (ADIABATIC) UPDRAFT / second part :
809! ---   FIND THE REST OF THE LIFTED PARCEL TEMPERATURES
810! ---   &
811! ---   COMPUTE THE PRECIPITATION EFFICIENCIES AND THE
812! ---   FRACTION OF PRECIPITATION FALLING OUTSIDE OF CLOUD
813! ---   &
814! ---   FIND THE LEVEL OF NEUTRAL BUOYANCY
815! -------------------------------------------------------------------
816
817    IF (iflag_con==3) THEN
818      CALL cv3_undilute2(nloc, ncum, nd, icb, icbs, nk, &              !na->nd
819                         tnk, qnk, gznk, hnk, t, q, qs, gz, &
820                         p, h, tv, lv, lf, pbase, buoybase, plcl, &
821                         inb, tp, tvp, clw, hp, ep, sigp, buoy, &
822                         frac)
823
824    END IF
825
826    IF (iflag_con==4) THEN
827      CALL cv_undilute2(nloc, ncum, nd, icb, nk, &
828                        tnk, qnk, gznk, t, q, qs, gz, &
829                        p, dph, h, tv, lv, &
830                        inb, inbis, tp, tvp, clw, hp, ep, sigp, frac)
831    END IF
832
833! -------------------------------------------------------------------
834! --- MIXING(1)   (if iflag_mix .ge. 1)
835! -------------------------------------------------------------------
836    IF (iflag_con==3) THEN
837      IF ((iflag_ice_thermo==1) .AND. (iflag_mix/=0)) THEN
838        WRITE (*, *) ' iflag_ice_thermo==1 requires iflag_mix==0', ' but iflag_mix=', iflag_mix, &
839          '. Might as well stop here.'
840        STOP
841      END IF
842      IF (iflag_mix>=1) THEN
843        CALL zilch(supmax, nloc*klev)
844        CALL cv3p_mixing(nloc, ncum, nd, nd, ntra, icb, nk, inb, &           ! na->nd
845                         ph, t, q, qs, u, v, tra, h, lv, qnk, &
846                         unk, vnk, hp, tv, tvp, ep, clw, sig, &
847                         ment, qent, hent, uent, vent, nent, &
848                         sigij, elij, supmax, ments, qents, traent)
849! print*, 'cv3p_mixing-> supmax ', (supmax(1,k), k=1,nd)
850
851      ELSE
852        CALL zilch(supmax, nloc*klev)
853      END IF
854    END IF
855! -------------------------------------------------------------------
856! --- CLOSURE
857! -------------------------------------------------------------------
858
859
860    IF (iflag_con==3) THEN
861      IF (iflag_clos==0) THEN
862        CALL cv3_closure(nloc, ncum, nd, icb, inb, &           ! na->nd
863                         pbase, p, ph, tv, buoy, &
864                         sig, w0, cape, m, iflag)
865      END IF
866
867      ok_inhib = iflag_mix == 2
868
869      IF (iflag_clos==1) THEN
870        PRINT *, ' pas d appel cv3p_closure'
871! c        CALL cv3p_closure(nloc,ncum,nd,icb,inb              ! na->nd
872! c    :                       ,pbase,plcl,p,ph,tv,tvp,buoy
873! c    :                       ,supmax
874! c    o                       ,sig,w0,ptop2,cape,cin,m)
875      END IF
876      IF (iflag_clos==2) THEN
877        CALL cv3p1_closure(nloc, ncum, nd, icb, inb, &         ! na->nd
878                           pbase, plcl, p, ph, tv, tvp, buoy, &
879                           supmax, ok_inhib, Ale, Alp, &
880                           sig, w0, ptop2, cape, cin, m, iflag, coef_clos, &
881                           Plim1, plim2, asupmax, supmax0, &
882                           asupmaxmin, cbmf, plfc, wbeff)
883
884        if (prt_level >= 10) &
885             PRINT *, 'cv3p1_closure-> plfc,wbeff ', plfc(1), wbeff(1)
886      END IF
887    END IF ! iflag_con.eq.3
888
889    IF (iflag_con==4) THEN
890      CALL cv_closure(nloc, ncum, nd, nk, icb, &
891                         tv, tvp, p, ph, dph, plcl, cpn, &
892                         iflag, cbmf)
893    END IF
894
895! print *,'cv_closure-> cape ',cape(1)
896
897! -------------------------------------------------------------------
898! --- MIXING(2)
899! -------------------------------------------------------------------
900
901    IF (iflag_con==3) THEN
902      IF (iflag_mix==0) THEN
903        CALL cv3_mixing(nloc, ncum, nd, nd, ntra, icb, nk, inb, &             ! na->nd
904                        ph, t, q, qs, u, v, tra, h, lv, lf, frac, qnk, &
905                        unk, vnk, hp, tv, tvp, ep, clw, m, sig, &
906                        ment, qent, uent, vent, nent, sigij, elij, ments, qents, traent)
907        CALL zilch(hent, nloc*klev*klev)
908      ELSE
909        CALL cv3_mixscale(nloc, ncum, nd, ment, m)
910        IF (debut) THEN
911          PRINT *, ' cv3_mixscale-> '
912        END IF !(debut) THEN
913      END IF
914    END IF
915
916    IF (iflag_con==4) THEN
917      CALL cv_mixing(nloc, ncum, nd, icb, nk, inb, inbis, &
918                     ph, t, q, qs, u, v, h, lv, qnk, &
919                     hp, tv, tvp, ep, clw, cbmf, &
920                     m, ment, qent, uent, vent, nent, sigij, elij)
921    END IF                                                                                         
922
923    IF (debut) THEN
924      PRINT *, ' cv_mixing ->'
925    END IF !(debut) THEN
926! do i = 1,klev
927! print*,'cv_mixing-> i,ment ',i,(ment(1,i,j),j=1,klev)
928! enddo
929
930! -------------------------------------------------------------------
931! --- UNSATURATED (PRECIPITATING) DOWNDRAFTS
932! -------------------------------------------------------------------
933    IF (iflag_con==3) THEN
934      IF (debut) THEN
935        PRINT *, ' cva_driver -> cv3_unsat '
936      END IF !(debut) THEN
937
938      CALL cv3_unsat(nloc, ncum, nd, nd, ntra, icb, inb, iflag, &              ! na->nd
939                     t_wake, q_wake, qs_wake, gz, u, v, tra, p, ph, &
940                     th_wake, tv_wake, lv_wake, lf_wake, cpn_wake, &
941                     ep, sigp, clw, &
942                     m, ment, elij, delt, plcl, coef_clos, &
943                     mp, qp, up, vp, trap, wt, water, evap, fondue, ice, &
944                     faci, b, sigd, &
945                     wdtrainA, wdtrainM)                                       ! RomP
946    END IF
947
948    IF (iflag_con==4) THEN
949      CALL cv_unsat(nloc, ncum, nd, inb, t, q, qs, gz, u, v, p, ph, &
950                     h, lv, ep, sigp, clw, m, ment, elij, &
951                     iflag, mp, qp, up, vp, wt, water, evap)
952    END IF
953
954    IF (debut) THEN
955      PRINT *, 'cv_unsat-> '
956    END IF !(debut) THEN
957
958! print *,'cv_unsat-> mp ',mp
959! print *,'cv_unsat-> water ',water
960! -------------------------------------------------------------------
961! --- YIELD
962! (tendencies, precipitation, variables of interface with other
963! processes, etc)
964! -------------------------------------------------------------------
965
966    IF (iflag_con==3) THEN
967
968      CALL cv3_yield(nloc, ncum, nd, nd, ntra, ok_conserv_q, &                      ! na->nd
969                     icb, inb, delt, &
970                     t, q, t_wake, q_wake, s_wake, u, v, tra, &
971                     gz, p, ph, h, hp, lv, lf, cpn, th, th_wake, &
972                     ep, clw, m, tp, mp, qp, up, vp, trap, &
973                     wt, water, ice, evap, fondue, faci, b, sigd, &
974                     ment, qent, hent, iflag_mix, uent, vent, &
975                     nent, elij, traent, sig, &
976                     tv, tvp, wghti, &
977                     iflag, precip, vprecip, ft, fq, fu, fv, ftra, &
978                     cbmf, upwd, dnwd, dnwd0, ma, mip, &
979                     tls, tps, qcondc, wd, &
980                     ftd, fqd)
981    END IF
982
983    IF (debut) THEN
984      PRINT *, ' cv3_yield -> fqd(1) = ', fqd(1, 1)
985    END IF !(debut) THEN
986
987    IF (iflag_con==4) THEN
988      CALL cv_yield(nloc, ncum, nd, nk, icb, inb, delt, &
989                     t, q, u, v, &
990                     gz, p, ph, h, hp, lv, cpn, &
991                     ep, clw, frac, m, mp, qp, up, vp, &
992                     wt, water, evap, &
993                     ment, qent, uent, vent, nent, elij, &
994                     tv, tvp, &
995                     iflag, wd, qprime, tprime, &
996                     precip, cbmf, ft, fq, fu, fv, ma, qcondc)
997    END IF
998
999!AC!
1000!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
1001!--- passive tracers
1002!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
1003
1004    IF (iflag_con==3) THEN
1005!RomP >>>
1006      CALL cv3_tracer(nloc, len, ncum, nd, nd, &
1007                     ment, sigij, da, phi, phi2, d1a, dam, &
1008                     ep, vprecip, elij, clw, epmlmMm, eplaMm, &
1009                     icb, inb)
1010!RomP <<<
1011    END IF
1012
1013!AC!
1014
1015! ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
1016! --- UNCOMPRESS THE FIELDS
1017! ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
1018
1019
1020    IF (iflag_con==3) THEN
1021      CALL cv3a_uncompress(nloc, len, ncum, nd, ntra, idcum, &
1022                           iflag, icb, inb, &
1023                           precip, cbmf, plcl, plfc, wbeff, sig, w0, ptop2, &
1024                           ft, fq, fu, fv, ftra, &
1025                           sigd, ma, mip, vprecip, upwd, dnwd, dnwd0, &
1026                           qcondc, wd, cape, cin, &
1027                           tvp, &
1028                           ftd, fqd, &
1029                           Plim1, plim2, asupmax, supmax0, &
1030                           asupmaxmin, &
1031                           da, phi, mp, phi2, d1a, dam, sigij, &         ! RomP
1032                           clw, elij, evap, ep, epmlmMm, eplaMm, &       ! RomP
1033                           wdtrainA, wdtrainM, &                         ! RomP
1034                           iflag1, kbas1, ktop1, &
1035                           precip1, cbmf1, plcl1, plfc1, wbeff1, sig1, w01, ptop21, &
1036                           ft1, fq1, fu1, fv1, ftra1, &
1037                           sigd1, ma1, mip1, vprecip1, upwd1, dnwd1, dnwd01, &
1038                           qcondc1, wd1, cape1, cin1, &
1039                           tvp1, &
1040                           ftd1, fqd1, &
1041                           Plim11, plim21, asupmax1, supmax01, &
1042                           asupmaxmin1, &
1043                           da1, phi1, mp1, phi21, d1a1, dam1, sigij1,  & ! RomP
1044                           clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, & ! RomP
1045                           wdtrainA1, wdtrainM1)                         ! RomP
1046    END IF
1047
1048    IF (iflag_con==4) THEN
1049      CALL cv_uncompress(nloc, len, ncum, nd, idcum, &
1050                           iflag, &
1051                           precip, cbmf, &
1052                           ft, fq, fu, fv, &
1053                           ma, qcondc, &
1054                           iflag1, &
1055                           precip1,cbmf1, &
1056                           ft1, fq1, fu1, fv1, &
1057                           ma1, qcondc1)
1058    END IF
1059
1060  END IF ! ncum>0
1061
1062  IF (debut) THEN
1063    PRINT *, ' cv_compress -> '
1064    debut = .FALSE.
1065  END IF  !(debut) THEN
1066
1067
1068  RETURN
1069END SUBROUTINE cva_driver
Note: See TracBrowser for help on using the repository browser.