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