source: LMDZ6/branches/Amaury_dev/libf/phylmdiso/concvl.F90 @ 5136

Last change on this file since 5136 was 5117, checked in by abarral, 2 months ago

rename modules properly lmdz_*
move some unused files to obsolete/
(lint) uppercase fortran keywords

  • 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_qcomp, 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, detrain, &
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 lmdz_print_control, 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_qcomp, 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, detrain
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!  Local
218!  ----
219  REAL, DIMENSION(klon,klev)                    :: em_p
220  REAL, DIMENSION(klon,klev+1)                  :: em_ph
221  REAL                                          :: em_sig1feed ! sigma at lower bound of feeding layer
222  REAL                                          :: em_sig2feed ! sigma at upper bound of feeding layer
223  REAL, DIMENSION(klev)                         :: em_wght ! weight density determining the feeding mixture
224  REAL, DIMENSION(klon,klev+1)                  :: Vprecipi                       !jyg
225!on enleve le save
226! SAVE em_sig1feed,em_sig2feed,em_wght
227
228  REAL, DIMENSION(klon)                         :: rflag
229  REAL, DIMENSION(klon)                         :: plim1, plim2
230  REAL, DIMENSION(klon)                         :: ptop2
231  REAL, DIMENSION(klon,klev)                    :: asupmax
232  REAL, DIMENSION(klon)                         :: supmax0, asupmaxmin
233  REAL                                          :: zx_t, zdelta, zx_qs, zcor
234
235!   INTEGER iflag_mix
236!   SAVE iflag_mix
237  INTEGER                                       :: noff, minorig
238  INTEGER                                       :: i,j, k, itra
239  REAL, DIMENSION(klon,klev)                    :: qs, qs_wake
240!LF          SAVE cbmf
241!IM/JYG      REAL, SAVE, ALLOCATABLE :: cbmf(:)
242!!!$OMP THREADPRIVATE(cbmf)!
243  REAL, DIMENSION(klon)                         :: cbmflast
244#ifdef ISO
245REAL, DIMENSION(ntraciso,klon,klev+1)                  :: xtVprecipi
246  INTEGER                                       :: ixt
247#endif
248
249
250! Variables supplementaires liees au bilan d'energie
251! Real paire(klon)
252!LF      Real ql(klon,klev)
253! Save paire
254!LF      Save ql
255!LF      Real t1(klon,klev),q1(klon,klev)
256!LF      Save t1,q1
257! Data paire /1./
258  REAL, SAVE, ALLOCATABLE :: ql(:, :), q1(:, :), t1(:, :)
259!$OMP THREADPRIVATE(ql, q1, t1)
260        ! pas besoin d'isos ici
261
262! Variables liees au bilan d'energie et d'enthAlpi
263  REAL ztsol(klon)
264  REAL        h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot, &
265              h_qs_tot, qw_tot, ql_tot, qs_tot, ec_tot
266  SAVE        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!$OMP THREADPRIVATE(h_vcol_tot, h_dair_tot, h_qw_tot, h_ql_tot)
269!$OMP THREADPRIVATE(h_qs_tot, qw_tot, ql_tot, qs_tot , ec_tot)
270  REAL        d_h_vcol, d_h_dair, d_qt, d_qw, d_ql, d_qs, d_ec
271  REAL        d_h_vcol_phy
272  REAL        fs_bound, fq_bound
273  SAVE        d_h_vcol_phy
274!$OMP THREADPRIVATE(d_h_vcol_phy)
275  REAL        zero_v(klon)
276  CHARACTER *15 ztit
277  INTEGER     ip_ebil ! PRINT level for energy conserv. diag.
278  SAVE        ip_ebil
279  DATA        ip_ebil/2/
280!$OMP THREADPRIVATE(ip_ebil)
281  INTEGER     if_ebil ! level for energy conserv. dignostics
282  SAVE        if_ebil
283  DATA        if_ebil/2/
284!$OMP THREADPRIVATE(if_ebil)
285!+jld ec_conser
286  REAL d_t_ec(klon, klev) ! tendance du a la conersion Ec -> E thermique
287  REAL zrcpd
288!-jld ec_conser
289!LF
290  INTEGER nloc
291  LOGICAL, SAVE            :: first = .TRUE.
292!$OMP THREADPRIVATE(first)
293  INTEGER, SAVE            :: itap, igout
294!$OMP THREADPRIVATE(itap, igout)
295
296
297  include "YOMCST.h"
298  include "YOMCST2.h"
299  include "YOETHF.h"
300  include "FCTTRE.h"
301!jyg<
302  include "conema3.h"
303!>jyg
304
305  IF (first) THEN
306! Allocate some variables LF 04/2008
307
308!IM/JYG allocate(cbmf(klon))
309    ALLOCATE (ql(klon,klev))
310    ALLOCATE (t1(klon,klev))
311    ALLOCATE (q1(klon,klev))
312
313    convoccur(:) = 0.
314
315    itap = 0
316    igout = klon/2 + 1/klon
317  END IF
318! Incrementer le compteur de la physique
319  itap = itap + 1
320
321! Copy T into Tconv
322  DO k = 1, klev
323    DO i = 1, klon
324      Tconv(i, k) = t(i, k)
325    END DO
326  END DO
327
328  IF (if_ebil>=1) THEN
329    DO i = 1, klon
330      ztsol(i) = t(i, 1)
331      zero_v(i) = 0.
332      DO k = 1, klev
333        ql(i, k) = 0.
334      END DO
335    END DO
336  END IF
337
338! ym
339  snow(:) = 0
340#ifdef ISO
341      xtsnow(:,:)=0
342#endif
343
344  IF (first) THEN
345    first = .FALSE.
346
347! ===========================================================================
348! READ IN PARAMETERS FOR THE CLOSURE AND THE MIXING DISTRIBUTION
349! ===========================================================================
350
351    IF (iflag_con==3) THEN
352!      CALL cv3_inicp()
353      CALL cv3_inip()
354    END IF
355
356! ===========================================================================
357! READ IN PARAMETERS FOR CONVECTIVE INHIBITION BY TROPOS. DRYNESS
358! ===========================================================================
359
360! c$$$         open (56,file='supcrit.data')
361! c$$$         read (56,*) Supcrit1, Supcrit2
362! c$$$         close (56)
363
364    IF (prt_level>=10) WRITE (lunout, *) 'supcrit1, supcrit2', supcrit1, supcrit2
365
366! ===========================================================================
367! Initialisation pour les bilans d'eau et d'energie
368! ===========================================================================
369    IF (if_ebil>=1) d_h_vcol_phy = 0.
370
371    DO i = 1, klon
372      cbmf(i) = 0.
373!!          plcl(i) = 0.
374      sigd(i) = 0.
375    END DO
376  END IF !(first)
377
378! Initialisation a chaque pas de temps
379  plfc(:) = 0.
380  wbeff(:) = 100.
381  plcl(:) = 0.
382
383  DO k = 1, klev + 1
384    DO i = 1, klon
385      em_ph(i, k) = paprs(i, k)/100.0
386      pmflxr(i, k) = 0.
387      pmflxs(i, k) = 0.
388    END DO
389  END DO
390
391  DO k = 1, klev
392    DO i = 1, klon
393      em_p(i, k) = pplay(i, k)/100.0
394    END DO
395  END DO
396
397
398! Feeding layer
399
400  em_sig1feed = 1.
401!jyg<
402!  em_sig2feed = 0.97
403  em_sig2feed = cvl_sig2feed
404!>jyg
405! em_sig2feed = 0.8
406! Relative Weight densities
407  DO k = 1, klev
408    em_wght(k) = 1.
409  END DO
410!CRtest: couche alim des tehrmiques ponderee par a*
411! DO i = 1, klon
412! do k=1,lalim_conv(i)
413! em_wght(k)=wght_th(i,k)
414! PRINT*,'em_wght=',em_wght(k),wght_th(i,k)
415! END DO
416! END DO
417
418  IF (iflag_con==4) THEN
419    DO k = 1, klev
420      DO i = 1, klon
421        zx_t = t(i, k)
422        zdelta = max(0., sign(1.,rtt-zx_t))
423        zx_qs = min(0.5, r2es*foeew(zx_t,zdelta)/em_p(i,k)/100.0)
424        zcor = 1./(1.-retv*zx_qs)
425        qs(i, k) = zx_qs*zcor
426      END DO
427      DO i = 1, klon
428        zx_t = t_wake(i, k)
429        zdelta = max(0., sign(1.,rtt-zx_t))
430        zx_qs = min(0.5, r2es*foeew(zx_t,zdelta)/em_p(i,k)/100.0)
431        zcor = 1./(1.-retv*zx_qs)
432        qs_wake(i, k) = zx_qs*zcor
433      END DO
434    END DO
435  ELSE ! iflag_con=3 (modif de puristes qui fait la diffce pour la convergence numerique)
436    DO k = 1, klev
437      DO i = 1, klon
438        zx_t = t(i, k)
439        zdelta = max(0., sign(1.,rtt-zx_t))
440        zx_qs = r2es*foeew(zx_t, zdelta)/em_p(i, k)/100.0
441        zx_qs = min(0.5, zx_qs)
442        zcor = 1./(1.-retv*zx_qs)
443        zx_qs = zx_qs*zcor
444        qs(i, k) = zx_qs
445      END DO
446      DO i = 1, klon
447        zx_t = t_wake(i, k)
448        zdelta = max(0., sign(1.,rtt-zx_t))
449        zx_qs = r2es*foeew(zx_t, zdelta)/em_p(i, k)/100.0
450        zx_qs = min(0.5, zx_qs)
451        zcor = 1./(1.-retv*zx_qs)
452        zx_qs = zx_qs*zcor
453        qs_wake(i, k) = zx_qs
454      END DO
455    END DO
456  END IF ! iflag_con
457
458! ------------------------------------------------------------------
459
460! Main driver for convection:
461!                   iflag_con=3 -> nvlle version de KE (JYG)
462!                   iflag_con = 30  -> equivAlent to convect3
463!                   iflag_con = 4  -> equivAlent to convect1/2
464
465
466  IF (iflag_con==30) THEN
467
468 
469#ifdef ISO         
470#ifdef ISOVERIF
471       do k = 1, klev
472        do i = 1, klon               
473         do ixt=1,ntraciso         
474             CALL iso_verif_noNaN(xt(ixt,i,k),'concvl 394')
475         enddo
476        enddo !do i = 1, klon
477       enddo !do k = 1, klev       
478       IF (iso_eau.gt.0) THEN
479       do k = 1, klev
480        do i = 1, klon   
481          CALL iso_verif_egalite_choix(xt(iso_eau,i,k),q(i,k), &
482                  'concvl 174',errmax,errmaxrel)
483        enddo !do i = 1, klon
484       enddo !do k = 1, klev   
485       endif !if (iso_eau.gt.0) THEN
486       IF (iso_HDO.gt.0) THEN
487       do k = 1, klev
488        do i = 1, klon         
489         IF (q(i,k).gt.ridicule) THEN
490          CALL iso_verif_aberrant(xt(iso_hdo,i,k)/q(i,k),'concvl 175')
491         endif ! if (q(i,k).gt.ridicule) THEN
492        enddo
493       enddo   
494       endif !if (iso_eau.gt.0) THEN
495#ifdef ISOTRAC
496        do k = 1, klev
497        do i = 1, klon   
498           CALL iso_verif_traceur(xt(1,i,k),'concvl 218')
499        enddo
500        enddo
501#endif       
502       WRITE(*,*) 'concvl 170: avant appel cv_driver'
503#endif
504        ! ISOVERIF ! end verif       
505#endif
506
507! print *, '-> cv_driver'      !jyg
508    CALL cv_driver(klon, klev, klevp1, ntra, iflag_con, &
509                   t, q, qs, u, v, tra, &
510                   em_p, em_ph, iflag, &
511                   d_t, d_q, d_u, d_v, d_tra, rain, &
512                   Vprecip, cbmf, sig1, w01, & !jyg
513                   kbas, ktop, &
514                   dtime, Ma, upwd, dnwd, dnwdbis, qcondc, wd, cape, &
515                   da, phi, mp, phii, d1a, dam, sij, clw, elij, &       !RomP
516                   evap, ep, epmlmMm, eplaMm, &                         !RomP
517                   wdtrainA, wdtrainM, &                                !RomP
518                   epmax_diag & ! epmax_cape
519#ifdef ISO
520                   , xt,d_xt &
521                   , xtrain,xtVprecip &
522                   , xtevap,xtclw,xtwdtrainA &
523#ifdef DIAGISO
524                , qlp,xtlp,qvp,xtvp & ! juste diagnostique
525                , fq_detrainement,fq_ddft,fq_fluxmasse,fq_evapprecip &
526                , fxt_detrainement,fxt_ddft,fxt_fluxmasse,fxt_evapprecip &
527                , f_detrainement,q_detrainement,xt_detrainement &
528#endif     
529#endif
530                    )
531!           print *, 'cv_driver ->'      !jyg
532
533#ifdef ISO
534      ! verif
535#ifdef ISOVERIF
536       WRITE(*,*) 'concvl 463: après appel cv_driver'
537       do k = 1, klev
538        do i = 1, klon
539        IF (iso_eau.gt.0) THEN
540            CALL iso_verif_egalite(xt(iso_eau,i,k),q(i,k),'concvl 203')
541            CALL iso_verif_egalite(xt_wake(iso_eau,i,k),q_wake(i,k),'concvl 204')
542            CALL iso_verif_egalite(d_xt(iso_eau,i,k),d_q(i,k), &
543                  'concvl 452')
544         endif !if (iso_eau.gt.0) THEN
545#ifdef DIAGISO
546         do ixt=1,ntraciso
547            CALL iso_verif_noNaN(xt(ixt,i,k),'concvl 460')
548            CALL iso_verif_noNaN(xtlp(ixt,i,k),'concvl 295')
549            CALL iso_verif_noNaN(xtvp(ixt,i,k),'concvl 260')
550          enddo
551#endif                 
552        enddo
553       enddo       
554#ifdef ISOTRAC
555           CALL iso_verif_traceur_vect(xt,klon,klev,'concvl 218')
556           CALL iso_verif_trac_masse_vect(d_xt,klon,klev, &
557                 'concvl 464',errmax,errmaxrel)
558#endif           
559#endif
560       ! end verif       
561#endif
562
563    DO i = 1, klon
564      cbmf(i) = Ma(i, kbas(i))
565    END DO
566
567!RL
568    wght(:, :) = 0.
569    DO i = 1, klon
570      wght(i, 1) = 1.
571    END DO
572!RL
573
574  ELSE
575
576!LF   necessary for gathered fields
577    nloc = klon
578#ifdef ISOVERIF
579        WRITE(*,*) 'concvl 581: juste avant appel de cva_driver'
580#endif
581    CALL cva_driver(klon, klev, klev+1, ntra, nloc, k_upper_cv, &
582                    iflag_con, iflag_mix, iflag_ice_thermo, &
583                    iflag_clos, ok_conserv_q, dtime, cvl_comp_threshold, &
584                    t, q, qs, t_wake, q_wake, qs_wake, s_wake, u, v, tra, &
585                    em_p, em_ph, &
586                    Ale, Alp, omega, &
587                    em_sig1feed, em_sig2feed, em_wght, &
588                    iflag, d_t, d_q, d_qcomp, d_u, d_v, d_tra, rain, kbas, ktop, &
589                    cbmf, plcl, plfc, wbeff, sig1, w01, ptop2, sigd, &
590                    Ma, mip, Vprecip, Vprecipi, upwd, dnwd, dnwdbis, qcondc, wd, &
591                    cape, cin, tvp, &
592                    dd_t, dd_q, plim1, plim2, asupmax, supmax0, &
593                    asupmaxmin, lalim_conv, &
594!AC!+!RomP+jyg
595!!                   da,phi,mp,phii,d1a,dam,sij,clw,elij, &               ! RomP
596!!                   evap,ep,epmlmMm,eplaMm,                              ! RomP
597                    da, phi, mp, phii, d1a, dam, sij, wght, &           ! RomP+RL
598                    qta, clw, elij, evap, ep, epmlmMm, eplaMm, &        ! RomP+RL
599                    wdtrainA, wdtrainS, wdtrainM, qtc, sigt, detrain, &
600                    tau_cld_cv, coefw_cld_cv, &                         ! RomP,AJ
601!AC!+!RomP+jyg
602                    epmax_diag & ! epmax_cape
603#ifdef ISO
604                   ,xt,xt_wake,d_xt, xtrain &
605                   ,xtvprecip,xtvprecipi &
606                   ,xtclw,dd_xt,xtevap,xtwdtrainA &
607#ifdef DIAGISO     
608                , qlp,xtlp,qvp,xtvp &
609                , fq_detrainement,fq_ddft,fq_fluxmasse,fq_evapprecip &
610                , fxt_detrainement,fxt_ddft,fxt_fluxmasse,fxt_evapprecip &
611                , f_detrainement,q_detrainement,xt_detrainement &
612#endif     
613#endif
614                )
615  END IF
616! ------------------------------------------------------------------
617  IF (prt_level>=10) WRITE (lunout, *) ' cva_driver -> cbmf,plcl,plfc,wbeff, d_t, d_q ', &
618                                         cbmf(1), plcl(1), plfc(1), wbeff(1), d_t(1,1), d_q(1,1)
619
620  DO i = 1, klon
621    rain(i) = rain(i)/86400.
622    rflag(i) = iflag(i)
623#ifdef ISO
624       do ixt = 1, ntraciso
625        xtrain(ixt,i) = xtrain(ixt,i)/86400.
626       enddo
627#endif
628  END DO
629
630  DO k = 1, klev
631    DO i = 1, klon
632      d_t(i, k) = dtime*d_t(i, k)
633      d_q(i, k) = dtime*d_q(i, k)
634      d_u(i, k) = dtime*d_u(i, k)
635      d_v(i, k) = dtime*d_v(i, k)
636#ifdef ISO
637           do ixt = 1, ntraciso
638            d_xt(ixt,i,k) = dtime*d_xt(ixt,i,k)
639           enddo
640#endif
641    END DO
642  END DO
643
644             
645#ifdef ISO
646#ifdef ISOVERIF     
647!        k=1
648!        i=  538
649        WRITE(*,*) 'concvl 640'
650!        WRITE(*,*) 'q(i,k),d_q(i,k)=', q(i,k),d_q(i,k)
651!        WRITE(*,*) 'xt(iso_HDO,i,k),d_xt(iso_HDO,i,k)=', &
652!     &          xt(iso_HDO,i,k),d_xt(iso_HDO,i,k)
653  DO k = 1, klev
654    DO i = 1, klon
655           IF (iso_HDO.gt.0) THEN
656             IF (q(i,k).gt.ridicule) THEN
657                 CALL iso_verif_aberrant_encadre((xt(iso_HDO,i,k) &
658              +d_xt(iso_HDO,i,k))/(q(i,k)+d_q(i,k)),'concvl 250')
659             endif !if (q_seri(i,k).gt.ridicule) THEN
660          endif !if (iso_HDO.gt.0) THEN
661           IF (iso_eau.gt.0) THEN
662             CALL iso_verif_egalite_choix(d_xt(iso_eau,i,k), &
663                d_q(i,k),'concvl 530',errmax*dtime,errmaxrel)
664          endif !if (iso_HDO.gt.0) THEN
665#ifdef ISOTRAC
666           CALL iso_verif_traceur_justmass(d_xt(1,i,k),'concvl 316')
667#endif 
668    END DO
669  END DO         
670#endif           
671#endif
672
673#ifdef ISO
674      IF ((iso_eau.gt.0).AND.(bidouille_anti_divergence)) THEN
675        do k=1,klev   
676        do i=1,klon
677            d_xt(iso_eau,i,k)=d_q(i,k)
678        enddo !do i=1,klon
679        enddo !do k=1,klev               
680#ifdef ISOTRAC 
681        CALL iso_verif_traceur_jbid_vect(d_xt, &
682                  klon,klev)
683#endif         
684      endif !if ((iso_eau.gt.0).AND.(bidouille_anti_divergence)) THEN
685#endif
686
687  IF (iflag_con==3) THEN
688    DO i = 1,klon
689      IF (wbeff(i) > 100. .OR. wbeff(i) == 0 .OR. iflag(i) > 3) THEN
690        wbeff(i) = 0.
691        convoccur(i) = 0. 
692      ELSE
693        convoccur(i) = 1.
694      ENDIF
695    ENDDO
696  ENDIF
697
698  IF (iflag_con==30) THEN
699    DO itra = 1, ntra
700      DO k = 1, klev
701        DO i = 1, klon
702!RL!            d_tra(i,k,itra) =dtime*d_tra(i,k,itra)
703          d_tra(i, k, itra) = 0.
704        END DO
705      END DO
706    END DO
707  END IF
708
709!!AC!
710  IF (iflag_con==3) THEN
711    DO itra = 1, ntra
712      DO k = 1, klev
713        DO i = 1, klon
714!RL!            d_tra(i,k,itra) =dtime*d_tra(i,k,itra)
715          d_tra(i, k, itra) = 0.
716        END DO
717      END DO
718    END DO
719  END IF
720!!AC!
721
722  DO k = 1, klev
723    DO i = 1, klon
724      t1(i, k) = t(i, k) + d_t(i, k)
725      q1(i, k) = q(i, k) + d_q(i, k)
726! juste diag: pas besoin d'isos ici
727    END DO
728  END DO
729
730!                                                     !jyg
731  IF (iflag_con == 30 .OR. iflag_ice_thermo ==0) THEN
732! --Separation neige/pluie (pour diagnostics)         !jyg
733    DO k = 1, klev                                    !jyg
734      DO i = 1, klon                                  !jyg
735        IF (t1(i,k)<rtt) THEN                         !jyg
736          pmflxs(i, k) = Vprecip(i, k)                !jyg
737        ELSE                                          !jyg
738          pmflxr(i, k) = Vprecip(i, k)                !jyg
739        END IF                                        !jyg
740      END DO                                          !jyg
741    END DO                                            !jyg
742  ELSE
743    DO k = 1, klev                                    !jyg
744      DO i = 1, klon                                  !jyg
745        pmflxs(i, k) = Vprecipi(i, k)                 !jyg
746        pmflxr(i, k) = Vprecip(i, k)-Vprecipi(i, k)   !jyg
747      END DO                                          !jyg
748    END DO                                            !jyg
749  ENDIF
750
751! c      IF (if_ebil.ge.2) THEN
752! c        ztit='after convect'
753! c        CALL diagetpq(paire,ztit,ip_ebil,2,2,dtime
754! c     e      , t1,q1,ql,qs,u,v,paprs,pplay
755! c     s      , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
756! c         CALL diagphy(paire,ztit,ip_ebil
757! c     e      , zero_v, zero_v, zero_v, zero_v, zero_v
758! c     e      , zero_v, rain, zero_v, ztsol
759! c     e      , d_h_vcol, d_qt, d_ec
760! c     s      , fs_bound, fq_bound )
761! c      END IF
762
763
764! les traceurs ne sont pas mis dans cette version de convect4:
765  IF (iflag_con==4) THEN
766    DO itra = 1, ntra
767      DO k = 1, klev
768        DO i = 1, klon
769          d_tra(i, k, itra) = 0.
770        END DO
771      END DO
772    END DO
773  END IF
774! PRINT*, 'concvl->: dd_t,dd_q ',dd_t(1,1),dd_q(1,1)
775
776  DO k = 1, klev
777    DO i = 1, klon
778      dtvpdt1(i, k) = 0.
779      dtvpdq1(i, k) = 0.
780    END DO
781  END DO
782  DO i = 1, klon
783    dplcldt(i) = 0.
784    dplcldr(i) = 0.
785  END DO
786
787  IF (prt_level>=20) THEN
788    DO k = 1, klev
789! PRINT*,'physiq apres_add_con i k it d_u d_v d_t d_q qdl0',igout, &
790!         k,itap,d_u_con(igout,k) ,d_v_con(igout,k), d_t_con(igout,k), &
791!         d_q_con(igout,k),dql0(igout,k)
792! PRINT*,'phys apres_add_con itap Ma cin ALE ALP wak t q undi t q', &
793!         itap,Ma(igout,k),cin(igout),ALE(igout), ALP(igout), &
794!         t_wake(igout,k),q_wake(igout,k),t_undi(igout,k),q_undi(igout,k)
795! PRINT*,'phy apres_add_con itap CON rain snow EMA wk1 wk2 Vpp mip', &
796!         itap,rain_con(igout),snow_con(igout),ema_work1(igout,k), &
797!         ema_work2(igout,k),Vprecip(igout,k), mip(igout,k)
798! PRINT*,'phy apres_add_con itap upwd dnwd dnwd0 cape tvp Tconv ', &
799!         itap,upwd(igout,k),dnwd(igout,k),dnwd0(igout,k),cape(igout), &
800!         tvp(igout,k),Tconv(igout,k)
801! PRINT*,'phy apres_add_con itap dtvpdt dtvdq dplcl dplcldr qcondc', &
802!         itap,dtvpdt1(igout,k),dtvpdq1(igout,k),dplcldt(igout), &
803!         dplcldr(igout),qcondc(igout,k)
804! PRINT*,'phy apres_add_con itap wd pmflxr Kpmflxr Kp1 Kpmflxs Kp1', &
805!         itap,wd(igout),pmflxr(igout,k),pmflxr(igout,k+1),pmflxs(igout,k), &
806!         pmflxs(igout,k+1)
807! PRINT*,'phy apres_add_con itap da phi mp ftd fqd lalim wgth', &
808!         itap,da(igout,k),phi(igout,k,k),mp(igout,k),ftd(igout,k), &
809!         fqd(igout,k),lalim_conv(igout),wght_th(igout,k)
810    END DO
811  END IF !(prt_level.EQ.20) THEN
812
813
814END SUBROUTINE concvl
815
Note: See TracBrowser for help on using the repository browser.