source: LMDZ4/trunk/libf/phylmd/cv_driver.F @ 588

Last change on this file since 588 was 559, checked in by lmdzadmin, 20 years ago

Initialisations diverses YM
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 22.7 KB
Line 
1!
2! $Header$
3!
4      SUBROUTINE cv_driver(len,nd,ndp1,ntra,iflag_con,
5     &                   t1,q1,qs1,u1,v1,tra1,
6     &                   p1,ph1,iflag1,ft1,fq1,fu1,fv1,ftra1,
7     &                   precip1,
8     &                   cbmf1,sig1,w01,
9     &                   delt,Ma1,upwd1,dnwd1,dnwd01,qcondc1,wd1,cape1)
10C
11      implicit none
12C
13C.............................START PROLOGUE............................
14C
15C PARAMETERS:
16C      Name            Type         Usage            Description
17C   ----------      ----------     -------  ----------------------------
18C
19C      len           Integer        Input        first (i) dimension
20C      nd            Integer        Input        vertical (k) dimension
21C      ndp1          Integer        Input        nd + 1
22C      ntra          Integer        Input        number of tracors
23C      iflag_con     Integer        Input        version of convect (3/4)
24C      t1            Real           Input        temperature
25C      q1            Real           Input        specific hum
26C      qs1           Real           Input        sat specific hum
27C      u1            Real           Input        u-wind
28C      v1            Real           Input        v-wind
29C      tra1          Real           Input        tracors
30C      p1            Real           Input        full level pressure
31C      ph1           Real           Input        half level pressure
32C      iflag1        Integer        Output       flag for Emanuel conditions
33C      ft1           Real           Output       temp tend
34C      fq1           Real           Output       spec hum tend
35C      fu1           Real           Output       u-wind tend
36C      fv1           Real           Output       v-wind tend
37C      ftra1         Real           Output       tracor tend
38C      precip1       Real           Output       precipitation
39C      cbmf1         Real           Output       cloud base mass flux
40C      sig1          Real           In/Out       section adiabatic updraft
41C      w01           Real           In/Out       vertical velocity within adiab updraft
42C      delt          Real           Input        time step
43C      Ma1           Real           Output       mass flux adiabatic updraft
44C      upwd1         Real           Output       total upward mass flux (adiab+mixed)
45C      dnwd1         Real           Output       saturated downward mass flux (mixed)
46C      dnwd01        Real           Output       unsaturated downward mass flux
47C      qcondc1       Real           Output       in-cld mixing ratio of condensed water
48C      wd1           Real           Output       downdraft velocity scale for sfc fluxes
49C      cape1         Real           Output       CAPE
50C
51C S. Bony, Mar 2002:
52C       * Several modules corresponding to different physical processes
53C       * Several versions of convect may be used:
54C               - iflag_con=3: version lmd  (previously named convect3)
55C               - iflag_con=4: version 4.3b (vect. version, previously convect1/2)
56C   + tard:     - iflag_con=5: version lmd with ice (previously named convectg)
57C S. Bony, Oct 2002:
58C       * Vectorization of convect3 (ie version lmd)
59C
60C..............................END PROLOGUE.............................
61c
62c
63#include "dimensions.h"
64#include "dimphy.h"
65
66      integer len
67      integer nd
68      integer ndp1
69      integer noff
70      integer iflag_con
71      integer ntra
72      real t1(len,nd)
73      real q1(len,nd)
74      real qs1(len,nd)
75      real u1(len,nd)
76      real v1(len,nd)
77      real p1(len,nd)
78      real ph1(len,ndp1)
79      integer iflag1(len)
80      real ft1(len,nd)
81      real fq1(len,nd)
82      real fu1(len,nd)
83      real fv1(len,nd)
84      real precip1(len)
85      real cbmf1(len)
86      real Ma1(len,nd)
87      real upwd1(len,nd)
88      real dnwd1(len,nd)
89      real dnwd01(len,nd)
90
91      real qcondc1(len,nd)     ! cld
92      real wd1(len)            ! gust
93      real cape1(len)     
94
95      real tra1(len,nd,ntra)
96      real ftra1(len,nd,ntra)
97
98      real delt
99
100!-------------------------------------------------------------------
101! --- ARGUMENTS
102!-------------------------------------------------------------------
103! --- On input:
104!
105!  t:   Array of absolute temperature (K) of dimension ND, with first
106!       index corresponding to lowest model level. Note that this array
107!       will be altered by the subroutine if dry convective adjustment
108!       occurs and if IPBL is not equal to 0.
109!
110!  q:   Array of specific humidity (gm/gm) of dimension ND, with first
111!       index corresponding to lowest model level. Must be defined
112!       at same grid levels as T. Note that this array will be altered
113!       if dry convective adjustment occurs and if IPBL is not equal to 0.
114!
115!  qs:  Array of saturation specific humidity of dimension ND, with first
116!       index corresponding to lowest model level. Must be defined
117!       at same grid levels as T. Note that this array will be altered
118!       if dry convective adjustment occurs and if IPBL is not equal to 0.
119!
120!  u:   Array of zonal wind velocity (m/s) of dimension ND, witth first
121!       index corresponding with the lowest model level. Defined at
122!       same levels as T. Note that this array will be altered if
123!       dry convective adjustment occurs and if IPBL is not equal to 0.
124!
125!  v:   Same as u but for meridional velocity.
126!
127!  tra: Array of passive tracer mixing ratio, of dimensions (ND,NTRA),
128!       where NTRA is the number of different tracers. If no
129!       convective tracer transport is needed, define a dummy
130!       input array of dimension (ND,1). Tracers are defined at
131!       same vertical levels as T. Note that this array will be altered
132!       if dry convective adjustment occurs and if IPBL is not equal to 0.
133!
134!  p:   Array of pressure (mb) of dimension ND, with first
135!       index corresponding to lowest model level. Must be defined
136!       at same grid levels as T.
137!
138!  ph:  Array of pressure (mb) of dimension ND+1, with first index
139!       corresponding to lowest level. These pressures are defined at
140!       levels intermediate between those of P, T, Q and QS. The first
141!       value of PH should be greater than (i.e. at a lower level than)
142!       the first value of the array P.
143!
144!  nl:  The maximum number of levels to which convection can penetrate, plus 1.
145!       NL MUST be less than or equal to ND-1.
146!
147!  delt: The model time step (sec) between calls to CONVECT
148!
149!----------------------------------------------------------------------------
150! ---   On Output:
151!
152!  iflag: An output integer whose value denotes the following:
153!       VALUE   INTERPRETATION
154!       -----   --------------
155!         0     Moist convection occurs.
156!         1     Moist convection occurs, but a CFL condition
157!               on the subsidence warming is violated. This
158!               does not cause the scheme to terminate.
159!         2     Moist convection, but no precip because ep(inb) lt 0.0001
160!         3     No moist convection because new cbmf is 0 and old cbmf is 0.
161!         4     No moist convection; atmosphere is not
162!               unstable
163!         6     No moist convection because ihmin le minorig.
164!         7     No moist convection because unreasonable
165!               parcel level temperature or specific humidity.
166!         8     No moist convection: lifted condensation
167!               level is above the 200 mb level.
168!         9     No moist convection: cloud base is higher
169!               then the level NL-1.
170!
171!  ft:   Array of temperature tendency (K/s) of dimension ND, defined at same
172!        grid levels as T, Q, QS and P.
173!
174!  fq:   Array of specific humidity tendencies ((gm/gm)/s) of dimension ND,
175!        defined at same grid levels as T, Q, QS and P.
176!
177!  fu:   Array of forcing of zonal velocity (m/s^2) of dimension ND,
178!        defined at same grid levels as T.
179!
180!  fv:   Same as FU, but for forcing of meridional velocity.
181!
182!  ftra: Array of forcing of tracer content, in tracer mixing ratio per
183!        second, defined at same levels as T. Dimensioned (ND,NTRA).
184!
185!  precip: Scalar convective precipitation rate (mm/day).
186!
187!  wd:   A convective downdraft velocity scale. For use in surface
188!        flux parameterizations. See convect.ps file for details.
189!
190!  tprime: A convective downdraft temperature perturbation scale (K).
191!          For use in surface flux parameterizations. See convect.ps
192!          file for details.
193!
194!  qprime: A convective downdraft specific humidity
195!          perturbation scale (gm/gm).
196!          For use in surface flux parameterizations. See convect.ps
197!          file for details.
198!
199!  cbmf: The cloud base mass flux ((kg/m**2)/s). THIS SCALAR VALUE MUST
200!        BE STORED BY THE CALLING PROGRAM AND RETURNED TO CONVECT AT
201!        ITS NEXT CALL. That is, the value of CBMF must be "remembered"
202!        by the calling program between calls to CONVECT.
203!
204!  det:   Array of detrainment mass flux of dimension ND.
205!
206!-------------------------------------------------------------------
207c
208c  Local arrays
209c
210
211      integer i,k,n,il,j
212      integer icbmax
213      integer nk1(klon)
214      integer icb1(klon)
215      integer icbs1(klon)
216
217      real plcl1(klon)
218      real tnk1(klon)
219      real qnk1(klon)
220      real gznk1(klon)
221      real pnk1(klon)
222      real qsnk1(klon)
223      real pbase1(klon)
224      real buoybase1(klon)
225
226      real lv1(klon,klev)
227      real cpn1(klon,klev)
228      real tv1(klon,klev)
229      real gz1(klon,klev)
230      real hm1(klon,klev)
231      real h1(klon,klev)
232      real tp1(klon,klev)
233      real tvp1(klon,klev)
234      real clw1(klon,klev)
235      real sig1(klon,klev)
236      real w01(klon,klev)
237      real th1(klon,klev)
238c
239      integer ncum
240c
241c (local) compressed fields:
242c
243      integer nloc
244      parameter (nloc=klon) ! pour l'instant
245
246      integer idcum(nloc)
247      integer iflag(nloc),nk(nloc),icb(nloc)
248      integer nent(nloc,klev)
249      integer icbs(nloc)
250      integer inb(nloc), inbis(nloc)
251
252      real cbmf(nloc),plcl(nloc),tnk(nloc),qnk(nloc),gznk(nloc)
253      real t(nloc,klev),q(nloc,klev),qs(nloc,klev)
254      real u(nloc,klev),v(nloc,klev)
255      real gz(nloc,klev),h(nloc,klev),lv(nloc,klev),cpn(nloc,klev)
256      real p(nloc,klev),ph(nloc,klev+1),tv(nloc,klev),tp(nloc,klev)
257      real clw(nloc,klev)
258      real dph(nloc,klev)
259      real pbase(nloc), buoybase(nloc), th(nloc,klev)
260      real tvp(nloc,klev)
261      real sig(nloc,klev), w0(nloc,klev)
262      real hp(nloc,klev), ep(nloc,klev), sigp(nloc,klev)
263      real frac(nloc), buoy(nloc,klev)
264      real cape(nloc)
265      real m(nloc,klev), ment(nloc,klev,klev), qent(nloc,klev,klev)
266      real uent(nloc,klev,klev), vent(nloc,klev,klev)
267      real ments(nloc,klev,klev), qents(nloc,klev,klev)
268      real sij(nloc,klev,klev), elij(nloc,klev,klev)
269      real mp(nloc,klev), qp(nloc,klev), up(nloc,klev), vp(nloc,klev)
270      real wt(nloc,klev), water(nloc,klev), evap(nloc,klev)
271      real b(nloc,klev), ft(nloc,klev), fq(nloc,klev)
272      real fu(nloc,klev), fv(nloc,klev)
273      real upwd(nloc,klev), dnwd(nloc,klev), dnwd0(nloc,klev)
274      real Ma(nloc,klev), mike(nloc,klev), tls(nloc,klev)
275      real tps(nloc,klev), qprime(nloc), tprime(nloc)
276      real precip(nloc)
277      real tra(nloc,klev,ntra), trap(nloc,klev,ntra)
278      real ftra(nloc,klev,ntra), traent(nloc,klev,klev,ntra)
279      real qcondc(nloc,klev)  ! cld
280      real wd(nloc)           ! gust
281
282!-------------------------------------------------------------------
283! --- SET CONSTANTS AND PARAMETERS
284!-------------------------------------------------------------------
285
286c -- set simulation flags:
287c   (common cvflag)
288
289       CALL cv_flag
290
291c -- set thermodynamical constants:
292c       (common cvthermo)
293
294       CALL cv_thermo(iflag_con)
295
296c -- set convect parameters
297c
298c       includes microphysical parameters and parameters that
299c       control the rate of approach to quasi-equilibrium)
300c       (common cvparam)
301
302      if (iflag_con.eq.3) then
303       CALL cv3_param(nd,delt)
304      endif
305
306      if (iflag_con.eq.4) then
307       CALL cv_param(nd)
308      endif
309
310!---------------------------------------------------------------------
311! --- INITIALIZE OUTPUT ARRAYS AND PARAMETERS
312!---------------------------------------------------------------------
313
314      do 20 k=1,nd
315        do 10 i=1,len
316         ft1(i,k)=0.0
317         fq1(i,k)=0.0
318         fu1(i,k)=0.0
319         fv1(i,k)=0.0
320         tvp1(i,k)=0.0
321         tp1(i,k)=0.0
322         clw1(i,k)=0.0
323cym
324         clw(i,k)=0.0   
325         gz1(i,k) = 0.
326
327         Ma1(i,k)=0.0
328         upwd1(i,k)=0.0
329         dnwd1(i,k)=0.0
330         dnwd01(i,k)=0.0
331         qcondc1(i,k)=0.0
332 10     continue
333 20   continue
334
335      do 30 j=1,ntra
336       do 31 k=1,nd
337        do 32 i=1,len
338         ftra1(i,k,j)=0.0
339 32     continue   
340 31    continue   
341 30   continue   
342
343      do 60 i=1,len
344        precip1(i)=0.0
345        iflag1(i)=0
346        wd1(i)=0.0
347        cape1(i)=0.0
348 60   continue
349
350      if (iflag_con.eq.3) then
351        do il=1,len
352         sig1(il,nd)=sig1(il,nd)+1.
353         sig1(il,nd)=amin1(sig1(il,nd),12.1)
354        enddo
355      endif
356
357!--------------------------------------------------------------------
358! --- CALCULATE ARRAYS OF GEOPOTENTIAL, HEAT CAPACITY & STATIC ENERGY
359!--------------------------------------------------------------------
360
361      if (iflag_con.eq.3) then
362       CALL cv3_prelim(len,nd,ndp1,t1,q1,p1,ph1            ! nd->na
363     o               ,lv1,cpn1,tv1,gz1,h1,hm1,th1)
364      endif
365
366      if (iflag_con.eq.4) then
367       CALL cv_prelim(len,nd,ndp1,t1,q1,p1,ph1
368     o               ,lv1,cpn1,tv1,gz1,h1,hm1)
369      endif
370
371!--------------------------------------------------------------------
372! --- CONVECTIVE FEED
373!--------------------------------------------------------------------
374
375      if (iflag_con.eq.3) then
376       CALL cv3_feed(len,nd,t1,q1,qs1,p1,ph1,hm1,gz1           ! nd->na
377     o         ,nk1,icb1,icbmax,iflag1,tnk1,qnk1,gznk1,plcl1)
378      endif
379
380      if (iflag_con.eq.4) then
381       CALL cv_feed(len,nd,t1,q1,qs1,p1,hm1,gz1
382     o         ,nk1,icb1,icbmax,iflag1,tnk1,qnk1,gznk1,plcl1)
383      endif
384
385!--------------------------------------------------------------------
386! --- UNDILUTE (ADIABATIC) UPDRAFT / 1st part
387! (up through ICB for convect4, up through ICB+1 for convect3)
388!     Calculates the lifted parcel virtual temperature at nk, the
389!     actual temperature, and the adiabatic liquid water content.
390!--------------------------------------------------------------------
391
392      if (iflag_con.eq.3) then
393       CALL cv3_undilute1(len,nd,t1,q1,qs1,gz1,plcl1,p1,nk1,icb1  ! nd->na
394     o                        ,tp1,tvp1,clw1,icbs1)
395      endif
396
397      if (iflag_con.eq.4) then
398       CALL cv_undilute1(len,nd,t1,q1,qs1,gz1,p1,nk1,icb1,icbmax
399     :                        ,tp1,tvp1,clw1)
400      endif
401
402!-------------------------------------------------------------------
403! --- TRIGGERING
404!-------------------------------------------------------------------
405
406      if (iflag_con.eq.3) then
407       CALL cv3_trigger(len,nd,icb1,plcl1,p1,th1,tv1,tvp1      ! nd->na
408     o                 ,pbase1,buoybase1,iflag1,sig1,w01)
409      endif
410
411      if (iflag_con.eq.4) then
412       CALL cv_trigger(len,nd,icb1,cbmf1,tv1,tvp1,iflag1)
413      endif
414
415!=====================================================================
416! --- IF THIS POINT IS REACHED, MOIST CONVECTIVE ADJUSTMENT IS NECESSARY
417!=====================================================================
418
419      ncum=0
420      do 400 i=1,len
421        if(iflag1(i).eq.0)then
422           ncum=ncum+1
423           idcum(ncum)=i
424        endif
425 400  continue
426
427c       print*,'klon, ncum = ',len,ncum
428
429      IF (ncum.gt.0) THEN
430
431!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
432! --- COMPRESS THE FIELDS
433!               (-> vectorization over convective gridpoints)
434!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
435
436      if (iflag_con.eq.3) then
437       CALL cv3_compress( len,nloc,ncum,nd,ntra
438     :    ,iflag1,nk1,icb1,icbs1
439     :    ,plcl1,tnk1,qnk1,gznk1,pbase1,buoybase1
440     :    ,t1,q1,qs1,u1,v1,gz1,th1
441     :    ,tra1
442     :    ,h1,lv1,cpn1,p1,ph1,tv1,tp1,tvp1,clw1
443     :    ,sig1,w01
444     o    ,iflag,nk,icb,icbs
445     o    ,plcl,tnk,qnk,gznk,pbase,buoybase
446     o    ,t,q,qs,u,v,gz,th
447     o    ,tra
448     o    ,h,lv,cpn,p,ph,tv,tp,tvp,clw
449     o    ,sig,w0  )
450      endif
451
452      if (iflag_con.eq.4) then
453       CALL cv_compress( len,nloc,ncum,nd
454     :    ,iflag1,nk1,icb1
455     :    ,cbmf1,plcl1,tnk1,qnk1,gznk1
456     :    ,t1,q1,qs1,u1,v1,gz1
457     :    ,h1,lv1,cpn1,p1,ph1,tv1,tp1,tvp1,clw1
458     o    ,iflag,nk,icb
459     o    ,cbmf,plcl,tnk,qnk,gznk
460     o    ,t,q,qs,u,v,gz,h,lv,cpn,p,ph,tv,tp,tvp,clw
461     o    ,dph )
462      endif
463
464!-------------------------------------------------------------------
465! --- UNDILUTE (ADIABATIC) UPDRAFT / second part :
466! ---   FIND THE REST OF THE LIFTED PARCEL TEMPERATURES
467! ---   &
468! ---   COMPUTE THE PRECIPITATION EFFICIENCIES AND THE
469! ---   FRACTION OF PRECIPITATION FALLING OUTSIDE OF CLOUD
470! ---   &
471! ---   FIND THE LEVEL OF NEUTRAL BUOYANCY
472!-------------------------------------------------------------------
473
474      if (iflag_con.eq.3) then
475       CALL cv3_undilute2(nloc,ncum,nd,icb,icbs,nk        !na->nd
476     :                        ,tnk,qnk,gznk,t,q,qs,gz
477     :                        ,p,h,tv,lv,pbase,buoybase,plcl
478     o                        ,inb,tp,tvp,clw,hp,ep,sigp,buoy)
479      endif
480
481      if (iflag_con.eq.4) then
482       CALL cv_undilute2(nloc,ncum,nd,icb,nk
483     :                        ,tnk,qnk,gznk,t,q,qs,gz
484     :                        ,p,dph,h,tv,lv
485     o             ,inb,inbis,tp,tvp,clw,hp,ep,sigp,frac)
486      endif
487
488!-------------------------------------------------------------------
489! --- CLOSURE
490!-------------------------------------------------------------------
491
492      if (iflag_con.eq.3) then
493       CALL cv3_closure(nloc,ncum,nd,icb,inb              ! na->nd
494     :                       ,pbase,p,ph,tv,buoy
495     o                       ,sig,w0,cape,m)
496      endif
497
498      if (iflag_con.eq.4) then
499       CALL cv_closure(nloc,ncum,nd,nk,icb
500     :                ,tv,tvp,p,ph,dph,plcl,cpn
501     o                ,iflag,cbmf)
502      endif
503
504!-------------------------------------------------------------------
505! --- MIXING
506!-------------------------------------------------------------------
507
508      if (iflag_con.eq.3) then
509       CALL cv3_mixing(nloc,ncum,nd,nd,ntra,icb,nk,inb    ! na->nd
510     :                     ,ph,t,q,qs,u,v,tra,h,lv,qnk
511     :                     ,hp,tv,tvp,ep,clw,m,sig
512     o ,ment,qent,uent,vent,sij,elij,ments,qents,traent)
513      endif
514
515      if (iflag_con.eq.4) then
516       CALL cv_mixing(nloc,ncum,nd,icb,nk,inb,inbis
517     :                     ,ph,t,q,qs,u,v,h,lv,qnk
518     :                     ,hp,tv,tvp,ep,clw,cbmf
519     o                     ,m,ment,qent,uent,vent,nent,sij,elij)
520      endif
521
522!-------------------------------------------------------------------
523! --- UNSATURATED (PRECIPITATING) DOWNDRAFTS
524!-------------------------------------------------------------------
525
526      if (iflag_con.eq.3) then
527       CALL cv3_unsat(nloc,ncum,nd,nd,ntra,icb,inb    ! na->nd
528     :               ,t,q,qs,gz,u,v,tra,p,ph
529     :               ,th,tv,lv,cpn,ep,sigp,clw
530     :               ,m,ment,elij,delt,plcl
531     o          ,mp,qp,up,vp,trap,wt,water,evap,b)
532      endif
533
534      if (iflag_con.eq.4) then
535       CALL cv_unsat(nloc,ncum,nd,inb,t,q,qs,gz,u,v,p,ph
536     :                   ,h,lv,ep,sigp,clw,m,ment,elij
537     o                   ,iflag,mp,qp,up,vp,wt,water,evap)
538      endif
539
540!-------------------------------------------------------------------
541! --- YIELD
542!     (tendencies, precipitation, variables of interface with other
543!      processes, etc)
544!-------------------------------------------------------------------
545
546      if (iflag_con.eq.3) then
547       CALL cv3_yield(nloc,ncum,nd,nd,ntra            ! na->nd
548     :                     ,icb,inb,delt
549     :                     ,t,q,u,v,tra,gz,p,ph,h,hp,lv,cpn,th
550     :                     ,ep,clw,m,tp,mp,qp,up,vp,trap
551     :                     ,wt,water,evap,b
552     :                     ,ment,qent,uent,vent,nent,elij,traent,sig
553     :                     ,tv,tvp
554     o                     ,iflag,precip,ft,fq,fu,fv,ftra
555     o                     ,upwd,dnwd,dnwd0,ma,mike,tls,tps,qcondc,wd)
556      endif
557
558      if (iflag_con.eq.4) then
559       CALL cv_yield(nloc,ncum,nd,nk,icb,inb,delt
560     :              ,t,q,u,v,gz,p,ph,h,hp,lv,cpn
561     :              ,ep,clw,frac,m,mp,qp,up,vp
562     :              ,wt,water,evap
563     :              ,ment,qent,uent,vent,nent,elij
564     :              ,tv,tvp
565     o              ,iflag,wd,qprime,tprime
566     o              ,precip,cbmf,ft,fq,fu,fv,Ma,qcondc)
567      endif
568
569!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
570! --- UNCOMPRESS THE FIELDS
571!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
572
573
574      if (iflag_con.eq.3) then
575       CALL cv3_uncompress(nloc,len,ncum,nd,ntra,idcum
576     :          ,iflag
577     :          ,precip,sig,w0
578     :          ,ft,fq,fu,fv,ftra
579     :          ,Ma,upwd,dnwd,dnwd0,qcondc,wd,cape
580     o          ,iflag1
581     o          ,precip1,sig1,w01
582     o          ,ft1,fq1,fu1,fv1,ftra1
583     o          ,Ma1,upwd1,dnwd1,dnwd01,qcondc1,wd1,cape1 )
584      endif
585
586      if (iflag_con.eq.4) then
587       CALL cv_uncompress(nloc,len,ncum,nd,idcum
588     :          ,iflag
589     :          ,precip,cbmf
590     :          ,ft,fq,fu,fv
591     :          ,Ma,qcondc           
592     o          ,iflag1
593     o          ,precip1,cbmf1
594     o          ,ft1,fq1,fu1,fv1
595     o          ,Ma1,qcondc1 )           
596      endif
597
598      ENDIF ! ncum>0
599
6009999  continue
601
602      return
603      end
604
605!==================================================================
606      SUBROUTINE cv_flag
607      implicit none
608
609#include "cvflag.h"
610
611c -- si .TRUE., on rend la gravite plus explicite et eventuellement
612c differente de 10.0 dans convect3:
613      cvflag_grav = .FALSE.
614
615      return
616      end
617
618!==================================================================
619      SUBROUTINE cv_thermo(iflag_con)
620          implicit none
621
622c-------------------------------------------------------------
623c Set thermodynamical constants for convectL
624c-------------------------------------------------------------
625
626#include "YOMCST.h"
627#include "cvthermo.h"
628
629      integer iflag_con
630
631
632c original set from convect:
633      if (iflag_con.eq.4) then
634       cpd=1005.7
635       cpv=1870.0
636       cl=4190.0
637       rrv=461.5
638       rrd=287.04
639       lv0=2.501E6
640       g=9.8
641       t0=273.15
642       grav=g
643      endif
644
645c constants consistent with LMDZ:
646      if (iflag_con.eq.3) then
647       cpd = RCPD
648       cpv = RCPV
649       cl  = RCW
650       rrv = RV
651       rrd = RD
652       lv0 = RLVTT
653       g   = RG     ! not used in convect3
654c ori      t0  = RTT
655       t0  = 273.15 ! convect3 (RTT=273.16)
656       grav= 10.    ! implicitely or explicitely used in convect3
657      endif
658
659      rowl=1000.0 !(a quelle variable de YOMCST cela correspond-il?)
660
661      clmcpv=cl-cpv
662      clmcpd=cl-cpd
663      cpdmcp=cpd-cpv
664      cpvmcpd=cpv-cpd
665      cpvmcl=cl-cpv ! for convect3
666      eps=rrd/rrv
667      epsi=1.0/eps
668      epsim1=epsi-1.0
669c      ginv=1.0/g
670      ginv=1.0/grav
671      hrd=0.5*rrd
672
673      return
674      end
675
Note: See TracBrowser for help on using the repository browser.