source: LMDZ5/trunk/libf/phylmd/cva_driver.F @ 1893

Last change on this file since 1893 was 1849, checked in by Ehouarn Millour, 11 years ago

Including the thermodynamics of ice in the convection scheme (iactive only if iflag_ice_thermo==1).
CR+JYG

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