Context Navigation

phytrac_mod.f90 @ 5447

Last change on this file since 5447 was 5447, checked in by jyg, 11 hours ago
output the convective wet deposit of tracers
Property copyright set to Name of program: LMDZ Creation date: 1984 Version: LMDZ5 License: CeCILL version 2 Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539 See the license file in the root directory Property svn:keywords set to `Id`
File size: 45.5 KB

Line
1	!$Id: phytrac_mod.f90 5447 2024-12-22 16:48:23Z jyg $
2	MODULE phytrac_mod
3	!=================================================================================
4	! Interface between the LMDZ physical package and tracer computation.
5	! Chemistry modules (INCA, Reprobus or the more specific traclmdz routine)
6	! are called from phytrac.
7	!
8	!======================================================================
9	! Auteur(s) FH
10	! Objet: Moniteur general des tendances traceurs
11	!
12	! iflag_vdf_trac : Options for activating transport by vertical diffusion :
13	! 1. notmal
14	! 0. emission is injected in the first layer only, without diffusion
15	! -1 no emission & no diffusion
16	! Modification 2013/07/22 : transformed into a module to pass tendencies to
17	! physics outputs. Additional keys for controling activation of sub processes.
18	! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
19	! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
20	!=================================================================================
21
22	!
23	! Tracer tendencies, for outputs
24	!-------------------------------
25	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_cl ! Td couche limite/traceur
26	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_dec !RomP
27	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_cv ! Td convection/traceur
28	! RomP >>>
29	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_insc
30	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_bcscav
31	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_evapls
32	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_ls
33	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_trsp
34	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_sscav
35	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_sat
36	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_uscav
37	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: flux_tr_wet ! tracer wet deposit (surface) jyg
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,flux_tr_wet,qPr,qDi)
51	!$OMP THREADPRIVATE(d_tr_insc,d_tr_bcscav,d_tr_evapls,d_tr_ls,qPrls)
52	!$OMP THREADPRIVATE(d_tr_cl,d_tr_dry,flux_tr_dry,d_tr_dec,d_tr_cv)
53
54	CONTAINS
55
56	SUBROUTINE phytrac_init()
57
58	USE dimphy
59	USE infotrac_phy, ONLY: nbtr, type_trac
60	USE tracco2i_mod, ONLY: tracco2i_init
61
62	IMPLICIT NONE
63
64	ALLOCATE(d_tr_cl(klon,klev,nbtr),d_tr_dry(klon,nbtr))
65	ALLOCATE(flux_tr_dry(klon,nbtr),d_tr_dec(klon,klev,nbtr),d_tr_cv(klon,klev,nbtr))
66	ALLOCATE(d_tr_insc(klon,klev,nbtr),d_tr_bcscav(klon,klev,nbtr))
67	ALLOCATE(d_tr_evapls(klon,klev,nbtr),d_tr_ls(klon,klev,nbtr))
68	ALLOCATE(qPrls(klon,nbtr),d_tr_trsp(klon,klev,nbtr))
69	ALLOCATE(d_tr_sscav(klon,klev,nbtr),d_tr_sat(klon,klev,nbtr))
70	ALLOCATE(d_tr_uscav(klon,klev,nbtr),qPr(klon,klev,nbtr),qDi(klon,klev,nbtr))
71	ALLOCATE(flux_tr_wet(klon,nbtr))
72	ALLOCATE(qPa(klon,klev,nbtr),qMel(klon,klev,nbtr))
73	ALLOCATE(qTrdi(klon,klev,nbtr),dtrcvMA(klon,klev,nbtr))
74	ALLOCATE(d_tr_th(klon,klev,nbtr))
75	ALLOCATE(d_tr_lessi_impa(klon,klev,nbtr),d_tr_lessi_nucl(klon,klev,nbtr))
76
77	!===============================================================================
78	! -- Do specific treatment according to chemestry model or local LMDZ tracers
79	!
80	!===============================================================================
81	! -- CO2 interactif --
82	IF (ANY(type_trac == ['co2i','inco'])) CALL tracco2i_init()
83
84	! -- type_trac == 'co2i' ! PC
85	! -- CO2 interactif --
86	! -- source is updated with FF and BB emissions
87	! -- and net fluxes from ocean and orchidee
88	! -- sign convention : positive into the atmosphere
89
90	END SUBROUTINE phytrac_init
91
92	SUBROUTINE phytrac( &
93	nstep, julien, gmtime, debutphy, &
94	lafin, pdtphys, u, v, t_seri, &
95	paprs, pplay, pmfu, pmfd, &
96	pen_u, pde_u, pen_d, pde_d, &
97	cdragh, coefh, fm_therm, entr_therm, &
98	yu1, yv1, ftsol, pctsrf, &
99	ustar, u10m, v10m, &
100	wstar, ale_bl, ale_wake, &
101	xlat, xlon, &
102	frac_impa,frac_nucl,beta_fisrt,beta_v1, &
103	presnivs, pphis, pphi, albsol, &
104	sh, ch, rh, cldfra, rneb, &
105	diafra, cldliq, itop_con, ibas_con, &
106	pmflxr, pmflxs, prfl, psfl, &
107	da, phi, mp, upwd, &
108	phi2, d1a, dam, sij, wght_cvfd, & ! RomP +RL
109	wdtrainA, wdtrainM, sigd, clw, elij, & ! RomP
110	evap, ep, epmlmMm, eplaMm, & ! RomP
111	dnwd, aerosol_couple, flxmass_w, &
112	tau_aero, piz_aero, cg_aero, ccm, &
113	rfname, &
114	d_tr_dyn, & ! RomP
115	tr_seri, init_source)
116	!
117	!======================================================================
118	! Auteur(s) FH
119	! Objet: Moniteur general des tendances traceurs
120	! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
121	! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
122	!======================================================================
123
124	USE ioipsl
125	USE phys_cal_mod, only : hour
126	USE dimphy
127	USE infotrac_phy, ONLY: nbtr, nqCO2, type_trac, conv_flg, pbl_flg
128	USE strings_mod, ONLY: int2str
129	USE mod_grid_phy_lmdz
130	USE mod_phys_lmdz_para
131	USE iophy
132	USE traclmdz_mod
133	USE tracinca_mod
134	USE tracreprobus_mod
135	USE indice_sol_mod
136	USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root
137	USE print_control_mod, ONLY: lunout
138	USE aero_mod, ONLY : naero_grp
139	USE lmdz_thermcell_dq, ONLY : thermcell_dq
140
141	USE tracco2i_mod
142
143	USE traccoag_mod
144	USE phys_local_var_mod, ONLY: mdw
145	USE phys_local_var_mod, ONLY: budg_dep_dry_ocs, budg_dep_wet_ocs
146	USE phys_local_var_mod, ONLY: budg_dep_dry_so2, budg_dep_wet_so2
147	USE phys_local_var_mod, ONLY: budg_dep_dry_h2so4, budg_dep_wet_h2so4
148	USE phys_local_var_mod, ONLY: budg_dep_dry_part, budg_dep_wet_part
149	USE infotrac_phy, ONLY: nbtr_sulgas, id_OCS_strat, id_SO2_strat, id_H2SO4_strat
150	USE strataer_nuc_mod, ONLY : tracstrataer_init
151	USE aerophys
152
153	USE lmdz_cppkeys_wrapper, ONLY: CPPKEY_STRATAER
154
155	USE yomcst_mod_h
156	USE clesphys_mod_h
157	IMPLICIT NONE
158
159
160	!==========================================================================
161	! -- ARGUMENT DESCRIPTION --
162	!==========================================================================
163
164	! Input arguments
165	!----------------
166	!Configuration grille,temps:
167	INTEGER,INTENT(IN) :: nstep ! Appel physique
168	INTEGER,INTENT(IN) :: julien ! Jour julien
169	REAL,INTENT(IN) :: gmtime ! Heure courante
170	REAL,INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
171	LOGICAL,INTENT(IN) :: debutphy ! le flag de l'initialisation de la physique
172	LOGICAL,INTENT(IN) :: lafin ! le flag de la fin de la physique
173
174	REAL,DIMENSION(klon),INTENT(IN) :: xlat ! latitudes pour chaque point
175	REAL,DIMENSION(klon),INTENT(IN) :: xlon ! longitudes pour chaque point
176	!
177	!Physique:
178	!--------
179	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
180	REAL,DIMENSION(klon,klev),INTENT(IN) :: u ! variable not used
181	REAL,DIMENSION(klon,klev),INTENT(IN) :: v ! variable not used
182	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
183	REAL,DIMENSION(klon,klev),INTENT(IN) :: rh ! humidite relative
184	REAL,DIMENSION(klon,klev),INTENT(IN) :: ch ! eau liquide (+ glace si le traceur existe)
185	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
186	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
187	REAL,DIMENSION(klon,klev),INTENT(IN) :: pphi ! geopotentiel
188	REAL,DIMENSION(klon),INTENT(IN) :: pphis
189	REAL,DIMENSION(klev),INTENT(IN) :: presnivs
190	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldliq ! eau condensee totale
191	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldfra ! fraction nuageuse (tous les nuages)
192	REAL,DIMENSION(klon,klev),INTENT(IN) :: diafra ! fraction nuageuse (convection ou stratus artificiels)
193	REAL,DIMENSION(klon,klev),INTENT(IN) :: rneb ! fraction nuageuse (grande echelle)
194	!
195	REAL :: ql_incl ! contenu en eau liquide nuageuse dans le nuage ! ql_incl=oliq/rneb
196	REAL,DIMENSION(klon,klev),INTENT(IN) :: beta_fisrt ! taux de conversion de l'eau cond (de fisrtilp)
197	REAL,DIMENSION(klon,klev),INTENT(out) :: beta_v1 ! -- (originale version)
198
199	!
200	INTEGER,DIMENSION(klon),INTENT(IN) :: itop_con
201	INTEGER,DIMENSION(klon),INTENT(IN) :: ibas_con
202	REAL,DIMENSION(klon),INTENT(IN) :: albsol ! albedo surface
203	!
204	!Dynamique
205	!--------
206	REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: d_tr_dyn
207	!
208	!Convection:
209	!----------
210	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu ! flux de masse dans le panache montant
211	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd ! flux de masse dans le panache descendant
212	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
213	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
214	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
215	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
216
217	!...Tiedke
218	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
219	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
220
221	LOGICAL,INTENT(IN) :: aerosol_couple
222	REAL,DIMENSION(klon,klev),INTENT(IN) :: flxmass_w
223	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: tau_aero
224	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: piz_aero
225	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: cg_aero
226	CHARACTER(len=4),DIMENSION(naero_grp),INTENT(IN) :: rfname
227	REAL,DIMENSION(klon,klev,2),INTENT(IN) :: ccm
228	!... K.Emanuel
229	REAL,DIMENSION(klon,klev),INTENT(IN) :: da
230	REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
231	! RomP >>>
232	REAL,DIMENSION(klon,klev),INTENT(IN) :: d1a,dam
233	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2
234	!
235	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainA
236	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainM
237	REAL,DIMENSION(klon),INTENT(IN) :: sigd
238	! ---- RomP flux entraine, detraine et precipitant kerry Emanuel
239	REAL,DIMENSION(klon,klev),INTENT(IN) :: evap
240	REAL,DIMENSION(klon,klev),INTENT(IN) :: ep
241	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij
242	REAL,DIMENSION(klon,klev),INTENT(IN) :: wght_cvfd !RL
243	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij
244	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm
245	REAL,DIMENSION(klon,klev),INTENT(IN) :: eplaMm
246	REAL,DIMENSION(klon,klev),INTENT(IN) :: clw
247	! RomP <<<
248	!
249	REAL,DIMENSION(klon,klev),INTENT(IN) :: mp
250	REAL,DIMENSION(klon,klev),INTENT(IN) :: upwd ! saturated updraft mass flux
251	REAL,DIMENSION(klon,klev),INTENT(IN) :: dnwd ! saturated downdraft mass flux
252	!
253	!Thermiques:
254	!----------
255	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: fm_therm
256	REAL,DIMENSION(klon,klev),INTENT(INOUT) :: entr_therm
257	!
258	!Couche limite:
259	!--------------
260	!
261	REAL,DIMENSION(:),INTENT(IN) :: cdragh ! (klon) coeff drag pour T et Q
262	REAL,DIMENSION(:,:),INTENT(IN) :: coefh ! (klon,klev) coeff melange CL (m**2/s)
263	REAL,DIMENSION(:),INTENT(IN) :: ustar,u10m,v10m ! (klon) u* & vent a 10m (m/s)
264	REAL,DIMENSION(:),INTENT(IN) :: wstar,ale_bl,ale_wake ! (klon) w* and Avail. Lifting Ener.
265	REAL,DIMENSION(:),INTENT(IN) :: yu1 ! (klon) vents au premier niveau
266	REAL,DIMENSION(:),INTENT(IN) :: yv1 ! (klon) vents au premier niveau
267	!
268	!Lessivage:
269	!----------
270	!
271	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrAA
272	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrENV
273	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: coefcoli
274	LOGICAL, DIMENSION(:), ALLOCATABLE, SAVE :: flag_cvltr
275	!$OMP THREADPRIVATE(ccntrAA,ccntrENV,coefcoli,flag_cvltr)
276	REAL, DIMENSION(klon,klev) :: ccntrAA_3d
277	REAL, DIMENSION(klon,klev) :: ccntrENV_3d
278	REAL, DIMENSION(klon,klev) :: coefcoli_3d
279	!
280	! pour le ON-LINE
281	!
282	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
283	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
284
285	! Arguments necessaires pour les sources et puits de traceur:
286	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
287	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
288
289	REAL,DIMENSION(klon) :: v_dep_dry !dry deposition velocity of stratospheric sulfate in m/s
290	! Output argument
291	!----------------
292	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
293	REAL,DIMENSION(klon,klev) :: sourceBE
294	REAL,DIMENSION(klon,nbtr), INTENT(IN) :: init_source
295
296	!=======================================================================================
297	! -- LOCAL VARIABLES --
298	!=======================================================================================
299
300	INTEGER :: i, k, it
301	INTEGER :: nsplit
302
303	!Sources et Reservoirs de traceurs (ex:Radon):
304	!--------------------------------------------
305	!
306	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source ! a voir lorsque le flux de surface est prescrit
307	!$OMP THREADPRIVATE(source)
308
309	!
310	!Entrees/Sorties:
311	!---------------
312	INTEGER :: iiq, ierr
313	INTEGER :: nhori, nvert
314	REAL :: zsto, zout, zjulian
315	INTEGER,SAVE :: nid_tra ! pointe vers le fichier histrac.nc
316	!$OMP THREADPRIVATE(nid_tra)
317	REAL,DIMENSION(klon) :: zx_tmp_fi2d ! variable temporaire grille physique
318	INTEGER :: itau_w ! pas de temps ecriture = nstep + itau_phy
319	LOGICAL,PARAMETER :: ok_sync=.TRUE.
320	!
321	! Nature du traceur
322	!------------------
323	LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol ! aerosol(it) = true => aerosol => lessivage
324	!$OMP THREADPRIVATE(aerosol) ! aerosol(it) = false => gaz
325	REAL,DIMENSION(klon,klev) :: delp ! epaisseur de couche (Pa)
326	!
327	! Tendances de traceurs (Td) et flux de traceurs:
328	!------------------------
329	REAL,DIMENSION(klon,klev) :: d_tr ! Td dans l'atmosphere
330	REAL,DIMENSION(klon,klev) :: Mint
331	REAL,DIMENSION(klon,klev,nbtr) :: zmfd1a
332	REAL,DIMENSION(klon,klev,nbtr) :: zmfdam
333	REAL,DIMENSION(klon,klev,nbtr) :: zmfphi2
334	! Physique
335	!----------
336	REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
337	REAL,DIMENSION(klon,klev) :: zmasse ! densite atmospherique Kg/m2
338	REAL,DIMENSION(klon,klev) :: ztra_th
339	!PhH
340	REAL,DIMENSION(klon,klev) :: zrho
341	REAL,DIMENSION(klon,klev) :: zdz
342	REAL :: evaplsc,dx,beta ! variable pour lessivage Genthon
343	REAL,DIMENSION(klon) :: his_dh ! ---
344	! in-cloud scav variables
345	REAL :: ql_incloud_ref ! ref value of in-cloud condensed water content
346
347	!Controles:
348	!---------
349	INTEGER,SAVE :: iflag_vdf_trac,iflag_con_trac,iflag_the_trac
350	INTEGER,SAVE :: iflag_con_trac_omp, iflag_vdf_trac_omp,iflag_the_trac_omp
351	!$OMP THREADPRIVATE(iflag_vdf_trac,iflag_con_trac,iflag_the_trac)
352
353	LOGICAL,SAVE :: lessivage
354	!$OMP THREADPRIVATE(lessivage)
355
356	!RomP >>>
357	INTEGER,SAVE :: iflag_lscav_omp,iflag_lscav
358	REAL, SAVE :: ccntrAA_in,ccntrAA_omp
359	REAL, SAVE :: ccntrENV_in,ccntrENV_omp
360	REAL, SAVE :: coefcoli_in,coefcoli_omp
361
362	LOGICAL,SAVE :: convscav_omp,convscav
363	!$OMP THREADPRIVATE(iflag_lscav)
364	!$OMP THREADPRIVATE(ccntrAA_in,ccntrENV_in,coefcoli_in)
365	!$OMP THREADPRIVATE(convscav)
366	!RomP <<<
367	!######################################################################
368	! -- INITIALIZATION --
369	!######################################################################
370
371	DO k=1,klev
372	DO i=1,klon
373	sourceBE(i,k)=0.
374	Mint(i,k)=0.
375	zrho(i,k)=0.
376	zdz(i,k)=0.
377	ENDDO
378	ENDDO
379
380	DO it=1, nbtr
381	DO k=1,klev
382	DO i=1,klon
383	d_tr_insc(i,k,it)=0.
384	d_tr_bcscav(i,k,it)=0.
385	d_tr_evapls(i,k,it)=0.
386	d_tr_ls(i,k,it)=0.
387	d_tr_cv(i,k,it)=0.
388	d_tr_cl(i,k,it)=0.
389	d_tr_trsp(i,k,it)=0.
390	d_tr_sscav(i,k,it)=0.
391	d_tr_sat(i,k,it)=0.
392	d_tr_uscav(i,k,it)=0.
393	d_tr_lessi_impa(i,k,it)=0.
394	d_tr_lessi_nucl(i,k,it)=0.
395	qDi(i,k,it)=0.
396	qPr(i,k,it)=0.
397	qPa(i,k,it)=0.
398	qMel(i,k,it)=0.
399	qTrdi(i,k,it)=0.
400	dtrcvMA(i,k,it)=0.
401	zmfd1a(i,k,it)=0.
402	zmfdam(i,k,it)=0.
403	zmfphi2(i,k,it)=0.
404	ENDDO
405	ENDDO
406	ENDDO
407
408	DO it=1, nbtr
409	DO i=1,klon
410	d_tr_dry(i,it)=0.
411	flux_tr_dry(i,it)=0.
412	flux_tr_wet(i,it)=0.
413	ENDDO
414	ENDDO
415
416	DO k = 1, klev
417	DO i = 1, klon
418	delp(i,k) = paprs(i,k)-paprs(i,k+1)
419	ENDDO
420	ENDDO
421
422	IF (debutphy) THEN
423	!!jyg
424	!$OMP BARRIER
425	ecrit_tra=86400. ! frequence de stokage en dur
426	! obsolete car remplace par des ecritures dans phys_output_write
427	!RomP >>>
428	!
429	!Config Key = convscav
430	!Config Desc = Convective scavenging switch: 0=off, 1=on.
431	!Config Def = .FALSE.
432	!Config Help =
433	!
434	!$OMP MASTER
435	convscav_omp=.FALSE.
436	call getin('convscav', convscav_omp)
437	iflag_vdf_trac_omp=1
438	call getin('iflag_vdf_trac', iflag_vdf_trac_omp)
439	iflag_con_trac_omp=1
440	call getin('iflag_con_trac', iflag_con_trac_omp)
441	iflag_the_trac_omp=1
442	call getin('iflag_the_trac', iflag_the_trac_omp)
443	!$OMP END MASTER
444	!$OMP BARRIER
445	convscav=convscav_omp
446	iflag_vdf_trac=iflag_vdf_trac_omp
447	iflag_con_trac=iflag_con_trac_omp
448	iflag_the_trac=iflag_the_trac_omp
449	write(lunout,*) 'phytrac passage dans routine conv avec lessivage', convscav
450	!
451	!Config Key = iflag_lscav
452	!Config Desc = Large scale scavenging parametrization: 0=none, 1=old(Genthon92),
453	! 2=1+PHeinrich, 3=Reddy_Boucher2004, 4=3+RPilon.
454	!Config Def = 1
455	!Config Help =
456	!
457	!$OMP MASTER
458	iflag_lscav_omp=1
459	call getin('iflag_lscav', iflag_lscav_omp)
460	ccntrAA_omp=1
461	ccntrENV_omp=1.
462	coefcoli_omp=0.001
463	call getin('ccntrAA', ccntrAA_omp)
464	call getin('ccntrENV', ccntrENV_omp)
465	call getin('coefcoli', coefcoli_omp)
466	!$OMP END MASTER
467	!$OMP BARRIER
468	iflag_lscav=iflag_lscav_omp
469	ccntrAA_in=ccntrAA_omp
470	ccntrENV_in=ccntrENV_omp
471	coefcoli_in=coefcoli_omp
472	!
473	SELECT CASE(iflag_lscav)
474	CASE(0)
475	WRITE(lunout,*) 'Large scale scavenging: none'
476	CASE(1)
477	WRITE(lunout,*) 'Large scale scavenging: C. Genthon, Tellus(1992), 44B, 371-389'
478	CASE(2)
479	WRITE(lunout,*) 'Large scale scavenging: C. Genthon, modified P. Heinrich'
480	CASE(3)
481	WRITE(lunout,*) 'Large scale scavenging: M. Shekkar Reddy and O. Boucher, JGR(2004), 109, D14202'
482	CASE(4)
483	WRITE(lunout,*) 'Large scale scavenging: Reddy and Boucher, modified R. Pilon'
484	END SELECT
485	!RomP <<<
486	WRITE(,) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
487	ALLOCATE( source(klon,nbtr), stat=ierr)
488	IF (ierr /= 0) CALL abort_physic('phytrac', 'pb in allocation 1',1)
489
490	ALLOCATE( aerosol(nbtr), stat=ierr)
491	IF (ierr /= 0) CALL abort_physic('phytrac', 'pb in allocation 2',1)
492
493
494	! Initialize module for specific tracers
495	IF (type_trac == 'inca') THEN
496	source(:,:)=init_source(:,:)
497	CALL tracinca_init(aerosol,lessivage)
498	ELSE IF (type_trac == 'repr') THEN
499	source(:,:)=0.
500	ELSE IF (type_trac == 'co2i') THEN
501	source(:,:)=0.
502	lessivage = .FALSE.
503	aerosol(:) = .FALSE.
504	pbl_flg(:) = 1
505	iflag_the_trac= 1
506	iflag_vdf_trac= 1
507	iflag_con_trac= 1
508	ELSE IF (type_trac == 'inco') THEN
509	source(:,1:nqCO2) = 0. ! from CO2i ModThL
510	source(:,nqCO2+1:nbtr)=init_source(:,:) ! from INCA ModThL
511	aerosol(1:nqCO2) = .FALSE. ! from CO2i ModThL
512	CALL tracinca_init(aerosol(nqCO2+1:nbtr),lessivage) ! from INCA ModThL
513	pbl_flg(1:nqCO2) = 1 ! From CO2i ModThL
514	iflag_the_trac = 1 ! From CO2i
515	iflag_vdf_trac = 1 ! From CO2i
516	iflag_con_trac = 1 ! From CO2i
517	ELSE IF (type_trac == 'coag' .AND. CPPKEY_STRATAER) THEN
518	source(:,:)=0.
519	CALL tracstrataer_init(aerosol,lessivage) ! init aerosols and lessivage param
520	ELSE IF (type_trac == 'lmdz') THEN
521	CALL traclmdz_init(pctsrf,xlat,xlon,ftsol,tr_seri,t_seri,pplay,sh,pdtphys,aerosol,lessivage)
522	ENDIF
523
524	!
525	!--initialising coefficients for scavenging in the case of NP
526	!
527	ALLOCATE(flag_cvltr(nbtr))
528	IF (iflag_con.EQ.3) THEN
529	!
530	ALLOCATE(ccntrAA(nbtr))
531	ALLOCATE(ccntrENV(nbtr))
532	ALLOCATE(coefcoli(nbtr))
533	!
534	DO it=1, nbtr
535	IF (type_trac == 'repr') THEN
536	flag_cvltr(it)=.FALSE.
537	ELSE IF (type_trac == 'inca') THEN
538	! IF ((it.EQ.id_Rn222) .OR. ((it.GE.id_SO2) .AND. (it.LE.id_NH3)) ) THEN
539	! !--gas-phase species
540	! flag_cvltr(it)=.FALSE.
541	!
542	! ELSEIF ( (it.GE.id_CIDUSTM) .AND. (it.LE.id_AIN) ) THEN
543	! !--insoluble aerosol species
544	! flag_cvltr(it)=.TRUE.
545	! ccntrAA(it)=0.7
546	! ccntrENV(it)=0.7
547	! coefcoli(it)=0.001
548	! ELSEIF ( (it.EQ.id_Pb210) .OR. ((it.GE.id_CSSSM) .AND. (it.LE.id_SSN))) THEN
549	! !--soluble aerosol species
550	! flag_cvltr(it)=.TRUE.
551	! ccntrAA(it)=0.9
552	! ccntrENV(it)=0.9
553	! coefcoli(it)=0.001
554	! ELSE
555	! WRITE(lunout,*) 'pb it=', it
556	! CALL abort_physic('phytrac','pb it scavenging',1)
557	! ENDIF
558	!--test OB
559	!--for now we do not scavenge in cvltr
560	flag_cvltr(it)=.FALSE.
561	ELSE IF (type_trac == 'co2i') THEN
562	!--co2 tracers are not scavenged
563	flag_cvltr(it)=.FALSE.
564	ELSE IF (type_trac == 'inco') THEN ! Add ThL
565	flag_cvltr(it)=.FALSE.
566	ELSE IF (type_trac == 'coag' .AND. CPPKEY_STRATAER) THEN
567	IF (convscav.and.aerosol(it)) THEN
568	flag_cvltr(it)=.TRUE.
569	ccntrAA(it) =ccntrAA_in
570	ccntrENV(it)=ccntrENV_in
571	coefcoli(it)=coefcoli_in
572	ELSE
573	flag_cvltr(it)=.FALSE.
574	ENDIF
575	ELSE IF (type_trac == 'lmdz') THEN
576	IF (convscav.and.aerosol(it)) THEN
577	flag_cvltr(it)=.TRUE.
578	ccntrAA(it) =ccntrAA_in !--a modifier par JYG a lire depuis fichier
579	ccntrENV(it)=ccntrENV_in
580	coefcoli(it)=coefcoli_in
581	ELSE
582	flag_cvltr(it)=.FALSE.
583	ENDIF
584	ENDIF
585	ENDDO
586	!
587	ELSE ! iflag_con .ne. 3
588	flag_cvltr(:) = .FALSE.
589	ENDIF
590	!
591	! print out all tracer flags
592	!
593	WRITE(lunout,*) 'print out all tracer flags'
594	WRITE(lunout,*) 'type_trac =', type_trac
595	WRITE(lunout,*) 'config_inca =', config_inca
596	WRITE(lunout,*) 'iflag_con_trac =', iflag_con_trac
597	WRITE(lunout,*) 'iflag_con =', iflag_con
598	WRITE(lunout,*) 'convscav =', convscav
599	WRITE(lunout,*) 'iflag_lscav =', iflag_lscav
600	WRITE(lunout,*) 'aerosol =', aerosol
601	WRITE(lunout,*) 'iflag_the_trac =', iflag_the_trac
602	WRITE(lunout,*) 'iflag_thermals =', iflag_thermals
603	WRITE(lunout,*) 'iflag_vdf_trac =', iflag_vdf_trac
604	WRITE(lunout,*) 'pbl_flg =', pbl_flg
605	WRITE(lunout,*) 'lessivage =', lessivage
606	write(lunout,*) 'flag_cvltr = ', flag_cvltr
607
608	IF (lessivage .AND. ANY(type_trac == ['inca','inco'])) &
609	CALL abort_physic('phytrac', 'lessivage=T config_inca=inca impossible',1)
610	!
611	ENDIF ! of IF (debutphy)
612	!############################################ END INITIALIZATION #######
613
614	DO k=1,klev
615	DO i=1,klon
616	zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
617	ENDDO
618	ENDDO
619	!
620	IF (id_be .GT. 0) THEN
621	DO k=1,klev
622	DO i=1,klon
623	sourceBE(i,k)=srcbe(i,k) !RomP -> pour sortie histrac
624	ENDDO
625	ENDDO
626	ENDIF
627
628	!===============================================================================
629	! -- Do specific treatment according to chemestry model or local LMDZ tracers
630	!
631	!===============================================================================
632	IF (type_trac == 'inca') THEN
633	! -- CHIMIE INCA config_inca = aero or chem --
634	! Appel fait en fin de phytrac pour avoir les emissions modifiees par
635	! la couche limite et la convection avant le calcul de la chimie
636
637	ELSE IF (type_trac == 'repr') THEN
638	! -- CHIMIE REPROBUS --
639	CALL tracreprobus(pdtphys, gmtime, debutphy, julien, &
640	presnivs, xlat, xlon, pphis, pphi, &
641	t_seri, pplay, paprs, sh , &
642	tr_seri)
643
644	ELSE IF (type_trac == 'co2i') THEN
645	! -- CO2 interactif --
646	! -- source is updated with FF and BB emissions
647	! -- and net fluxes from ocean and orchidee
648	! -- sign convention : positive into the atmosphere
649
650	CALL tracco2i(pdtphys, debutphy, &
651	xlat, xlon, pphis, pphi, &
652	t_seri, pplay, paprs, tr_seri, source)
653	ELSE IF (type_trac == 'inco') THEN ! Add ThL
654	CALL tracco2i(pdtphys, debutphy, &
655	xlat, xlon, pphis, pphi, &
656	t_seri, pplay, paprs, tr_seri, source)
657
658	ELSE IF (type_trac == 'coag' .AND. CPPKEY_STRATAER) THEN
659	! --STRATOSPHERIC AER IN THE STRAT --
660	CALL traccoag(pdtphys, gmtime, debutphy, julien, &
661	presnivs, xlat, xlon, pphis, pphi, &
662	t_seri, pplay, paprs, sh, rh , &
663	tr_seri)
664	ELSE IF (type_trac == 'lmdz') THEN
665	! -- Traitement des traceurs avec traclmdz
666	CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
667	cdragh, coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,iflag_vdf_trac>=0,sh, &
668	rh, pphi, ustar, wstar, ale_bl, ale_wake, u10m, v10m, &
669	tr_seri, source, d_tr_cl,d_tr_dec, zmasse) !RomP
670	ENDIF
671	!======================================================================
672	! -- Calcul de l'effet de la convection --
673	!======================================================================
674
675	IF (iflag_con_trac==1) THEN
676
677	DO it=1, nbtr
678	IF ( conv_flg(it) == 0 ) CYCLE
679	IF (iflag_con.LT.2) THEN
680	!--pas de transport convectif
681	d_tr_cv(:,:,it)=0.
682
683	ELSE IF (iflag_con.EQ.2) THEN
684	!--ancien transport convectif de Tiedtke
685
686	CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
687	pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
688	ELSE
689	!--nouveau transport convectif de Emanuel
690
691	IF (flag_cvltr(it)) THEN
692	!--nouveau transport convectif de Emanuel avec lessivage convectif
693	!
694	!
695	ccntrAA_3d(:,:) =ccntrAA(it)
696	ccntrENV_3d(:,:)=ccntrENV(it)
697	coefcoli_3d(:,:)=coefcoli(it)
698
699	!--beware this interface is a bit weird because it is called for each tracer
700	!--with the full array tr_seri even if only item it is processed
701
702	CALL cvltr_scav(pdtphys, da, phi,phi2,d1a,dam, mp,ep, &
703	sigd,sij,wght_cvfd,clw,elij,epmlmMm,eplaMm, &
704	pmflxr,pmflxs,evap,t_seri,wdtrainA,wdtrainM, &
705	paprs,it,tr_seri,upwd,dnwd,itop_con,ibas_con, &
706	ccntrAA_3d,ccntrENV_3d,coefcoli_3d, &
707	d_tr_cv,d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,flux_tr_wet, &
708	qDi,qPr, &
709	qPa,qMel,qTrdi,dtrcvMA,Mint, &
710	zmfd1a,zmfphi2,zmfdam)
711
712
713	ELSE !---flag_cvltr(it).EQ.FALSE
714	!--nouveau transport convectif de Emanuel mais pas de lessivage convectif
715
716	!--beware this interface is a bit weird because it is called for each tracer
717	!--with the full array tr_seri even if only item it is processed
718	!
719	CALL cvltr_noscav(it,pdtphys, da, phi,mp,wght_cvfd,paprs,pplay, & !jyg
720	tr_seri,upwd,dnwd,d_tr_cv) !jyg
721
722	ENDIF
723
724	ENDIF !--iflag
725
726	!--on ajoute les tendances
727
728	DO k = 1, klev
729	DO i = 1, klon
730	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
731	ENDDO
732	ENDDO
733
734	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//TRIM(int2str(it)))
735
736	ENDDO ! nbtr
737
738	IF (CPPKEY_STRATAER) THEN
739	IF (type_trac=='coag') THEN
740	! initialize wet deposition flux of sulfur
741	budg_dep_wet_ocs(:)=0.0
742	budg_dep_wet_so2(:)=0.0
743	budg_dep_wet_h2so4(:)=0.0
744	budg_dep_wet_part(:)=0.0
745	! compute wet deposition flux of sulfur (sum over gases and particles)
746	! and convert to kg(S)/m2/s
747	DO i = 1, klon
748	DO k = 1, klev
749	DO it = 1, nbtr
750	!do not include SO2 because most of it comes trom the troposphere
751	IF (it==id_OCS_strat) THEN
752	budg_dep_wet_ocs(i)=budg_dep_wet_ocs(i)+d_tr_cv(i,k,it)*(mSatom/mOCSmol) &
753	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
754	ELSEIF (it==id_SO2_strat) THEN
755	budg_dep_wet_so2(i)=budg_dep_wet_so2(i)+d_tr_cv(i,k,it)*(mSatom/mSO2mol) &
756	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
757	ELSEIF (it==id_H2SO4_strat) THEN
758	budg_dep_wet_h2so4(i)=budg_dep_wet_h2so4(i)+d_tr_cv(i,k,it)*(mSatom/mH2SO4mol) &
759	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
760	ELSEIF (it.GT.nbtr_sulgas) THEN
761	budg_dep_wet_part(i)=budg_dep_wet_part(i)+d_tr_cv(i,k,it)*(mSatom/mH2SO4mol) &
762	& dens_aer_dry4./3.RPI(mdw(it-nbtr_sulgas)/2.)**3 &
763	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
764	ENDIF
765	ENDDO
766	ENDDO
767	ENDDO
768	ENDIF
769	END IF
770
771	ENDIF ! convection
772
773	!======================================================================
774	! -- Calcul de l'effet des thermiques --
775	!======================================================================
776
777	DO it=1,nbtr
778	DO k=1,klev
779	DO i=1,klon
780	d_tr_th(i,k,it)=0.
781	tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
782	! the next safeguard causes some problem for stratospheric aerosol tracers (particle number)
783	! and there is little justification for it so it is commented out (4 December 2017) by OB
784	! if reinstated please keep the ifndef CPP_StratAer
785	!#ifndef CPP_StratAer
786	! tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
787	!#endif
788	ENDDO
789	ENDDO
790	ENDDO
791
792	IF (iflag_thermals.GT.0.AND.iflag_the_trac>0) THEN
793
794	DO it=1, nbtr
795
796	CALL thermcell_dq(klon,klev,1,pdtphys,fm_therm,entr_therm, &
797	zmasse,tr_seri(1:klon,1:klev,it), &
798	d_tr_th(1:klon,1:klev,it),ztra_th,0 )
799
800	DO k=1,klev
801	DO i=1,klon
802	d_tr_th(i,k,it)=pdtphys*d_tr_th(i,k,it)
803	tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr_th(i,k,it),0.)
804	ENDDO
805	ENDDO
806
807	ENDDO ! it
808
809	ENDIF ! Thermiques
810
811	!======================================================================
812	! -- Calcul de l'effet de la couche limite --
813	!======================================================================
814
815	IF (iflag_vdf_trac==1) THEN
816
817	! Injection during BL mixing
818	!
819	IF (CPPKEY_STRATAER) THEN
820	IF (type_trac=='coag') THEN
821
822	! initialize dry deposition flux of sulfur
823	budg_dep_dry_ocs(:)=0.0
824	budg_dep_dry_so2(:)=0.0
825	budg_dep_dry_h2so4(:)=0.0
826	budg_dep_dry_part(:)=0.0
827
828	! compute dry deposition velocity as function of surface type (numbers
829	! from IPSL note 23, 2002)
830	v_dep_dry(:) = pctsrf(:,is_ter) * 2.5e-3 &
831	& + pctsrf(:,is_oce) * 0.5e-3 &
832	& + pctsrf(:,is_lic) * 2.5e-3 &
833	& + pctsrf(:,is_sic) * 2.5e-3
834
835	! compute surface dry deposition flux
836	zrho(:,1)=pplay(:,1)/t_seri(:,1)/RD
837
838	DO it=1, nbtr
839	source(:,it) = - v_dep_dry(:) * tr_seri(:,1,it) * zrho(:,1)
840	ENDDO
841
842	ENDIF
843	END IF
844
845	DO it=1, nbtr
846	!
847	IF (pbl_flg(it) /= 0) THEN
848	!
849	CALL cltrac(pdtphys, coefh,t_seri, &
850	tr_seri(:,:,it), source(:,it), &
851	paprs, pplay, delp, &
852	d_tr_cl(:,:,it),d_tr_dry(:,it),flux_tr_dry(:,it))
853	!
854	tr_seri(:,:,it)=tr_seri(:,:,it)+d_tr_cl(:,:,it)
855	!
856	IF (CPPKEY_STRATAER) THEN
857	IF (type_trac=='coag') THEN
858	! compute dry deposition flux of sulfur (sum over gases and particles)
859	IF (it==id_OCS_strat) THEN
860	budg_dep_dry_ocs(:)=budg_dep_dry_ocs(:)-source(:,it)*(mSatom/mOCSmol)
861	ELSEIF (it==id_SO2_strat) THEN
862	budg_dep_dry_so2(:)=budg_dep_dry_so2(:)-source(:,it)*(mSatom/mSO2mol)
863	ELSEIF (it==id_H2SO4_strat) THEN
864	budg_dep_dry_h2so4(:)=budg_dep_dry_h2so4(:)-source(:,it)*(mSatom/mH2SO4mol)
865	ELSEIF (it.GT.nbtr_sulgas) THEN
866	budg_dep_dry_part(:)=budg_dep_dry_part(:)-source(:,it)(mSatom/mH2SO4mol)dens_aer_dry &
867	& 4./3.RPI(mdw(it-nbtr_sulgas)/2.)*3
868	ENDIF
869	ENDIF
870	END IF
871	!
872	ENDIF
873	!
874	ENDDO
875	!
876	ELSE IF (iflag_vdf_trac==0) THEN
877	!
878	! Injection of source in the first model layer
879	!
880	DO it=1,nbtr
881	d_tr_cl(:,1,it)=source(:,it)RG/delp(:,1)pdtphys
882	tr_seri(:,1,it)=tr_seri(:,1,it)+d_tr_cl(:,1,it)
883	ENDDO
884	d_tr_cl(:,2:klev,1:nbtr)=0.
885	!
886	ELSE IF (iflag_vdf_trac==-1) THEN
887	!
888	! Nothing happens
889	d_tr_cl=0.
890	!
891	ELSE
892	!
893	CALL abort_physic('iflag_vdf_trac', 'cas non prevu',1)
894	!
895	ENDIF ! couche limite
896
897	!======================================================================
898	! Calcul de l'effet de la precipitation grande echelle
899	! POUR INCA le lessivage est fait directement dans INCA
900	!======================================================================
901
902	IF (lessivage) THEN
903
904	ql_incloud_ref = 10.e-4
905	ql_incloud_ref = 5.e-4
906
907
908	! calcul du contenu en eau liquide au sein du nuage
909	ql_incl = ql_incloud_ref
910	! choix du lessivage
911	!
912	IF (iflag_lscav .EQ. 3 .OR. iflag_lscav .EQ. 4) THEN
913	! ******** Olivier Boucher version (3) possibly with modified ql_incl (4)
914	!
915	DO it = 1, nbtr
916
917	IF (aerosol(it)) THEN
918	! incloud scavenging and removal by large scale rain ! orig : ql_incl was replaced by 0.5e-3 kg/kg
919	! the value 0.5e-3 kg/kg is from Giorgi and Chameides (1986), JGR
920	! Liu (2001) proposed to use 1.5e-3 kg/kg
921
922	CALL lsc_scav(pdtphys,it,iflag_lscav,aerosol,ql_incl,prfl,psfl,rneb,beta_fisrt, &
923	beta_v1,pplay,paprs,t_seri,tr_seri,d_tr_insc,d_tr_bcscav,d_tr_evapls,qPrls)
924
925	!large scale scavenging tendency
926	DO k = 1, klev
927	DO i = 1, klon
928	d_tr_ls(i,k,it)=d_tr_insc(i,k,it)+d_tr_bcscav(i,k,it)+d_tr_evapls(i,k,it)
929	tr_seri(i,k,it)=tr_seri(i,k,it)+d_tr_ls(i,k,it)
930	ENDDO
931	ENDDO
932	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'lsc scav it = '//TRIM(int2str(it)))
933	ENDIF
934
935	ENDDO !tr
936
937	IF (CPPKEY_STRATAER) THEN
938	IF (type_trac=='coag') THEN
939	! compute wet deposition flux of sulfur (sum over gases and
940	! particles) and convert to kg(S)/m2/s
941	! adding contribution of d_tr_ls to d_tr_cv (above)
942	DO i = 1, klon
943	DO k = 1, klev
944	DO it = 1, nbtr
945	IF (it==id_OCS_strat) THEN
946	budg_dep_wet_ocs(i)=budg_dep_wet_ocs(i)+d_tr_ls(i,k,it)*(mSatom/mOCSmol) &
947	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
948	ELSEIF (it==id_SO2_strat) THEN
949	budg_dep_wet_so2(i)=budg_dep_wet_so2(i)+d_tr_ls(i,k,it)*(mSatom/mSO2mol) &
950	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
951	ELSEIF (it==id_H2SO4_strat) THEN
952	budg_dep_wet_h2so4(i)=budg_dep_wet_h2so4(i)+d_tr_ls(i,k,it)*(mSatom/mH2SO4mol) &
953	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
954	ELSEIF (it.GT.nbtr_sulgas) THEN
955	budg_dep_wet_part(i)=budg_dep_wet_part(i)+d_tr_ls(i,k,it)*(mSatom/mH2SO4mol) &
956	& dens_aer_dry4./3.RPI(mdw(it-nbtr_sulgas)/2.)**3 &
957	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
958	ENDIF
959	ENDDO
960	ENDDO
961	ENDDO
962	ENDIF
963	END IF
964
965	ELSE IF (iflag_lscav .EQ. 2) THEN ! frac_impa, frac_nucl
966	! ********* modified old version
967
968	d_tr_lessi_nucl(:,:,:) = 0.
969	d_tr_lessi_impa(:,:,:) = 0.
970	flestottr(:,:,:) = 0.
971	! Tendance des aerosols nuclees et impactes
972	DO it = 1, nbtr
973	IF (aerosol(it)) THEN
974	his_dh(:)=0.
975	DO k = 1, klev
976	DO i = 1, klon
977	!PhH
978	zrho(i,k)=pplay(i,k)/t_seri(i,k)/RD
979	zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
980	!
981	ENDDO
982	ENDDO
983
984	DO k=klev-1, 1, -1
985	DO i=1, klon
986	! d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.)
987	dx=d_tr_ls(i,k,it)
988	his_dh(i)=his_dh(i)-dxzrho(i,k)zdz(i,k)/pdtphys ! kg/m2/s
989	evaplsc = prfl(i,k) - prfl(i,k+1) + psfl(i,k) - psfl(i,k+1)
990	! Evaporation Partielle -> Liberation Partielle 0.5*evap
991	IF ( evaplsc .LT.0..and.abs(prfl(i,k+1)+psfl(i,k+1)).gt.1.e-10) THEN
992	evaplsc = (-evaplsc)/(prfl(i,k+1)+psfl(i,k+1))
993	! evaplsc est donc positif, his_dh(i) est positif
994	!--------------
995	d_tr_evapls(i,k,it)=0.5evaplsc(d_tr_lessi_nucl(i,k+1,it) &
996	+d_tr_lessi_impa(i,k+1,it))
997	!------------- d_tr_evapls(i,k,it)=-0.5evaplsc(d_tr_lsc(i,k+1,it))
998	beta=0.5*evaplsc
999	if ((prfl(i,k)+psfl(i,k)).lt.1.e-10) THEN
1000	beta=1.0*evaplsc
1001	endif
1002	dx=betahis_dh(i)/zrho(i,k)/zdz(i,k)pdtphys
1003	his_dh(i)=(1.-beta)*his_dh(i) ! tracer from
1004	d_tr_evapls(i,k,it)=dx
1005	ENDIF
1006	d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.) &
1007	+d_tr_evapls(i,k,it)
1008
1009	!--------------
1010	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
1011	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
1012	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
1013	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
1014	!
1015	! Flux lessivage total
1016	flestottr(i,k,it) = flestottr(i,k,it) - &
1017	( d_tr_lessi_nucl(i,k,it) + &
1018	d_tr_lessi_impa(i,k,it) ) * &
1019	( paprs(i,k)-paprs(i,k+1) ) / &
1020	(RG * pdtphys)
1021	!! Mise a jour des traceurs due a l'impaction,nucleation
1022	! tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
1023	!! calcul de la tendance liee au lessivage stratiforme
1024	! d_tr_ls(i,k,it)=tr_seri(i,k,it)*&
1025	! (1.-1./(frac_impa(i,k)*frac_nucl(i,k)))
1026	!--------------
1027	ENDDO
1028	ENDDO
1029	ENDIF
1030	ENDDO
1031	! ********* end modified old version
1032
1033	ELSE IF (iflag_lscav .EQ. 1) THEN ! frac_impa, frac_nucl
1034	! ********* old version
1035
1036	d_tr_lessi_nucl(:,:,:) = 0.
1037	d_tr_lessi_impa(:,:,:) = 0.
1038	flestottr(:,:,:) = 0.
1039	!=========================
1040	! LESSIVAGE LARGE SCALE :
1041	!=========================
1042
1043	! Tendance des aerosols nuclees et impactes
1044	! -----------------------------------------
1045	DO it = 1, nbtr
1046	IF (aerosol(it)) THEN
1047	DO k = 1, klev
1048	DO i = 1, klon
1049	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
1050	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
1051	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
1052	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
1053
1054	!
1055	! Flux lessivage total
1056	! ------------------------------------------------------------
1057	flestottr(i,k,it) = flestottr(i,k,it) - &
1058	( d_tr_lessi_nucl(i,k,it) + &
1059	d_tr_lessi_impa(i,k,it) ) * &
1060	( paprs(i,k)-paprs(i,k+1) ) / &
1061	(RG * pdtphys)
1062	!
1063	! Mise a jour des traceurs due a l'impaction,nucleation
1064	! ----------------------------------------------------------------------
1065	tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
1066	ENDDO
1067	ENDDO
1068	ENDIF
1069	ENDDO
1070
1071	! ********* end old version
1072	ENDIF ! iflag_lscav .EQ. 1, 2, 3 or 4
1073	!
1074	ENDIF ! lessivage
1075
1076
1077	! -- CHIMIE INCA config_inca = aero or chem --
1078	IF (ANY(type_trac == ['inca','inco'])) THEN ! ModThL
1079
1080	CALL tracinca(&
1081	nstep, julien, gmtime, lafin, &
1082	pdtphys, t_seri, paprs, pplay, &
1083	pmfu, upwd, ftsol, pctsrf, pphis, &
1084	pphi, albsol, sh, ch, rh, &
1085	cldfra, rneb, diafra, cldliq, &
1086	itop_con, ibas_con, pmflxr, pmflxs, &
1087	prfl, psfl, aerosol_couple, flxmass_w, &
1088	tau_aero, piz_aero, cg_aero, ccm, &
1089	rfname, &
1090	tr_seri(:,:,1+nqCO2:nbtr), source(:,1+nqCO2:nbtr)) ! ModThL
1091	ENDIF
1092
1093	END SUBROUTINE phytrac
1094
1095	END MODULE

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format