LMDZ
cvltr_scav.F90
Go to the documentation of this file.
1 !
2 ! $Id $
3 !
4 SUBROUTINE cvltr_scav(pdtime, da, phi,phi2,d1a,dam, mpIN,epIN, &
5  sigd,sij,wght_cvfd,clw,elij,epmlmmm,eplamm, &
6  pmflxrin,pmflxsin,ev,te,wdtraina,wdtrainm, &
7  paprs,it,tr,upd,dnd,inb,icb, &
8  ccntraa_3d,ccntrenv_3d,coefcoli_3d, &
9  dtrcv,trsptd,dtrsscav,dtrsat,dtruscav,qdi,qpr, &
10  qpa,qmel,qtrdi,dtrcvma,mint, &
11  zmfd1a,zmfphi2,zmfdam)
12  !
13  USE ioipsl
14  USE dimphy
15  USE infotrac_phy, ONLY : nbtr,tname
16  IMPLICIT NONE
17  !=====================================================================
18  ! Objet : convection des traceurs / KE
19  ! Auteurs: M-A Filiberti and J-Y Grandpeix
20  ! modifiee par R Pilon : lessivage des traceurs / KE
21  !=====================================================================
22 
23  include "YOMCST.h"
24  include "YOECUMF.h"
25  include "conema3.h"
26  include "chem.h"
27 
28  ! Entree
29  REAL,INTENT(IN) :: pdtime
30  REAL,DIMENSION(klon,klev),INTENT(IN) :: da
31  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi
32  ! RomP
33  REAL,DIMENSION(klon,klev),INTENT(IN) :: d1a,dam ! matrices pour simplifier
34  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2 ! l'ecriture des tendances
35  !
36  REAL,DIMENSION(klon,klev),INTENT(IN) :: mpIN
37  REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression aux 1/2 couches (bas en haut)
38  INTEGER,INTENT(IN) :: it ! numero du traceur
39  REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: tr ! q de traceur (bas en haut)
40  REAL,DIMENSION(klon,klev),INTENT(IN) :: upd ! saturated updraft mass flux
41  REAL,DIMENSION(klon,klev),INTENT(IN) :: dnd ! saturated downdraft mass flux
42  !
43  REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainA ! masses precipitantes de l'asc adiab
44  REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainM ! masses precipitantes des melanges
45  !JE REAL,DIMENSION(klon,klev),INTENT(IN) :: pmflxrIN ! vprecip: eau
46  REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxrIN ! vprecip: eau
47  !JE REAL,DIMENSION(klon,klev),INTENT(IN) :: pmflxsIN ! vprecip: neige
48  REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxsIN ! vprecip: neige
49  REAL,DIMENSION(klon,klev),INTENT(IN) :: ev ! evaporation cv30_routine
50  REAL,DIMENSION(klon,klev),INTENT(IN) :: epIN
51  REAL,DIMENSION(klon,klev),INTENT(IN) :: te
52  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij ! fraction dair de lenv
53  REAL,DIMENSION(klon,klev),INTENT(IN) :: wght_cvfd ! weights of the layers feeding convection
54  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij ! contenu en eau condensée spécifique/conc deau condensée massique
55  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm ! eau condensee precipitee dans mel masse dair sat
56  REAL,DIMENSION(klon,klev),INTENT(IN) :: eplaMm ! eau condensee precipitee dans aa masse dair sat
57 
58  REAL,DIMENSION(klon,klev),INTENT(IN) :: clw ! contenu en eau condensée dans lasc adiab
59  REAL,DIMENSION(klon),INTENT(IN) :: sigd
60  INTEGER,DIMENSION(klon),INTENT(IN) :: icb,inb
61  !
62  REAL,DIMENSION(klon,klev),INTENT(IN) :: ccntrAA_3d
63  REAL,DIMENSION(klon,klev),INTENT(IN) :: ccntrENV_3d
64  REAL,DIMENSION(klon,klev),INTENT(IN) :: coefcoli_3d
65  !
66  ! Sortie
67  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrcv ! tendance totale (bas en haut)
68  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrcvMA ! M-A Filiberti
69  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: trsptd ! tendance du transport
70  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrSscav ! tendance du lessivage courant sat
71  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrsat ! tendance trsp+sat scav
72  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrUscav ! tendance du lessivage courant unsat
73  !
74  ! Variables locales
75  INTEGER :: i,j,k
76  REAL,DIMENSION(klon,klev) :: dxpres ! difference de pression entre niveau (j+1) et (j)
77  REAL :: pdtimeRG ! pas de temps * gravite
78  REAL,DIMENSION(klon,nbtr) :: qfeed ! tracer concentration feeding convection
79  ! variables pour les courants satures
80  REAL,DIMENSION(klon,klev,klev) :: zmd
81  REAL,DIMENSION(klon,klev,klev) :: za
82  REAL,DIMENSION(klon,klev,nbtr) :: zmfd,zmfa
83  REAL,DIMENSION(klon,klev,nbtr) :: zmfp,zmfu
84 
85  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfd1a
86  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfdam
87  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfphi2
88 
89  ! RomP ! les variables sont nettoyees des valeurs aberrantes
90  REAL,DIMENSION(klon,klev) :: Pa, Pm ! pluie AA et mélanges, var temporaire
91  REAL,DIMENSION(klon,klev) :: pmflxs,pmflxr ! pmflxrIN,pmflxsIN sans valeur aberante
92  REAL,DIMENSION(klon,klev) :: mp ! flux de masse
93  REAL,DIMENSION(klon,klev) :: ep ! fraction d'eau convertie en precipitation
94  REAL,DIMENSION(klon,klev) :: evap ! evaporation : variable temporaire
95  REAL,DIMENSION(klon,klev) :: rho !environmental density
96 
97  REAL,DIMENSION(klon,klev) :: kappa ! denominateur du au calcul de la matrice
98  ! pour obtenir qd et qp
99  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qTrdi ! traceurs descente air insature transport MA
100  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qDi ! traceurs descente insaturees
101  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qPr ! traceurs colonne precipitante
102  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qPa ! traceurs dans les precip issues lasc. adiab.
103  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qMel ! traceurs dans les precip issues des melanges
104  REAL,DIMENSION(klon,klev,nbtr) :: qMeltmp ! variable temporaire
105  REAL,DIMENSION(klon,klev,nbtr) :: qpmMint
106  REAL,DIMENSION(klon,klev),INTENT(OUT) :: Mint
107  ! tendances
108  REAL :: tdcvMA ! terme de transport de traceur (schema Marie Angele)
109  REAL :: trsptrac ! terme de transport de traceur par l'air
110  REAL :: scavtrac ! terme de lessivage courant sature
111  REAL :: uscavtrac ! terme de lessivage courant insature
112  ! impaction
113 !!! Correction apres discussion Romain P. / Olivier B.
114 !!! REAL,PARAMETER :: rdrop=2.5e-3 ! rayon des gouttes d'eau
115  REAL,PARAMETER :: rdrop=1.e-3 ! rayon des gouttes d'eau
116 !!!
117  REAL,DIMENSION(klon,klev) :: imp ! coefficient d'impaction
118  !
119  LOGICAL,DIMENSION(klon,klev) :: NO_precip
120  ! var tmp tests
121  REAL :: conserv
122  real :: conservMA
123 
124 !jyg<
125 !! ! ======================================================
126 !! ! calcul de l'impaction
127 !! ! ======================================================
128 !!
129 !! ! impaction sur la surface de la colonne de la descente insaturee
130 !! ! On prend la moyenne des precip entre le niveau i+1 et i
131 !! ! I=3/4* (P(1+1)+P(i))/2 / (sigd*r*rho_l)
132 !! ! 1000kg/m3= densite de l'eau
133 !! ! 0.75e-3 = 3/4 /1000
134 !! ! Par la suite, I est tout le temps multiplie par sig_d pour avoir l'impaction sur la surface de la maille
135 !!!! ! on le neglige ici pour simplifier le code
136 !!
137 !! DO j=1,klev-1
138 !! DO i=1,klon
139 !! imp(i,j) = coefcoli_3d(i,j)*0.75e-3/rdrop *&
140 !! 0.5*(pmflxr(i,j+1)+pmflxs(i,j+1)+pmflxr(i,j)+pmflxs(i,j))
141 !! ENDDO
142 !! ENDDO
144  !
145  ! initialisation pour flux de traceurs, td et autre
146  !
147  trsptrac = 0.
148  scavtrac = 0.
149  uscavtrac = 0.
150  qfeed(:,it) = 0. !RL
151  DO j=1,klev
152  DO i=1,klon
153  zmfd(i,j,it)=0.
154  zmfa(i,j,it)=0.
155  zmfu(i,j,it)=0.
156  zmfp(i,j,it)=0.
157  zmfphi2(i,j,it)=0.
158  zmfd1a(i,j,it)=0.
159  zmfdam(i,j,it)=0.
160  qdi(i,j,it)=0.
161  qpr(i,j,it)=0.
162  qpa(i,j,it)=0.
163  qmel(i,j,it)=0.
164  qmeltmp(i,j,it)=0.
165  qtrdi(i,j,it)=0.
166  kappa(i,j)=0.
167  trsptd(i,j,it)=0.
168  dtrsat(i,j,it)=0.
169  dtrsscav(i,j,it)=0.
170  dtruscav(i,j,it)=0.
171  dtrcv(i,j,it)=0.
172  dtrcvma(i,j,it)=0.
173  evap(i,j)=0.
174  dxpres(i,j)=0.
175  qpmmint(i,j,it)=0.
176  mint(i,j)=0.
177  END DO
178  END DO
179 
180  ! suppression des valeurs très faibles (~1e-320)
181  ! multiplication de levaporation pour lavoir par unite de temps
182  ! et par unite de surface de la maille
183  ! -> cv30_unsat : evap : masse evaporee/s/(m2 de la descente)
184  DO j=1,klev
185  DO i=1,klon
186  IF(ev(i,j).lt.1.e-16) THEN
187  evap(i,j)=0.
188  ELSE
189  evap(i,j)=ev(i,j)*sigd(i)
190  ENDIF
191  END DO
192  END DO
193 
194  DO j=1,klev
195  DO i=1,klon
196  IF(j.LT.klev) THEN
197  IF(epin(i,j).LT.1.e-32) THEN
198  ep(i,j)=0.
199  ELSE
200  ep(i,j)=epin(i,j)
201  ENDIF
202  ELSE
203  ep(i,j)=epmax
204  ENDIF
205  IF(mpin(i,j).LT.1.e-32) THEN
206  mp(i,j)=0.
207  ELSE
208  mp(i,j)=mpin(i,j)
209  ENDIF
210  IF(pmflxsin(i,j).LT.1.e-32) THEN
211  pmflxs(i,j)=0.
212  ELSE
213  pmflxs(i,j)=pmflxsin(i,j)
214  ENDIF
215  IF(pmflxrin(i,j).LT.1.e-32) THEN
216  pmflxr(i,j)=0.
217  ELSE
218  pmflxr(i,j)=pmflxrin(i,j)
219  ENDIF
220  IF(wdtraina(i,j).LT.1.e-32) THEN
221  pa(i,j)=0.
222  ELSE
223  pa(i,j)=wdtraina(i,j)
224  ENDIF
225  IF(wdtrainm(i,j).LT.1.e-32) THEN
226  pm(i,j)=0.
227  ELSE
228  pm(i,j)=wdtrainm(i,j)
229  ENDIF
230  END DO
231  END DO
232 
233  !==========================================
234  DO j = klev-1,1,-1
235  DO i = 1,klon
236  no_precip(i,j) = (pmflxr(i,j+1)+pmflxs(i,j+1)).LT.1.e-10&
237  .AND.pa(i,j).LT.1.e-10.AND.pm(i,j).LT.1.e-10
238  END DO
239  END DO
240 
241 !jyg<
242  ! ======================================================
243  ! calcul de l'impaction
244  ! ======================================================
245 
246  ! impaction sur la surface de la colonne de la descente insaturee
247  ! On prend la moyenne des precip entre le niveau i+1 et i
248  ! I=3/4* (P(1+1)+P(i))/2 / (sigd*r*rho_l)
249  ! 1000kg/m3= densite de l'eau
250  ! 0.75e-3 = 3/4 /1000
251  ! Par la suite, I est tout le temps multiplie par sig_d pour avoir l'impaction sur la surface de la maille
252  ! on le neglige ici pour simplifier le code
253 
254  DO j=1,klev-1
255  DO i=1,klon
256  imp(i,j) = coefcoli_3d(i,j)*0.75e-3/rdrop *&
257  0.5*(pmflxr(i,j+1)+pmflxs(i,j+1)+pmflxr(i,j)+pmflxs(i,j))
258  ENDDO
259  ENDDO
261  ! =========================================
262  ! calcul des tendances liees au downdraft
263  ! =========================================
264  !cdir collapse
265  DO k=1,klev
266  DO j=1,klev
267  DO i=1,klon
268  zmd(i,j,k)=0.
269  za(i,j,k)=0.
270  END DO
271  END DO
272  END DO
273  ! calcul de la matrice d echange
274  ! matrice de distribution de la masse entrainee en k
275  ! commmentaire RomP : mp > 0
276  DO k=1,klev-1
277  DO i=1,klon
278  zmd(i,k,k)=max(0.,mp(i,k)-mp(i,k+1)) ! ~ mk(k)
279  END DO
280  END DO
281  DO k=2,klev
282  DO j=k-1,1,-1
283  DO i=1,klon
284  IF(mp(i,j+1).GT.1.e-10) THEN
285  zmd(i,j,k)=zmd(i,j+1,k)*min(1.,mp(i,j)/mp(i,j+1)) !det ~ mk(j)=mk(j+1)*mp(i,j)/mp(i,j+1)
286  ENDIF
287  END DO
288  END DO
289  END DO
290  DO k=1,klev
291  DO j=1,klev-1
292  DO i=1,klon
293  za(i,j,k)=max(0.,zmd(i,j+1,k)-zmd(i,j,k))
294  END DO
295  END DO
296  END DO
297 !!!!! quantite de traceur dans la descente d'air insaturee : 4 juin 2012
298  DO k=1,klev
299  DO j=1,klev-1
300  DO i=1,klon
301  IF(mp(i,j+1).GT.1.e-10) THEN
302  qtrdi(i,j+1,it)=qtrdi(i,j+1,it)+(zmd(i,j+1,k)/mp(i,j+1))*tr(i,k,it)
303  ELSE
304  qtrdi(i,j,it)=0.!tr(i,j,it)
305  ENDIF
306  ENDDO
307  ENDDO
308  ENDDO
309 !!!!!
310  !
311  ! rajout du terme lie a l ascendance induite
312  !
313  DO j=2,klev
314  DO i=1,klon
315  za(i,j,j-1)=za(i,j,j-1)+mp(i,j)
316  END DO
317  END DO
318  !
319  ! tendance courants insatures ! sans lessivage ancien schema
320  !
321  DO k=1,klev
322  DO j=1,klev
323  DO i=1,klon
324  zmfd(i,j,it)=zmfd(i,j,it)+za(i,j,k)*(tr(i,k,it)-tr(i,j,it))
325  END DO
326  END DO
327  END DO
328  !
329  ! =========================================
330  ! calcul des tendances liees aux courants satures j <-> z ; k <-> z'
331  ! =========================================
332  !RL
333  ! Feeding concentrations
334  DO j=1,klev
335  DO i=1,klon
336  qfeed(i,it)=qfeed(i,it)+wght_cvfd(i,j)*tr(i,j,it)
337  END DO
338  END DO
339  !RL
340  !
341  DO j=1,klev
342  DO i=1,klon
343  !RL
344  !! zmfa(i,j,it)=da(i,j)*(tr(i,1,it)-tr(i,j,it)) ! da
345  zmfa(i,j,it)=da(i,j)*(qfeed(i,it)-tr(i,j,it)) ! da
346  !RL
347  END DO
348  END DO
349  !
350  DO k=1,klev
351  DO j=1,klev
352  DO i=1,klon
353  zmfp(i,j,it)=zmfp(i,j,it)+phi(i,j,k)*(tr(i,k,it)-tr(i,j,it)) ! phi
354  END DO
355  END DO
356  END DO
357  ! RomP ajout des matrices liees au lessivage
358  DO j=1,klev
359  DO i=1,klon
360  zmfd1a(i,j,it)=d1a(i,j)*tr(i,1,it) ! da1
361  zmfdam(i,j,it)=dam(i,j)*tr(i,1,it) ! dam
362  END DO
363  END DO
364  DO k=1,klev
365  DO j=1,klev
366  DO i=1,klon
367  zmfphi2(i,j,it)=zmfphi2(i,j,it)+phi2(i,j,k)*tr(i,k,it) ! psi
368  END DO
369  END DO
370  END DO
371  DO j=1,klev-1
372  DO i=1,klon
373  zmfu(i,j,it)=max(0.,upd(i,j+1)+dnd(i,j+1))*(tr(i,j+1,it)-tr(i,j,it))
374  END DO
375  END DO
376  DO j=2,klev
377  DO i=1,klon
378  zmfu(i,j,it)=zmfu(i,j,it)+min(0.,upd(i,j)+dnd(i,j))*(tr(i,j,it)-tr(i,j-1,it))
379  END DO
380  END DO
381  ! ===================================================
382  ! calcul des tendances liees aux courants insatures
383  ! ===================================================
384  ! pression
385  DO k=1, klev
386  DO i=1, klon
387  dxpres(i,k)=paprs(i,k)-paprs(i,k+1)
388  ENDDO
389  ENDDO
390  pdtimerg=pdtime*rg
391 
392  ! q_pa et q_pm traceurs issues des courants satures se retrouvant dans les precipitations
393  DO j=1,klev
394  DO i=1,klon
395  IF(j.GE.icb(i).AND.j.LE.inb(i)) THEN
396  IF(clw(i,j).GT.1.e-16) THEN
397  !JE qPa(i,j,it)=ccntrAA_coef*tr(i,1,it)/clw(i,j)
398  qpa(i,j,it)=ccntraa_3d(i,j)*tr(i,1,it)/clw(i,j)
399  ELSE
400  qpa(i,j,it)=0.
401  ENDIF
402  ENDIF
403  END DO
404  END DO
405 
406  ! calcul de q_pm en 2 parties :
407  ! 1) calcul de sa valeur pour un niveau z' donne
408  ! 2) integration sur la verticale sur z'
409  DO j=1,klev
410  DO k=1,j-1
411  DO i=1,klon
412  IF(k.GE.icb(i).AND.k.LE.inb(i).AND.&
413  j.LE.inb(i)) THEN
414  IF(elij(i,k,j).GT.1.e-16) THEN
415  !JE qMeltmp(i,j,it)=((1-ep(i,k))*ccntrAA_coef*tr(i,1,it)&
416  !JE *(1.-sij(i,k,j)) +ccntrENV_coef&
417  !JE *tr(i,k,it)*sij(i,k,j)) / elij(i,k,j)
418  qmeltmp(i,j,it)=((1-ep(i,k))*ccntraa_3d(i,k)*tr(i,1,it)&
419  *(1.-sij(i,k,j)) +ccntrenv_3d(i,k)&
420  *tr(i,k,it)*sij(i,k,j)) / elij(i,k,j)
421  ELSE
422  qmeltmp(i,j,it)=0.
423  ENDIF
424  qpmmint(i,j,it)=qpmmint(i,j,it) + qmeltmp(i,j,it)*epmlmmm(i,j,k)
425  mint(i,j)=mint(i,j) + epmlmmm(i,j,k)
426  ENDIF ! end if dans nuage
427  END DO
428  END DO
429  END DO
430 
431  DO j=1,klev
432  DO i=1,klon
433  IF(mint(i,j).GT.1.e-16) THEN
434  qmel(i,j,it)=qpmmint(i,j,it)/mint(i,j)
435  ELSE
436  qmel(i,j,it)=0.
437  ENDIF
438  END DO
439  END DO
440 
441  ! calcul de q_d et q_p traceurs de la descente precipitante
442  DO j=klev-1,1,-1
443  DO i=1,klon
444  IF(mp(i,j+1).GT.mp(i,j).AND.mp(i,j+1).GT.1.e-10) THEN ! detrainement
445  kappa(i,j)=((pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))*&
446  (-mp(i,j+1)-imp(i,j)/rg*dxpres(i,j))&
447  + (imp(i,j)/rg*dxpres(i,j))*(evap(i,j)/rg*dxpres(i,j)))
448 
449  ELSEIF(mp(i,j).GT.mp(i,j+1).AND.mp(i,j).GT.1.e-10) then! entrainement
450  IF(j.eq.1) THEN
451  kappa(i,j)=((pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))*&
452  (-mp(i,2)-imp(i,j)/rg*dxpres(i,j))&
453  + (imp(i,j)/rg*dxpres(i,j))*(evap(i,j)/rg*dxpres(i,j)))
454  ELSE
455  kappa(i,j)=((pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))*&
456  (-mp(i,j)-imp(i,j)/rg*dxpres(i,j))&
457  + (imp(i,j)/rg*dxpres(i,j))*(evap(i,j)/rg*dxpres(i,j)))
458  ENDIF
459  ELSE
460  kappa(i,j)=1.
461  ENDIF
462  ENDDO
463  ENDDO
464 
465  DO j=klev-1,1,-1
466  DO i=1,klon
467  IF (abs(kappa(i,j)).LT.1.e-25) THEN !si denominateur nul (il peut y avoir des mp!=0)
468  kappa(i,j)=1.
469  IF(j.eq.1) THEN
470  qdi(i,j,it)=qdi(i,j+1,it) !orig tr(i,j,it) ! mp(1)=0 donc tout vient de la couche supérieure
471  ELSEIF(mp(i,j+1).GT.mp(i,j).AND.mp(i,j+1).GT.1.e-10) THEN
472  qdi(i,j,it)=qdi(i,j+1,it)
473  ELSEIF(mp(i,j).GT.mp(i,j+1).AND.mp(i,j).GT.1.e-10) then! entrainement
474  qdi(i,j,it)=(-mp(i,j+1)*(qdi(i,j+1,it)-tr(i,j,it))-mp(i,j)*tr(i,j,it))/(-mp(i,j))
475  ELSE ! si mp (i)=0 et mp(j+1)=0
476  qdi(i,j,it)=tr(i,j,it) ! orig 0.
477  ENDIF
478 
479  IF(no_precip(i,j)) THEN
480  qpr(i,j,it)=0.
481  ELSE
482  qpr(i,j,it)=((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
483  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it)&
484  +imp(i,j)/rg*dxpres(i,j)*qdi(i,j,it))/&
485  (pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))
486  ENDIF
487  ELSE ! denominateur non nul
488  kappa(i,j)=1./kappa(i,j)
489  ! calcul de qd et qp
490  !!jyg (20130119) correction pour le sommet du nuage
491  !! if(j.GE.inb(i)) THEN !au-dessus du nuage, sommet inclu
492  if(j.GT.inb(i)) THEN !au-dessus du nuage
493  qdi(i,j,it)=tr(i,j,it) ! pas de descente => environnement = descente insaturee
494  qpr(i,j,it)=0.
495 
496  ! vvv premiere couche du modele ou mp(1)=0 ! det tout le temps vvv
497  ELSEIF(j.eq.1) THEN
498  if(mp(i,2).GT.1.e-10) THEN !mp(2) non nul -> detrainement (car mp(1) = 0) !ent pas possible
499  if(no_precip(i,j)) THEN !pas de precip en (i)
500  qdi(i,j,it)=qdi(i,j+1,it)
501  qpr(i,j,it)=0.
502  ELSE
503  qdi(i,j,it)=kappa(i,j)*(&
504  (-evap(i,j)/rg*dxpres(i,j))*((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
505  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it)) +&
506  (pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))*&
507  (-mp(i,j+1)*qdi(i,j+1,it)))
508 
509  qpr(i,j,it)=kappa(i,j)*(&
510  (-mp(i,j+1)-imp(i,j)/rg*dxpres(i,j))*&
511  ((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
512  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it))&
513  +(-mp(i,j+1)*qdi(i,j+1,it)) * (imp(i,j)/rg*dxpres(i,j)))
514  ENDIF
515 
516  ELSE !mp(2) nul -> plus de descente insaturee -> pluie agit sur environnement
517  qdi(i,j,it)=tr(i,j,it) ! orig 0.
518  if(no_precip(i,j)) THEN
519  qpr(i,j,it)=0.
520  ELSE
521  qpr(i,j,it)=((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
522  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it)&
523  +imp(i,j)/rg*dxpres(i,j)*tr(i,j,it))/&
524  (pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))
525  ENDIF
526 
527  ENDIF !mp(2) nul ou non
528 
529  ! vvv (j!=1.AND.j.LT.inb(i)) en-dessous du sommet nuage vvv
530  ELSE
531  !------------------------------------------------------------- detrainement
532  if(mp(i,j+1).GT.mp(i,j).AND.mp(i,j+1).GT.1.e-10) THEN !mp(i,j).GT.1.e-10) THEN
533  if(no_precip(i,j)) THEN
534  qdi(i,j,it)=qdi(i,j+1,it)
535  qpr(i,j,it)=0.
536  ELSE
537  qdi(i,j,it)=kappa(i,j)*(&
538  (-evap(i,j)/rg*dxpres(i,j))*((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
539  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it)) +&
540  (pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))*&
541  (-mp(i,j+1)*qdi(i,j+1,it)))
542  !
543  qpr(i,j,it)=kappa(i,j)*(&
544  (-mp(i,j+1)-imp(i,j)/rg*dxpres(i,j))*&
545  ((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
546  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it))&
547  +(-mp(i,j+1)*qdi(i,j+1,it)) * (imp(i,j)/rg*dxpres(i,j)))
548  ENDIF !precip
549  !------------------------------------------------------------- entrainement
550  ELSEIF(mp(i,j).GT.mp(i,j+1).AND.mp(i,j).GT.1.e-10) THEN
551  if(no_precip(i,j)) THEN
552  qdi(i,j,it)=(-mp(i,j+1)*(qdi(i,j+1,it)-tr(i,j,it))-mp(i,j)*tr(i,j,it))/(-mp(i,j))
553  qpr(i,j,it)=0.
554  ELSE
555  qdi(i,j,it)=kappa(i,j)*(&
556  (-evap(i,j)/rg*dxpres(i,j))*((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
557  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it)) +&
558  (pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))*&
559  (-mp(i,j+1)*(qdi(i,j+1,it)-tr(i,j,it))-mp(i,j)*tr(i,j,it)))
560  !
561  qpr(i,j,it)=kappa(i,j)*(&
562  (-mp(i,j)-imp(i,j)/rg*dxpres(i,j))*&
563  ((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
564  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it))&
565  +(-mp(i,j+1)*(qdi(i,j+1,it)-tr(i,j,it))-mp(i,j)*tr(i,j,it))*&
566  (imp(i,j)/rg*dxpres(i,j)))
567  ENDIF !precip
568  !------------------------------------------------------------- ENDIF ! ent/det
569  ELSE !mp nul
570  qdi(i,j,it)=tr(i,j,it) ! orig 0.
571  if(no_precip(i,j)) THEN
572  qpr(i,j,it)=0.
573  ELSE
574  qpr(i,j,it)=((pmflxr(i,j+1)+pmflxs(i,j+1))*qpr(i,j+1,it)+&
575  pa(i,j)*qpa(i,j,it)+pm(i,j)*qmel(i,j,it)&
576  +imp(i,j)/rg*dxpres(i,j)*tr(i,j,it))/&
577  (pmflxr(i,j+1)+pmflxs(i,j+1)+pa(i,j)+pm(i,j))
578  ENDIF
579  ENDIF ! mp nul ou non
580  ENDIF ! condition sur j
581  ENDIF ! kappa
582  ENDDO
583  ENDDO
584 
585  !! print test descente insaturee
586  ! DO j=klev,1,-1
587  ! DO i=1,klon
588  ! if(it.eq.3) THEN
589  ! write(*,'(I2,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12)') j,&
590  !! 'zmfdam',zmfdam(i,j,it),'zmfpsi',zmfphi2(i,j,it),&
591  ! 'zmfdam+zmfpsi',zmfdam(i,j,it)+zmfphi2(i,j,it),'qpmMint',qpmMint(i,j,it),&
592  ! 'Pm',Pm(i,j),'Mint',Mint(i,j),&
593  !! 'zmfa',zmfa(i,j,it),'zmfp',zmfp(i,j,it),&
594  ! 'zmfdam',zmfdam(i,j,it),'zmfpsi',zmfphi2(i,j,it),'zmfd1a',zmfd1a(i,j,it)
595  !! 'Pa',Pa(i,j),'eplaMm',eplaMm(i,j)
596  !! 'zmfd1a=da1*qa',zmfd1a(i,j,it),'Pa*qPa',wdtrainA(i,j)*qPa(i,j,it),'da1',d1a(i,j)
597  ! ENDIF
598  ! ENDDO
599  ! ENDDO
600 
601 
602  ! ===================================================
603  ! calcul final des tendances
604  ! ===================================================
605 
606  DO k=klev-1,1,-1
607  DO i=1, klon
608  ! transport
609  tdcvma=zmfd(i,k,it)+zmfu(i,k,it)+zmfa(i,k,it)+zmfp(i,k,it) ! double comptage des downdraft insatures
610  trsptrac=zmfu(i,k,it)+zmfa(i,k,it)+zmfp(i,k,it)
611  ! lessivage courants satures
612  !JE scavtrac=-ccntrAA_coef*zmfd1a(i,k,it)&
613  !JE -zmfphi2(i,k,it)*ccntrENV_coef&
614  !JE -zmfdam(i,k,it)*ccntrAA_coef
615  scavtrac=-ccntraa_3d(i,k)*zmfd1a(i,k,it)&
616  -zmfphi2(i,k,it)*ccntrenv_3d(i,k)&
617  -zmfdam(i,k,it)*ccntraa_3d(i,k)
618  ! lessivage courants insatures
619  if(k.LE.inb(i).AND.k.GT.1) THEN ! tendances dans le nuage
620  !------------------------------------------------------------- detrainement
621  if(mp(i,k+1).GT.mp(i,k).AND.mp(i,k+1).GT.1.e-10) THEN
622  uscavtrac= (-mp(i,k)+mp(i,k+1))*(qdi(i,k,it)-tr(i,k,it))&
623  + mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))
624  !
625  ! if(it.eq.3) write(*,'(I2,1X,a,5X,e20.12,82X,a,e20.12)')k,' det incloud',&
626  ! (-mp(i,k)+mp(i,k+1))*(qDi(i,k,it)-tr(i,k,it))*pdtimeRG/dxpres(i,k)+&
627  ! mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))*pdtimeRG/dxpres(i,k),&
628  ! 'mp',mp(i,k)
629  !------------------------------------------------------------- entrainement
630  ELSEIF(mp(i,k).GT.mp(i,k+1).AND.mp(i,k).GT.1.e-10) THEN
631  uscavtrac= mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))
632  !
633  ! if(it.eq.3) write(*,'(I2,1X,a,5X,e20.12,82X,a,e20.12)')k,' ent incloud',uscavtrac*pdtimeRG/dxpres(i,k), 'mp',mp(i,k)
634  !=!------------------------------------------------------------- end ent/det
635  ELSE ! mp(i,k+1)=0. et mp(i,k)=0. pluie directement sur l environnement
636 
637  if(no_precip(i,k)) THEN
638  uscavtrac=0.
639  ! if(it.eq.3) write(*,'(I2,1X,a,e20.12,82X,a,e20.12)')k,' no P ent incloud',uscavtrac*pdtimeRG/dxpres(i,k), 'mp',mp(i,k)
640  ELSE
641  uscavtrac=-imp(i,k)*tr(i,k,it)*dxpres(i,k)/rg+evap(i,k)*qpr(i,k,it)*dxpres(i,k)/rg
642  ! if(it.eq.3) write(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac*pdtimeRG/dxpres(i,k), 'mp',mp(i,k)
643  !!JE adds
644  ! if(it.eq.3) write(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'imp',imp(i,k)
645  ! if(it.eq.3) write(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'tr',tr(i,k,it)
646  ! if(it.eq.3) write(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'evap',evap(i,k)
647  ! if(it.eq.3) write(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'qPr',qPr(i,k,it)
648  ! if(it.eq.3) write(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'dxpres',dxpres(i,k)
649  !!Je end
650 
651  ENDIF
652  ENDIF ! mp/det/ent
653  !------------------------------------------------------------- premiere couche
654  ELSEIF(k.eq.1) THEN ! mp(1)=0.
655  if(mp(i,2).GT.1.e-10) THEN !detrainement
656  uscavtrac= (-0.+mp(i,2))*(qdi(i,k,it)-tr(i,k,it)) !&
657  ! + mp(i,2)*(0.-tr(i,k,it))
658  !
659  ! if(it.eq.3) write(*,'(I2,1X,a,e20.12,84X,a,e20.12)')k,' 1 det',&
660  ! (-0.+mp(i,2))*(qDi(i,k,it)-tr(i,k,it))*pdtimeRG/dxpres(i,k)+&
661  ! mp(i,2)*(0.-tr(i,k,it))*pdtimeRG/dxpres(i,k),&
662  ! 'mp',mp(i,k)
663  ELSE ! mp(2) = 0 = mp(1) pas de descente insaturee, rien ne se passe s'il ne pleut pas, sinon pluie->env
664  if(no_precip(i,1)) THEN
665  uscavtrac=0.
666  ELSE
667  uscavtrac=-imp(i,k)*tr(i,k,it)*dxpres(i,k)/rg+evap(i,k)*qpr(i,k,it)*dxpres(i,k)/rg
668  ENDIF
669  ! if(it.eq.3) write(*,'(I2,1X,a,2X,e20.12,82X,a,e20.12)')k,'1 P env incloud',uscavtrac*pdtimeRG/dxpres(i,k), 'mp',mp(i,k)
670  ENDIF
671 
672  ELSE ! k > INB au-dessus du nuage
673  uscavtrac=0.
674  ENDIF
675 
676  ! ===== tendances finales ======
677  trsptd(i,k,it)=trsptrac*pdtimerg/dxpres(i,k) ! td transport sans eau dans courants satures
678  dtrsscav(i,k,it)=scavtrac*pdtimerg/dxpres(i,k) ! td du lessivage dans courants satures
679  dtruscav(i,k,it)=uscavtrac*pdtimerg/dxpres(i,k) ! td courant insat
680  dtrsat(i,k,it)=(trsptrac+scavtrac)*pdtimerg/dxpres(i,k) ! td courant sat
681  dtrcv(i,k,it)=(trsptrac+scavtrac+uscavtrac)*pdtimerg/dxpres(i,k)!dtrsat(i,k,it)+dtrUscav(i,k,it) td conv
682 !!!!!!
683  dtrcvma(i,k,it)=tdcvma*pdtimerg/dxpres(i,k) ! MA tendance convection
684  ENDDO
685  ENDDO
686 
687  ! test de conservation du traceur
688  !print*,"_____________________________________________________________"
689  !print*," "
690  ! conserv=0.
691  ! conservMA=0.
692  ! DO k= klev-1,1,-1
693  ! DO i=1, klon
694  ! conserv=conserv+dtrcv(i,k,it)* &
695  ! (paprs(i,k)-paprs(i,k+1))/RG
696  ! conservMA=conservMA+dtrcvMA(i,k,it)* &
697  ! (paprs(i,k)-paprs(i,k+1))/RG
698  !
699  ! if(it.eq.3) write(*,'(I2,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12)') k,&
700  ! 'MA td ',dtrcvMA(i,k,it)*dxpres(i,k)/RG,&
701  ! ' td',dtrcv(i,k,it)*dxpres(i,k)/RG,' conservMA ',conservMA,'conserv ',conserv
702  !!
703  ! ENDDO
704  ! ENDDO
705  ! if(it.eq.3) print *,'it',it,'conserv ',conserv,'conservMA ',conservMA
706 
707 END SUBROUTINE cvltr_scav
infotrac_phy::nbtr
integer, save nbtr
Definition: infotrac_phy.F90:20
dimphy::klon
integer, save klon
Definition: dimphy.F90:3
infotrac_phy
Definition: infotrac_phy.F90:4
dimphy::klev
integer, save klev
Definition: dimphy.F90:7
it
!$Id Turb_fcg_gcssold get_uvd it
Definition: 1D_interp_cases.h:16
paprs
!$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 &zphi geo500!IM on interpole a chaque pas de temps le paprs
Definition: calcul_STDlev.h:144
infotrac_phy::tname
character(len=20), dimension(:), allocatable, save tname
Definition: infotrac_phy.F90:27
epmax
!$Id sig2feed!common comconema2 iflag_cvl_sigd common comconema1 epmax
Definition: conema3.h:15
cvltr_scav
subroutine cvltr_scav(pdtime, da, phi, phi2, d1a, dam, mpIN, epIN, sigd, sij, wght_cvfd, clw, elij, epmlmMm, eplaMm, pmflxrIN, pmflxsIN, ev, te, wdtrainA, wdtrainM, paprs, it, tr, upd, dnd, inb, icb, ccntrAA_3d, ccntrENV_3d, coefcoli_3d, dtrcv, trsptd, dtrSscav, dtrsat, dtrUscav, qDi, qPr, qPa, qMel, qTrdi, dtrcvMA, Mint, zmfd1a, zmfphi2, zmfdam)
Definition: cvltr_scav.F90:12
dimphy
Definition: dimphy.F90:1
sigd
!$Id!Parameters for nlm real sigd
Definition: cv30param.h:5
rg
real rg
Definition: comcstphy.h:1