LMDZ
conlmd.F90
Go to the documentation of this file.
1 
2 ! $Id: conlmd.F90 2346 2015-08-21 15:13:46Z emillour $
3 
4 SUBROUTINE conlmd(dtime, paprs, pplay, t, q, conv_q, d_t, d_q, rain, snow, &
5  ibas, itop)
6  USE dimphy
7  IMPLICIT NONE
8  ! ======================================================================
9  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818
10  ! Objet: Schema de convection utilis'e dans le modele du LMD
11  ! Ajustement humide (Manabe) + Ajustement convectif (Kuo)
12  ! ======================================================================
13  include "YOMCST.h"
14  include "YOETHF.h"
15 
16  ! Arguments:
17 
18  REAL dtime ! pas d'integration (s)
19  REAL paprs(klon, klev+1) ! pression inter-couche (Pa)
20  REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
21  REAL t(klon, klev) ! temperature (K)
22  REAL q(klon, klev) ! humidite specifique (kg/kg)
23  REAL conv_q(klon, klev) ! taux de convergence humidite (g/g/s)
24 
25  REAL d_t(klon, klev) ! incrementation temperature
26  REAL d_q(klon, klev) ! incrementation humidite
27  REAL rain(klon) ! pluies (mm/s)
28  REAL snow(klon) ! neige (mm/s)
29  INTEGER ibas(klon) ! niveau du bas
30  INTEGER itop(klon) ! niveau du haut
31 
32  LOGICAL usekuo ! utiliser convection profonde (schema Kuo)
33  parameter(usekuo=.true.)
34 
35  REAL d_t_bis(klon, klev)
36  REAL d_q_bis(klon, klev)
37  REAL rain_bis(klon)
38  REAL snow_bis(klon)
39  INTEGER ibas_bis(klon)
40  INTEGER itop_bis(klon)
41  REAL d_ql(klon, klev), d_ql_bis(klon, klev)
42  REAL rneb(klon, klev), rneb_bis(klon, klev)
43 
44  INTEGER i, k
45  REAL zlvdcp, zlsdcp, zdelta, zz, za, zb
46 
47  ! cc CALL fiajh ! ancienne version de Convection Manabe
48  CALL conman & ! nouvelle version de Convection
49  ! Manabe
50  (dtime, paprs, pplay, t, q, d_t, d_q, d_ql, rneb, rain, snow, ibas, itop)
51 
52  IF (usekuo) THEN
53  ! cc CALL fiajc ! ancienne version de Convection Kuo
54  CALL conkuo & ! nouvelle version de Convection
55  ! Kuo
56  (dtime, paprs, pplay, t, q, conv_q, d_t_bis, d_q_bis, d_ql_bis, &
57  rneb_bis, rain_bis, snow_bis, ibas_bis, itop_bis)
58  DO k = 1, klev
59  DO i = 1, klon
60  d_t(i, k) = d_t(i, k) + d_t_bis(i, k)
61  d_q(i, k) = d_q(i, k) + d_q_bis(i, k)
62  d_ql(i, k) = d_ql(i, k) + d_ql_bis(i, k)
63  END DO
64  END DO
65  DO i = 1, klon
66  rain(i) = rain(i) + rain_bis(i)
67  snow(i) = snow(i) + snow_bis(i)
68  ibas(i) = min(ibas(i), ibas_bis(i))
69  itop(i) = max(itop(i), itop_bis(i))
70  END DO
71  END IF
72 
73  ! L'eau liquide convective est dispersee dans l'air:
74 
75  DO k = 1, klev
76  DO i = 1, klon
77  zlvdcp = rlvtt/rcpd/(1.0+rvtmp2*q(i,k))
78  zlsdcp = rlstt/rcpd/(1.0+rvtmp2*q(i,k))
79  zdelta = max(0., sign(1.,rtt-t(i,k)))
80  zz = d_ql(i, k) ! re-evap. de l'eau liquide
81  zb = max(0.0, zz)
82  za = -max(0.0, zz)*(zlvdcp*(1.-zdelta)+zlsdcp*zdelta)
83  d_t(i, k) = d_t(i, k) + za
84  d_q(i, k) = d_q(i, k) + zb
85  END DO
86  END DO
87 
88  RETURN
89 END SUBROUTINE conlmd
90 SUBROUTINE conman(dtime, paprs, pplay, t, q, d_t, d_q, d_ql, rneb, rain, &
91  snow, ibas, itop)
92  USE dimphy
93  IMPLICIT NONE
94  ! ======================================================================
95  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19970324
96  ! Objet: ajustement humide convectif avec la possibilite de faire
97  ! l'ajustement sur une fraction de la maille.
98  ! Methode: On impose une distribution uniforme pour la vapeur d'eau
99  ! au sein d'une maille. On applique la procedure d'ajustement
100  ! successivement a la totalite, 75%, 50%, 25% et 5% de la maille
101  ! jusqu'a ce que l'ajustement a lieu. J'espere que ceci augmente
102  ! les activites convectives et corrige le biais "trop froid et sec"
103  ! du modele.
104  ! ======================================================================
105  include "YOMCST.h"
106 
107  REAL dtime ! pas d'integration (s)
108  REAL t(klon, klev) ! temperature (K)
109  REAL q(klon, klev) ! humidite specifique (kg/kg)
110  REAL paprs(klon, klev+1) ! pression inter-couche (Pa)
111  REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
112 
113  REAL d_t(klon, klev) ! incrementation temperature
114  REAL d_q(klon, klev) ! incrementation humidite
115  REAL d_ql(klon, klev) ! incrementation eau liquide
116  REAL rneb(klon, klev) ! nebulosite
117  REAL rain(klon) ! pluies (mm/s)
118  REAL snow(klon) ! neige (mm/s)
119  INTEGER ibas(klon) ! niveau du bas
120  INTEGER itop(klon) ! niveau du haut
121 
122  LOGICAL afaire(klon) ! .TRUE. implique l'ajustement
123  LOGICAL accompli(klon) ! .TRUE. si l'ajustement est effectif
124 
125  INTEGER nb ! nombre de sous-fractions a considere
126  parameter(nb=1)
127  ! cc PARAMETER (nb=3)
128 
129  REAL ratqs ! largeur de la distribution pour vapeur d'eau
130  parameter(ratqs=0.05)
131 
132  REAL w_q(klon, klev)
133  REAL w_d_t(klon, klev), w_d_q(klon, klev), w_d_ql(klon, klev)
134  REAL w_rneb(klon, klev)
135  REAL w_rain(klon), w_snow(klon)
136  INTEGER w_ibas(klon), w_itop(klon)
137  REAL zq1, zq2
138  INTEGER i, k, n
139 
140  REAL t_coup
141  parameter(t_coup=234.0)
142  REAL zdp1, zdp2
143  REAL zqs1, zqs2, zdqs1, zdqs2
144  REAL zgamdz
145  REAL zflo ! flotabilite
146  REAL zsat ! sur-saturation
147  REAL zdelta, zcor, zcvm5
148  LOGICAL imprim
149 
150  INTEGER ncpt
151  SAVE ncpt
152  !$OMP THREADPRIVATE(ncpt)
153  REAL frac(nb) ! valeur de la maille fractionnelle
154  SAVE frac
155  !$OMP THREADPRIVATE(frac)
156  INTEGER opt_cld(nb) ! option pour le modele nuageux
157  SAVE opt_cld
158  !$OMP THREADPRIVATE(opt_cld)
159  LOGICAL appel1er
160  SAVE appel1er
161  !$OMP THREADPRIVATE(appel1er)
162 
163  ! Fonctions thermodynamiques:
164 
165  include "YOETHF.h"
166  include "FCTTRE.h"
167 
168  DATA frac/1.0/
169  DATA opt_cld/4/
170  ! cc DATA frac / 1.0, 0.50, 0.25/
171  ! cc DATA opt_cld / 4, 4, 4/
172 
173  DATA appel1er/.true./
174  DATA ncpt/0/
175 
176  IF (appel1er) THEN
177  print *, 'conman, nb:', nb
178  print *, 'conman, frac:', frac
179  print *, 'conman, opt_cld:', opt_cld
180  appel1er = .false.
181  END IF
182 
183  ! Initialiser les sorties a zero:
184 
185  DO k = 1, klev
186  DO i = 1, klon
187  d_t(i, k) = 0.0
188  d_q(i, k) = 0.0
189  d_ql(i, k) = 0.0
190  rneb(i, k) = 0.0
191  END DO
192  END DO
193  DO i = 1, klon
194  ibas(i) = klev
195  itop(i) = 1
196  rain(i) = 0.0
197  snow(i) = 0.0
198  END DO
199 
200  ! S'il n'y a pas d'instabilite conditionnelle,
201  ! pas la penne de se fatiguer:
202 
203  DO i = 1, klon
204  afaire(i) = .false.
205  END DO
206  DO k = 1, klev - 1
207  DO i = 1, klon
208  IF (thermcep) THEN
209  zdelta = max(0., sign(1.,rtt-t(i,k)))
210  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
211  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*q(i,k))
212  zqs1 = r2es*foeew(t(i,k), zdelta)/pplay(i, k)
213  zqs1 = min(0.5, zqs1)
214  zcor = 1./(1.-retv*zqs1)
215  zqs1 = zqs1*zcor
216  zdqs1 = foede(t(i,k), zdelta, zcvm5, zqs1, zcor)
217 
218  zdelta = max(0., sign(1.,rtt-t(i,k+1)))
219  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
220  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*q(i,k+1))
221  zqs2 = r2es*foeew(t(i,k+1), zdelta)/pplay(i, k+1)
222  zqs2 = min(0.5, zqs2)
223  zcor = 1./(1.-retv*zqs2)
224  zqs2 = zqs2*zcor
225  zdqs2 = foede(t(i,k+1), zdelta, zcvm5, zqs2, zcor)
226  ELSE
227  IF (t(i,k)<t_coup) THEN
228  zqs1 = qsats(t(i,k))/pplay(i, k)
229  zdqs1 = dqsats(t(i,k), zqs1)
230 
231  zqs2 = qsats(t(i,k+1))/pplay(i, k+1)
232  zdqs2 = dqsats(t(i,k+1), zqs2)
233  ELSE
234  zqs1 = qsatl(t(i,k))/pplay(i, k)
235  zdqs1 = dqsatl(t(i,k), zqs1)
236 
237  zqs2 = qsatl(t(i,k+1))/pplay(i, k+1)
238  zdqs2 = dqsatl(t(i,k+1), zqs2)
239  END IF
240  END IF
241  zdp1 = paprs(i, k) - paprs(i, k+1)
242  zdp2 = paprs(i, k+1) - paprs(i, k+2)
243  zgamdz = -(pplay(i,k)-pplay(i,k+1))/paprs(i, k+1)/rcpd*(rd*(t(i, &
244  k)*zdp1+t(i,k+1)*zdp2)/(zdp1+zdp2)+rlvtt*(zqs1*zdp1+zqs2*zdp2)/(zdp1+ &
245  zdp2))/(1.0+(zdqs1*zdp1+zdqs2*zdp2)/(zdp1+zdp2))
246  zflo = t(i, k) + zgamdz - t(i, k+1)
247  zsat = (q(i,k)-zqs1)*zdp1 + (q(i,k+1)-zqs2)*zdp2
248  IF (zflo>0.0) afaire(i) = .true.
249  ! erreur IF (zflo.GT.0.0 .AND. zsat.GT.0.0) afaire(i) = .TRUE.
250  END DO
251  END DO
252 
253  imprim = mod(ncpt, 48) == 0
254  DO n = 1, nb
255 
256  DO k = 1, klev
257  DO i = 1, klon
258  IF (afaire(i)) THEN
259  zq1 = q(i, k)*(1.0-ratqs)
260  zq2 = q(i, k)*(1.0+ratqs)
261  w_q(i, k) = zq2 - frac(n)/2.0*(zq2-zq1)
262  END IF
263  END DO
264  END DO
265 
266  CALL conmanv(dtime, paprs, pplay, t, w_q, afaire, opt_cld(n), w_d_t, &
267  w_d_q, w_d_ql, w_rneb, w_rain, w_snow, w_ibas, w_itop, accompli, &
268  imprim)
269  DO k = 1, klev
270  DO i = 1, klon
271  IF (afaire(i) .AND. accompli(i)) THEN
272  d_t(i, k) = w_d_t(i, k)*frac(n)
273  d_q(i, k) = w_d_q(i, k)*frac(n)
274  d_ql(i, k) = w_d_ql(i, k)*frac(n)
275  IF (nint(w_rneb(i,k))==1) rneb(i, k) = frac(n)
276  END IF
277  END DO
278  END DO
279  DO i = 1, klon
280  IF (afaire(i) .AND. accompli(i)) THEN
281  rain(i) = w_rain(i)*frac(n)
282  snow(i) = w_snow(i)*frac(n)
283  ibas(i) = min(ibas(i), w_ibas(i))
284  itop(i) = max(itop(i), w_itop(i))
285  END IF
286  END DO
287  DO i = 1, klon
288  IF (afaire(i) .AND. accompli(i)) afaire(i) = .false.
289  END DO
290 
291  END DO
292 
293  ncpt = ncpt + 1
294 
295  RETURN
296 END SUBROUTINE conman
297 SUBROUTINE conmanv(dtime, paprs, pplay, t, q, afaire, opt_cld, d_t, d_q, &
298  d_ql, rneb, rain, snow, ibas, itop, accompli, imprim)
299  USE dimphy
300  IMPLICIT NONE
301  ! ======================================================================
302  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818
303  ! Objet: ajustement humide (convection proposee par Manabe).
304  ! Pour une colonne verticale, il peut avoir plusieurs blocs
305  ! necessitant l'ajustement. ibas est le bas du plus bas bloc
306  ! et itop est le haut du plus haut bloc
307  ! ======================================================================
308  include "YOMCST.h"
309 
310  ! Arguments:
311 
312  REAL dtime ! pas d'integration (s)
313  REAL t(klon, klev) ! temperature (K)
314  REAL q(klon, klev) ! humidite specifique (kg/kg)
315  REAL paprs(klon, klev+1) ! pression inter-couche (Pa)
316  REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
317  INTEGER opt_cld ! comment traiter l'eau liquide
318  LOGICAL afaire(klon) ! .TRUE. si le point est a faire (Input)
319  LOGICAL imprim ! .T. pour imprimer quelques diagnostiques
320 
321  REAL d_t(klon, klev) ! incrementation temperature
322  REAL d_q(klon, klev) ! incrementation humidite
323  REAL d_ql(klon, klev) ! incrementation eau liquide
324  REAL rneb(klon, klev) ! nebulosite
325  REAL rain(klon) ! pluies (mm/s)
326  REAL snow(klon) ! neige (mm/s)
327  INTEGER ibas(klon) ! niveau du bas
328  INTEGER itop(klon) ! niveau du haut
329  LOGICAL accompli(klon) ! .TRUE. si l'ajustement a eu lieu (Output)
330 
331  ! Quelques options:
332 
333  LOGICAL new_top ! re-calculer sommet quand re-ajustement est fait
334  parameter(new_top=.false.)
335  LOGICAL evap_prec ! evaporation de pluie au-dessous de convection
336  parameter(evap_prec=.true.)
337  REAL coef_eva
338  parameter(coef_eva=1.0e-05)
339  REAL t_coup
340  parameter(t_coup=234.0)
341  REAL seuil_vap
342  parameter(seuil_vap=1.0e-10)
343  LOGICAL old_tau ! implique precip nulle, si vrai.
344  parameter(old_tau=.false.)
345  REAL toliq(klon) ! rapport entre l'eau nuageuse et l'eau precipitante
346  REAL dpmin, tomax !Epaisseur faible, rapport eau liquide plus grande
347  parameter(dpmin=0.15, tomax=0.97)
348  REAL dpmax, tomin !Epaisseur grande, rapport eau liquide plus faible
349  parameter(dpmax=0.30, tomin=0.05)
350  REAL deep_sig, deep_to ! au dela de deep_sig, utiliser deep_to
351  parameter(deep_sig=0.50, deep_to=0.05)
352  LOGICAL exigent ! implique un calcul supplementaire pour Qs
353  parameter(exigent=.false.)
354 
355  INTEGER kbase
356  parameter(kbase=0)
357 
358  ! Variables locales:
359 
360  INTEGER nexpo
361  INTEGER i, k, k1min, k1max, k2min, k2max, is
362  REAL zgamdz(klon, klev-1)
363  REAL zt(klon, klev), zq(klon, klev)
364  REAL zqs(klon, klev), zdqs(klon, klev)
365  REAL zqmqsdp(klon, klev)
366  REAL ztnew(klon, klev), zqnew(klon, klev)
367  REAL zcond(klon), zvapo(klon), zrapp(klon)
368  REAL zrfl(klon), zrfln, zqev, zqevt
369  REAL zsat(klon) ! sur-saturation
370  REAL zflo(klon) ! flotabilite
371  REAL za(klon), zb(klon), zc(klon)
372  INTEGER k1(klon), k2(klon)
373  REAL zdelta, zcor, zcvm5
374  REAL delp(klon, klev)
375  LOGICAL possible(klon), todo(klon), etendre(klon)
376  LOGICAL aller(klon), todobis(klon)
377  REAL zalfa
378  INTEGER nbtodo, nbdone
379 
380  ! Fonctions thermodynamiques:
381 
382  include "YOETHF.h"
383  include "FCTTRE.h"
384 
385  DO k = 1, klev
386  DO i = 1, klon
387  delp(i, k) = paprs(i, k) - paprs(i, k+1)
388  END DO
389  END DO
390 
391  ! Initialiser les sorties a zero
392 
393  DO k = 1, klev
394  DO i = 1, klon
395  d_t(i, k) = 0.0
396  d_q(i, k) = 0.0
397  d_ql(i, k) = 0.0
398  rneb(i, k) = 0.0
399  END DO
400  END DO
401  DO i = 1, klon
402  ibas(i) = klev
403  itop(i) = 1
404  rain(i) = 0.0
405  snow(i) = 0.0
406  accompli(i) = .false.
407  END DO
408 
409  ! Preparations
410 
411  DO k = 1, klev
412  DO i = 1, klon
413  IF (afaire(i)) THEN
414  zt(i, k) = t(i, k)
415  zq(i, k) = q(i, k)
416 
417  ! Calculer Qs et L/Cp*dQs/dT
418 
419  IF (thermcep) THEN
420  zdelta = max(0., sign(1.,rtt-zt(i,k)))
421  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
422  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*zq(i,k))
423  zqs(i, k) = r2es*foeew(zt(i,k), zdelta)/pplay(i, k)
424  zqs(i, k) = min(0.5, zqs(i,k))
425  zcor = 1./(1.-retv*zqs(i,k))
426  zqs(i, k) = zqs(i, k)*zcor
427  zdqs(i, k) = foede(zt(i,k), zdelta, zcvm5, zqs(i,k), zcor)
428  ELSE
429  IF (zt(i,k)<t_coup) THEN
430  zqs(i, k) = qsats(zt(i,k))/pplay(i, k)
431  zdqs(i, k) = dqsats(zt(i,k), zqs(i,k))
432  ELSE
433  zqs(i, k) = qsatl(zt(i,k))/pplay(i, k)
434  zdqs(i, k) = dqsatl(zt(i,k), zqs(i,k))
435  END IF
436  END IF
437 
438  ! Calculer (q-qs)*dp
439  zqmqsdp(i, k) = (zq(i,k)-zqs(i,k))*delp(i, k)
440  END IF
441  END DO
442  END DO
443 
444  ! -----zgama is the moist convective lapse rate (-dT/dz).
445  ! -----zgamdz(*,k) est la difference minimale autorisee des temperatures
446  ! -----entre deux couches (k et k+1), c.a.d. si T(k+1)-T(k) est inferieur
447  ! -----a zgamdz(*,k), alors ces 2 couches sont instables conditionnellement
448 
449  DO k = 1, klev - 1
450  DO i = 1, klon
451  IF (afaire(i)) THEN
452  zgamdz(i, k) = -(pplay(i,k)-pplay(i,k+1))/paprs(i, k+1)/rcpd*(rd*(zt( &
453  i,k)*delp(i,k)+zt(i,k+1)*delp(i,k+1))/(delp(i,k)+delp(i, &
454  k+1))+rlvtt*(zqs(i,k)*delp(i,k)+zqs(i,k+1)*delp(i,k+1))/(delp(i, &
455  k)+delp(i,k+1)))/(1.0+(zdqs(i,k)*delp(i,k)+zdqs(i,k+1)*delp(i, &
456  k+1))/(delp(i,k)+delp(i,k+1)))
457  END IF
458  END DO
459  END DO
460 
461  ! On cherche la presence simultanee d'instabilite conditionnelle
462  ! et de sur-saturation. Sinon, pas la penne de se fatiguer:
463 
464  DO i = 1, klon
465  possible(i) = .false.
466  END DO
467  DO k = 2, klev
468  DO i = 1, klon
469  IF (afaire(i)) THEN
470  zflo(i) = zt(i, k-1) + zgamdz(i, k-1) - zt(i, k)
471  zsat(i) = zqmqsdp(i, k) + zqmqsdp(i, k-1)
472  IF (zflo(i)>0.0 .AND. zsat(i)>0.0) possible(i) = .true.
473  END IF
474  END DO
475  END DO
476 
477  DO i = 1, klon
478  IF (possible(i)) THEN
479  k1(i) = kbase
480  k2(i) = k1(i) + 1
481  END IF
482  END DO
483 
484 810 CONTINUE ! chercher le bas de la colonne a ajuster
485 
486  k2min = klev
487  DO i = 1, klon
488  todo(i) = .false.
489  aller(i) = .true.
490  IF (possible(i)) k2min = min(k2min, k2(i))
491  END DO
492  IF (k2min==klev) GO TO 860
493  DO k = k2min, klev - 1
494  DO i = 1, klon
495  IF (possible(i) .AND. k>=k2(i) .AND. aller(i)) THEN
496  zflo(i) = zt(i, k) + zgamdz(i, k) - zt(i, k+1)
497  zsat(i) = zqmqsdp(i, k) + zqmqsdp(i, k+1)
498  IF (zflo(i)>0.0 .AND. zsat(i)>0.0) THEN
499  k1(i) = k
500  k2(i) = k + 1
501  todo(i) = .true.
502  aller(i) = .false.
503  END IF
504  END IF
505  END DO
506  END DO
507  DO i = 1, klon
508  IF (possible(i) .AND. aller(i)) THEN
509  todo(i) = .false.
510  k1(i) = klev
511  k2(i) = klev
512  END IF
513  END DO
514 
515  ! CC DO i = 1, klon
516  ! CC IF (possible(i)) THEN
517  ! CC 811 k2(i) = k2(i) + 1
518  ! CC IF (k2(i) .GT. klev) THEN
519  ! CC todo(i) = .FALSE.
520  ! CC GOTO 812
521  ! CC ENDIF
522  ! CC k = k2(i)
523  ! CC zflo(i) = zt(i,k-1) + zgamdz(i,k-1) - zt(i,k)
524  ! CC zsat(i) = zqmqsdp(i,k) + zqmqsdp(i,k-1)
525  ! CC IF (zflo(i).LE.0.0 .OR. zsat(i).LE.0.0) GOTO 811
526  ! CC k1(i) = k2(i) - 1
527  ! CC todo(i) = .TRUE.
528  ! CC ENDIF
529  ! CC 812 CONTINUE
530  ! CC ENDDO
531 
532 820 CONTINUE ! chercher le haut de la colonne
533 
534  k2min = klev
535  DO i = 1, klon
536  aller(i) = .true.
537  IF (todo(i)) k2min = min(k2min, k2(i))
538  END DO
539  IF (k2min<klev) THEN
540  DO k = k2min, klev
541  DO i = 1, klon
542  IF (todo(i) .AND. k>k2(i) .AND. aller(i)) THEN
543  zsat(i) = zsat(i) + zqmqsdp(i, k)
544  zflo(i) = zt(i, k-1) + zgamdz(i, k-1) - zt(i, k)
545  IF (zflo(i)<=0.0 .OR. zsat(i)<=0.0) THEN
546  aller(i) = .false.
547  ELSE
548  k2(i) = k
549  END IF
550  END IF
551  END DO
552  END DO
553  ! error is = 0
554  ! error DO i = 1, klon
555  ! error IF(todo(i).AND.aller(i)) THEN
556  ! error is = is + 1
557  ! error todo(i) = .FALSE.
558  ! error k2(i) = klev
559  ! error ENDIF
560  ! error ENDDO
561  ! error IF (is.GT.0) THEN
562  ! error PRINT*, "Bizard. je pourrais continuer mais j arrete"
563  ! error CALL abort
564  ! error ENDIF
565  END IF
566 
567  ! CC DO i = 1, klon
568  ! CC IF (todo(i)) THEN
569  ! CC 821 CONTINUE
570  ! CC IF (k2(i) .EQ. klev) GOTO 822
571  ! CC k = k2(i) + 1
572  ! CC zsat(i) = zsat(i) + zqmqsdp(i,k)
573  ! CC zflo(i) = zt(i,k-1) + zgamdz(i,k-1) - zt(i,k)
574  ! CC IF (zflo(i).LE.0.0 .OR. zsat(i).LE.0.0) GOTO 822
575  ! CC k2(i) = k
576  ! CC GOTO 821
577  ! CC ENDIF
578  ! CC 822 CONTINUE
579  ! CC ENDDO
580 
581 830 CONTINUE ! faire l'ajustement en sachant k1 et k2
582 
583  is = 0
584  DO i = 1, klon
585  IF (todo(i)) THEN
586  IF (k2(i)<=k1(i)) is = is + 1
587  END IF
588  END DO
589  IF (is>0) THEN
590  print *, 'Impossible: k1 trop grand ou k2 trop petit'
591  print *, 'is=', is
592  CALL abort
593  END IF
594 
595  k1min = klev
596  k1max = 1
597  k2max = 1
598  DO i = 1, klon
599  IF (todo(i)) THEN
600  k1min = min(k1min, k1(i))
601  k1max = max(k1max, k1(i))
602  k2max = max(k2max, k2(i))
603  END IF
604  END DO
605 
606  DO i = 1, klon
607  IF (todo(i)) THEN
608  k = k1(i)
609  za(i) = 0.
610  zb(i) = (rcpd*(1.+zdqs(i,k))*(zt(i,k)-za(i))-rlvtt*(zqs(i,k)-zq(i, &
611  k)))*delp(i, k)
612  zc(i) = delp(i, k)*rcpd*(1.+zdqs(i,k))
613  END IF
614  END DO
615 
616  DO k = k1min, k2max
617  DO i = 1, klon
618  IF (todo(i) .AND. k>=(k1(i)+1) .AND. k<=k2(i)) THEN
619  za(i) = za(i) + zgamdz(i, k-1)
620  zb(i) = zb(i) + (rcpd*(1.+zdqs(i,k))*(zt(i,k)-za(i))-rlvtt*(zqs(i, &
621  k)-zq(i,k)))*delp(i, k)
622  zc(i) = zc(i) + delp(i, k)*rcpd*(1.+zdqs(i,k))
623  END IF
624  END DO
625  END DO
626 
627  DO i = 1, klon
628  IF (todo(i)) THEN
629  k = k1(i)
630  ztnew(i, k) = zb(i)/zc(i)
631  zqnew(i, k) = zqs(i, k) + (ztnew(i,k)-zt(i,k))*rcpd/rlvtt*zdqs(i, k)
632  END IF
633  END DO
634 
635  DO k = k1min, k2max
636  DO i = 1, klon
637  IF (todo(i) .AND. k>=(k1(i)+1) .AND. k<=k2(i)) THEN
638  ztnew(i, k) = ztnew(i, k-1) + zgamdz(i, k-1)
639  zqnew(i, k) = zqs(i, k) + (ztnew(i,k)-zt(i,k))*rcpd/rlvtt*zdqs(i, k)
640  END IF
641  END DO
642  END DO
643 
644  ! Quantite de condensation produite pendant l'ajustement:
645 
646  DO i = 1, klon
647  zcond(i) = 0.0
648  END DO
649  DO k = k1min, k2max
650  DO i = 1, klon
651  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i)) THEN
652  rneb(i, k) = 1.0
653  zcond(i) = zcond(i) + (zq(i,k)-zqnew(i,k))*delp(i, k)/rg
654  END IF
655  END DO
656  END DO
657 
658  ! Si condensation negative, effort completement perdu:
659 
660  DO i = 1, klon
661  IF (todo(i) .AND. zcond(i)<=0.) todo(i) = .false.
662  END DO
663 
664  ! L'ajustement a ete accompli, meme les calculs accessoires
665  ! ne sont pas encore faits:
666 
667  DO i = 1, klon
668  IF (todo(i)) accompli(i) = .true.
669  END DO
670 
671  ! =====
672  ! Une fois que la condensation a lieu, on doit construire un
673  ! "modele nuageux" pour partager la condensation entre l'eau
674  ! liquide nuageuse et la precipitation (leur rapport toliq
675  ! est calcule selon l'epaisseur nuageuse). Je suppose que
676  ! toliq=tomax quand l'epaisseur nuageuse est inferieure a dpmin,
677  ! et que toliq=tomin quand l'epaisseur depasse dpmax (interpolation
678  ! lineaire entre dpmin et dpmax).
679  ! =====
680  DO i = 1, klon
681  IF (todo(i)) THEN
682  toliq(i) = tomax - ((paprs(i,k1(i))-paprs(i,k2(i)+1))/paprs(i,1)-dpmin) &
683  *(tomax-tomin)/(dpmax-dpmin)
684  toliq(i) = max(tomin, min(tomax,toliq(i)))
685  IF (pplay(i,k2(i))/paprs(i,1)<=deep_sig) toliq(i) = deep_to
686  IF (old_tau) toliq(i) = 1.0
687  END IF
688  END DO
689  ! =====
690  ! On doit aussi determiner la distribution verticale de
691  ! l'eau nuageuse. Plusieurs options sont proposees:
692 
693  ! (0) La condensation precipite integralement (toliq ne sera
694  ! pas utilise).
695  ! (1) L'eau liquide est distribuee entre k1 et k2 et proportionnelle
696  ! a la vapeur d'eau locale.
697  ! (2) Elle est distribuee entre k1 et k2 avec une valeur constante.
698  ! (3) Elle est seulement distribuee aux couches ou la vapeur d'eau
699  ! est effectivement diminuee pendant le processus d'ajustement.
700  ! (4) Elle est en fonction (lineaire ou exponentielle) de la
701  ! distance (epaisseur en pression) avec le niveau k1 (la couche
702  ! k1 n'aura donc pas d'eau liquide).
703  ! =====
704 
705  IF (opt_cld==0) THEN
706 
707  DO i = 1, klon
708  IF (todo(i)) zrfl(i) = zcond(i)/dtime
709  END DO
710 
711  ELSE IF (opt_cld==1) THEN
712 
713  DO i = 1, klon
714  IF (todo(i)) zvapo(i) = 0.0 ! quantite integrale de vapeur d'eau
715  END DO
716  DO k = k1min, k2max
717  DO i = 1, klon
718  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i)) zvapo(i) = zvapo(i) + &
719  zqnew(i, k)*delp(i, k)/rg
720  END DO
721  END DO
722  DO i = 1, klon
723  IF (todo(i)) THEN
724  zrapp(i) = toliq(i)*zcond(i)/zvapo(i)
725  zrapp(i) = max(0., min(1.,zrapp(i)))
726  zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
727  END IF
728  END DO
729  DO k = k1min, k2max
730  DO i = 1, klon
731  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i)) THEN
732  d_ql(i, k) = d_ql(i, k) + zrapp(i)*zqnew(i, k)
733  END IF
734  END DO
735  END DO
736 
737  ELSE IF (opt_cld==2) THEN
738 
739  DO i = 1, klon
740  IF (todo(i)) zvapo(i) = 0.0 ! quantite integrale de masse
741  END DO
742  DO k = k1min, k2max
743  DO i = 1, klon
744  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i)) zvapo(i) = zvapo(i) + &
745  delp(i, k)/rg
746  END DO
747  END DO
748  DO k = k1min, k2max
749  DO i = 1, klon
750  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i)) THEN
751  d_ql(i, k) = d_ql(i, k) + toliq(i)*zcond(i)/zvapo(i)
752  END IF
753  END DO
754  END DO
755  DO i = 1, klon
756  IF (todo(i)) zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
757  END DO
758 
759  ELSE IF (opt_cld==3) THEN
760 
761  DO i = 1, klon
762  IF (todo(i)) zvapo(i) = 0.0 ! quantite de l'eau strictement condensee
763  END DO
764  DO k = k1min, k2max
765  DO i = 1, klon
766  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i)) zvapo(i) = zvapo(i) + &
767  max(0.0, zq(i,k)-zqnew(i,k))*delp(i, k)/rg
768  END DO
769  END DO
770  DO k = k1min, k2max
771  DO i = 1, klon
772  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i) .AND. zvapo(i)>0.0) d_ql(i, &
773  k) = d_ql(i, k) + toliq(i)*zcond(i)/zvapo(i)*max(0.0, zq(i,k)-zqnew &
774  (i,k))
775  END DO
776  END DO
777  DO i = 1, klon
778  IF (todo(i)) zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
779  END DO
780 
781  ELSE IF (opt_cld==4) THEN
782 
783  nexpo = 3
784  ! cc nexpo = 1 ! distribution lineaire
785 
786  DO i = 1, klon
787  IF (todo(i)) zvapo(i) = 0.0 ! quantite integrale de masse
788  END DO ! (avec ponderation)
789  DO k = k1min, k2max
790  DO i = 1, klon
791  IF (todo(i) .AND. k>=(k1(i)+1) .AND. k<=k2(i)) zvapo(i) = zvapo(i) + &
792  delp(i, k)/rg*(pplay(i,k1(i))-pplay(i,k))**nexpo
793  END DO
794  END DO
795  DO k = k1min, k2max
796  DO i = 1, klon
797  IF (todo(i) .AND. k>=(k1(i)+1) .AND. k<=k2(i)) d_ql(i, k) = d_ql(i, &
798  k) + toliq(i)*zcond(i)/zvapo(i)*(pplay(i,k1(i))-pplay(i,k))**nexpo
799  END DO
800  END DO
801  DO i = 1, klon
802  IF (todo(i)) zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
803  END DO
804 
805  ELSE ! valeur non-prevue pour opt_cld
806 
807  print *, 'opt_cld est faux:', opt_cld
808  CALL abort
809 
810  END IF ! fin de opt_cld
811 
812  ! L'eau precipitante peut etre evaporee:
813 
814  zalfa = 0.05
815  IF (evap_prec .AND. (k1max>=2)) THEN
816  DO k = k1max - 1, 1, -1
817  DO i = 1, klon
818  IF (todo(i) .AND. k<k1(i) .AND. zrfl(i)>0.0) THEN
819  zqev = max(0.0, (zqs(i,k)-zq(i,k))*zalfa)
820  zqevt = coef_eva*(1.0-zq(i,k)/zqs(i,k))*sqrt(zrfl(i))*delp(i, k)/ &
821  pplay(i, k)*zt(i, k)*rd/rg
822  zqevt = max(0.0, min(zqevt,zrfl(i)))*rg*dtime/delp(i, k)
823  zqev = min(zqev, zqevt)
824  zrfln = zrfl(i) - zqev*(delp(i,k))/rg/dtime
825  zq(i, k) = zq(i, k) - (zrfln-zrfl(i))*(rg/(delp(i,k)))*dtime
826  zt(i, k) = zt(i, k) + (zrfln-zrfl(i))*(rg/(delp(i, &
827  k)))*dtime*rlvtt/rcpd/(1.0+rvtmp2*zq(i,k))
828  zrfl(i) = zrfln
829  END IF
830  END DO
831  END DO
832  END IF
833 
834  ! La temperature de la premiere couche determine la pluie ou la neige:
835 
836  DO i = 1, klon
837  IF (todo(i)) THEN
838  IF (zt(i,1)>rtt) THEN
839  rain(i) = rain(i) + zrfl(i)
840  ELSE
841  snow(i) = snow(i) + zrfl(i)
842  END IF
843  END IF
844  END DO
845 
846  ! Mise a jour de la temperature et de l'humidite
847 
848  DO k = k1min, k2max
849  DO i = 1, klon
850  IF (todo(i) .AND. k>=k1(i) .AND. k<=k2(i)) THEN
851  zt(i, k) = ztnew(i, k)
852  zq(i, k) = zqnew(i, k)
853  END IF
854  END DO
855  END DO
856 
857  ! Re-calculer certaines variables pour etendre et re-ajuster la colonne
858 
859  IF (exigent) THEN
860  DO k = 1, klev
861  DO i = 1, klon
862  IF (todo(i)) THEN
863  IF (thermcep) THEN
864  zdelta = max(0., sign(1.,rtt-zt(i,k)))
865  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
866  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*zq(i,k))
867  zqs(i, k) = r2es*foeew(zt(i,k), zdelta)/pplay(i, k)
868  zqs(i, k) = min(0.5, zqs(i,k))
869  zcor = 1./(1.-retv*zqs(i,k))
870  zqs(i, k) = zqs(i, k)*zcor
871  zdqs(i, k) = foede(zt(i,k), zdelta, zcvm5, zqs(i,k), zcor)
872  ELSE
873  IF (zt(i,k)<t_coup) THEN
874  zqs(i, k) = qsats(zt(i,k))/pplay(i, k)
875  zdqs(i, k) = dqsats(zt(i,k), zqs(i,k))
876  ELSE
877  zqs(i, k) = qsatl(zt(i,k))/pplay(i, k)
878  zdqs(i, k) = dqsatl(zt(i,k), zqs(i,k))
879  END IF
880  END IF
881  END IF
882  END DO
883  END DO
884  END IF
885 
886  IF (exigent) THEN
887  DO k = 1, klev - 1
888  DO i = 1, klon
889  IF (todo(i)) THEN
890  zgamdz(i, k) = -(pplay(i,k)-pplay(i,k+1))/paprs(i, k+1)/rcpd*(rd*( &
891  zt(i,k)*delp(i,k)+zt(i,k+1)*delp(i,k+1))/(delp(i,k)+delp(i, &
892  k+1))+rlvtt*(zqs(i,k)*delp(i,k)+zqs(i,k+1)*delp(i,k+1))/(delp(i, &
893  k)+delp(i,k+1)))/(1.0+(zdqs(i,k)*delp(i,k)+zdqs(i,k+1)*delp(i, &
894  k+1))/(delp(i,k)+delp(i,k+1)))
895  END IF
896  END DO
897  END DO
898  END IF
899 
900  ! Puisque l'humidite a ete modifiee, on re-fait (q-qs)*dp
901 
902  DO k = 1, klev
903  DO i = 1, klon
904  IF (todo(i)) THEN
905  zqmqsdp(i, k) = (zq(i,k)-zqs(i,k))*delp(i, k)
906  END IF
907  END DO
908  END DO
909 
910  ! Verifier si l'on peut etendre le bas de la colonne
911 
912  DO i = 1, klon
913  etendre(i) = .false.
914  END DO
915 
916  k1max = 1
917  DO i = 1, klon
918  IF (todo(i) .AND. k1(i)>(kbase+1)) THEN
919  k = k1(i)
920  zflo(i) = zt(i, k-1) + zgamdz(i, k-1) - zt(i, k)
921  zsat(i) = zqmqsdp(i, k) + zqmqsdp(i, k-1)
922  ! sc voici l'ancienne ligne:
923  ! sc IF (zflo(i).LE.0.0 .OR. zsat(i).LE.0.0) THEN
924  ! sc sylvain: il faut RESPECTER les 2 criteres:
925  IF (zflo(i)>0.0 .AND. zsat(i)>0.0) THEN
926  etendre(i) = .true.
927  k1(i) = k1(i) - 1
928  k1max = max(k1max, k1(i))
929  aller(i) = .true.
930  END IF
931  END IF
932  END DO
933 
934  IF (k1max>(kbase+1)) THEN
935  DO k = k1max, kbase + 1, -1
936  DO i = 1, klon
937  IF (etendre(i) .AND. k<k1(i) .AND. aller(i)) THEN
938  zsat(i) = zsat(i) + zqmqsdp(i, k)
939  zflo(i) = zt(i, k) + zgamdz(i, k) - zt(i, k+1)
940  IF (zsat(i)<=0.0 .OR. zflo(i)<=0.0) THEN
941  aller(i) = .false.
942  ELSE
943  k1(i) = k
944  END IF
945  END IF
946  END DO
947  END DO
948  DO i = 1, klon
949  IF (etendre(i) .AND. aller(i)) THEN
950  k1(i) = 1
951  END IF
952  END DO
953  END IF
954 
955  ! CC DO i = 1, klon
956  ! CC IF (etendre(i)) THEN
957  ! CC 840 k = k1(i)
958  ! CC IF (k.GT.1) THEN
959  ! CC zsat(i) = zsat(i) + zqmqsdp(i,k-1)
960  ! CC zflo(i) = zt(i,k-1) + zgamdz(i,k-1) - zt(i,k)
961  ! CC IF (zflo(i).GT.0.0 .AND. zsat(i).GT.0.0) THEN
962  ! CC k1(i) = k - 1
963  ! CC GOTO 840
964  ! CC ENDIF
965  ! CC ENDIF
966  ! CC ENDIF
967  ! CC ENDDO
968 
969  DO i = 1, klon
970  todobis(i) = todo(i)
971  todo(i) = .false.
972  END DO
973  is = 0
974  DO i = 1, klon
975  IF (etendre(i)) THEN
976  todo(i) = .true.
977  is = is + 1
978  END IF
979  END DO
980  IF (is>0) THEN
981  IF (new_top) THEN
982  GO TO 820 ! chercher de nouveau le sommet k2
983  ELSE
984  GO TO 830 ! supposer que le sommet est celui deja trouve
985  END IF
986  END IF
987 
988  DO i = 1, klon
989  possible(i) = .false.
990  END DO
991  is = 0
992  DO i = 1, klon
993  IF (todobis(i) .AND. k2(i)<klev) THEN
994  is = is + 1
995  possible(i) = .true.
996  END IF
997  END DO
998  IF (is>0) GO TO 810 !on cherche en haut d'autres blocks
999  ! a ajuster a partir du sommet de la colonne precedente
1000 
1001 860 CONTINUE ! Calculer les tendances et diagnostiques
1002  ! cc print*, "Apres 860"
1003 
1004  DO k = 1, klev
1005  DO i = 1, klon
1006  IF (accompli(i)) THEN
1007  d_t(i, k) = zt(i, k) - t(i, k)
1008  zq(i, k) = max(zq(i,k), seuil_vap)
1009  d_q(i, k) = zq(i, k) - q(i, k)
1010  END IF
1011  END DO
1012  END DO
1013 
1014  DO i = 1, klon
1015  IF (accompli(i)) THEN
1016  DO k = 1, klev
1017  IF (rneb(i,k)>0.0) THEN
1018  ibas(i) = k
1019  GO TO 807
1020  END IF
1021  END DO
1022 807 CONTINUE
1023  DO k = klev, 1, -1
1024  IF (rneb(i,k)>0.0) THEN
1025  itop(i) = k
1026  GO TO 808
1027  END IF
1028  END DO
1029 808 CONTINUE
1030  END IF
1031  END DO
1032 
1033  IF (imprim) THEN
1034  nbtodo = 0
1035  nbdone = 0
1036  DO i = 1, klon
1037  IF (afaire(i)) nbtodo = nbtodo + 1
1038  IF (accompli(i)) nbdone = nbdone + 1
1039  END DO
1040  print *, 'nbTodo, nbDone=', nbtodo, nbdone
1041  END IF
1042 
1043  RETURN
1044 END SUBROUTINE conmanv
1045 SUBROUTINE conkuo(dtime, paprs, pplay, t, q, conv_q, d_t, d_q, d_ql, rneb, &
1046  rain, snow, ibas, itop)
1048  IMPLICIT NONE
1049  ! ======================================================================
1050  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818
1051  ! Objet: Schema de convection de type Kuo (1965).
1052  ! Cette version du code peut calculer le niveau de depart
1053  ! N.B. version vectorielle (le 6 oct. 1997)
1054  ! ======================================================================
1055  include "YOMCST.h"
1056 
1057  ! Arguments:
1058 
1059  REAL dtime ! intervalle du temps (s)
1060  REAL paprs(klon, klev+1) ! pression a inter-couche (Pa)
1061  REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
1062  REAL t(klon, klev) ! temperature (K)
1063  REAL q(klon, klev) ! humidite specifique
1064  REAL conv_q(klon, klev) ! taux de convergence humidite (g/g/s)
1065 
1066  REAL d_t(klon, klev) ! incrementation temperature
1067  REAL d_q(klon, klev) ! incrementation humidite
1068  REAL d_ql(klon, klev) ! incrementation eau liquide
1069  REAL rneb(klon, klev) ! nebulosite
1070  REAL rain(klon) ! pluies (mm/s)
1071  REAL snow(klon) ! neige (mm/s)
1072  INTEGER itop(klon) ! niveau du sommet
1073  INTEGER ibas(klon) ! niveau du bas
1074 
1075  LOGICAL ldcum(klon) ! convection existe
1076  LOGICAL todo(klon)
1077 
1078  ! Quelsques options:
1079 
1080  LOGICAL calcfcl ! calculer le niveau de convection libre
1081  parameter(calcfcl=.true.)
1082  INTEGER ldepar ! niveau fixe de convection libre
1083  parameter(ldepar=4)
1084  INTEGER opt_cld ! comment traiter l'eau liquide
1085  parameter(opt_cld=4) ! valeur possible: 0, 1, 2, 3 ou 4
1086  LOGICAL evap_prec ! evaporation de pluie au-dessous de convection
1087  parameter(evap_prec=.true.)
1088  REAL coef_eva
1089  parameter(coef_eva=1.0e-05)
1090  LOGICAL new_deh ! nouvelle facon de calculer dH
1091  parameter(new_deh=.false.)
1092  REAL t_coup
1093  parameter(t_coup=234.0)
1094  LOGICAL old_tau ! implique precipitation nulle
1095  parameter(old_tau=.false.)
1096  REAL toliq(klon) ! rapport entre l'eau nuageuse et l'eau precipitante
1097  REAL dpmin, tomax !Epaisseur faible, rapport eau liquide plus grande
1098  parameter(dpmin=0.15, tomax=0.97)
1099  REAL dpmax, tomin !Epaisseur grande, rapport eau liquide plus faible
1100  parameter(dpmax=0.30, tomin=0.05)
1101  REAL deep_sig, deep_to ! au dela de deep_sig, utiliser deep_to
1102  parameter(deep_sig=0.50, deep_to=0.05)
1103 
1104  ! Variables locales:
1105 
1106  INTEGER nexpo
1107  LOGICAL nuage(klon)
1108  INTEGER i, k, kbmin, kbmax, khmax
1109  REAL ztotal(klon, klev), zdeh(klon, klev)
1110  REAL zgz(klon, klev)
1111  REAL zqs(klon, klev)
1112  REAL zdqs(klon, klev)
1113  REAL ztemp(klon, klev)
1114  REAL zpres(klon, klev)
1115  REAL zconv(klon) ! convergence d'humidite
1116  REAL zvirt(klon) ! convergence virtuelle d'humidite
1117  REAL zfrac(klon) ! fraction convective
1118  INTEGER kb(klon), kh(klon)
1119 
1120  REAL zcond(klon), zvapo(klon), zrapp(klon)
1121  REAL zrfl(klon), zrfln, zqev, zqevt
1122  REAL zdelta, zcvm5, zcor
1123  REAL zvar
1124 
1125  LOGICAL appel1er
1126  SAVE appel1er
1127  !$OMP THREADPRIVATE(appel1er)
1128 
1129  ! Fonctions thermodynamiques
1130 
1131  include "YOETHF.h"
1132  include "FCTTRE.h"
1133 
1134  DATA appel1er/.true./
1135 
1136  IF (appel1er) THEN
1137  print *, 'conkuo, calcfcl:', calcfcl
1138  IF (.NOT. calcfcl) print *, 'conkuo, ldepar:', ldepar
1139  print *, 'conkuo, opt_cld:', opt_cld
1140  print *, 'conkuo, evap_prec:', evap_prec
1141  print *, 'conkuo, new_deh:', new_deh
1142  appel1er = .false.
1143  END IF
1144 
1145  ! Initialiser les sorties a zero
1146 
1147  DO k = 1, klev
1148  DO i = 1, klon
1149  d_q(i, k) = 0.0
1150  d_t(i, k) = 0.0
1151  d_ql(i, k) = 0.0
1152  rneb(i, k) = 0.0
1153  END DO
1154  END DO
1155  DO i = 1, klon
1156  rain(i) = 0.0
1157  snow(i) = 0.0
1158  ibas(i) = 0
1159  itop(i) = 0
1160  END DO
1161 
1162  ! Calculer la vapeur d'eau saturante Qs et sa derive L/Cp * dQs/dT
1163 
1164  DO k = 1, klev
1165  DO i = 1, klon
1166  IF (thermcep) THEN
1167  zdelta = max(0., sign(1.,rtt-t(i,k)))
1168  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
1169  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*q(i,k))
1170  zqs(i, k) = r2es*foeew(t(i,k), zdelta)/pplay(i, k)
1171  zqs(i, k) = min(0.5, zqs(i,k))
1172  zcor = 1./(1.-retv*zqs(i,k))
1173  zqs(i, k) = zqs(i, k)*zcor
1174  zdqs(i, k) = foede(t(i,k), zdelta, zcvm5, zqs(i,k), zcor)
1175  ELSE
1176  IF (t(i,k)<t_coup) THEN
1177  zqs(i, k) = qsats(t(i,k))/pplay(i, k)
1178  zdqs(i, k) = dqsats(t(i,k), zqs(i,k))
1179  ELSE
1180  zqs(i, k) = qsatl(t(i,k))/pplay(i, k)
1181  zdqs(i, k) = dqsatl(t(i,k), zqs(i,k))
1182  END IF
1183  END IF
1184  END DO
1185  END DO
1186 
1187  ! Calculer gz (energie potentielle)
1188 
1189  DO i = 1, klon
1190  zgz(i, 1) = rd*t(i, 1)/(0.5*(paprs(i,1)+pplay(i, &
1191  1)))*(paprs(i,1)-pplay(i,1))
1192  END DO
1193  DO k = 2, klev
1194  DO i = 1, klon
1195  zgz(i, k) = zgz(i, k-1) + rd*0.5*(t(i,k-1)+t(i,k))/paprs(i, k)*(pplay(i &
1196  ,k-1)-pplay(i,k))
1197  END DO
1198  END DO
1199 
1200  ! Calculer l'energie statique humide saturee (Cp*T + gz + L*Qs)
1201 
1202  DO k = 1, klev
1203  DO i = 1, klon
1204  ztotal(i, k) = rcpd*t(i, k) + rlvtt*zqs(i, k) + zgz(i, k)
1205  END DO
1206  END DO
1207 
1208  ! Determiner le niveau de depart et calculer la difference de
1209  ! l'energie statique humide saturee (ztotal) entre la couche
1210  ! de depart et chaque couche au-dessus.
1211 
1212  IF (calcfcl) THEN
1213  DO k = 1, klev
1214  DO i = 1, klon
1215  zpres(i, k) = pplay(i, k)
1216  ztemp(i, k) = t(i, k)
1217  END DO
1218  END DO
1219  CALL kuofcl(ztemp, q, zgz, zpres, ldcum, kb)
1220  DO i = 1, klon
1221  IF (ldcum(i)) THEN
1222  k = kb(i)
1223  IF (new_deh) THEN
1224  zdeh(i, k) = ztotal(i, k-1) - ztotal(i, k)
1225  ELSE
1226  zdeh(i, k) = rcpd*(t(i,k-1)-t(i,k)) - rd*0.5*(t(i,k-1)+t(i,k))/ &
1227  paprs(i, k)*(pplay(i,k-1)-pplay(i,k)) + &
1228  rlvtt*(zqs(i,k-1)-zqs(i,k))
1229  END IF
1230  zdeh(i, k) = zdeh(i, k)*0.5
1231  END IF
1232  END DO
1233  DO k = 1, klev
1234  DO i = 1, klon
1235  IF (ldcum(i) .AND. k>=(kb(i)+1)) THEN
1236  IF (new_deh) THEN
1237  zdeh(i, k) = zdeh(i, k-1) + (ztotal(i,k-1)-ztotal(i,k))
1238  ELSE
1239  zdeh(i, k) = zdeh(i, k-1) + rcpd*(t(i,k-1)-t(i,k)) - &
1240  rd*0.5*(t(i,k-1)+t(i,k))/paprs(i, k)* &
1241  (pplay(i,k-1)-pplay(i,k)) + rlvtt*(zqs(i,k-1)-zqs(i,k))
1242  END IF
1243  END IF
1244  END DO
1245  END DO
1246  ELSE
1247  DO i = 1, klon
1248  k = ldepar
1249  kb(i) = ldepar
1250  ldcum(i) = .true.
1251  IF (new_deh) THEN
1252  zdeh(i, k) = ztotal(i, k-1) - ztotal(i, k)
1253  ELSE
1254  zdeh(i, k) = rcpd*(t(i,k-1)-t(i,k)) - rd*0.5*(t(i,k-1)+t(i,k))/paprs( &
1255  i, k)*(pplay(i,k-1)-pplay(i,k)) + rlvtt*(zqs(i,k-1)-zqs(i,k))
1256  END IF
1257  zdeh(i, k) = zdeh(i, k)*0.5
1258  END DO
1259  DO k = ldepar + 1, klev
1260  DO i = 1, klon
1261  IF (new_deh) THEN
1262  zdeh(i, k) = zdeh(i, k-1) + (ztotal(i,k-1)-ztotal(i,k))
1263  ELSE
1264  zdeh(i, k) = zdeh(i, k-1) + rcpd*(t(i,k-1)-t(i,k)) - &
1265  rd*0.5*(t(i,k-1)+t(i,k))/paprs(i, k)*(pplay(i,k-1)-pplay(i,k)) + &
1266  rlvtt*(zqs(i,k-1)-zqs(i,k))
1267  END IF
1268  END DO
1269  END DO
1270  END IF
1271 
1272  ! -----Chercher le sommet du nuage
1273  ! -----Calculer la convergence de l'humidite (en kg/m**2 a un facteur
1274  ! -----psolpa/RG pres) du bas jusqu'au sommet du nuage.
1275  ! -----Calculer la convergence virtuelle pour que toute la maille
1276  ! -----deviennt nuageuse (du bas jusqu'au sommet du nuage)
1277 
1278  DO i = 1, klon
1279  nuage(i) = .true.
1280  zconv(i) = 0.0
1281  zvirt(i) = 0.0
1282  kh(i) = -999
1283  END DO
1284  DO k = 1, klev
1285  DO i = 1, klon
1286  IF (k>=kb(i) .AND. ldcum(i)) THEN
1287  nuage(i) = nuage(i) .AND. zdeh(i, k) > 0.0
1288  IF (nuage(i)) THEN
1289  kh(i) = k
1290  zconv(i) = zconv(i) + conv_q(i, k)*dtime*(paprs(i,k)-paprs(i,k+1))
1291  zvirt(i) = zvirt(i) + (zdeh(i,k)/rlvtt+zqs(i,k)-q(i,k))*(paprs(i,k) &
1292  -paprs(i,k+1))
1293  END IF
1294  END IF
1295  END DO
1296  END DO
1297 
1298  DO i = 1, klon
1299  todo(i) = ldcum(i) .AND. kh(i) > kb(i) .AND. zconv(i) > 0.0
1300  END DO
1301 
1302  kbmin = klev
1303  kbmax = 0
1304  khmax = 0
1305  DO i = 1, klon
1306  IF (todo(i)) THEN
1307  kbmin = min(kbmin, kb(i))
1308  kbmax = max(kbmax, kb(i))
1309  khmax = max(khmax, kh(i))
1310  END IF
1311  END DO
1312 
1313  ! -----Calculer la surface couverte par le nuage
1314 
1315  DO i = 1, klon
1316  IF (todo(i)) THEN
1317  zfrac(i) = max(0.0, min(zconv(i)/zvirt(i),1.0))
1318  END IF
1319  END DO
1320 
1321  ! -----Calculs essentiels:
1322 
1323  DO i = 1, klon
1324  IF (todo(i)) THEN
1325  zcond(i) = 0.0
1326  END IF
1327  END DO
1328  DO k = kbmin, khmax
1329  DO i = 1, klon
1330  IF (todo(i) .AND. k>=kb(i) .AND. k<=kh(i)) THEN
1331  zvar = zdeh(i, k)/(1.+zdqs(i,k))
1332  d_t(i, k) = zvar*zfrac(i)/rcpd
1333  d_q(i, k) = (zvar*zdqs(i,k)/rlvtt+zqs(i,k)-q(i,k))*zfrac(i) - &
1334  conv_q(i, k)*dtime
1335  zcond(i) = zcond(i) - d_q(i, k)*(paprs(i,k)-paprs(i,k+1))/rg
1336  rneb(i, k) = zfrac(i)
1337  END IF
1338  END DO
1339  END DO
1340 
1341  DO i = 1, klon
1342  IF (todo(i) .AND. zcond(i)<0.0) THEN
1343  print *, 'WARNING: cond. negative (Kuo) ', i, kb(i), kh(i), zcond(i)
1344  zcond(i) = 0.0
1345  DO k = kb(i), kh(i)
1346  d_t(i, k) = 0.0
1347  d_q(i, k) = 0.0
1348  END DO
1349  todo(i) = .false. ! effort totalement perdu
1350  END IF
1351  END DO
1352 
1353  ! =====
1354  ! Une fois que la condensation a lieu, on doit construire un
1355  ! "modele nuageux" pour partager la condensation entre l'eau
1356  ! liquide nuageuse et la precipitation (leur rapport toliq
1357  ! est calcule selon l'epaisseur nuageuse). Je suppose que
1358  ! toliq=tomax quand l'epaisseur nuageuse est inferieure a dpmin,
1359  ! et que toliq=tomin quand l'epaisseur depasse dpmax (interpolation
1360  ! lineaire entre dpmin et dpmax).
1361  ! =====
1362  DO i = 1, klon
1363  IF (todo(i)) THEN
1364  toliq(i) = tomax - ((paprs(i,kb(i))-paprs(i,kh(i)+1))/paprs(i,1)-dpmin) &
1365  *(tomax-tomin)/(dpmax-dpmin)
1366  toliq(i) = max(tomin, min(tomax,toliq(i)))
1367  IF (pplay(i,kh(i))/paprs(i,1)<=deep_sig) toliq(i) = deep_to
1368  IF (old_tau) toliq(i) = 1.0
1369  END IF
1370  END DO
1371  ! =====
1372  ! On doit aussi determiner la distribution verticale de
1373  ! l'eau nuageuse. Plusieurs options sont proposees:
1374 
1375  ! (0) La condensation precipite integralement (toliq ne sera
1376  ! pas utilise).
1377  ! (1) L'eau liquide est distribuee entre k1 et k2 et proportionnelle
1378  ! a la vapeur d'eau locale.
1379  ! (2) Elle est distribuee entre k1 et k2 avec une valeur constante.
1380  ! (3) Elle est seulement distribuee aux couches ou la vapeur d'eau
1381  ! est effectivement diminuee pendant le processus d'ajustement.
1382  ! (4) Elle est en fonction (lineaire ou exponentielle) de la
1383  ! distance (epaisseur en pression) avec le niveau k1 (la couche
1384  ! k1 n'aura donc pas d'eau liquide).
1385  ! =====
1386 
1387  IF (opt_cld==0) THEN
1388 
1389  DO i = 1, klon
1390  IF (todo(i)) zrfl(i) = zcond(i)/dtime
1391  END DO
1392 
1393  ELSE IF (opt_cld==1) THEN
1394 
1395  DO i = 1, klon
1396  IF (todo(i)) zvapo(i) = 0.0 ! quantite integrale de vapeur d'eau
1397  END DO
1398  DO k = kbmin, khmax
1399  DO i = 1, klon
1400  IF (todo(i) .AND. k>=kb(i) .AND. k<=kh(i)) THEN
1401  zvapo(i) = zvapo(i) + (q(i,k)+d_q(i,k))*(paprs(i,k)-paprs(i,k+1))/ &
1402  rg
1403  END IF
1404  END DO
1405  END DO
1406  DO i = 1, klon
1407  IF (todo(i)) THEN
1408  zrapp(i) = toliq(i)*zcond(i)/zvapo(i)
1409  zrapp(i) = max(0., min(1.,zrapp(i)))
1410  END IF
1411  END DO
1412  DO k = kbmin, khmax
1413  DO i = 1, klon
1414  IF (todo(i) .AND. k>=kb(i) .AND. k<=kh(i)) THEN
1415  d_ql(i, k) = zrapp(i)*(q(i,k)+d_q(i,k))
1416  END IF
1417  END DO
1418  END DO
1419  DO i = 1, klon
1420  IF (todo(i)) THEN
1421  zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
1422  END IF
1423  END DO
1424 
1425  ELSE IF (opt_cld==2) THEN
1426 
1427  DO i = 1, klon
1428  IF (todo(i)) zvapo(i) = 0.0 ! quantite integrale de masse
1429  END DO
1430  DO k = kbmin, khmax
1431  DO i = 1, klon
1432  IF (todo(i) .AND. k>=kb(i) .AND. k<=kh(i)) THEN
1433  zvapo(i) = zvapo(i) + (paprs(i,k)-paprs(i,k+1))/rg
1434  END IF
1435  END DO
1436  END DO
1437  DO k = kbmin, khmax
1438  DO i = 1, klon
1439  IF (todo(i) .AND. k>=kb(i) .AND. k<=kh(i)) THEN
1440  d_ql(i, k) = toliq(i)*zcond(i)/zvapo(i)
1441  END IF
1442  END DO
1443  END DO
1444  DO i = 1, klon
1445  IF (todo(i)) THEN
1446  zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
1447  END IF
1448  END DO
1449 
1450  ELSE IF (opt_cld==3) THEN
1451 
1452  DO i = 1, klon
1453  IF (todo(i)) THEN
1454  zvapo(i) = 0.0 ! quantite de l'eau strictement condensee
1455  END IF
1456  END DO
1457  DO k = kbmin, khmax
1458  DO i = 1, klon
1459  IF (todo(i) .AND. k>=kb(i) .AND. k<=kh(i)) THEN
1460  zvapo(i) = zvapo(i) + max(0.0, -d_q(i,k))*(paprs(i,k)-paprs(i,k+1)) &
1461  /rg
1462  END IF
1463  END DO
1464  END DO
1465  DO k = kbmin, khmax
1466  DO i = 1, klon
1467  IF (todo(i) .AND. k>=kb(i) .AND. k<=kh(i) .AND. zvapo(i)>0.0) THEN
1468  d_ql(i, k) = d_ql(i, k) + toliq(i)*zcond(i)/zvapo(i)*max(0.0, -d_q( &
1469  i,k))
1470  END IF
1471  END DO
1472  END DO
1473  DO i = 1, klon
1474  IF (todo(i)) THEN
1475  zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
1476  END IF
1477  END DO
1478 
1479  ELSE IF (opt_cld==4) THEN
1480 
1481  nexpo = 3
1482  ! cc nexpo = 1 ! distribution lineaire
1483 
1484  DO i = 1, klon
1485  IF (todo(i)) THEN
1486  zvapo(i) = 0.0 ! quantite integrale de masse (avec ponderation)
1487  END IF
1488  END DO
1489  DO k = kbmin, khmax
1490  DO i = 1, klon
1491  IF (todo(i) .AND. k>=(kb(i)+1) .AND. k<=kh(i)) THEN
1492  zvapo(i) = zvapo(i) + (paprs(i,k)-paprs(i,k+1))/rg*(pplay(i,kb(i))- &
1493  pplay(i,k))**nexpo
1494  END IF
1495  END DO
1496  END DO
1497  DO k = kbmin, khmax
1498  DO i = 1, klon
1499  IF (todo(i) .AND. k>=(kb(i)+1) .AND. k<=kh(i)) THEN
1500  d_ql(i, k) = d_ql(i, k) + toliq(i)*zcond(i)/zvapo(i)*(pplay(i,kb(i) &
1501  )-pplay(i,k))**nexpo
1502  END IF
1503  END DO
1504  END DO
1505  DO i = 1, klon
1506  IF (todo(i)) THEN
1507  zrfl(i) = (1.0-toliq(i))*zcond(i)/dtime
1508  END IF
1509  END DO
1510 
1511  ELSE ! valeur non-prevue pour opt_cld
1512 
1513  print *, 'opt_cld est faux:', opt_cld
1514  CALL abort
1515 
1516  END IF ! fin de opt_cld
1517 
1518  ! L'eau precipitante peut etre re-evaporee:
1519 
1520  IF (evap_prec .AND. kbmax>=2) THEN
1521  DO k = kbmax, 1, -1
1522  DO i = 1, klon
1523  IF (todo(i) .AND. k<=(kb(i)-1) .AND. zrfl(i)>0.0) THEN
1524  zqev = max(0.0, (zqs(i,k)-q(i,k))*zfrac(i))
1525  zqevt = coef_eva*(1.0-q(i,k)/zqs(i,k))*sqrt(zrfl(i))* &
1526  (paprs(i,k)-paprs(i,k+1))/pplay(i, k)*t(i, k)*rd/rg
1527  zqevt = max(0.0, min(zqevt,zrfl(i)))*rg*dtime/ &
1528  (paprs(i,k)-paprs(i,k+1))
1529  zqev = min(zqev, zqevt)
1530  zrfln = zrfl(i) - zqev*(paprs(i,k)-paprs(i,k+1))/rg/dtime
1531  d_q(i, k) = -(zrfln-zrfl(i))*(rg/(paprs(i,k)-paprs(i,k+1)))*dtime
1532  d_t(i, k) = (zrfln-zrfl(i))*(rg/(paprs(i,k)-paprs(i, &
1533  k+1)))*dtime*rlvtt/rcpd
1534  zrfl(i) = zrfln
1535  END IF
1536  END DO
1537  END DO
1538  END IF
1539 
1540  ! La temperature de la premiere couche determine la pluie ou la neige:
1541 
1542  DO i = 1, klon
1543  IF (todo(i)) THEN
1544  IF (t(i,1)>rtt) THEN
1545  rain(i) = rain(i) + zrfl(i)
1546  ELSE
1547  snow(i) = snow(i) + zrfl(i)
1548  END IF
1549  END IF
1550  END DO
1551 
1552  RETURN
1553 END SUBROUTINE conkuo
1554 SUBROUTINE kuofcl(pt, pq, pg, pp, ldcum, kcbot)
1556  IMPLICIT NONE
1557  ! ======================================================================
1558  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19940927
1559  ! adaptation du code de Tiedtke du ECMWF
1560  ! Objet: calculer le niveau de convection libre
1561  ! (FCL: Free Convection Level)
1562  ! ======================================================================
1563  ! Arguments:
1564  ! pt---input-R- temperature (K)
1565  ! pq---input-R- vapeur d'eau (kg/kg)
1566  ! pg---input-R- geopotentiel (g*z ou z est en metre)
1567  ! pp---input-R- pression (Pa)
1568 
1569  ! LDCUM---output-L- Y-t-il la convection
1570  ! kcbot---output-I- Niveau du bas de la convection
1571  ! ======================================================================
1572  include "YOMCST.h"
1573  include "YOETHF.h"
1574 
1575  REAL pt(klon, klev), pq(klon, klev), pg(klon, klev), pp(klon, klev)
1576  INTEGER kcbot(klon)
1577  LOGICAL ldcum(klon)
1578 
1579  REAL ztu(klon, klev), zqu(klon, klev), zlu(klon, klev)
1580  REAL zqold(klon), zbuo
1581  INTEGER is, i, k
1582 
1583  ! klab=1: on est sous le nuage convectif
1584  ! klab=2: le bas du nuage convectif
1585  ! klab=0: autres couches
1586  INTEGER klab(klon, klev)
1587 
1588  ! quand lflag=.true., on est sous le nuage, il faut donc appliquer
1589  ! le processus d'elevation.
1590  LOGICAL lflag(klon)
1591 
1592  DO k = 1, klev
1593  DO i = 1, klon
1594  ztu(i, k) = pt(i, k)
1595  zqu(i, k) = pq(i, k)
1596  zlu(i, k) = 0.0
1597  klab(i, k) = 0
1598  END DO
1599  END DO
1600  ! ----------------------------------------------------------------------
1601  DO i = 1, klon
1602  klab(i, 1) = 1
1603  kcbot(i) = 2
1604  ldcum(i) = .false.
1605  END DO
1606 
1607  DO k = 2, klev - 1
1608 
1609  is = 0
1610  DO i = 1, klon
1611  IF (klab(i,k-1)==1) is = is + 1
1612  lflag(i) = .false.
1613  IF (klab(i,k-1)==1) lflag(i) = .true.
1614  END DO
1615  IF (is==0) GO TO 290
1616 
1617  ! on eleve le parcel d'air selon l'adiabatique sec
1618 
1619  DO i = 1, klon
1620  IF (lflag(i)) THEN
1621  zqu(i, k) = zqu(i, k-1)
1622  ztu(i, k) = ztu(i, k-1) + (pg(i,k-1)-pg(i,k))/rcpd
1623  zbuo = ztu(i, k)*(1.+retv*zqu(i,k)) - pt(i, k)*(1.+retv*pq(i,k)) + &
1624  0.5
1625  IF (zbuo>0.) klab(i, k) = 1
1626  zqold(i) = zqu(i, k)
1627  END IF
1628  END DO
1629 
1630  ! on calcule la condensation eventuelle
1631 
1632  CALL adjtq(pp(1,k), ztu(1,k), zqu(1,k), lflag, 1)
1633 
1634  ! s'il y a la condensation et la "buoyancy" force est positive
1635  ! c'est bien le bas de la tour de convection
1636 
1637  DO i = 1, klon
1638  IF (lflag(i) .AND. zqu(i,k)/=zqold(i)) THEN
1639  klab(i, k) = 2
1640  zlu(i, k) = zlu(i, k) + zqold(i) - zqu(i, k)
1641  zbuo = ztu(i, k)*(1.+retv*zqu(i,k)) - pt(i, k)*(1.+retv*pq(i,k)) + &
1642  0.5
1643  IF (zbuo>0.) THEN
1644  kcbot(i) = k
1645  ldcum(i) = .true.
1646  END IF
1647  END IF
1648  END DO
1649 
1650 290 END DO
1651 
1652  RETURN
1653 END SUBROUTINE kuofcl
1654 SUBROUTINE adjtq(pp, pt, pq, ldflag, kcall)
1656  IMPLICIT NONE
1657  ! ======================================================================
1658  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19940927
1659  ! adaptation du code de Tiedtke du ECMWF
1660  ! Objet: ajustement entre T et Q
1661  ! ======================================================================
1662  ! Arguments:
1663  ! pp---input-R- pression (Pa)
1664  ! pt---input/output-R- temperature (K)
1665  ! pq---input/output-R- vapeur d'eau (kg/kg)
1666  ! ======================================================================
1667  ! TO PRODUCE T,Q AND L VALUES FOR CLOUD ASCENT
1668 
1669  ! NOTE: INPUT PARAMETER KCALL DEFINES CALCULATION AS
1670  ! KCALL=0 ENV. T AND QS IN*CUINI*
1671  ! KCALL=1 CONDENSATION IN UPDRAFTS (E.G. CUBASE, CUASC)
1672  ! KCALL=2 EVAPORATION IN DOWNDRAFTS (E.G. CUDLFS,CUDDRAF)
1673 
1674  include "YOMCST.h"
1675 
1676  REAL pt(klon), pq(klon), pp(klon)
1677  LOGICAL ldflag(klon)
1678  INTEGER kcall
1679 
1680  REAL t_coup
1681  parameter(t_coup=234.0)
1682 
1683  REAL zcond(klon), zcond1
1684  REAL zdelta, zcvm5, zldcp, zqsat, zcor, zdqsat
1685  INTEGER is, i
1686  include "YOETHF.h"
1687  include "FCTTRE.h"
1688 
1689  DO i = 1, klon
1690  zcond(i) = 0.0
1691  END DO
1692 
1693  DO i = 1, klon
1694  IF (ldflag(i)) THEN
1695  zdelta = max(0., sign(1.,rtt-pt(i)))
1696  zldcp = rlvtt*(1.-zdelta) + zdelta*rlstt
1697  zldcp = zldcp/rcpd/(1.0+rvtmp2*pq(i))
1698  IF (thermcep) THEN
1699  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
1700  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*pq(i))
1701  zqsat = r2es*foeew(pt(i), zdelta)/pp(i)
1702  zqsat = min(0.5, zqsat)
1703  zcor = 1./(1.-retv*zqsat)
1704  zqsat = zqsat*zcor
1705  zdqsat = foede(pt(i), zdelta, zcvm5, zqsat, zcor)
1706  ELSE
1707  IF (pt(i)<t_coup) THEN
1708  zqsat = qsats(pt(i))/pp(i)
1709  zdqsat = dqsats(pt(i), zqsat)
1710  ELSE
1711  zqsat = qsatl(pt(i))/pp(i)
1712  zdqsat = dqsatl(pt(i), zqsat)
1713  END IF
1714  END IF
1715  zcond(i) = (pq(i)-zqsat)/(1.+zdqsat)
1716  IF (kcall==1) zcond(i) = max(zcond(i), 0.)
1717  IF (kcall==2) zcond(i) = min(zcond(i), 0.)
1718  pt(i) = pt(i) + zldcp*zcond(i)
1719  pq(i) = pq(i) - zcond(i)
1720  END IF
1721  END DO
1722 
1723  is = 0
1724  DO i = 1, klon
1725  IF (zcond(i)/=0.) is = is + 1
1726  END DO
1727  IF (is==0) GO TO 230
1728 
1729  DO i = 1, klon
1730  IF (ldflag(i) .AND. zcond(i)/=0.) THEN
1731  zdelta = max(0., sign(1.,rtt-pt(i)))
1732  zldcp = rlvtt*(1.-zdelta) + zdelta*rlstt
1733  zldcp = zldcp/rcpd/(1.0+rvtmp2*pq(i))
1734  IF (thermcep) THEN
1735  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
1736  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*pq(i))
1737  zqsat = r2es*foeew(pt(i), zdelta)/pp(i)
1738  zqsat = min(0.5, zqsat)
1739  zcor = 1./(1.-retv*zqsat)
1740  zqsat = zqsat*zcor
1741  zdqsat = foede(pt(i), zdelta, zcvm5, zqsat, zcor)
1742  ELSE
1743  IF (pt(i)<t_coup) THEN
1744  zqsat = qsats(pt(i))/pp(i)
1745  zdqsat = dqsats(pt(i), zqsat)
1746  ELSE
1747  zqsat = qsatl(pt(i))/pp(i)
1748  zdqsat = dqsatl(pt(i), zqsat)
1749  END IF
1750  END IF
1751  zcond1 = (pq(i)-zqsat)/(1.+zdqsat)
1752  pt(i) = pt(i) + zldcp*zcond1
1753  pq(i) = pq(i) - zcond1
1754  END IF
1755  END DO
1756 
1757 230 CONTINUE
1758  RETURN
1759 END SUBROUTINE adjtq
1760 SUBROUTINE fiajh(dtime, paprs, pplay, t, q, d_t, d_q, d_ql, rneb, rain, snow, &
1761  ibas, itop)
1763  IMPLICIT NONE
1764 
1765  ! Ajustement humide (Schema de convection de Manabe)
1766  ! .
1767  include "YOMCST.h"
1768 
1769  ! Arguments:
1770 
1771  REAL dtime ! intervalle du temps (s)
1772  REAL t(klon, klev) ! temperature (K)
1773  REAL q(klon, klev) ! humidite specifique (kg/kg)
1774  REAL paprs(klon, klev+1) ! pression a inter-couche (Pa)
1775  REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
1776 
1777  REAL d_t(klon, klev) ! incrementation pour la temperature
1778  REAL d_q(klon, klev) ! incrementation pour vapeur d'eau
1779  REAL d_ql(klon, klev) ! incrementation pour l'eau liquide
1780  REAL rneb(klon, klev) ! fraction nuageuse
1781 
1782  REAL rain(klon) ! variable non utilisee
1783  REAL snow(klon) ! variable non utilisee
1784  INTEGER ibas(klon) ! variable non utilisee
1785  INTEGER itop(klon) ! variable non utilisee
1786 
1787  REAL t_coup
1788  parameter(t_coup=234.0)
1789  REAL seuil_vap
1790  parameter(seuil_vap=1.0e-10)
1791 
1792  ! Variables locales:
1793 
1794  INTEGER i, k
1795  INTEGER k1, k1p, k2, k2p
1796  LOGICAL itest(klon)
1797  REAL delta_q(klon, klev)
1798  REAL cp_new_t(klev)
1799  REAL cp_delta_t(klev)
1800  REAL new_qb(klev)
1801  REAL v_cptj(klev), v_cptjk1, v_ssig
1802  REAL v_cptt(klon, klev), v_p, v_t
1803  REAL v_qs(klon, klev), v_qsd(klon, klev)
1804  REAL zq1(klon), zq2(klon)
1805  REAL gamcpdz(klon, 2:klev)
1806  REAL zdp, zdpm
1807 
1808  REAL zsat ! sur-saturation
1809  REAL zflo ! flotabilite
1810 
1811  REAL local_q(klon, klev), local_t(klon, klev)
1812 
1813  REAL zdelta, zcor, zcvm5
1814 
1815  include "YOETHF.h"
1816  include "FCTTRE.h"
1817 
1818  DO k = 1, klev
1819  DO i = 1, klon
1820  local_q(i, k) = q(i, k)
1821  local_t(i, k) = t(i, k)
1822  rneb(i, k) = 0.0
1823  d_ql(i, k) = 0.0
1824  d_t(i, k) = 0.0
1825  d_q(i, k) = 0.0
1826  END DO
1827  END DO
1828  DO i = 1, klon
1829  rain(i) = 0.0
1830  snow(i) = 0.0
1831  ibas(i) = 0
1832  itop(i) = 0
1833  END DO
1834 
1835  ! Calculer v_qs et v_qsd:
1836 
1837  DO k = 1, klev
1838  DO i = 1, klon
1839  v_cptt(i, k) = rcpd*local_t(i, k)
1840  v_t = local_t(i, k)
1841  v_p = pplay(i, k)
1842 
1843  IF (thermcep) THEN
1844  zdelta = max(0., sign(1.,rtt-v_t))
1845  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
1846  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*local_q(i,k))
1847  v_qs(i, k) = r2es*foeew(v_t, zdelta)/v_p
1848  v_qs(i, k) = min(0.5, v_qs(i,k))
1849  zcor = 1./(1.-retv*v_qs(i,k))
1850  v_qs(i, k) = v_qs(i, k)*zcor
1851  v_qsd(i, k) = foede(v_t, zdelta, zcvm5, v_qs(i,k), zcor)
1852  ELSE
1853  IF (v_t<t_coup) THEN
1854  v_qs(i, k) = qsats(v_t)/v_p
1855  v_qsd(i, k) = dqsats(v_t, v_qs(i,k))
1856  ELSE
1857  v_qs(i, k) = qsatl(v_t)/v_p
1858  v_qsd(i, k) = dqsatl(v_t, v_qs(i,k))
1859  END IF
1860  END IF
1861  END DO
1862  END DO
1863 
1864  ! Calculer Gamma * Cp * dz: (gamm est le gradient critique)
1865 
1866  DO k = 2, klev
1867  DO i = 1, klon
1868  zdp = paprs(i, k) - paprs(i, k+1)
1869  zdpm = paprs(i, k-1) - paprs(i, k)
1870  gamcpdz(i, k) = ((rd/rcpd/(zdpm+zdp)*(v_cptt(i,k-1)*zdpm+ &
1871  v_cptt(i,k)*zdp)+rlvtt/(zdpm+zdp)*(v_qs(i,k-1)*zdpm+ &
1872  v_qs(i,k)*zdp))*(pplay(i,k-1)-pplay(i,k))/paprs(i,k))/(1.0+(v_qsd(i, &
1873  k-1)*zdpm+v_qsd(i,k)*zdp)/(zdpm+zdp))
1874  END DO
1875  END DO
1876 
1877  ! ------------------------------------ modification des profils instables
1878  DO i = 1, klon
1879  itest(i) = .false.
1880 
1881  k1 = 0
1882  k2 = 1
1883 
1884 810 CONTINUE ! chercher k1, le bas de la colonne
1885  k2 = k2 + 1
1886  IF (k2>klev) GO TO 9999
1887  zflo = v_cptt(i, k2-1) - v_cptt(i, k2) - gamcpdz(i, k2)
1888  zsat = (local_q(i,k2-1)-v_qs(i,k2-1))*(paprs(i,k2-1)-paprs(i,k2)) + &
1889  (local_q(i,k2)-v_qs(i,k2))*(paprs(i,k2)-paprs(i,k2+1))
1890  IF (zflo<=0.0 .OR. zsat<=0.0) GO TO 810
1891  k1 = k2 - 1
1892  itest(i) = .true.
1893 
1894 820 CONTINUE ! chercher k2, le haut de la colonne
1895  IF (k2==klev) GO TO 821
1896  k2p = k2 + 1
1897  zsat = zsat + (paprs(i,k2p)-paprs(i,k2p+1))*(local_q(i,k2p)-v_qs(i,k2p))
1898  zflo = v_cptt(i, k2p-1) - v_cptt(i, k2p) - gamcpdz(i, k2p)
1899  IF (zflo<=0.0 .OR. zsat<=0.0) GO TO 821
1900  k2 = k2p
1901  GO TO 820
1902 821 CONTINUE
1903 
1904  ! ------------------------------------------------------ ajustement local
1905 830 CONTINUE ! ajustement proprement dit
1906  v_cptj(k1) = 0.0
1907  zdp = paprs(i, k1) - paprs(i, k1+1)
1908  v_cptjk1 = ((1.0+v_qsd(i,k1))*(v_cptt(i,k1)+v_cptj(k1))+rlvtt*(local_q(i, &
1909  k1)-v_qs(i,k1)))*zdp
1910  v_ssig = zdp*(1.0+v_qsd(i,k1))
1911 
1912  k1p = k1 + 1
1913  DO k = k1p, k2
1914  zdp = paprs(i, k) - paprs(i, k+1)
1915  v_cptj(k) = v_cptj(k-1) + gamcpdz(i, k)
1916  v_cptjk1 = v_cptjk1 + zdp*((1.0+v_qsd(i,k))*(v_cptt(i, &
1917  k)+v_cptj(k))+rlvtt*(local_q(i,k)-v_qs(i,k)))
1918  v_ssig = v_ssig + zdp*(1.0+v_qsd(i,k))
1919  END DO
1920 
1921  DO k = k1, k2
1922  cp_new_t(k) = v_cptjk1/v_ssig - v_cptj(k)
1923  cp_delta_t(k) = cp_new_t(k) - v_cptt(i, k)
1924  new_qb(k) = v_qs(i, k) + v_qsd(i, k)*cp_delta_t(k)/rlvtt
1925  local_q(i, k) = new_qb(k)
1926  local_t(i, k) = cp_new_t(k)/rcpd
1927  END DO
1928 
1929  ! --------------------------------------------------- sondage vers le bas
1930  ! -- on redefinit les variables prognostiques dans
1931  ! -- la colonne qui vient d'etre ajustee
1932 
1933  DO k = k1, k2
1934  v_cptt(i, k) = rcpd*local_t(i, k)
1935  v_t = local_t(i, k)
1936  v_p = pplay(i, k)
1937 
1938  IF (thermcep) THEN
1939  zdelta = max(0., sign(1.,rtt-v_t))
1940  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
1941  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*local_q(i,k))
1942  v_qs(i, k) = r2es*foeew(v_t, zdelta)/v_p
1943  v_qs(i, k) = min(0.5, v_qs(i,k))
1944  zcor = 1./(1.-retv*v_qs(i,k))
1945  v_qs(i, k) = v_qs(i, k)*zcor
1946  v_qsd(i, k) = foede(v_t, zdelta, zcvm5, v_qs(i,k), zcor)
1947  ELSE
1948  IF (v_t<t_coup) THEN
1949  v_qs(i, k) = qsats(v_t)/v_p
1950  v_qsd(i, k) = dqsats(v_t, v_qs(i,k))
1951  ELSE
1952  v_qs(i, k) = qsatl(v_t)/v_p
1953  v_qsd(i, k) = dqsatl(v_t, v_qs(i,k))
1954  END IF
1955  END IF
1956  END DO
1957  DO k = 2, klev
1958  zdpm = paprs(i, k-1) - paprs(i, k)
1959  zdp = paprs(i, k) - paprs(i, k+1)
1960  gamcpdz(i, k) = ((rd/rcpd/(zdpm+zdp)*(v_cptt(i,k-1)*zdpm+ &
1961  v_cptt(i,k)*zdp)+rlvtt/(zdpm+zdp)*(v_qs(i,k-1)*zdpm+ &
1962  v_qs(i,k)*zdp))*(pplay(i,k-1)-pplay(i,k))/paprs(i,k))/(1.0+(v_qsd(i, &
1963  k-1)*zdpm+v_qsd(i,k)*zdp)/(zdpm+zdp))
1964  END DO
1965 
1966  ! Verifier si l'on peut etendre la colonne vers le bas
1967 
1968  IF (k1==1) GO TO 841 ! extension echouee
1969  zflo = v_cptt(i, k1-1) - v_cptt(i, k1) - gamcpdz(i, k1)
1970  zsat = (local_q(i,k1-1)-v_qs(i,k1-1))*(paprs(i,k1-1)-paprs(i,k1)) + &
1971  (local_q(i,k1)-v_qs(i,k1))*(paprs(i,k1)-paprs(i,k1+1))
1972  IF (zflo<=0.0 .OR. zsat<=0.0) GO TO 841 ! extension echouee
1973 
1974 840 CONTINUE
1975  k1 = k1 - 1
1976  IF (k1==1) GO TO 830 ! GOTO 820 (a tester, Z.X.Li, mars 1995)
1977  zsat = zsat + (local_q(i,k1-1)-v_qs(i,k1-1))*(paprs(i,k1-1)-paprs(i,k1))
1978  zflo = v_cptt(i, k1-1) - v_cptt(i, k1) - gamcpdz(i, k1)
1979  IF (zflo>0.0 .AND. zsat>0.0) THEN
1980  GO TO 840
1981  ELSE
1982  GO TO 830 ! GOTO 820 (a tester, Z.X.Li, mars 1995)
1983  END IF
1984 841 CONTINUE
1985 
1986  GO TO 810 ! chercher d'autres blocks en haut
1987 
1988 9999 END DO ! boucle sur tous les points
1989  ! -----------------------------------------------------------------------
1990 
1991  ! Determiner la fraction nuageuse (hypothese: la nebulosite a lieu
1992  ! a l'endroit ou la vapeur d'eau est diminuee par l'ajustement):
1993 
1994  DO k = 1, klev
1995  DO i = 1, klon
1996  IF (itest(i)) THEN
1997  delta_q(i, k) = local_q(i, k) - q(i, k)
1998  IF (delta_q(i,k)<0.) rneb(i, k) = 1.0
1999  END IF
2000  END DO
2001  END DO
2002 
2003  ! Distribuer l'eau condensee en eau liquide nuageuse (hypothese:
2004  ! l'eau liquide est distribuee aux endroits ou la vapeur d'eau
2005  ! diminue et d'une maniere proportionnelle a cet diminution):
2006 
2007  DO i = 1, klon
2008  IF (itest(i)) THEN
2009  zq1(i) = 0.0
2010  zq2(i) = 0.0
2011  END IF
2012  END DO
2013  DO k = 1, klev
2014  DO i = 1, klon
2015  IF (itest(i)) THEN
2016  zdp = paprs(i, k) - paprs(i, k+1)
2017  zq1(i) = zq1(i) - delta_q(i, k)*zdp
2018  zq2(i) = zq2(i) - min(0.0, delta_q(i,k))*zdp
2019  END IF
2020  END DO
2021  END DO
2022  DO k = 1, klev
2023  DO i = 1, klon
2024  IF (itest(i)) THEN
2025  IF (zq2(i)/=0.0) d_ql(i, k) = -min(0.0, delta_q(i,k))*zq1(i)/zq2(i)
2026  END IF
2027  END DO
2028  END DO
2029 
2030  DO k = 1, klev
2031  DO i = 1, klon
2032  local_q(i, k) = max(local_q(i,k), seuil_vap)
2033  END DO
2034  END DO
2035 
2036  DO k = 1, klev
2037  DO i = 1, klon
2038  d_t(i, k) = local_t(i, k) - t(i, k)
2039  d_q(i, k) = local_q(i, k) - q(i, k)
2040  END DO
2041  END DO
2042 
2043  RETURN
2044 END SUBROUTINE fiajh
2045 SUBROUTINE fiajc(dtime, paprs, pplay, t, q, conv_q, d_t, d_q, d_ql, rneb, &
2046  rain, snow, ibas, itop)
2048  IMPLICIT NONE
2049 
2050  include "YOMCST.h"
2051 
2052  ! Options:
2053 
2054  INTEGER plb ! niveau de depart pour la convection
2055  parameter(plb=4)
2056 
2057  ! Mystere: cette option n'est pas innocente pour les resultats !
2058  ! Qui peut resoudre ce mystere ? (Z.X.Li mars 1995)
2059  LOGICAL vector ! calcul vectorise
2060  parameter(vector=.false.)
2061 
2062  REAL t_coup
2063  parameter(t_coup=234.0)
2064 
2065  ! Arguments:
2066 
2067  REAL q(klon, klev) ! humidite specifique (kg/kg)
2068  REAL t(klon, klev) ! temperature (K)
2069  REAL paprs(klon, klev+1) ! pression a inter-couche (Pa)
2070  REAL pplay(klon, klev) ! pression au milieu de couche (Pa)
2071  REAL dtime ! intervalle du temps (s)
2072  REAL conv_q(klon, klev) ! taux de convergence de l'humidite
2073  REAL rneb(klon, klev) ! fraction nuageuse
2074  REAL d_q(klon, klev) ! incrementaion pour la vapeur d'eau
2075  REAL d_ql(klon, klev) ! incrementation pour l'eau liquide
2076  REAL d_t(klon, klev) ! incrementation pour la temperature
2077  REAL rain(klon) ! variable non-utilisee
2078  REAL snow(klon) ! variable non-utilisee
2079  INTEGER itop(klon) ! variable non-utilisee
2080  INTEGER ibas(klon) ! variable non-utilisee
2081 
2082  INTEGER kh(klon), i, k
2083  LOGICAL nuage(klon), test(klon, klev)
2084  REAL zconv(klon), zdeh(klon, klev), zvirt(klon)
2085  REAL zdqs(klon, klev), zqs(klon, klev)
2086  REAL ztt, zvar, zfrac(klon)
2087  REAL zq1(klon), zq2(klon)
2088  REAL zdelta, zcor, zcvm5
2089 
2090  include "YOETHF.h"
2091  include "FCTTRE.h"
2092 
2093  ! Initialiser les sorties:
2094 
2095  DO k = 1, klev
2096  DO i = 1, klon
2097  rneb(i, k) = 0.0
2098  d_ql(i, k) = 0.0
2099  d_t(i, k) = 0.0
2100  d_q(i, k) = 0.0
2101  END DO
2102  END DO
2103  DO i = 1, klon
2104  itop(i) = 0
2105  ibas(i) = 0
2106  rain(i) = 0.0
2107  snow(i) = 0.0
2108  END DO
2109 
2110  ! Calculer Qs et L/Cp * dQs/dT:
2111 
2112  DO k = 1, klev
2113  DO i = 1, klon
2114  ztt = t(i, k)
2115  IF (thermcep) THEN
2116  zdelta = max(0., sign(1.,rtt-ztt))
2117  zcvm5 = r5les*rlvtt*(1.-zdelta) + zdelta*r5ies*rlstt
2118  zcvm5 = zcvm5/rcpd/(1.0+rvtmp2*q(i,k))
2119  zqs(i, k) = r2es*foeew(ztt, zdelta)/pplay(i, k)
2120  zqs(i, k) = min(0.5, zqs(i,k))
2121  zcor = 1./(1.-retv*zqs(i,k))
2122  zqs(i, k) = zqs(i, k)*zcor
2123  zdqs(i, k) = foede(ztt, zdelta, zcvm5, zqs(i,k), zcor)
2124  ELSE
2125  IF (ztt<t_coup) THEN
2126  zqs(i, k) = qsats(ztt)/pplay(i, k)
2127  zdqs(i, k) = dqsats(ztt, zqs(i,k))
2128  ELSE
2129  zqs(i, k) = qsatl(ztt)/pplay(i, k)
2130  zdqs(i, k) = dqsatl(ztt, zqs(i,k))
2131  END IF
2132  END IF
2133  END DO
2134  END DO
2135 
2136  ! Determiner la difference de l'energie totale saturee:
2137 
2138  DO i = 1, klon
2139  k = plb
2140  zdeh(i, k) = rcpd*(t(i,k-1)-t(i,k)) - rd*0.5*(t(i,k-1)+t(i,k))/paprs(i, k &
2141  )*(pplay(i,k-1)-pplay(i,k)) + rlvtt*(zqs(i,k-1)-zqs(i,k))
2142  zdeh(i, k) = zdeh(i, k)*0.5 ! on prend la moitie
2143  END DO
2144  DO k = plb + 1, klev
2145  DO i = 1, klon
2146  zdeh(i, k) = zdeh(i, k-1) + rcpd*(t(i,k-1)-t(i,k)) - &
2147  rd*0.5*(t(i,k-1)+t(i,k))/paprs(i, k)*(pplay(i,k-1)-pplay(i,k)) + &
2148  rlvtt*(zqs(i,k-1)-zqs(i,k))
2149  END DO
2150  END DO
2151 
2152  ! Determiner le sommet du nuage selon l'instabilite
2153  ! Calculer les convergences d'humidite (reelle et virtuelle)
2154 
2155  DO i = 1, klon
2156  nuage(i) = .true.
2157  zconv(i) = 0.0
2158  zvirt(i) = 0.0
2159  kh(i) = -999
2160  END DO
2161  DO k = plb, klev
2162  DO i = 1, klon
2163  nuage(i) = nuage(i) .AND. zdeh(i, k) > 0.0
2164  IF (nuage(i)) THEN
2165  kh(i) = k
2166  zconv(i) = zconv(i) + conv_q(i, k)*dtime*(paprs(i,k)-paprs(i,k+1))
2167  zvirt(i) = zvirt(i) + (zdeh(i,k)/rlvtt+zqs(i,k)-q(i,k))*(paprs(i,k)- &
2168  paprs(i,k+1))
2169  END IF
2170  END DO
2171  END DO
2172 
2173  IF (vector) THEN
2174 
2175 
2176  DO k = plb, klev
2177  DO i = 1, klon
2178  IF (k<=kh(i) .AND. kh(i)>plb .AND. zconv(i)>0.0) THEN
2179  test(i, k) = .true.
2180  zfrac(i) = max(0.0, min(zconv(i)/zvirt(i),1.0))
2181  ELSE
2182  test(i, k) = .false.
2183  END IF
2184  END DO
2185  END DO
2186 
2187  DO k = plb, klev
2188  DO i = 1, klon
2189  IF (test(i,k)) THEN
2190  zvar = zdeh(i, k)/(1.0+zdqs(i,k))
2191  d_q(i, k) = (zvar*zdqs(i,k)/rlvtt+zqs(i,k)-q(i,k))*zfrac(i) - &
2192  conv_q(i, k)*dtime
2193  d_t(i, k) = zvar*zfrac(i)/rcpd
2194  END IF
2195  END DO
2196  END DO
2197 
2198  DO i = 1, klon
2199  zq1(i) = 0.0
2200  zq2(i) = 0.0
2201  END DO
2202  DO k = plb, klev
2203  DO i = 1, klon
2204  IF (test(i,k)) THEN
2205  IF (d_q(i,k)<0.0) rneb(i, k) = zfrac(i)
2206  zq1(i) = zq1(i) - d_q(i, k)*(paprs(i,k)-paprs(i,k+1))
2207  zq2(i) = zq2(i) - min(0.0, d_q(i,k))*(paprs(i,k)-paprs(i,k+1))
2208  END IF
2209  END DO
2210  END DO
2211 
2212  DO k = plb, klev
2213  DO i = 1, klon
2214  IF (test(i,k)) THEN
2215  IF (zq2(i)/=0.) d_ql(i, k) = -min(0.0, d_q(i,k))*zq1(i)/zq2(i)
2216  END IF
2217  END DO
2218  END DO
2219 
2220  ELSE ! (.NOT. vector)
2221 
2222  DO i = 1, klon
2223  IF (kh(i)>plb .AND. zconv(i)>0.0) THEN
2224  ! cc IF (kh(i).LE.plb) GOTO 999 ! il n'y a pas d'instabilite
2225  ! cc IF (zconv(i).LE.0.0) GOTO 999 ! convergence insuffisante
2226  zfrac(i) = max(0.0, min(zconv(i)/zvirt(i),1.0))
2227  DO k = plb, kh(i)
2228  zvar = zdeh(i, k)/(1.0+zdqs(i,k))
2229  d_q(i, k) = (zvar*zdqs(i,k)/rlvtt+zqs(i,k)-q(i,k))*zfrac(i) - &
2230  conv_q(i, k)*dtime
2231  d_t(i, k) = zvar*zfrac(i)/rcpd
2232  END DO
2233 
2234  zq1(i) = 0.0
2235  zq2(i) = 0.0
2236  DO k = plb, kh(i)
2237  IF (d_q(i,k)<0.0) rneb(i, k) = zfrac(i)
2238  zq1(i) = zq1(i) - d_q(i, k)*(paprs(i,k)-paprs(i,k+1))
2239  zq2(i) = zq2(i) - min(0.0, d_q(i,k))*(paprs(i,k)-paprs(i,k+1))
2240  END DO
2241  DO k = plb, kh(i)
2242  IF (zq2(i)/=0.) d_ql(i, k) = -min(0.0, d_q(i,k))*zq1(i)/zq2(i)
2243  END DO
2244  END IF
2245  END DO
2246 
2247  END IF ! fin de teste sur vector
2248 
2249  RETURN
2250 END SUBROUTINE fiajc
subroutine conman(dtime, paprs, pplay, t, q, d_t, d_q, d_ql, rneb, rain, snow, ibas, itop)
Definition: conlmd.F90:92
subroutine conlmd(dtime, paprs, pplay, t, q, conv_q, d_t, d_q, rain, snow, ibas, itop)
Definition: conlmd.F90:6
subroutine adjtq(pp, pt, pq, ldflag, kcall)
Definition: conlmd.F90:1655
integer, save klon
Definition: dimphy.F90:3
integer, save klev
Definition: dimphy.F90:7
!$Id itapm1 ENDIF!IM on interpole les champs sur les niveaux STD de pression!IM a chaque pas de temps de la physique c!positionnement de l argument logique a false c!pour ne pas recalculer deux fois la meme chose!c!a cet effet un appel a plevel_new a ete deplace c!a la fin de la serie d appels c!la boucle DO nlevSTD a ete internalisee c!dans d ou la creation de cette routine c c!CALL false
Definition: calcul_STDlev.h:26
subroutine fiajh(dtime, paprs, pplay, t, q, d_t, d_q, d_ql, rneb, rain, snow, ibas, itop)
Definition: conlmd.F90:1762
!$Header!integer nvarmx parameter(nfmx=10, imx=200, jmx=150, lmx=200, nvarmx=1000) real xd(imx
!$Id itapm1 ENDIF!IM on interpole les champs sur les niveaux STD de pression!IM a chaque pas de temps de la physique c!positionnement de l argument logique a false c!pour ne pas recalculer deux fois la meme chose!c!a cet effet un appel a plevel_new a ete deplace c!a la fin de la serie d appels c!la boucle DO nlevSTD a ete internalisee c!dans d ou la creation de cette routine c c!CALL ulevSTD CALL &zphi philevSTD CALL &zx_rh rhlevSTD!DO klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon du jour ou toutes les read_climoz CALL true
subroutine conkuo(dtime, paprs, pplay, t, q, conv_q, d_t, d_q, d_ql, rneb, rain, snow, ibas, itop)
Definition: conlmd.F90:1047
subroutine fiajc(dtime, paprs, pplay, t, q, conv_q, d_t, d_q, d_ql, rneb, rain, snow, ibas, itop)
Definition: conlmd.F90:2047
subroutine conmanv(dtime, paprs, pplay, t, q, afaire, opt_cld, d_t, d_q, d_ql, rneb, rain, snow, ibas, itop, accompli, imprim)
Definition: conlmd.F90:299
Definition: dimphy.F90:1
subroutine kuofcl(pt, pq, pg, pp, ldcum, kcbot)
Definition: conlmd.F90:1555
real rg
Definition: comcstphy.h:1