source: LMDZ6/branches/contrails/libf/phylmdiso/concvl.F90 @ 5418

Last change on this file since 5418 was 5304, checked in by abarral, 7 weeks ago

Turn YOMCST2.h.h into module

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