source: LMDZ6/branches/Amaury_dev/libf/phylmdiso/lmdz_concvl.F90 @ 5209

Last change on this file since 5209 was 5160, checked in by abarral, 7 weeks ago

Put .h into modules

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