LMDZ
cvltr.F90
Go to the documentation of this file.
1 !
2 ! $Id $
3 !
4 SUBROUTINE cvltr(pdtime, da, phi,phi2,d1a,dam, mpIN,epIN, &
5 !! sigd,sij,clw,elij,epmlmMm,eplaMm, & !RL
6  sigd,sij,wght_cvfd,clw,elij,epmlmmm,eplamm, & !RL
7  pmflxrin,pmflxsin,ev,te,wdtraina,wdtrainm, &
8  paprs,it,tr,upd,dnd,inb,icb, &
9  dtrcv,trsptd,dtrsscav,dtrsat,dtruscav,qdi,qpr, &
10  qpa,qmel,qtrdi,dtrcvma,mint, &
11  zmfd1a,zmfphi2,zmfdam)
12  USE ioipsl
13  USE dimphy
14  USE infotrac_phy, ONLY : nbtr,tname
15  IMPLICIT NONE
16 !=====================================================================
17 ! Objet : convection des traceurs / KE
18 ! Auteurs: M-A Filiberti and J-Y Grandpeix
19 ! modifiee par R Pilon : lessivage des traceurs / KE
20 !=====================================================================
21 
22  include "YOMCST.h"
23  include "YOECUMF.h"
24  include "conema3.h"
25 
26 ! Entree
27  REAL,INTENT(IN) :: pdtime
28  REAL,DIMENSION(klon,klev),INTENT(IN) :: da
29  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi
30 ! RomP
31  REAL,DIMENSION(klon,klev),INTENT(IN) :: d1a,dam ! matrices pour simplifier
32  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2 ! l'ecriture des tendances
33 !
34  REAL,DIMENSION(klon,klev),INTENT(IN) :: mpIN
35  REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression aux 1/2 couches (bas en haut)
36 ! REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression aux 1/2 couches (bas en haut)
37  REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: tr ! q de traceur (bas en haut)
38  INTEGER,INTENT(IN) :: it
39  REAL,DIMENSION(klon,klev),INTENT(IN) :: upd ! saturated updraft mass flux
40  REAL,DIMENSION(klon,klev),INTENT(IN) :: dnd ! saturated downdraft mass flux
41 !
42  REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainA ! masses precipitantes de l'asc adiab
43  REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainM ! masses precipitantes des melanges
44  REAL,DIMENSION(klon,klev),INTENT(IN) :: pmflxrIN ! vprecip: eau
45  REAL,DIMENSION(klon,klev),INTENT(IN) :: pmflxsIN ! vprecip: neige
46  REAL,DIMENSION(klon,klev),INTENT(IN) :: ev ! evaporation cv30_routine
47  REAL,DIMENSION(klon,klev),INTENT(IN) :: epIN
48  REAL,DIMENSION(klon,klev),INTENT(IN) :: te
49  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij ! fraction dair de lenv
50  REAL,DIMENSION(klon,klev),INTENT(IN) :: wght_cvfd ! weights of the layers feeding convection
51  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij ! contenu en eau condensée spécifique/conc deau condensée massique
52  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm ! eau condensee precipitee dans mel masse dair sat
53  REAL,DIMENSION(klon,klev),INTENT(IN) :: eplaMm ! eau condensee precipitee dans aa masse dair sat
54 
55  REAL,DIMENSION(klon,klev),INTENT(IN) :: clw ! contenu en eau condensée dans lasc adiab
56  REAL,DIMENSION(klon),INTENT(IN) :: sigd
57  INTEGER,DIMENSION(klon),INTENT(IN) :: icb,inb
58 ! Sortie
59  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrcv ! tendance totale (bas en haut)
60  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrcvMA ! M-A Filiberti
61  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: trsptd ! tendance du transport
62  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrSscav ! tendance du lessivage courant sat
63  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrsat ! tendance trsp+sat scav
64  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrUscav ! tendance du lessivage courant unsat
65 !
66 ! Variables locales
67  INTEGER :: i,j,k
68  REAL,DIMENSION(klon,klev) :: dxpres ! difference de pression entre niveau (j+1) et (j)
69  REAL :: pdtimeRG ! pas de temps * gravite
70 ! variables pour les courants satures
71  REAL,DIMENSION(klon,klev,klev) :: zmd
72  REAL,DIMENSION(klon,klev,klev) :: za
73  REAL,DIMENSION(klon,klev,nbtr) :: zmfd,zmfa
74  REAL,DIMENSION(klon,klev,nbtr) :: zmfp,zmfu
75  REAL,DIMENSION(klon,nbtr) :: qfeed ! tracer concentration feeding convection
76 
77  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfd1a
78  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfdam
79  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfphi2
80 
81 ! RomP ! les variables sont nettoyees des valeurs aberrantes
82  REAL,DIMENSION(klon,klev) :: Pa, Pm ! pluie AA et mélanges, var temporaire
83  REAL,DIMENSION(klon,klev) :: pmflxs,pmflxr ! pmflxrIN,pmflxsIN sans valeur aberante
84  REAL,DIMENSION(klon,klev) :: mp ! flux de masse
85  REAL,DIMENSION(klon,klev) :: ep ! fraction d'eau convertie en precipitation
86  REAL,DIMENSION(klon,klev) :: evap ! evaporation : variable temporaire
87  REAL,DIMENSION(klon,klev) :: rho !environmental density
88 
89  REAL,DIMENSION(klon,klev) :: kappa ! denominateur du au calcul de la matrice
90  ! pour obtenir qd et qp
91  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qTrdi ! traceurs descente air insature transport MA
92  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qDi ! traceurs descente insaturees
93  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qPr ! traceurs colonne precipitante
94  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qPa ! traceurs dans les precip issues lasc. adiab.
95  REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qMel ! traceurs dans les precip issues des melanges
96  REAL,DIMENSION(klon,klev,nbtr) :: qMeltmp ! variable temporaire
97  REAL,DIMENSION(klon,klev,nbtr) :: qpmMint
98  REAL,DIMENSION(klon,klev),INTENT(OUT) :: Mint
99 ! tendances
100  REAL :: tdcvMA ! terme de transport de traceur (schema Marie Angele)
101  REAL :: trsptrac ! terme de transport de traceur par l'air
102  REAL :: scavtrac ! terme de lessivage courant sature
103  REAL :: uscavtrac ! terme de lessivage courant insature
104 ! impaction
105 !!! Correction apres discussion Romain P. / Olivier B.
106 !!! REAL,PARAMETER :: rdrop=2.5e-3 ! rayon des gouttes d'eau
107  REAL,PARAMETER :: rdrop=1.e-3 ! rayon des gouttes d'eau
108 !!!
109  REAL,DIMENSION(klon,klev) :: imp ! coefficient d'impaction
110 ! parametres lessivage
111  REAL :: ccntrAA_coef ! \alpha_a : fract aerosols de l'AA convertis en CCN
112  REAL :: ccntrENV_coef ! \beta_m : fract aerosols de l'env convertis en CCN
113  REAL :: coefcoli ! coefficient de collision des gouttes sur les aerosols
114 !
115  LOGICAL,DIMENSION(klon,klev) :: NO_precip
116 ! LOGICAL :: scavON
117 ! var tmp tests
118  REAL :: conserv
119  real :: conservMA
120 
121 ! coefficient lessivage
122  ccntraa_coef = 0.
123  ccntrenv_coef = 0.
124  coefcoli = 0.
125 
126 !$OMP MASTER
127  call getin('ccntrAA_coef',ccntraa_coef)
128  call getin('ccntrENV_coef',ccntrenv_coef)
129  call getin('coefcoli',coefcoli)
130 !$OMP END MASTER
131 !$OMP BARRIER
132  print*,'cvltr coef lessivage convectif', ccntraa_coef,ccntrenv_coef,coefcoli
133 
134 ! scavON=.TRUE.
135 ! if(scavON) then
136 ! ccntrAA_coef = 1.
137 ! ccntrENV_coef = 1.
138 ! coefcoli = 1.
139 ! else
140 ! ccntrAA_coef = 0.
141 ! ccntrENV_coef = 0.
142 ! coefcoli = 0.
143 ! endif
144 
145 ! ======================================================
146 ! calcul de l'impaction
147 ! ======================================================
148 !initialisation
149  do j=1,klev
150  do i=1,klon
151  imp(i,j)=0.
152  enddo
153  enddo
154 ! impaction sur la surface de la colonne de la descente insaturee
155 ! On prend la moyenne des precip entre le niveau i+1 et i
156 ! I=3/4* (P(1+1)+P(i))/2 / (sigd*r*rho_l)
157 ! 1000kg/m3= densité de l'eau
158 ! 0.75e-3 = 3/4 /1000
159 ! Par la suite, I est tout le temps multiplié par sig_d pour avoir l'impaction sur la surface de la maille
160 ! on le néglige ici pour simplifier le code
161  do j=1,klev-1
162  do i=1,klon
163  imp(i,j) = coefcoli*0.75e-3/rdrop *&
164  0.5*(pmflxr(i,j+1)+pmflxs(i,j+1)+pmflxr(i,j)+pmflxs(i,j))
165 ! rho(i,j)=pplay(i,j)/(rd*te(i,j))
166  enddo
167  enddo
168 !
169 ! initialisation pour flux de traceurs, td et autre
170  trsptrac = 0.
171  scavtrac = 0.
172  uscavtrac = 0.
173  qfeed(:,it) = 0. !RL
174  DO j=1,klev
175  DO i=1,klon
176  zmfd(i,j,it)=0.
177  zmfa(i,j,it)=0.
178  zmfu(i,j,it)=0.
179  zmfp(i,j,it)=0.
180  zmfphi2(i,j,it)=0.
181  zmfd1a(i,j,it)=0.
182  zmfdam(i,j,it)=0.
183  qdi(i,j,it)=0.
184  qpr(i,j,it)=0.
185  qpa(i,j,it)=0.
186  qmel(i,j,it)=0.
187  qmeltmp(i,j,it)=0.
188  qtrdi(i,j,it)=0.
189  kappa(i,j)=0.
190  trsptd(i,j,it)=0.
191  dtrsat(i,j,it)=0.
192  dtrsscav(i,j,it)=0.
193  dtruscav(i,j,it)=0.
194  dtrcv(i,j,it)=0.
195  dtrcvma(i,j,it)=0.
196  evap(i,j)=0.
197  dxpres(i,j)=0.
198  qpmmint(i,j,it)=0.
199  mint(i,j)=0.
200  END DO
201  END DO
202 
203 ! suppression des valeurs très faibles (~1e-320)
204 ! multiplication de levaporation pour lavoir par unite de temps
205 ! et par unite de surface de la maille
206 ! -> cv30_unsat : evap : masse evaporee/s/(m2 de la descente)
207  DO j=1,klev
208  DO i=1,klon
209  if(ev(i,j).lt.1.e-16) then
210  evap(i,j)=0.
211  else
212  evap(i,j)=ev(i,j)*sigd(i)
213  endif
214  END DO
215  END DO
216 
217  DO j=1,klev
218  DO i=1,klon
219  if(j.lt.klev) then
220  if(epin(i,j).lt.1.e-32) then
221  ep(i,j)=0.
222  else
223  ep(i,j)=epin(i,j)
224  endif
225  else
226  ep(i,j)=epmax
227  endif
228  if(mpin(i,j).lt.1.e-32) then
229  mp(i,j)=0.
230  else
231  mp(i,j)=mpin(i,j)
232  endif
233  if(pmflxsin(i,j).lt.1.e-32) then
234  pmflxs(i,j)=0.
235  else
236  pmflxs(i,j)=pmflxsin(i,j)
237  endif
238  if(pmflxrin(i,j).lt.1.e-32) then
239  pmflxr(i,j)=0.
240  else
241  pmflxr(i,j)=pmflxrin(i,j)
242  endif
243  if(wdtraina(i,j).lt.1.e-32) then
244  pa(i,j)=0.
245  else
246  pa(i,j)=wdtraina(i,j)
247  endif
248  if(wdtrainm(i,j).lt.1.e-32) then
249  pm(i,j)=0.
250  else
251  pm(i,j)=wdtrainm(i,j)
252  endif
253  END DO
254  END DO
255 
256 !==========================================
257  DO j = klev-1,1,-1
258  DO i = 1,klon
259  no_precip(i,j) = (pmflxr(i,j+1)+pmflxs(i,j+1)).lt.1.e-10&
260  .and.pa(i,j).lt.1.e-10.and.pm(i,j).lt.1.e-10
261  END DO
262  END DO
263 
264 ! =========================================
265 ! calcul des tendances liees au downdraft
266 ! =========================================
267 !cdir collapse
268  DO k=1,klev
269  DO j=1,klev
270  DO i=1,klon
271  zmd(i,j,k)=0.
272  za(i,j,k)=0.
273  END DO
274  END DO
275  END DO
276 ! calcul de la matrice d echange
277 ! matrice de distribution de la masse entrainee en k
278 ! commmentaire RomP : mp > 0
279  DO k=1,klev-1
280  DO i=1,klon
281  zmd(i,k,k)=max(0.,mp(i,k)-mp(i,k+1)) ! ~ mk(k)
282  END DO
283  END DO
284  DO k=2,klev
285  DO j=k-1,1,-1
286  DO i=1,klon
287  if(mp(i,j+1).gt.1.e-10) then
288  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)
289  ENDif
290  END DO
291  END DO
292  END DO
293  DO k=1,klev
294  DO j=1,klev-1
295  DO i=1,klon
296  za(i,j,k)=max(0.,zmd(i,j+1,k)-zmd(i,j,k))
297  END DO
298  END DO
299  END DO
300 !!!!! quantite de traceur dans la descente d'air insaturee : 4 juin 2012
301  DO k=1,klev
302  DO j=1,klev-1
303  DO i=1,klon
304  if(mp(i,j+1).gt.1.e-10) then
305  qtrdi(i,j+1,it)=qtrdi(i,j+1,it)+(zmd(i,j+1,k)/mp(i,j+1))*tr(i,k,it)
306  else
307  qtrdi(i,j,it)=0.!tr(i,j,it)
308  endif
309  ENDDO
310  ENDDO
311  ENDDO
312 !!!!!
313 !
314 ! rajout du terme lie a l ascendance induite
315 !
316  DO j=2,klev
317  DO i=1,klon
318  za(i,j,j-1)=za(i,j,j-1)+mp(i,j)
319  END DO
320  END DO
321 !
322 ! tendance courants insatures ! sans lessivage ancien schema
323 !
324  DO k=1,klev
325  DO j=1,klev
326  DO i=1,klon
327  zmfd(i,j,it)=zmfd(i,j,it)+za(i,j,k)*(tr(i,k,it)-tr(i,j,it))
328  END DO
329  END DO
330  END DO
331 !
332 ! =========================================
333 ! calcul des tendances liees aux courants satures j <-> z ; k <-> z'
334 ! =========================================
335 !
336 !RL
337 ! Feeding concentrations
338  DO j=1,klev
339  DO i=1,klon
340  qfeed(i,it)=qfeed(i,it)+wght_cvfd(i,j)*tr(i,j,it)
341  END DO
342  END DO
343 !RL
344 !
345  DO j=1,klev
346  DO i=1,klon
347 !RL
348 !! zmfa(i,j,it)=da(i,j)*(tr(i,1,it)-tr(i,j,it)) ! da
349  zmfa(i,j,it)=da(i,j)*(qfeed(i,it)-tr(i,j,it)) ! da
350 !RL
351  END DO
352  END DO
353 !
354  DO k=1,klev
355  DO j=1,klev
356  DO i=1,klon
357  zmfp(i,j,it)=zmfp(i,j,it)+phi(i,j,k)*(tr(i,k,it)-tr(i,j,it)) ! phi
358  END DO
359  END DO
360  END DO
361 ! RomP ajout des matrices liees au lessivage
362  DO j=1,klev
363  DO i=1,klon
364  zmfd1a(i,j,it)=d1a(i,j)*tr(i,1,it) ! da1
365  zmfdam(i,j,it)=dam(i,j)*tr(i,1,it) ! dam
366  END DO
367  END DO
368  DO k=1,klev
369  DO j=1,klev
370  DO i=1,klon
371  zmfphi2(i,j,it)=zmfphi2(i,j,it)+phi2(i,j,k)*tr(i,k,it) ! psi
372  END DO
373  END DO
374  END DO
375  DO j=1,klev-1
376  DO i=1,klon
377  zmfu(i,j,it)=max(0.,upd(i,j+1)+dnd(i,j+1))*(tr(i,j+1,it)-tr(i,j,it))
378  END DO
379  END DO
380  DO j=2,klev
381  DO i=1,klon
382  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))
383  END DO
384  END DO
385 ! ===================================================
386 ! calcul des tendances liees aux courants insatures
387 ! ===================================================
388 ! pression
389  DO k=1, klev
390  DO i=1, klon
391  dxpres(i,k)=paprs(i,k)-paprs(i,k+1)
392  ENDDO
393  ENDDO
394  pdtimerg=pdtime*rg
395 
396 ! q_pa et q_pm traceurs issues des courants satures se retrouvant dans les precipitations
397  DO j=1,klev
398  DO i=1,klon
399  if(j.ge.icb(i).and.j.le.inb(i)) then
400  if(clw(i,j).gt.1.e-16) then
401  qpa(i,j,it)=ccntraa_coef*tr(i,1,it)/clw(i,j)
402  else
403  qpa(i,j,it)=0.
404  endif
405  endif
406  END DO
407  END DO
408 
409 ! calcul de q_pm en 2 parties :
410 ! 1) calcul de sa valeur pour un niveau z' donne
411 ! 2) integration sur la verticale sur z'
412  DO j=1,klev
413  DO k=1,j-1
414  DO i=1,klon
415  if(k.ge.icb(i).and.k.le.inb(i).and.&
416  j.le.inb(i)) then
417  if(elij(i,k,j).gt.1.e-16) then
418  qmeltmp(i,j,it)=((1-ep(i,k))*ccntraa_coef*tr(i,1,it)&
419  *(1.-sij(i,k,j)) +ccntrenv_coef&
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  scavtrac=-ccntraa_coef*zmfd1a(i,k,it)&
613  -zmfphi2(i,k,it)*ccntrenv_coef&
614  -zmfdam(i,k,it)*ccntraa_coef
615 ! lessivage courants insatures
616  if(k.le.inb(i).and.k.gt.1) then ! tendances dans le nuage
617 !------------------------------------------------------------- detrainement
618  if(mp(i,k+1).gt.mp(i,k).and.mp(i,k+1).gt.1.e-10) then
619  uscavtrac= (-mp(i,k)+mp(i,k+1))*(qdi(i,k,it)-tr(i,k,it))&
620  + mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))
621 !
622 ! if(it.eq.3) write(*,'(I2,1X,a,5X,e20.12,82X,a,e20.12)')k,' det incloud',&
623 ! (-mp(i,k)+mp(i,k+1))*(qDi(i,k,it)-tr(i,k,it))*pdtimeRG/dxpres(i,k)+&
624 ! mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))*pdtimeRG/dxpres(i,k),&
625 ! 'mp',mp(i,k)
626 !------------------------------------------------------------- entrainement
627  elseif(mp(i,k).gt.mp(i,k+1).and.mp(i,k).gt.1.e-10) then
628  uscavtrac= mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))
629 !
630 ! 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)
631 !=!------------------------------------------------------------- end ent/det
632  else ! mp(i,k+1)=0. et mp(i,k)=0. pluie directement sur l environnement
633 
634  if(no_precip(i,k)) then
635  uscavtrac=0.
636 ! 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)
637  else
638  uscavtrac=-imp(i,k)*tr(i,k,it)*dxpres(i,k)/rg+evap(i,k)*qpr(i,k,it)*dxpres(i,k)/rg
639 ! 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)
640  endif
641  endif ! mp/det/ent
642 !------------------------------------------------------------- premiere couche
643  elseif(k.eq.1) then ! mp(1)=0.
644  if(mp(i,2).gt.1.e-10) then !detrainement
645  uscavtrac= (-0.+mp(i,2))*(qdi(i,k,it)-tr(i,k,it)) !&
646 ! + mp(i,2)*(0.-tr(i,k,it))
647 !
648 ! if(it.eq.3) write(*,'(I2,1X,a,e20.12,84X,a,e20.12)')k,' 1 det',&
649 ! (-0.+mp(i,2))*(qDi(i,k,it)-tr(i,k,it))*pdtimeRG/dxpres(i,k)+&
650 ! mp(i,2)*(0.-tr(i,k,it))*pdtimeRG/dxpres(i,k),&
651 ! 'mp',mp(i,k)
652  else ! mp(2) = 0 = mp(1) pas de descente insaturee, rien ne se passe s'il ne pleut pas, sinon pluie->env
653  if(no_precip(i,1)) then
654  uscavtrac=0.
655  else
656  uscavtrac=-imp(i,k)*tr(i,k,it)*dxpres(i,k)/rg+evap(i,k)*qpr(i,k,it)*dxpres(i,k)/rg
657  endif
658 ! 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)
659  endif
660 
661  else ! k > INB au-dessus du nuage
662  uscavtrac=0.
663  endif
664 
665 ! ===== tendances finales ======
666  trsptd(i,k,it)=trsptrac*pdtimerg/dxpres(i,k) ! td transport sans eau dans courants satures
667  dtrsscav(i,k,it)=scavtrac*pdtimerg/dxpres(i,k) ! td du lessivage dans courants satures
668  dtruscav(i,k,it)=uscavtrac*pdtimerg/dxpres(i,k) ! td courant insat
669  dtrsat(i,k,it)=(trsptrac+scavtrac)*pdtimerg/dxpres(i,k) ! td courant sat
670  dtrcv(i,k,it)=(trsptrac+scavtrac+uscavtrac)*pdtimerg/dxpres(i,k)!dtrsat(i,k,it)+dtrUscav(i,k,it) td conv
671 !!!!!!
672  dtrcvma(i,k,it)=tdcvma*pdtimerg/dxpres(i,k) ! MA tendance convection
673  ENDDO
674  ENDDO
675 
676 ! test de conservation du traceur
677 !print*,"_____________________________________________________________"
678 !print*," "
679 ! conserv=0.
680 ! conservMA=0.
681 ! DO k= klev-1,1,-1
682 ! DO i=1, klon
683 ! conserv=conserv+dtrcv(i,k,it)* &
684 ! (paprs(i,k)-paprs(i,k+1))/RG
685 ! conservMA=conservMA+dtrcvMA(i,k,it)* &
686 ! (paprs(i,k)-paprs(i,k+1))/RG
687 !
688 ! if(it.eq.3) write(*,'(I2,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12)') k,&
689 ! 'MA td ',dtrcvMA(i,k,it)*dxpres(i,k)/RG,&
690 ! ' td',dtrcv(i,k,it)*dxpres(i,k)/RG,' conservMA ',conservMA,'conserv ',conserv
691 !!
692 ! ENDDO
693 ! ENDDO
694 ! if(it.eq.3) print *,'it',it,'conserv ',conserv,'conservMA ',conservMA
695 
696 END SUBROUTINE cvltr
infotrac_phy::nbtr
integer, save nbtr
Definition: infotrac_phy.F90:20
cvltr
subroutine cvltr(pdtime, da, phi, phi2, d1a, dam, mpIN, epIN,
Definition: cvltr.F90:6
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
dimphy
Definition: dimphy.F90:1
sigd
!$Id!Parameters for nlm real sigd
Definition: cv30param.h:5
rg
real rg
Definition: comcstphy.h:1