source: LMDZ6/trunk/libf/phylmdiso/concvl.F90 @ 4143

Last change on this file since 4143 was 4143, checked in by dcugnet, 2 years ago
  • Some variables are renamed or replaced by direct equivalents:
    • iso_indnum -> tracers(:)%iso_iName
    • niso_possibles -> niso
    • iqiso -> iqIsoPha ; index_trac -> itZonIso
    • ok_iso_verif -> isoCheck
    • ntraceurs_zone -> nzone ; ntraciso -> ntiso
    • qperemin -> min_qparent ; masseqmin -> min_qmass ; ratiomin -> min_ratio
  • Some renamed variables are only aliased with the older name (using USE <module>, ONLY: <oldName> => <newName>) in routines where they are repeated many times.
  • Few hard-coded indexes are now computed (examples: ilic, iso, ivap, irneb, iq_vap, iq_liq, iso_H2O, iso_HDO, iso_HTO, iso_O17, iso_O18).
  • The IF(isoCheck) test is now embedded in the check_isotopes_seq and check_isotopes_loc routines (lighter calling).
  • Property svn:keywords set to Id
File size: 28.6 KB
Line 
1SUBROUTINE concvl(iflag_clos, &
2                  dtime, paprs, pplay, k_upper_cv, &
3                  t, q, t_wake, q_wake, s_wake, u, v, tra, ntra, &
4                  Ale, Alp, sig1, w01, &
5                  d_t, d_q, d_u, d_v, d_tra, &
6                  rain, snow, kbas, ktop, sigd, &
7                  cbmf, plcl, plfc, wbeff, convoccur, &
8                  upwd, dnwd, dnwdbis, &
9                  Ma, mip, Vprecip, &
10                  cape, cin, tvp, Tconv, iflag, &
11                  pbase, bbase, dtvpdt1, dtvpdq1, dplcldt, dplcldr, &
12                  qcondc, wd, pmflxr, pmflxs, &
13!RomP >>>
14!!     .             da,phi,mp,dd_t,dd_q,lalim_conv,wght_th)
15                  da, phi, mp, phii, d1a, dam, sij, qta, clw, elij, &! RomP
16                  dd_t, dd_q, lalim_conv, wght_th, &                 ! RomP
17                  evap, ep, epmlmMm, eplaMm, &                       ! RomP
18                  wdtrainA, wdtrainS, wdtrainM, wght, qtc, sigt, &
19                  tau_cld_cv, coefw_cld_cv, &                           ! RomP+RL, AJ
20!RomP <<<
21                  epmax_diag & ! epmax_cape
22#ifdef ISO
23     &       ,xt,xt_wake,d_xt,xtrain,xtsnow &
24     &       ,xtVprecip,xtVprecipi   &
25     &       ,xtclw,dd_xt,xtevap,xtwdtrainA &
26#ifdef DIAGISO
27     &       , qlp,xtlp,qvp,xtvp & ! juste diagnostique 
28     &       , fq_detrainement,fq_ddft,fq_fluxmasse,fq_evapprecip &
29     &       , fxt_detrainement,fxt_ddft,fxt_fluxmasse,fxt_evapprecip &
30     &       , f_detrainement,q_detrainement,xt_detrainement &
31#endif         
32#endif
33     &        ) ! **************************************************************
34! *
35! CONCVL                                                      *
36! *
37! *
38! written by   : Sandrine Bony-Lena, 17/05/2003, 11.16.04    *
39! modified by :                                               *
40! **************************************************************
41
42
43  USE dimphy
44  USE infotrac_phy, ONLY: nbtr
45#ifdef ISO
46  USE infotrac_phy, ONLY: ntraciso=>ntiso
47  USE isotopes_mod, ONLY: iso_eau, bidouille_anti_divergence, ridicule, &
48        iso_eau,iso_HDO
49#endif
50#ifdef ISOVERIF
51      USE isotopes_verif_mod, ONLY: errmax,errmaxrel, &
52        iso_verif_egalite_choix,iso_verif_aberrant,iso_verif_egalite, &
53        iso_verif_noNaN,iso_verif_aberrant_encadre
54#endif
55#ifdef ISOTRAC
56      USE isotrac_routines_mod, only: iso_verif_traceur_jbid_vect
57#ifdef ISOVERIF
58      USE isotopes_verif_mod, ONLY: iso_verif_traceur_vect, &
59&       iso_verif_trac_masse_vect, iso_verif_traceur,  &
60&       iso_verif_traceur_justmass
61#endif
62#endif
63  USE phys_local_var_mod, ONLY: omega
64  USE print_control_mod, ONLY: prt_level, lunout
65  IMPLICIT NONE
66! ======================================================================
67! Auteur(s): S. Bony-Lena (LMD/CNRS) date: ???
68! Objet: schema de convection de Emanuel (1991) interface
69! ======================================================================
70! Arguments:
71! dtime--input-R-pas d'integration (s)
72! s-------input-R-la vAleur "s" pour chaque couche
73! sigs----input-R-la vAleur "sigma" de chaque couche
74! sig-----input-R-la vAleur de "sigma" pour chaque niveau
75! psolpa--input-R-la pression au sol (en Pa)
76! pskapa--input-R-exponentiel kappa de psolpa
77! h-------input-R-enthAlpie potentielle (Cp*T/P**kappa)
78! q-------input-R-vapeur d'eau (en kg/kg)
79
80! work*: input et output: deux variables de travail,
81! on peut les mettre a 0 au debut
82! ALE--------input-R-energie disponible pour soulevement
83! ALP--------input-R-puissance disponible pour soulevement
84
85! d_h--------output-R-increment de l'enthAlpie potentielle (h)
86! d_q--------output-R-increment de la vapeur d'eau
87! rain-------output-R-la pluie (mm/s)
88! snow-------output-R-la neige (mm/s)
89! upwd-------output-R-saturated updraft mass flux (kg/m**2/s)
90! dnwd-------output-R-saturated downdraft mass flux (kg/m**2/s)
91! dnwd0------output-R-unsaturated downdraft mass flux (kg/m**2/s)
92! Ma---------output-R-adiabatic ascent mass flux (kg/m2/s)
93! mip--------output-R-mass flux shed by adiabatic ascent (kg/m2/s)
94! Vprecip----output-R-vertical profile of total precipitation (kg/m2/s)
95! Tconv------output-R-environment temperature seen by convective scheme (K)
96! Cape-------output-R-CAPE (J/kg)
97! Cin -------output-R-CIN  (J/kg)
98! Tvp--------output-R-Temperature virtuelle d'une parcelle soulevee
99! adiabatiquement a partir du niveau 1 (K)
100! deltapb----output-R-distance entre LCL et base de la colonne (<0 ; Pa)
101! Ice_flag---input-L-TRUE->prise en compte de la thermodynamique de la glace
102! dd_t-------output-R-increment de la temperature du aux descentes precipitantes
103! dd_q-------output-R-increment de la vapeur d'eau du aux desc precip
104! lalim_conv-
105! wght_th----
106! evap-------output-R
107! ep---------output-R
108! epmlmMm----output-R
109! eplaMm-----output-R
110! wdtrainA---output-R
111! wdtrainS---output-R
112! wdtrainM---output-R
113! wght-------output-R
114! ======================================================================
115
116
117  include "clesphys.h"
118
119  INTEGER, INTENT(IN)                           :: iflag_clos
120  REAL, INTENT(IN)                              :: dtime
121  REAL, DIMENSION(klon,klev),   INTENT(IN)      :: pplay
122  REAL, DIMENSION(klon,klev+1), INTENT(IN)      :: paprs
123  INTEGER,                      INTENT(IN)      :: k_upper_cv
124  REAL, DIMENSION(klon,klev),   INTENT(IN)      :: t, q, u, v
125  REAL, DIMENSION(klon,klev),   INTENT(IN)      :: t_wake, q_wake
126  REAL, DIMENSION(klon),        INTENT(IN)      :: s_wake
127  REAL, DIMENSION(klon,klev, nbtr),INTENT(IN)   :: tra
128  INTEGER,                      INTENT(IN)      :: ntra
129  REAL, DIMENSION(klon),        INTENT(IN)      :: Ale, Alp
130!CR:test: on passe lentr et alim_star des thermiques
131  INTEGER, DIMENSION(klon),     INTENT(IN)      :: lalim_conv
132  REAL, DIMENSION(klon,klev),   INTENT(IN)      :: wght_th
133#ifdef ISO
134  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(IN)    ::  xt
135  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(IN)    ::  xt_wake
136#endif
137
138  REAL, DIMENSION(klon,klev),   INTENT(INOUT)   :: sig1, w01
139
140  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: d_t, d_q, d_u, d_v
141  REAL, DIMENSION(klon,klev, nbtr),INTENT(OUT)  :: d_tra
142  REAL, DIMENSION(klon),        INTENT(OUT)     :: rain, snow
143
144  INTEGER, DIMENSION(klon),     INTENT(OUT)     :: kbas, ktop
145  REAL, DIMENSION(klon),        INTENT(OUT)     :: sigd
146  REAL, DIMENSION(klon),        INTENT(OUT)     :: cbmf, plcl, plfc, wbeff
147  REAL, DIMENSION(klon),        INTENT(OUT)     :: convoccur
148  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: upwd, dnwd, dnwdbis
149
150!!       REAL Ma(klon,klev), mip(klon,klev),Vprecip(klon,klev)                    !jyg
151  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: Ma, mip
152  REAL, DIMENSION(klon,klev+1), INTENT(OUT)     :: Vprecip                        !jyg
153  REAL, DIMENSION(klon),        INTENT(OUT)     :: cape, cin
154  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: tvp
155  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: Tconv
156  INTEGER, DIMENSION(klon),     INTENT(OUT)     :: iflag
157  REAL, DIMENSION(klon),        INTENT(OUT)     :: pbase, bbase
158  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: dtvpdt1, dtvpdq1
159  REAL, DIMENSION(klon),        INTENT(OUT)     :: dplcldt, dplcldr
160  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: qcondc
161  REAL, DIMENSION(klon),        INTENT(OUT)     :: wd
162  REAL, DIMENSION(klon,klev+1), INTENT(OUT)     :: pmflxr, pmflxs
163
164#ifdef ISO
165  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)    ::  d_xt
166  REAL, DIMENSION(ntraciso,klon),   INTENT(OUT)    ::  xtrain
167  REAL, DIMENSION(ntraciso,klon),   INTENT(OUT)    ::  xtsnow
168  REAL, DIMENSION(ntraciso,klon,klev+1),   INTENT(OUT)    ::  xtVprecip
169#endif
170
171
172  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: da, mp
173  REAL, DIMENSION(klon,klev,klev),INTENT(OUT)   :: phi
174! RomP >>>
175  REAL, DIMENSION(klon,klev,klev),INTENT(OUT)   :: phii
176  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: d1a, dam
177  REAL, DIMENSION(klon,klev,klev),INTENT(OUT)   :: sij, elij
178  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: qta
179  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: clw
180  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: dd_t, dd_q
181  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: evap, ep
182  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: eplaMm
183  REAL, DIMENSION(klon,klev,klev), INTENT(OUT)  :: epmlmMm
184  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: wdtrainA, wdtrainS, wdtrainM
185! RomP <<<
186  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: wght                       !RL
187  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: qtc
188  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: sigt
189  REAL,                         INTENT(OUT)     :: tau_cld_cv, coefw_cld_cv
190  REAL, DIMENSION(klon),        INTENT(OUT)     :: epmax_diag                ! epmax_cape
191
192#ifdef ISO
193  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: xtevap
194  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: xtwdtrainA
195  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: xtclw
196  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: dd_xt
197       ! juste diagnostique
198#ifdef DIAGISO
199  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: qlp
200  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: xtlp
201  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: qvp
202  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: xtvp
203  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: fq_detrainement
204  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: fq_fluxmasse
205  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: fq_ddft
206  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: fq_evapprecip
207  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: fxt_detrainement
208  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: fxt_fluxmasse
209  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: fxt_ddft
210  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: fxt_evapprecip
211  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: f_detrainement
212  REAL, DIMENSION(klon,klev),   INTENT(OUT)              :: q_detrainement
213  REAL, DIMENSION(ntraciso,klon,klev),   INTENT(OUT)     :: xt_detrainement
214#endif         
215#endif
216
217
218!
219!  Local
220!  ----
221  REAL, DIMENSION(klon,klev)                    :: em_p
222  REAL, DIMENSION(klon,klev+1)                  :: em_ph
223  REAL                                          :: em_sig1feed ! sigma at lower bound of feeding layer
224  REAL                                          :: em_sig2feed ! sigma at upper bound of feeding layer
225  REAL, DIMENSION(klev)                         :: em_wght ! weight density determining the feeding mixture
226  REAL, DIMENSION(klon,klev+1)                  :: Vprecipi                       !jyg
227!on enleve le save
228! SAVE em_sig1feed,em_sig2feed,em_wght
229
230  REAL, DIMENSION(klon)                         :: rflag
231  REAL, DIMENSION(klon)                         :: plim1, plim2
232  REAL, DIMENSION(klon)                         :: ptop2
233  REAL, DIMENSION(klon,klev)                    :: asupmax
234  REAL, DIMENSION(klon)                         :: supmax0, asupmaxmin
235  REAL                                          :: zx_t, zdelta, zx_qs, zcor
236!
237!   INTEGER iflag_mix
238!   SAVE iflag_mix
239  INTEGER                                       :: noff, minorig
240  INTEGER                                       :: i,j, k, itra
241  REAL, DIMENSION(klon,klev)                    :: qs, qs_wake
242!LF          SAVE cbmf
243!IM/JYG      REAL, SAVE, ALLOCATABLE :: cbmf(:)
244!!!$OMP THREADPRIVATE(cbmf)!
245  REAL, DIMENSION(klon)                         :: cbmflast
246#ifdef ISO
247REAL, DIMENSION(ntraciso,klon,klev+1)                  :: xtVprecipi
248  INTEGER                                       :: ixt
249#endif
250
251
252! Variables supplementaires liees au bilan d'energie
253! Real paire(klon)
254!LF      Real ql(klon,klev)
255! Save paire
256!LF      Save ql
257!LF      Real t1(klon,klev),q1(klon,klev)
258!LF      Save t1,q1
259! Data paire /1./
260  REAL, SAVE, ALLOCATABLE :: ql(:, :), q1(:, :), t1(:, :)
261!$OMP THREADPRIVATE(ql, q1, t1)
262        ! pas besoin d'isos ici
263
264! Variables liees au bilan d'energie et d'enthAlpi
265  REAL ztsol(klon)
266  REAL        h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot, &
267              h_qs_tot, qw_tot, ql_tot, qs_tot, ec_tot
268  SAVE        h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot, &
269              h_qs_tot, qw_tot, ql_tot, qs_tot, ec_tot
270!$OMP THREADPRIVATE(h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot)
271!$OMP THREADPRIVATE(h_qs_tot, qw_tot, ql_tot, qs_tot , ec_tot)
272  REAL        d_h_vcol, d_h_dair, d_qt, d_qw, d_ql, d_qs, d_ec
273  REAL        d_h_vcol_phy
274  REAL        fs_bound, fq_bound
275  SAVE        d_h_vcol_phy
276!$OMP THREADPRIVATE(d_h_vcol_phy)
277  REAL        zero_v(klon)
278  CHARACTER *15 ztit
279  INTEGER     ip_ebil ! PRINT level for energy conserv. diag.
280  SAVE        ip_ebil
281  DATA        ip_ebil/2/
282!$OMP THREADPRIVATE(ip_ebil)
283  INTEGER     if_ebil ! level for energy conserv. dignostics
284  SAVE        if_ebil
285  DATA        if_ebil/2/
286!$OMP THREADPRIVATE(if_ebil)
287!+jld ec_conser
288  REAL d_t_ec(klon, klev) ! tendance du a la conersion Ec -> E thermique
289  REAL zrcpd
290!-jld ec_conser
291!LF
292  INTEGER nloc
293  LOGICAL, SAVE            :: first = .TRUE.
294!$OMP THREADPRIVATE(first)
295  INTEGER, SAVE            :: itap, igout
296!$OMP THREADPRIVATE(itap, igout)
297
298
299  include "YOMCST.h"
300  include "YOMCST2.h"
301  include "YOETHF.h"
302  include "FCTTRE.h"
303!jyg<
304  include "conema3.h"
305!>jyg
306
307  IF (first) THEN
308! Allocate some variables LF 04/2008
309
310!IM/JYG allocate(cbmf(klon))
311    ALLOCATE (ql(klon,klev))
312    ALLOCATE (t1(klon,klev))
313    ALLOCATE (q1(klon,klev))
314!
315    convoccur(:) = 0.
316!
317    itap = 0
318    igout = klon/2 + 1/klon
319  END IF
320! Incrementer le compteur de la physique
321  itap = itap + 1
322
323! Copy T into Tconv
324  DO k = 1, klev
325    DO i = 1, klon
326      Tconv(i, k) = t(i, k)
327    END DO
328  END DO
329
330  IF (if_ebil>=1) THEN
331    DO i = 1, klon
332      ztsol(i) = t(i, 1)
333      zero_v(i) = 0.
334      DO k = 1, klev
335        ql(i, k) = 0.
336      END DO
337    END DO
338  END IF
339
340! ym
341  snow(:) = 0
342#ifdef ISO
343      xtsnow(:,:)=0
344#endif
345
346  IF (first) THEN
347    first = .FALSE.
348
349! ===========================================================================
350! READ IN PARAMETERS FOR THE CLOSURE AND THE MIXING DISTRIBUTION
351! ===========================================================================
352
353    IF (iflag_con==3) THEN
354!      CALL cv3_inicp()
355      CALL cv3_inip()
356    END IF
357
358! ===========================================================================
359! READ IN PARAMETERS FOR CONVECTIVE INHIBITION BY TROPOS. DRYNESS
360! ===========================================================================
361
362! c$$$         open (56,file='supcrit.data')
363! c$$$         read (56,*) Supcrit1, Supcrit2
364! c$$$         close (56)
365
366    IF (prt_level>=10) WRITE (lunout, *) 'supcrit1, supcrit2', supcrit1, supcrit2
367
368! ===========================================================================
369! Initialisation pour les bilans d'eau et d'energie
370! ===========================================================================
371    IF (if_ebil>=1) d_h_vcol_phy = 0.
372
373    DO i = 1, klon
374      cbmf(i) = 0.
375!!          plcl(i) = 0.
376      sigd(i) = 0.
377    END DO
378  END IF !(first)
379
380! Initialisation a chaque pas de temps
381  plfc(:) = 0.
382  wbeff(:) = 100.
383  plcl(:) = 0.
384
385  DO k = 1, klev + 1
386    DO i = 1, klon
387      em_ph(i, k) = paprs(i, k)/100.0
388      pmflxr(i, k) = 0.
389      pmflxs(i, k) = 0.
390    END DO
391  END DO
392
393  DO k = 1, klev
394    DO i = 1, klon
395      em_p(i, k) = pplay(i, k)/100.0
396    END DO
397  END DO
398
399
400! Feeding layer
401
402  em_sig1feed = 1.
403!jyg<
404!  em_sig2feed = 0.97
405  em_sig2feed = cvl_sig2feed
406!>jyg
407! em_sig2feed = 0.8
408! Relative Weight densities
409  DO k = 1, klev
410    em_wght(k) = 1.
411  END DO
412!CRtest: couche alim des tehrmiques ponderee par a*
413! DO i = 1, klon
414! do k=1,lalim_conv(i)
415! em_wght(k)=wght_th(i,k)
416! print*,'em_wght=',em_wght(k),wght_th(i,k)
417! end do
418! END DO
419
420  IF (iflag_con==4) THEN
421    DO k = 1, klev
422      DO i = 1, klon
423        zx_t = t(i, k)
424        zdelta = max(0., sign(1.,rtt-zx_t))
425        zx_qs = min(0.5, r2es*foeew(zx_t,zdelta)/em_p(i,k)/100.0)
426        zcor = 1./(1.-retv*zx_qs)
427        qs(i, k) = zx_qs*zcor
428      END DO
429      DO i = 1, klon
430        zx_t = t_wake(i, k)
431        zdelta = max(0., sign(1.,rtt-zx_t))
432        zx_qs = min(0.5, r2es*foeew(zx_t,zdelta)/em_p(i,k)/100.0)
433        zcor = 1./(1.-retv*zx_qs)
434        qs_wake(i, k) = zx_qs*zcor
435      END DO
436    END DO
437  ELSE ! iflag_con=3 (modif de puristes qui fait la diffce pour la convergence numerique)
438    DO k = 1, klev
439      DO i = 1, klon
440        zx_t = t(i, k)
441        zdelta = max(0., sign(1.,rtt-zx_t))
442        zx_qs = r2es*foeew(zx_t, zdelta)/em_p(i, k)/100.0
443        zx_qs = min(0.5, zx_qs)
444        zcor = 1./(1.-retv*zx_qs)
445        zx_qs = zx_qs*zcor
446        qs(i, k) = zx_qs
447      END DO
448      DO i = 1, klon
449        zx_t = t_wake(i, k)
450        zdelta = max(0., sign(1.,rtt-zx_t))
451        zx_qs = r2es*foeew(zx_t, zdelta)/em_p(i, k)/100.0
452        zx_qs = min(0.5, zx_qs)
453        zcor = 1./(1.-retv*zx_qs)
454        zx_qs = zx_qs*zcor
455        qs_wake(i, k) = zx_qs
456      END DO
457    END DO
458  END IF ! iflag_con
459
460! ------------------------------------------------------------------
461
462! Main driver for convection:
463!                   iflag_con=3 -> nvlle version de KE (JYG)
464!                   iflag_con = 30  -> equivAlent to convect3
465!                   iflag_con = 4  -> equivAlent to convect1/2
466
467
468  IF (iflag_con==30) THEN
469
470 
471#ifdef ISO         
472#ifdef ISOVERIF
473       do k = 1, klev
474        do i = 1, klon               
475         do ixt=1,ntraciso         
476             call iso_verif_noNaN(xt(ixt,i,k),'concvl 394')
477         enddo
478        enddo !do i = 1, klon
479       enddo !do k = 1, klev       
480       if (iso_eau.gt.0) then
481       do k = 1, klev
482        do i = 1, klon   
483          call iso_verif_egalite_choix(xt(iso_eau,i,k),q(i,k), &
484     &            'concvl 174',errmax,errmaxrel)
485        enddo !do i = 1, klon
486       enddo !do k = 1, klev   
487       endif !if (iso_eau.gt.0) then
488       if (iso_HDO.gt.0) then
489       do k = 1, klev
490        do i = 1, klon         
491         if (q(i,k).gt.ridicule) then
492          call iso_verif_aberrant(xt(iso_hdo,i,k)/q(i,k),'concvl 175')
493         endif ! if (q(i,k).gt.ridicule) then
494        enddo
495       enddo   
496       endif !if (iso_eau.gt.0) then
497#ifdef ISOTRAC
498        do k = 1, klev
499        do i = 1, klon   
500           call iso_verif_traceur(xt(1,i,k),'concvl 218')
501        enddo
502        enddo
503#endif       
504       write(*,*) 'concvl 170: avant appel cv_driver'
505#endif
506        ! ISOVERIF ! end verif       
507#endif
508
509! print *, '-> cv_driver'      !jyg
510    CALL cv_driver(klon, klev, klevp1, ntra, iflag_con, &
511                   t, q, qs, u, v, tra, &
512                   em_p, em_ph, iflag, &
513                   d_t, d_q, d_u, d_v, d_tra, rain, &
514                   Vprecip, cbmf, sig1, w01, & !jyg
515                   kbas, ktop, &
516                   dtime, Ma, upwd, dnwd, dnwdbis, qcondc, wd, cape, &
517                   da, phi, mp, phii, d1a, dam, sij, clw, elij, &       !RomP
518                   evap, ep, epmlmMm, eplaMm, &                         !RomP
519                   wdtrainA, wdtrainM, &                                !RomP
520                   epmax_diag & ! epmax_cape
521#ifdef ISO
522     &             , xt,d_xt &
523     &             , xtrain,xtVprecip &
524     &             , xtevap,xtclw,xtwdtrainA &
525#ifdef DIAGISO
526     &          , qlp,xtlp,qvp,xtvp & ! juste diagnostique 
527     &          , fq_detrainement,fq_ddft,fq_fluxmasse,fq_evapprecip &
528     &          , fxt_detrainement,fxt_ddft,fxt_fluxmasse,fxt_evapprecip &
529     &          , f_detrainement,q_detrainement,xt_detrainement &
530#endif     
531#endif
532     &              )
533!           print *, 'cv_driver ->'      !jyg
534
535#ifdef ISO
536      ! verif
537#ifdef ISOVERIF
538       write(*,*) 'concvl 463: après appel cv_driver'
539       do k = 1, klev
540        do i = 1, klon
541        if (iso_eau.gt.0) then
542            call iso_verif_egalite(xt(iso_eau,i,k),q(i,k),'concvl 203')
543            call iso_verif_egalite(xt_wake(iso_eau,i,k),q_wake(i,k),'concvl 204')
544            call iso_verif_egalite(d_xt(iso_eau,i,k),d_q(i,k), &
545     &            'concvl 452')
546         endif !if (iso_eau.gt.0) then
547#ifdef DIAGISO
548         do ixt=1,ntraciso
549            call iso_verif_noNaN(xt(ixt,i,k),'concvl 460')
550            call iso_verif_noNaN(xtlp(ixt,i,k),'concvl 295')
551            call iso_verif_noNaN(xtvp(ixt,i,k),'concvl 260')
552          enddo
553#endif                 
554        enddo
555       enddo       
556#ifdef ISOTRAC
557           call iso_verif_traceur_vect(xt,klon,klev,'concvl 218')
558           call iso_verif_trac_masse_vect(d_xt,klon,klev, &
559     &           'concvl 464',errmax,errmaxrel)
560#endif           
561#endif
562       ! end verif       
563#endif
564
565    DO i = 1, klon
566      cbmf(i) = Ma(i, kbas(i))
567    END DO
568
569!RL
570    wght(:, :) = 0.
571    DO i = 1, klon
572      wght(i, 1) = 1.
573    END DO
574!RL
575
576  ELSE
577
578!LF   necessary for gathered fields
579    nloc = klon
580#ifdef ISOVERIF
581        write(*,*) 'concvl 581: juste avant appel de cva_driver'
582#endif
583    CALL cva_driver(klon, klev, klev+1, ntra, nloc, k_upper_cv, &
584                    iflag_con, iflag_mix, iflag_ice_thermo, &
585                    iflag_clos, ok_conserv_q, dtime, cvl_comp_threshold, &
586                    t, q, qs, t_wake, q_wake, qs_wake, s_wake, u, v, tra, &
587                    em_p, em_ph, &
588                    Ale, Alp, omega, &
589                    em_sig1feed, em_sig2feed, em_wght, &
590                    iflag, d_t, d_q, d_u, d_v, d_tra, rain, kbas, ktop, &
591                    cbmf, plcl, plfc, wbeff, sig1, w01, ptop2, sigd, &
592                    Ma, mip, Vprecip, Vprecipi, upwd, dnwd, dnwdbis, qcondc, wd, &
593                    cape, cin, tvp, &
594                    dd_t, dd_q, plim1, plim2, asupmax, supmax0, &
595                    asupmaxmin, lalim_conv, &
596!AC!+!RomP+jyg
597!!                   da,phi,mp,phii,d1a,dam,sij,clw,elij, &               ! RomP
598!!                   evap,ep,epmlmMm,eplaMm,                              ! RomP
599                    da, phi, mp, phii, d1a, dam, sij, wght, &           ! RomP+RL
600                    qta, clw, elij, evap, ep, epmlmMm, eplaMm, &        ! RomP+RL
601                    wdtrainA, wdtrainS, wdtrainM, qtc, sigt, &
602                    tau_cld_cv, coefw_cld_cv, &                         ! RomP,AJ
603!AC!+!RomP+jyg
604                    epmax_diag & ! epmax_cape
605#ifdef ISO
606     &             ,xt,xt_wake,d_xt, xtrain &
607     &             ,xtvprecip,xtvprecipi &
608     &             ,xtclw,dd_xt,xtevap,xtwdtrainA &
609#ifdef DIAGISO     
610     &          , qlp,xtlp,qvp,xtvp &
611     &          , fq_detrainement,fq_ddft,fq_fluxmasse,fq_evapprecip &
612     &          , fxt_detrainement,fxt_ddft,fxt_fluxmasse,fxt_evapprecip &
613     &          , f_detrainement,q_detrainement,xt_detrainement &
614#endif     
615#endif
616     &          )     
617  END IF
618! ------------------------------------------------------------------
619  IF (prt_level>=10) WRITE (lunout, *) ' cva_driver -> cbmf,plcl,plfc,wbeff, d_t, d_q ', &
620                                         cbmf(1), plcl(1), plfc(1), wbeff(1), d_t(1,1), d_q(1,1)
621
622  DO i = 1, klon
623    rain(i) = rain(i)/86400.
624    rflag(i) = iflag(i)
625#ifdef ISO
626       do ixt = 1, ntraciso
627        xtrain(ixt,i) = xtrain(ixt,i)/86400.
628       enddo
629#endif
630  END DO
631
632  DO k = 1, klev
633    DO i = 1, klon
634      d_t(i, k) = dtime*d_t(i, k)
635      d_q(i, k) = dtime*d_q(i, k)
636      d_u(i, k) = dtime*d_u(i, k)
637      d_v(i, k) = dtime*d_v(i, k)
638#ifdef ISO
639           do ixt = 1, ntraciso
640            d_xt(ixt,i,k) = dtime*d_xt(ixt,i,k)
641           enddo
642#endif
643    END DO
644  END DO
645
646             
647#ifdef ISO
648#ifdef ISOVERIF     
649!        k=1
650!        i=  538
651        write(*,*) 'concvl 640'
652!        write(*,*) 'q(i,k),d_q(i,k)=', q(i,k),d_q(i,k)
653!        write(*,*) 'xt(iso_HDO,i,k),d_xt(iso_HDO,i,k)=', &
654!     &          xt(iso_HDO,i,k),d_xt(iso_HDO,i,k)
655  DO k = 1, klev
656    DO i = 1, klon
657           if (iso_HDO.gt.0) then
658             if (q(i,k).gt.ridicule) then
659                 call iso_verif_aberrant_encadre((xt(iso_HDO,i,k) &
660     &        +d_xt(iso_HDO,i,k))/(q(i,k)+d_q(i,k)),'concvl 250')
661             endif !if (q_seri(i,k).gt.ridicule) then
662          endif !if (iso_HDO.gt.0) then
663           if (iso_eau.gt.0) then
664             call iso_verif_egalite_choix(d_xt(iso_eau,i,k), &
665     &          d_q(i,k),'concvl 530',errmax*dtime,errmaxrel)
666          endif !if (iso_HDO.gt.0) then
667#ifdef ISOTRAC
668           call iso_verif_traceur_justmass(d_xt(1,i,k),'concvl 316')
669#endif 
670    END DO
671  END DO         
672#endif           
673#endif
674
675#ifdef ISO
676      if ((iso_eau.gt.0).and.(bidouille_anti_divergence)) then
677        do k=1,klev   
678        do i=1,klon
679            d_xt(iso_eau,i,k)=d_q(i,k)
680        enddo !do i=1,klon
681        enddo !do k=1,klev               
682#ifdef ISOTRAC 
683        call iso_verif_traceur_jbid_vect(d_xt, &
684     &            klon,klev)   
685#endif         
686      endif !if ((iso_eau.gt.0).and.(bidouille_anti_divergence)) then 
687#endif
688
689  IF (iflag_con==3) THEN
690    DO i = 1,klon
691      IF (wbeff(i) > 100. .OR. wbeff(i) == 0 .OR. iflag(i) > 3) THEN
692        wbeff(i) = 0.
693        convoccur(i) = 0. 
694      ELSE
695        convoccur(i) = 1.
696      ENDIF
697    ENDDO
698  ENDIF
699
700  IF (iflag_con==30) THEN
701    DO itra = 1, ntra
702      DO k = 1, klev
703        DO i = 1, klon
704!RL!            d_tra(i,k,itra) =dtime*d_tra(i,k,itra)
705          d_tra(i, k, itra) = 0.
706        END DO
707      END DO
708    END DO
709  END IF
710
711!!AC!
712  IF (iflag_con==3) THEN
713    DO itra = 1, ntra
714      DO k = 1, klev
715        DO i = 1, klon
716!RL!            d_tra(i,k,itra) =dtime*d_tra(i,k,itra)
717          d_tra(i, k, itra) = 0.
718        END DO
719      END DO
720    END DO
721  END IF
722!!AC!
723
724  DO k = 1, klev
725    DO i = 1, klon
726      t1(i, k) = t(i, k) + d_t(i, k)
727      q1(i, k) = q(i, k) + d_q(i, k)
728! juste diag: pas besoin d'isos ici
729    END DO
730  END DO
731
732!                                                     !jyg
733  IF (iflag_con == 30 .OR. iflag_ice_thermo ==0) THEN
734! --Separation neige/pluie (pour diagnostics)         !jyg
735    DO k = 1, klev                                    !jyg
736      DO i = 1, klon                                  !jyg
737        IF (t1(i,k)<rtt) THEN                         !jyg
738          pmflxs(i, k) = Vprecip(i, k)                !jyg
739        ELSE                                          !jyg
740          pmflxr(i, k) = Vprecip(i, k)                !jyg
741        END IF                                        !jyg
742      END DO                                          !jyg
743    END DO                                            !jyg
744  ELSE
745    DO k = 1, klev                                    !jyg
746      DO i = 1, klon                                  !jyg
747        pmflxs(i, k) = Vprecipi(i, k)                 !jyg
748        pmflxr(i, k) = Vprecip(i, k)-Vprecipi(i, k)   !jyg
749      END DO                                          !jyg
750    END DO                                            !jyg
751  ENDIF
752
753! c      IF (if_ebil.ge.2) THEN
754! c        ztit='after convect'
755! c        CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
756! c     e      , t1,q1,ql,qs,u,v,paprs,pplay
757! c     s      , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
758! c         call diagphy(paire,ztit,ip_ebil
759! c     e      , zero_v, zero_v, zero_v, zero_v, zero_v
760! c     e      , zero_v, rain, zero_v, ztsol
761! c     e      , d_h_vcol, d_qt, d_ec
762! c     s      , fs_bound, fq_bound )
763! c      END IF
764
765
766! les traceurs ne sont pas mis dans cette version de convect4:
767  IF (iflag_con==4) THEN
768    DO itra = 1, ntra
769      DO k = 1, klev
770        DO i = 1, klon
771          d_tra(i, k, itra) = 0.
772        END DO
773      END DO
774    END DO
775  END IF
776! print*, 'concvl->: dd_t,dd_q ',dd_t(1,1),dd_q(1,1)
777
778  DO k = 1, klev
779    DO i = 1, klon
780      dtvpdt1(i, k) = 0.
781      dtvpdq1(i, k) = 0.
782    END DO
783  END DO
784  DO i = 1, klon
785    dplcldt(i) = 0.
786    dplcldr(i) = 0.
787  END DO
788
789  IF (prt_level>=20) THEN
790    DO k = 1, klev
791! print*,'physiq apres_add_con i k it d_u d_v d_t d_q qdl0',igout, &
792!         k,itap,d_u_con(igout,k) ,d_v_con(igout,k), d_t_con(igout,k), &
793!         d_q_con(igout,k),dql0(igout,k)
794! print*,'phys apres_add_con itap Ma cin ALE ALP wak t q undi t q', &
795!         itap,Ma(igout,k),cin(igout),ALE(igout), ALP(igout), &
796!         t_wake(igout,k),q_wake(igout,k),t_undi(igout,k),q_undi(igout,k)
797! print*,'phy apres_add_con itap CON rain snow EMA wk1 wk2 Vpp mip', &
798!         itap,rain_con(igout),snow_con(igout),ema_work1(igout,k), &
799!         ema_work2(igout,k),Vprecip(igout,k), mip(igout,k)
800! print*,'phy apres_add_con itap upwd dnwd dnwd0 cape tvp Tconv ', &
801!         itap,upwd(igout,k),dnwd(igout,k),dnwd0(igout,k),cape(igout), &
802!         tvp(igout,k),Tconv(igout,k)
803! print*,'phy apres_add_con itap dtvpdt dtvdq dplcl dplcldr qcondc', &
804!         itap,dtvpdt1(igout,k),dtvpdq1(igout,k),dplcldt(igout), &
805!         dplcldr(igout),qcondc(igout,k)
806! print*,'phy apres_add_con itap wd pmflxr Kpmflxr Kp1 Kpmflxs Kp1', &
807!         itap,wd(igout),pmflxr(igout,k),pmflxr(igout,k+1),pmflxs(igout,k), &
808!         pmflxs(igout,k+1)
809! print*,'phy apres_add_con itap da phi mp ftd fqd lalim wgth', &
810!         itap,da(igout,k),phi(igout,k,k),mp(igout,k),ftd(igout,k), &
811!         fqd(igout,k),lalim_conv(igout),wght_th(igout,k)
812    END DO
813  END IF !(prt_level.EQ.20) THEN
814
815  RETURN
816END SUBROUTINE concvl
817
Note: See TracBrowser for help on using the repository browser.