source: LMDZ5/trunk/libf/phylmd/phytrac_mod.F90 @ 1816

Last change on this file since 1816 was 1813, checked in by idelkadi, 11 years ago
  1. transform phytrac into a module, in order to pass some variables

(tracer tendencies) to the standard physiq ouput codes.

  1. Correct a (big) bug in the call to phytrac.
  2. Include w*, and ALEs in the call to phytrac and traclmdz.

physiq.F

  • Bug correction in the call of phytrac from the physics u10m,v10m, ustar -> zu10m, zv10m, zustar

phytrac.F90 -> phytrac_mod.F90

  • Tranformation of routine phytrac into a module phytrac_mod, in order to tranfer the tracer tendencies from phytrac to

phys_output...

  • Inclusion of w*, Ale bl/wake in the call to phytrac and traclmdz

(to be used for dust emmission)

by respectively, vertical diffusion, thermal plumes and convection

  • desactivation of ini_histrac.h and write_histrac.h
  • USE phys_output_mod removed since it was creating a circular

dependency

between phytrac_mod and phys_output_mod.
So the automatic computation of ecrit_tra is desactivated

ini_histrac.h and write_histrac.h

Descactivated in phytrac but kept for backard compatibility
couchelimite -> iflag_vdf_trac>0

phys_output_ctrlout_mod.F90

New variables : o_dtr_vdf, o_dtr_the ... for output of tracer tendencies

phys_output_mod.F90

Default definition for these new output variables.

phys_output_write_F90.h

disapears, included directly in phys_output_write_mod.F90

phys_output_write_mod.F90

writing of the tracer tendencies

phys_state_var_mod.F90

New declaration (wstar)

traclmdz_mod.F90

  • Inclusion of w*, Ale bl/wake in the call to phytrac and traclmdz

(to be used for dust emmission)

File size: 28.7 KB
Line 
1!$Id $
2!$Id $
3MODULE phytrac_mod
4!=================================================================================
5! Interface between the LMDZ physical package and tracer computation.
6! Chemistry modules (INCA, Reprobus or the more specific traclmdz routine)
7! are called from phytrac.
8!
9!======================================================================
10! Auteur(s) FH
11! Objet: Moniteur general des tendances traceurs
12!
13! iflag_vdf_trac : Options for activating transport by vertical diffusion :
14!     1. notmal
15!     0. emission is injected in the first layer only, without diffusion
16!    -1  no emission & no diffusion
17! Modification 2013/07/22 : transformed into a module to pass tendencies to
18!     physics outputs. Additional keys for controling activation of sub processes.
19! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
20! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
21!=================================================================================
22
23!
24! Tracer tendencies, for outputs
25!-------------------------------
26  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_cl  ! Td couche limite/traceur
27  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_dec                            !RomP
28  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_cv  ! Td convection/traceur
29! RomP >>>
30  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_insc
31  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_bcscav
32  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_evapls
33  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_ls
34  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_trsp
35  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_sscav
36  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_sat
37  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_uscav
38  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qPr,qDi ! concentration tra dans pluie,air descente insaturee
39  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qPa,qMel
40  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qTrdi,dtrcvMA ! conc traceur descente air insaturee et td convective MA
41! RomP <<<
42  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_th  ! Td thermique
43  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_lessi_impa ! Td du lessivage par impaction
44  REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_lessi_nucl ! Td du lessivage par nucleation
45  REAL,DIMENSION(:,:),ALLOCATABLE,SAVE      :: qPrls      !jyg: concentration tra dans pluie LS a la surf.
46  REAL,DIMENSION(:,:),ALLOCATABLE,SAVE      :: d_tr_dry ! Td depot sec/traceur (1st layer),ALLOCATABLE,SAVE  jyg
47  REAL,DIMENSION(:,:),ALLOCATABLE,SAVE      :: flux_tr_dry ! depot sec/traceur (surface),ALLOCATABLE,SAVE    jyg
48
49!$OMP THREADPRIVATE(qPa,qMel,qTrdi,dtrcvMA,d_tr_th,d_tr_lessi_impa,d_tr_lessi_nucl)
50!$OMP THREADPRIVATE(d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qPr,qDi)
51!$OMP THREADPRIVATE(d_tr_insc,d_tr_bcscav,d_tr_evapls,d_tr_ls,qPrls)
52!$OMP THREADPRIVATE(d_tr,d_tr_cl,d_tr_dry,flux_tr_dry,d_tr_dec,d_tr_cv)
53
54
55CONTAINS
56
57SUBROUTINE phytrac(                            &
58     nstep,     julien,   gmtime,   debutphy,  &
59     lafin,     pdtphys,  u, v,     t_seri,    &
60     paprs,     pplay,    pmfu,     pmfd,      &
61     pen_u,     pde_u,    pen_d,    pde_d,     &
62     cdragh,    coefh,    fm_therm, entr_therm,&
63     yu1,       yv1,      ftsol,    pctsrf,    &
64     ustar,     u10m,      v10m,               &
65     wstar,     ale_bl,      ale_wake,         &
66     xlat,      xlon,                          &
67     frac_impa,frac_nucl,beta_fisrt,beta_v1,   &
68     presnivs,  pphis,    pphi,     albsol,    &
69     sh,        rh,       cldfra,   rneb,      &
70     diafra,    cldliq,   itop_con, ibas_con,  &
71     pmflxr,    pmflxs,   prfl,     psfl,      &
72     da,        phi,      mp,       upwd,      &
73     phi2,      d1a,      dam,      sij,       &   ! RomP
74     wdtrainA,  wdtrainM, sigd,     clw,elij,  &   ! RomP
75     evap,      ep,       epmlmMm,  eplaMm,    &   ! RomP
76     dnwd,      aerosol_couple,     flxmass_w, &
77     tau_aero,  piz_aero,  cg_aero, ccm,       &
78     rfname,                                   &
79     d_tr_dyn,                                 &   ! RomP
80     tr_seri)         
81!
82!======================================================================
83! Auteur(s) FH
84! Objet: Moniteur general des tendances traceurs
85! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
86! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
87!======================================================================
88
89  USE ioipsl
90  USE phys_cal_mod, only : hour
91  USE dimphy
92  USE infotrac
93  USE mod_grid_phy_lmdz
94  USE mod_phys_lmdz_para
95  USE comgeomphy
96  USE iophy
97  USE traclmdz_mod
98  USE tracinca_mod
99  USE tracreprobus_mod
100  USE control_mod
101  USE indice_sol_mod
102
103  IMPLICIT NONE
104
105  INCLUDE "YOMCST.h"
106  INCLUDE "dimensions.h"
107  INCLUDE "clesphys.h"
108  INCLUDE "temps.h"
109  INCLUDE "paramet.h"
110  INCLUDE "thermcell.h"
111  INCLUDE "iniprint.h"
112!==========================================================================
113!                   -- ARGUMENT DESCRIPTION --
114!==========================================================================
115
116! Input arguments
117!----------------
118!Configuration grille,temps:
119  INTEGER,INTENT(IN) :: nstep      ! Appel physique
120  INTEGER,INTENT(IN) :: julien     ! Jour julien
121  REAL,INTENT(IN)    :: gmtime     ! Heure courante
122  REAL,INTENT(IN)    :: pdtphys    ! Pas d'integration pour la physique (seconde)
123  LOGICAL,INTENT(IN) :: debutphy   ! le flag de l'initialisation de la physique
124  LOGICAL,INTENT(IN) :: lafin      ! le flag de la fin de la physique
125 
126  REAL,DIMENSION(klon),INTENT(IN) :: xlat    ! latitudes pour chaque point
127  REAL,DIMENSION(klon),INTENT(IN) :: xlon    ! longitudes pour chaque point
128!
129!Physique:
130!--------
131  REAL,DIMENSION(klon,klev),INTENT(IN)   :: t_seri  ! Temperature
132  REAL,DIMENSION(klon,klev),INTENT(IN)   :: u       ! variable not used
133  REAL,DIMENSION(klon,klev),INTENT(IN)   :: v       ! variable not used
134  REAL,DIMENSION(klon,klev),INTENT(IN)   :: sh      ! humidite specifique
135  REAL,DIMENSION(klon,klev),INTENT(IN)   :: rh      ! humidite relative
136  REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs   ! pression pour chaque inter-couche (en Pa)
137  REAL,DIMENSION(klon,klev),INTENT(IN)   :: pplay   ! pression pour le mileu de chaque couche (en Pa)
138  REAL,DIMENSION(klon,klev),INTENT(IN)   :: pphi    ! geopotentiel
139  REAL,DIMENSION(klon),INTENT(IN)        :: pphis
140  REAL,DIMENSION(klev),INTENT(IN)        :: presnivs
141  REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldliq  ! eau liquide nuageuse
142  REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldfra  ! fraction nuageuse (tous les nuages)
143  REAL,DIMENSION(klon,klev),INTENT(IN)   :: diafra  ! fraction nuageuse (convection ou stratus artificiels)
144  REAL,DIMENSION(klon,klev),INTENT(IN)   :: rneb    ! fraction nuageuse (grande echelle)
145!
146  REAL                                   :: ql_incl ! contenu en eau liquide nuageuse dans le nuage ! ql_incl=oliq/rneb
147  REAL,DIMENSION(klon,klev),INTENT(IN)   :: beta_fisrt ! taux de conversion de l'eau cond (de fisrtilp)
148  REAL,DIMENSION(klon,klev),INTENT(out)  :: beta_v1    ! -- (originale version)
149
150!
151  INTEGER,DIMENSION(klon),INTENT(IN)     :: itop_con
152  INTEGER,DIMENSION(klon),INTENT(IN)     :: ibas_con
153  REAL,DIMENSION(klon),INTENT(IN)        :: albsol  ! albedo surface
154!
155!Dynamique
156!--------
157  REAL,DIMENSION(klon,klev,nbtr),INTENT(IN)    :: d_tr_dyn
158!
159!Convection:
160!----------
161  REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu  ! flux de masse dans le panache montant
162  REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd  ! flux de masse dans le panache descendant
163  REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
164  REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
165  REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
166  REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
167
168!...Tiedke     
169  REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
170  REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
171
172  LOGICAL,INTENT(IN)                       :: aerosol_couple
173  REAL,DIMENSION(klon,klev),INTENT(IN)     :: flxmass_w
174  REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: tau_aero
175  REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: piz_aero
176  REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: cg_aero
177  CHARACTER(len=4),DIMENSION(9),INTENT(IN) :: rfname
178  REAL,DIMENSION(klon,klev,2),INTENT(IN)   :: ccm
179!... K.Emanuel
180  REAL,DIMENSION(klon,klev),INTENT(IN)     :: da
181  REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
182! RomP >>>
183  REAL,DIMENSION(klon,klev),INTENT(IN)      :: d1a,dam
184  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2
185!
186  REAL,DIMENSION(klon,klev),INTENT(IN)      :: wdtrainA
187  REAL,DIMENSION(klon,klev),INTENT(IN)      :: wdtrainM
188  REAL,DIMENSION(klon),INTENT(IN)           :: sigd
189! ---- RomP flux entraine, detraine et precipitant kerry Emanuel
190  REAL,DIMENSION(klon,klev),INTENT(IN)      :: evap
191  REAL,DIMENSION(klon,klev),INTENT(IN)      :: ep
192  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij
193  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij
194  REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm
195  REAL,DIMENSION(klon,klev),INTENT(IN)      :: eplaMm
196  REAL,DIMENSION(klon,klev),INTENT(IN)      :: clw
197! RomP <<<
198
199!
200  REAL,DIMENSION(klon,klev),INTENT(IN)     :: mp
201  REAL,DIMENSION(klon,klev),INTENT(IN)     :: upwd      ! saturated updraft mass flux
202  REAL,DIMENSION(klon,klev),INTENT(IN)     :: dnwd      ! saturated downdraft mass flux
203!
204!Thermiques:
205!----------
206  REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: fm_therm
207  REAL,DIMENSION(klon,klev),INTENT(IN)     :: entr_therm
208!
209!Couche limite:
210!--------------
211!
212!
213  REAL,DIMENSION(:),INTENT(IN)   :: cdragh          ! (klon) coeff drag pour T et Q
214  REAL,DIMENSION(:,:),INTENT(IN) :: coefh           ! (klon,klev) coeff melange CL (m**2/s)
215  REAL,DIMENSION(:),INTENT(IN)   :: ustar,u10m,v10m ! (klon) u* & vent a 10m (m/s)
216  REAL,DIMENSION(:),INTENT(IN)   :: wstar,ale_bl,ale_wake ! (klon) w* and Avail. Lifting Ener.
217  REAL,DIMENSION(:),INTENT(IN)   :: yu1             ! (klon) vents au premier niveau
218  REAL,DIMENSION(:),INTENT(IN)   :: yv1             ! (klon) vents au premier niveau
219
220!
221!Lessivage:
222!----------
223!
224! pour le ON-LINE
225!
226  REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
227  REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
228
229! Arguments necessaires pour les sources et puits de traceur:
230  REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol  ! Temperature du sol (surf)(Kelvin)
231  REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
232
233
234! Output argument
235!----------------
236  REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
237  REAL,DIMENSION(klon,klev)                    :: sourceBE
238!=======================================================================================
239!                        -- LOCAL VARIABLES --
240!=======================================================================================
241
242  INTEGER :: i, k, it
243  INTEGER :: nsplit
244
245!Sources et Reservoirs de traceurs (ex:Radon):
246!--------------------------------------------
247!
248  REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source  ! a voir lorsque le flux de surface est prescrit
249!$OMP THREADPRIVATE(source)
250
251!
252!Entrees/Sorties: (cf ini_histrac.h et write_histrac.h) 
253!---------------
254  INTEGER                   :: iiq, ierr
255  INTEGER                   :: nhori, nvert
256  REAL                      :: zsto, zout, zjulian
257  INTEGER,SAVE              :: nid_tra     ! pointe vers le fichier histrac.nc         
258!$OMP THREADPRIVATE(nid_tra)
259  REAL,DIMENSION(klon)      :: zx_tmp_fi2d ! variable temporaire grille physique
260  INTEGER                   :: itau_w      ! pas de temps ecriture = nstep + itau_phy
261  LOGICAL,PARAMETER         :: ok_sync=.TRUE.
262!$OMP THREADPRIVATE(chtratimestep)
263!
264! Nature du traceur
265!------------------
266  LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol  ! aerosol(it) = true  => aerosol => lessivage
267!$OMP THREADPRIVATE(aerosol)                        ! aerosol(it) = false => gaz
268  REAL,DIMENSION(klon,klev)             :: delp     ! epaisseur de couche (Pa)
269!
270! Tendances de traceurs (Td) et flux de traceurs:
271!------------------------
272  REAL,DIMENSION(klon,klev)      :: d_tr     ! Td dans l'atmosphere
273  REAL,DIMENSION(klon,klev)      :: Mint
274  REAL,DIMENSION(klon,klev,nbtr) :: zmfd1a
275  REAL,DIMENSION(klon,klev,nbtr) :: zmfdam
276  REAL,DIMENSION(klon,klev,nbtr) :: zmfphi2
277! Physique
278!----------
279  REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
280  REAL,DIMENSION(klon,klev)      :: zmasse    ! densité atmosphérique Kg/m2
281  REAL,DIMENSION(klon,klev)      :: ztra_th
282!PhH
283  REAL,DIMENSION(klon,klev)      :: zrho
284  REAL,DIMENSION(klon,klev)      :: zdz
285  REAL                           :: evaplsc,dx,beta ! variable pour lessivage Genthon
286  REAL,DIMENSION(klon)           :: his_dh          ! ---
287! in-cloud scav variables
288  REAL           :: ql_incloud_ref     ! ref value of in-cloud condensed water content
289 
290!Controles:
291!---------
292  INTEGER,SAVE :: iflag_vdf_trac,iflag_con_trac,iflag_the_trac
293  INTEGER :: iflag_con_trac_omp, iflag_vdf_trac_omp,iflag_the_trac_omp
294!$OMP THREADPRIVATE(iflag_lscav,convscav,iflag_the_trac)
295
296  LOGICAL,SAVE :: lessivage
297!$OMP THREADPRIVATE(lessivage)
298
299  CHARACTER(len=8),DIMENSION(nbtr) :: solsym
300!RomP >>>
301  INTEGER,SAVE  :: iflag_lscav_omp,iflag_lscav
302  LOGICAL,SAVE  :: convscav_omp,convscav
303!$OMP THREADPRIVATE(iflag_lscav,convscav)
304!RomP <<<
305!######################################################################
306!                    -- INITIALIZATION --
307!######################################################################
308IF (debutphy) THEN
309ALLOCATE(d_tr_cl(klon,klev,nbtr),d_tr_dry(klon,nbtr))
310ALLOCATE(flux_tr_dry(klon,nbtr),d_tr_dec(klon,klev,nbtr),d_tr_cv(klon,klev,nbtr))
311ALLOCATE(d_tr_insc(klon,klev,nbtr),d_tr_bcscav(klon,klev,nbtr))
312ALLOCATE(d_tr_evapls(klon,klev,nbtr),d_tr_ls(klon,klev,nbtr))
313ALLOCATE(qPrls(klon,nbtr),d_tr_trsp(klon,klev,nbtr))
314ALLOCATE(d_tr_sscav(klon,klev,nbtr),d_tr_sat(klon,klev,nbtr))
315ALLOCATE(d_tr_uscav(klon,klev,nbtr),qPr(klon,klev,nbtr),qDi(klon,klev,nbtr))
316ALLOCATE(qPa(klon,klev,nbtr),qMel(klon,klev,nbtr))
317ALLOCATE(qTrdi(klon,klev,nbtr),dtrcvMA(klon,klev,nbtr))
318ALLOCATE(d_tr_th(klon,klev,nbtr))
319ALLOCATE(d_tr_lessi_impa(klon,klev,nbtr),d_tr_lessi_nucl(klon,klev,nbtr))
320ENDIF
321
322  DO k=1,klev
323     DO i=1,klon
324      sourceBE(i,k)=0.
325      Mint(i,k)=0.
326      zrho(i,k)=0.
327      zdz(i,k)=0.
328     END DO
329  END DO
330
331  DO it=1, nbtr
332   DO k=1,klev
333    DO i=1,klon
334    d_tr_insc(i,k,it)=0.
335    d_tr_bcscav(i,k,it)=0.
336    d_tr_evapls(i,k,it)=0.
337    d_tr_ls(i,k,it)=0.
338    d_tr_cv(i,k,it)=0.
339    d_tr_cl(i,k,it)=0.
340    d_tr_trsp(i,k,it)=0.
341    d_tr_sscav(i,k,it)=0.
342    d_tr_sat(i,k,it)=0.
343    d_tr_uscav(i,k,it)=0.
344    d_tr_lessi_impa(i,k,it)=0.
345    d_tr_lessi_nucl(i,k,it)=0.
346    qDi(i,k,it)=0.
347    qPr(i,k,it)=0.
348    qPa(i,k,it)=0.
349    qMel(i,k,it)=0.
350    qTrdi(i,k,it)=0.
351    dtrcvMA(i,k,it)=0.
352    zmfd1a(i,k,it)=0.
353    zmfdam(i,k,it)=0.
354    zmfphi2(i,k,it)=0.
355    END DO
356   END DO
357  END DO
358  IF (debutphy) THEN
359!!jyg
360!$OMP BARRIER
361   ecrit_tra=86400. ! frequence de stokage en dur
362                    ! obsolete car remplace par des ecritures dans phys_output_write
363!RomP >>>
364!
365!Config Key  = convscav
366!Config Desc = Convective scavenging switch: 0=off, 1=on.
367!Config Def  = .false.
368!Config Help =
369!
370!$OMP MASTER
371  convscav_omp=.false.
372  call getin('convscav', convscav_omp)
373  iflag_vdf_trac_omp=1
374  call getin('iflag_vdf_trac', iflag_vdf_trac_omp)
375  iflag_con_trac_omp=1
376  call getin('iflag_con_trac', iflag_con_trac_omp)
377  iflag_the_trac_omp=1
378  call getin('iflag_the_trac', iflag_the_trac_omp)
379!$OMP END MASTER
380!$OMP BARRIER
381  convscav=convscav_omp
382  iflag_vdf_trac=iflag_vdf_trac_omp
383  iflag_con_trac=iflag_con_trac_omp
384  iflag_the_trac=iflag_the_trac_omp
385  print*,'phytrac passage dans routine conv avec lessivage', convscav
386!
387!Config Key  = iflag_lscav
388!Config Desc = Large scale scavenging parametrization: 0=none, 1=old(Genthon92),
389!              2=1+PHeinrich, 3=Reddy_Boucher2004, 4=3+RPilon.
390!Config Def  = 1
391!Config Help =
392!
393!$OMP MASTER
394  iflag_lscav_omp=1
395  call getin('iflag_lscav', iflag_lscav_omp)
396!$OMP END MASTER
397!$OMP BARRIER
398  iflag_lscav=iflag_lscav_omp
399!
400  SELECT CASE(iflag_lscav)
401  CASE(0)
402   PRINT*, 'Large scale scavenging: none'
403  CASE(1)
404   PRINT*, 'Large scale scavenging: C. Genthon, Tellus(1992), 44B, 371-389'
405  CASE(2)
406   PRINT*, 'Large scale scavenging: C. Genthon, modified P. Heinrich'
407  CASE(3)
408   PRINT*, 'Large scale scavenging: M. Shekkar Reddy and O. Boucher, JGR(2004), 109, D14202'
409  CASE(4)
410   PRINT*, 'Large scale scavenging: Reddy and Boucher, modified R. Pilon'
411  END SELECT
412!RomP <<<
413     WRITE(*,*) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
414     ALLOCATE( source(klon,nbtr), stat=ierr)
415     IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 1',1)
416     
417     ALLOCATE( aerosol(nbtr), stat=ierr)
418     IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 2',1)
419     
420
421     ! Initialize module for specific tracers
422     SELECT CASE(type_trac)
423     CASE('lmdz')
424        CALL traclmdz_init(pctsrf, xlat, xlon, ftsol, tr_seri, t_seri, pplay, sh, pdtphys, aerosol, lessivage)
425     CASE('inca')
426        source(:,:)=0.
427        CALL tracinca_init(aerosol,lessivage)
428     CASE('repr')
429        source(:,:)=0.
430     END SELECT
431!
432! Initialize diagnostic output
433! ----------------------------
434#ifdef CPP_IOIPSL
435!     INCLUDE "ini_histrac.h"
436#endif
437  END IF ! of IF (debutphy)
438!############################################ END INITIALIZATION #######
439
440  DO k=1,klev
441     DO i=1,klon
442        zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
443     END DO
444  END DO
445!
446  IF (id_be .GT. 0) THEN
447  DO k=1,klev
448     DO i=1,klon
449        sourceBE(i,k)=srcbe(i,k)       !RomP  -> pour sortie histrac
450     END DO
451  END DO
452  ENDIF
453
454!===============================================================================
455!    -- Do specific treatment according to chemestry model or local LMDZ tracers
456!     
457!===============================================================================
458  SELECT CASE(type_trac)
459  CASE('lmdz')
460     !    -- Traitement des traceurs avec traclmdz
461     CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
462          cdragh,  coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,iflag_vdf_trac>=0,sh, &
463           rh, pphi, ustar, wstar, ale_bl, ale_wake,  u10m, v10m, &
464           tr_seri, source, solsym, d_tr_cl,d_tr_dec, zmasse)               !RomP
465
466  CASE('inca')
467     !    -- CHIMIE INCA  config_inca = aero or chem --
468
469     CALL tracinca(&
470          nstep,    julien,   gmtime,         lafin,     &
471          pdtphys,  t_seri,   paprs,          pplay,     &
472          pmfu,     ftsol,    pctsrf,         pphis,     &
473          pphi,     albsol,   sh,             rh,        &
474          cldfra,   rneb,     diafra,         cldliq,    &
475          itop_con, ibas_con, pmflxr,         pmflxs,    &
476          prfl,     psfl,     aerosol_couple, flxmass_w, &
477          tau_aero, piz_aero, cg_aero,        ccm,       &
478          rfname,                                        &
479          tr_seri,  source,   solsym)     
480
481  CASE('repr')
482     !   -- CHIMIE REPROBUS --
483
484     CALL tracreprobus(pdtphys, gmtime, debutphy, julien, &
485          presnivs, xlat, xlon, pphis, pphi, &
486          t_seri, pplay, paprs, sh , &
487          tr_seri, solsym)
488     
489  END SELECT
490!======================================================================
491!       -- Calcul de l'effet de la convection --
492!======================================================================
493
494  IF (iflag_con_trac==1) THEN
495     DO it=1, nbtr
496        IF ( conv_flg(it) == 0 ) CYCLE
497        IF (iflag_con.LT.2) THEN
498           d_tr_cv(:,:,it)=0.
499        ELSE IF (iflag_con.EQ.2) THEN
500!..Tiedke
501           CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
502                pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
503! RomP >>>               
504        ELSE   
505!..K.Emanuel                  !RomP modif arg
506        if (convscav.and.aerosol(it)) then    ! lessivage convectif pour aerosol
507!
508          CALL cvltr(pdtphys, da, phi,phi2,d1a,dam, mp,ep,              &
509               sigd,sij,clw,elij,epmlmMm,eplaMm,                        &
510               pmflxr,pmflxs,evap,t_seri,wdtrainA,wdtrainM,             &   
511               paprs,it,tr_seri,upwd,dnwd,itop_con,ibas_con,            &
512               d_tr_cv,d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qDi,qPr,&
513               qPa,qMel,qTrdi,dtrcvMA,Mint,                             &
514               zmfd1a,zmfphi2,zmfdam)
515        else  !pas de lessivage convectif ou n'est pas un aerosol
516           CALL cvltrorig(it,pdtphys, da, phi,mp,paprs,pplay,tr_seri,&
517                    upwd,dnwd,d_tr_cv)
518        endif
519        END IF
520! RomP <<<
521
522        DO k = 1, klev
523           DO i = 1, klon       
524              tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
525           END DO
526        END DO
527
528        CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
529             
530     END DO ! nbtr
531  END IF ! convection
532
533!======================================================================
534!    -- Calcul de l'effet des thermiques --
535!======================================================================
536
537  DO it=1,nbtr
538     DO k=1,klev
539        DO i=1,klon
540           d_tr_th(i,k,it)=0.
541           tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
542           tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
543        END DO
544     END DO
545  END DO
546 
547  IF (iflag_thermals.GT.0.AND.iflag_the_trac>0) THEN   
548     nsplit=10
549     DO it=1, nbtr
550        DO isplit=1,nsplit
551
552           CALL dqthermcell(klon,klev,pdtphys/nsplit, &
553                fm_therm,entr_therm,zmasse, &
554                tr_seri(1:klon,1:klev,it),d_tr,ztra_th)
555
556           DO k=1,klev
557              DO i=1,klon
558                 d_tr(i,k)=pdtphys*d_tr(i,k)/nsplit
559                 d_tr_th(i,k,it)=d_tr_th(i,k,it)+d_tr(i,k)
560                 tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr(i,k),0.)
561              END DO
562           END DO
563        END DO ! nsplit
564     END DO ! it
565  END IF ! Thermiques
566
567!======================================================================
568!     -- Calcul de l'effet de la couche limite --
569!======================================================================
570
571  DO k = 1, klev
572     DO i = 1, klon
573        delp(i,k) = paprs(i,k)-paprs(i,k+1)
574     END DO
575  END DO
576
577  IF (iflag_vdf_trac==1) THEN
578     DO it=1, nbtr
579        print*,'trac pbl ',it,pbl_flg(it)
580        IF( pbl_flg(it) /= 0 ) THEN
581           CALL cltrac(pdtphys, coefh,t_seri,       &
582                tr_seri(:,:,it), source(:,it),      &
583                paprs, pplay, delp,                 &
584                d_tr_cl(:,:,it),d_tr_dry(:,it),flux_tr_dry(:,it))
585           tr_seri(:,:,it)=tr_seri(:,:,it)+d_tr_cl(:,:,it)
586        END IF
587     END DO
588  ELSE IF (iflag_vdf_trac==0) THEN
589!   Injection of source in the first model layer
590    DO it=1,nbtr
591       d_tr_cl(:,1,it)=source(:,it)*rg/delp(:,1)*pdtphys
592    ENDDO
593    tr_seri(:,1,it)=tr_seri(:,1,it)+d_tr_cl(:,1,it)
594    d_tr_cl(:,2:klev,1:nbtr)=0.
595  ELSE IF (iflag_vdf_trac==-1) THEN
596    d_tr_cl=0.
597  ELSE
598    CALL abort_gcm('iflag_vdf_trac', 'cas non prevu',1)
599  END IF ! couche limite
600
601
602
603!======================================================================
604!   Calcul de l'effet de la precipitation grande echelle
605!======================================================================
606  IF (lessivage) THEN
607
608   ql_incloud_ref = 10.e-4
609   ql_incloud_ref =  5.e-4
610
611
612! calcul du contenu en eau liquide au sein du nuage
613     ql_incl = ql_incloud_ref
614! choix du lessivage
615!
616  IF (iflag_lscav .EQ. 3 .OR. iflag_lscav .EQ. 4) THEN
617! ********  Olivier Boucher version (3) possibly with modified ql_incl (4)
618!
619   DO it = 1, nbtr
620!  incloud scavenging and removal by large scale rain ! orig : ql_incl was replaced by 0.5e-3 kg/kg
621! the value 0.5e-3 kg/kg is from Giorgi and Chameides (1986), JGR
622! Liu (2001) proposed to use 1.5e-3 kg/kg
623
624    CALL lsc_scav(pdtphys,it,iflag_lscav,ql_incl,prfl,psfl,rneb,beta_fisrt,  &
625                  beta_v1,pplay,paprs,t_seri,tr_seri,d_tr_insc,   &
626                  d_tr_bcscav,d_tr_evapls,qPrls)
627
628!large scale scavenging tendency
629   DO k = 1, klev
630    DO i = 1, klon
631    d_tr_ls(i,k,it)=d_tr_insc(i,k,it)+d_tr_bcscav(i,k,it)+d_tr_evapls(i,k,it)
632    tr_seri(i,k,it)=tr_seri(i,k,it)+d_tr_ls(i,k,it)
633    ENDDO
634   ENDDO
635     CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'lsc scav it = '//solsym(it))
636   END DO  !tr
637
638 ELSE IF (iflag_lscav .EQ. 2) THEN ! frac_impa, frac_nucl
639! *********   modified  old version
640
641     d_tr_lessi_nucl(:,:,:) = 0.
642     d_tr_lessi_impa(:,:,:) = 0.
643     flestottr(:,:,:) = 0.
644! Tendance des aerosols nuclees et impactes
645     DO it = 1, nbtr
646        IF (aerosol(it)) THEN
647        his_dh(:)=0.
648           DO k = 1, klev
649              DO i = 1, klon
650!PhH
651              zrho(i,k)=pplay(i,k)/t_seri(i,k)/RD
652              zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
653!
654              END DO
655           END DO
656
657          DO k=klev-1, 1, -1
658            DO i=1, klon
659!             d_tr_ls(i,k,it)=tr_seri(i,k,it)*(frac_impa(i,k)*frac_nucl(i,k)-1.)
660             dx=d_tr_ls(i,k,it)
661             his_dh(i)=his_dh(i)-dx*zrho(i,k)*zdz(i,k)/pdtphys !  kg/m2/s
662             evaplsc = prfl(i,k) - prfl(i,k+1) + psfl(i,k) - psfl(i,k+1)
663! Evaporation Partielle -> Liberation Partielle 0.5*evap
664            IF ( evaplsc .LT.0..and.abs(prfl(i,k+1)+psfl(i,k+1)).gt.1.e-10) THEN
665                evaplsc = (-evaplsc)/(prfl(i,k+1)+psfl(i,k+1))
666! evaplsc est donc positif, his_dh(i) est positif
667!-------------- 
668             d_tr_evapls(i,k,it)=0.5*evaplsc*(d_tr_lessi_nucl(i,k+1,it) &
669                                  +d_tr_lessi_impa(i,k+1,it))
670!-------------   d_tr_evapls(i,k,it)=-0.5*evaplsc*(d_tr_lsc(i,k+1,it))
671             beta=0.5*evaplsc
672              if ((prfl(i,k)+psfl(i,k)).lt.1.e-10) THEN
673               beta=1.0*evaplsc
674              endif
675            dx=beta*his_dh(i)/zrho(i,k)/zdz(i,k)*pdtphys
676            his_dh(i)=(1.-beta)*his_dh(i)   ! tracer from
677            d_tr_evapls(i,k,it)=dx
678            ENDIF
679            d_tr_ls(i,k,it)=tr_seri(i,k,it)*(frac_impa(i,k)*frac_nucl(i,k)-1.) &
680                            +d_tr_evapls(i,k,it)
681
682!--------------
683                 d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) +    &
684                      ( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
685                 d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) +    &
686                      ( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
687!
688! Flux lessivage total
689                 flestottr(i,k,it) = flestottr(i,k,it) -   &
690                      ( d_tr_lessi_nucl(i,k,it)   +        &
691                      d_tr_lessi_impa(i,k,it) ) *          &
692                      ( paprs(i,k)-paprs(i,k+1) ) /        &
693                      (RG * pdtphys)
694!! Mise a jour des traceurs due a l'impaction,nucleation
695!                 tr_seri(i,k,it)=tr_seri(i,k,it)*frac_impa(i,k)*frac_nucl(i,k)
696!!  calcul de la tendance liee au lessivage stratiforme
697!                 d_tr_ls(i,k,it)=tr_seri(i,k,it)*&
698!                                (1.-1./(frac_impa(i,k)*frac_nucl(i,k)))
699!--------------
700              END DO
701           END DO
702        END IF
703     END DO
704! *********   end modified old version
705
706 ELSE IF (iflag_lscav .EQ. 1) THEN ! frac_impa, frac_nucl
707! *********    old version
708
709     d_tr_lessi_nucl(:,:,:) = 0.
710     d_tr_lessi_impa(:,:,:) = 0.
711     flestottr(:,:,:) = 0.
712!=========================
713! LESSIVAGE LARGE SCALE :
714!=========================
715
716! Tendance des aerosols nuclees et impactes
717! -----------------------------------------
718     DO it = 1, nbtr
719        IF (aerosol(it)) THEN
720           DO k = 1, klev
721              DO i = 1, klon
722                 d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) +    &
723                      ( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
724                 d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) +    &
725                      ( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
726
727!
728! Flux lessivage total
729! ------------------------------------------------------------
730                 flestottr(i,k,it) = flestottr(i,k,it) -   &
731                      ( d_tr_lessi_nucl(i,k,it)   +        &
732                      d_tr_lessi_impa(i,k,it) ) *          &
733                      ( paprs(i,k)-paprs(i,k+1) ) /        &
734                      (RG * pdtphys)
735!
736! Mise a jour des traceurs due a l'impaction,nucleation
737! ----------------------------------------------------------------------
738                 tr_seri(i,k,it)=tr_seri(i,k,it)*frac_impa(i,k)*frac_nucl(i,k)
739              END DO
740           END DO
741        END IF
742     END DO
743     
744! *********   end old version
745  ENDIF  !  iflag_lscav . EQ. 1, 2, 3 or 4
746!
747  END IF !  lessivage
748
749!=============================================================
750!   Ecriture des sorties
751!=============================================================
752#ifdef CPP_IOIPSL
753! INCLUDE "write_histrac.h"
754#endif
755
756END SUBROUTINE phytrac
757
758END MODULE
Note: See TracBrowser for help on using the repository browser.