Context Navigation

phytrac_mod.F90 @ 4049

Last change on this file since 4049 was 2394, checked in by acozic, 9 years ago

to fit with inca we need to increase dimension of naero_grp in aero_mod.
When lmdz is running without inca dim 10 to 12 will be empty and there is no change in results

these modifications only impact the old physics coupled with inca

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

File size: 36.6 KB

Line
1	!$Id$
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 :: qPr,qDi ! concentration tra dans pluie,air descente insaturee
38	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qPa,qMel
39	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qTrdi,dtrcvMA ! conc traceur descente air insaturee et td convective MA
40	! RomP <<<
41	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_th ! Td thermique
42	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_lessi_impa ! Td du lessivage par impaction
43	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_lessi_nucl ! Td du lessivage par nucleation
44	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: qPrls !jyg: concentration tra dans pluie LS a la surf.
45	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: d_tr_dry ! Td depot sec/traceur (1st layer),ALLOCATABLE,SAVE jyg
46	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: flux_tr_dry ! depot sec/traceur (surface),ALLOCATABLE,SAVE jyg
47
48	!$OMP THREADPRIVATE(qPa,qMel,qTrdi,dtrcvMA,d_tr_th,d_tr_lessi_impa,d_tr_lessi_nucl)
49	!$OMP THREADPRIVATE(d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qPr,qDi)
50	!$OMP THREADPRIVATE(d_tr_insc,d_tr_bcscav,d_tr_evapls,d_tr_ls,qPrls)
51	!$OMP THREADPRIVATE(d_tr,d_tr_cl,d_tr_dry,flux_tr_dry,d_tr_dec,d_tr_cv)
52
53
54	CONTAINS
55
56	SUBROUTINE phytrac( &
57	nstep, julien, gmtime, debutphy, &
58	lafin, pdtphys, u, v, t_seri, &
59	paprs, pplay, pmfu, pmfd, &
60	pen_u, pde_u, pen_d, pde_d, &
61	cdragh, coefh, fm_therm, entr_therm, &
62	yu1, yv1, ftsol, pctsrf, &
63	ustar, u10m, v10m, &
64	wstar, ale_bl, ale_wake, &
65	xlat, xlon, &
66	frac_impa,frac_nucl,beta_fisrt,beta_v1, &
67	presnivs, pphis, pphi, albsol, &
68	sh, rh, cldfra, rneb, &
69	diafra, cldliq, itop_con, ibas_con, &
70	pmflxr, pmflxs, prfl, psfl, &
71	da, phi, mp, upwd, &
72	phi2, d1a, dam, sij, wght_cvfd, & ! RomP +RL
73	wdtrainA, wdtrainM, sigd, clw, elij, & ! RomP
74	evap, ep, epmlmMm, eplaMm, & ! RomP
75	dnwd, aerosol_couple, flxmass_w, &
76	tau_aero, piz_aero, cg_aero, ccm, &
77	rfname, &
78	d_tr_dyn, & ! RomP
79	tr_seri)
80	!
81	!======================================================================
82	! Auteur(s) FH
83	! Objet: Moniteur general des tendances traceurs
84	! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
85	! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
86	!======================================================================
87
88	USE ioipsl
89	USE phys_cal_mod, only : hour
90	USE dimphy
91	USE infotrac_phy, ONLY: nbtr, type_trac, conv_flg, solsym, pbl_flg
92	USE mod_grid_phy_lmdz
93	USE mod_phys_lmdz_para
94	USE iophy
95	USE traclmdz_mod
96	USE tracinca_mod
97	USE tracreprobus_mod
98	USE indice_sol_mod
99
100	USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root
101	USE print_control_mod, ONLY: lunout
102	USE aero_mod, ONLY : naero_grp
103
104	IMPLICIT NONE
105
106	INCLUDE "YOMCST.h"
107	INCLUDE "clesphys.h"
108	INCLUDE "thermcell.h"
109	!==========================================================================
110	! -- ARGUMENT DESCRIPTION --
111	!==========================================================================
112
113	! Input arguments
114	!----------------
115	!Configuration grille,temps:
116	INTEGER,INTENT(IN) :: nstep ! Appel physique
117	INTEGER,INTENT(IN) :: julien ! Jour julien
118	REAL,INTENT(IN) :: gmtime ! Heure courante
119	REAL,INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
120	LOGICAL,INTENT(IN) :: debutphy ! le flag de l'initialisation de la physique
121	LOGICAL,INTENT(IN) :: lafin ! le flag de la fin de la physique
122
123	REAL,DIMENSION(klon),INTENT(IN) :: xlat ! latitudes pour chaque point
124	REAL,DIMENSION(klon),INTENT(IN) :: xlon ! longitudes pour chaque point
125	!
126	!Physique:
127	!--------
128	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
129	REAL,DIMENSION(klon,klev),INTENT(IN) :: u ! variable not used
130	REAL,DIMENSION(klon,klev),INTENT(IN) :: v ! variable not used
131	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
132	REAL,DIMENSION(klon,klev),INTENT(IN) :: rh ! humidite relative
133	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
134	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
135	REAL,DIMENSION(klon,klev),INTENT(IN) :: pphi ! geopotentiel
136	REAL,DIMENSION(klon),INTENT(IN) :: pphis
137	REAL,DIMENSION(klev),INTENT(IN) :: presnivs
138	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldliq ! eau liquide nuageuse
139	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldfra ! fraction nuageuse (tous les nuages)
140	REAL,DIMENSION(klon,klev),INTENT(IN) :: diafra ! fraction nuageuse (convection ou stratus artificiels)
141	REAL,DIMENSION(klon,klev),INTENT(IN) :: rneb ! fraction nuageuse (grande echelle)
142	!
143	REAL :: ql_incl ! contenu en eau liquide nuageuse dans le nuage ! ql_incl=oliq/rneb
144	REAL,DIMENSION(klon,klev),INTENT(IN) :: beta_fisrt ! taux de conversion de l'eau cond (de fisrtilp)
145	REAL,DIMENSION(klon,klev),INTENT(out) :: beta_v1 ! -- (originale version)
146
147	!
148	INTEGER,DIMENSION(klon),INTENT(IN) :: itop_con
149	INTEGER,DIMENSION(klon),INTENT(IN) :: ibas_con
150	REAL,DIMENSION(klon),INTENT(IN) :: albsol ! albedo surface
151	!
152	!Dynamique
153	!--------
154	REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: d_tr_dyn
155	!
156	!Convection:
157	!----------
158	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu ! flux de masse dans le panache montant
159	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd ! flux de masse dans le panache descendant
160	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
161	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
162	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
163	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
164
165	!...Tiedke
166	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
167	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
168
169	LOGICAL,INTENT(IN) :: aerosol_couple
170	REAL,DIMENSION(klon,klev),INTENT(IN) :: flxmass_w
171	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: tau_aero
172	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: piz_aero
173	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: cg_aero
174	CHARACTER(len=4),DIMENSION(naero_grp),INTENT(IN) :: rfname
175	REAL,DIMENSION(klon,klev,2),INTENT(IN) :: ccm
176	!... K.Emanuel
177	REAL,DIMENSION(klon,klev),INTENT(IN) :: da
178	REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
179	! RomP >>>
180	REAL,DIMENSION(klon,klev),INTENT(IN) :: d1a,dam
181	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2
182	!
183	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainA
184	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainM
185	REAL,DIMENSION(klon),INTENT(IN) :: sigd
186	! ---- RomP flux entraine, detraine et precipitant kerry Emanuel
187	REAL,DIMENSION(klon,klev),INTENT(IN) :: evap
188	REAL,DIMENSION(klon,klev),INTENT(IN) :: ep
189	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij
190	REAL,DIMENSION(klon,klev),INTENT(IN) :: wght_cvfd !RL
191	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij
192	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm
193	REAL,DIMENSION(klon,klev),INTENT(IN) :: eplaMm
194	REAL,DIMENSION(klon,klev),INTENT(IN) :: clw
195	! RomP <<<
196
197	!
198	REAL,DIMENSION(klon,klev),INTENT(IN) :: mp
199	REAL,DIMENSION(klon,klev),INTENT(IN) :: upwd ! saturated updraft mass flux
200	REAL,DIMENSION(klon,klev),INTENT(IN) :: dnwd ! saturated downdraft mass flux
201	!
202	!Thermiques:
203	!----------
204	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: fm_therm
205	REAL,DIMENSION(klon,klev),INTENT(IN) :: entr_therm
206	!
207	!Couche limite:
208	!--------------
209	!
210	!
211	REAL,DIMENSION(:),INTENT(IN) :: cdragh ! (klon) coeff drag pour T et Q
212	REAL,DIMENSION(:,:),INTENT(IN) :: coefh ! (klon,klev) coeff melange CL (m**2/s)
213	REAL,DIMENSION(:),INTENT(IN) :: ustar,u10m,v10m ! (klon) u* & vent a 10m (m/s)
214	REAL,DIMENSION(:),INTENT(IN) :: wstar,ale_bl,ale_wake ! (klon) w* and Avail. Lifting Ener.
215	REAL,DIMENSION(:),INTENT(IN) :: yu1 ! (klon) vents au premier niveau
216	REAL,DIMENSION(:),INTENT(IN) :: yv1 ! (klon) vents au premier niveau
217
218	!
219	!Lessivage:
220	!----------
221	!
222	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrAA
223	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrENV
224	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: coefcoli
225	LOGICAL, DIMENSION(:), ALLOCATABLE, SAVE :: flag_cvltr
226	!$OMP THREADPRIVATE(ccntrAA,ccntrENV,coefcoli,flag_cvltr)
227	REAL, DIMENSION(klon,klev) :: ccntrAA_3d
228	REAL, DIMENSION(klon,klev) :: ccntrENV_3d
229	REAL, DIMENSION(klon,klev) :: coefcoli_3d
230	!
231	! pour le ON-LINE
232	!
233	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
234	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
235
236	! Arguments necessaires pour les sources et puits de traceur:
237	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
238	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
239
240
241	! Output argument
242	!----------------
243	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
244	REAL,DIMENSION(klon,klev) :: sourceBE
245	!=======================================================================================
246	! -- LOCAL VARIABLES --
247	!=======================================================================================
248
249	INTEGER :: i, k, it
250	INTEGER :: nsplit
251
252	!Sources et Reservoirs de traceurs (ex:Radon):
253	!--------------------------------------------
254	!
255	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source ! a voir lorsque le flux de surface est prescrit
256	!$OMP THREADPRIVATE(source)
257
258	!
259	!Entrees/Sorties: (cf ini_histrac.h et write_histrac.h)
260	!---------------
261	INTEGER :: iiq, ierr
262	INTEGER :: nhori, nvert
263	REAL :: zsto, zout, zjulian
264	INTEGER,SAVE :: nid_tra ! pointe vers le fichier histrac.nc
265	!$OMP THREADPRIVATE(nid_tra)
266	REAL,DIMENSION(klon) :: zx_tmp_fi2d ! variable temporaire grille physique
267	INTEGER :: itau_w ! pas de temps ecriture = nstep + itau_phy
268	LOGICAL,PARAMETER :: ok_sync=.TRUE.
269
270	!
271	! Nature du traceur
272	!------------------
273	LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol ! aerosol(it) = true => aerosol => lessivage
274	!$OMP THREADPRIVATE(aerosol) ! aerosol(it) = false => gaz
275	REAL,DIMENSION(klon,klev) :: delp ! epaisseur de couche (Pa)
276	!
277	! Tendances de traceurs (Td) et flux de traceurs:
278	!------------------------
279	REAL,DIMENSION(klon,klev) :: d_tr ! Td dans l'atmosphere
280	REAL,DIMENSION(klon,klev) :: Mint
281	REAL,DIMENSION(klon,klev,nbtr) :: zmfd1a
282	REAL,DIMENSION(klon,klev,nbtr) :: zmfdam
283	REAL,DIMENSION(klon,klev,nbtr) :: zmfphi2
284	! Physique
285	!----------
286	REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
287	REAL,DIMENSION(klon,klev) :: zmasse ! densité atmosphérique Kg/m2
288	REAL,DIMENSION(klon,klev) :: ztra_th
289	!PhH
290	REAL,DIMENSION(klon,klev) :: zrho
291	REAL,DIMENSION(klon,klev) :: zdz
292	REAL :: evaplsc,dx,beta ! variable pour lessivage Genthon
293	REAL,DIMENSION(klon) :: his_dh ! ---
294	! in-cloud scav variables
295	REAL :: ql_incloud_ref ! ref value of in-cloud condensed water content
296
297	!Controles:
298	!---------
299	INTEGER,SAVE :: iflag_vdf_trac,iflag_con_trac,iflag_the_trac
300	INTEGER,SAVE :: iflag_con_trac_omp, iflag_vdf_trac_omp,iflag_the_trac_omp
301	!$OMP THREADPRIVATE(iflag_vdf_trac,iflag_con_trac,iflag_the_trac)
302
303	LOGICAL,SAVE :: lessivage
304	!$OMP THREADPRIVATE(lessivage)
305
306	!RomP >>>
307	INTEGER,SAVE :: iflag_lscav_omp,iflag_lscav
308	REAL, SAVE :: ccntrAA_in,ccntrAA_omp
309	REAL, SAVE :: ccntrENV_in,ccntrENV_omp
310	REAL, SAVE :: coefcoli_in,coefcoli_omp
311
312	LOGICAL,SAVE :: convscav_omp,convscav
313	!$OMP THREADPRIVATE(iflag_lscav)
314	!$OMP THREADPRIVATE(ccntrAA_in,ccntrENV_in,coefcoli_in)
315	!$OMP THREADPRIVATE(convscav)
316	!RomP <<<
317	!######################################################################
318	! -- INITIALIZATION --
319	!######################################################################
320	IF (debutphy) THEN
321	ALLOCATE(d_tr_cl(klon,klev,nbtr),d_tr_dry(klon,nbtr))
322	ALLOCATE(flux_tr_dry(klon,nbtr),d_tr_dec(klon,klev,nbtr),d_tr_cv(klon,klev,nbtr))
323	ALLOCATE(d_tr_insc(klon,klev,nbtr),d_tr_bcscav(klon,klev,nbtr))
324	ALLOCATE(d_tr_evapls(klon,klev,nbtr),d_tr_ls(klon,klev,nbtr))
325	ALLOCATE(qPrls(klon,nbtr),d_tr_trsp(klon,klev,nbtr))
326	ALLOCATE(d_tr_sscav(klon,klev,nbtr),d_tr_sat(klon,klev,nbtr))
327	ALLOCATE(d_tr_uscav(klon,klev,nbtr),qPr(klon,klev,nbtr),qDi(klon,klev,nbtr))
328	ALLOCATE(qPa(klon,klev,nbtr),qMel(klon,klev,nbtr))
329	ALLOCATE(qTrdi(klon,klev,nbtr),dtrcvMA(klon,klev,nbtr))
330	ALLOCATE(d_tr_th(klon,klev,nbtr))
331	ALLOCATE(d_tr_lessi_impa(klon,klev,nbtr),d_tr_lessi_nucl(klon,klev,nbtr))
332	ENDIF
333
334	DO k=1,klev
335	DO i=1,klon
336	sourceBE(i,k)=0.
337	Mint(i,k)=0.
338	zrho(i,k)=0.
339	zdz(i,k)=0.
340	END DO
341	END DO
342
343	DO it=1, nbtr
344	DO k=1,klev
345	DO i=1,klon
346	d_tr_insc(i,k,it)=0.
347	d_tr_bcscav(i,k,it)=0.
348	d_tr_evapls(i,k,it)=0.
349	d_tr_ls(i,k,it)=0.
350	d_tr_cv(i,k,it)=0.
351	d_tr_cl(i,k,it)=0.
352	d_tr_trsp(i,k,it)=0.
353	d_tr_sscav(i,k,it)=0.
354	d_tr_sat(i,k,it)=0.
355	d_tr_uscav(i,k,it)=0.
356	d_tr_lessi_impa(i,k,it)=0.
357	d_tr_lessi_nucl(i,k,it)=0.
358	qDi(i,k,it)=0.
359	qPr(i,k,it)=0.
360	qPa(i,k,it)=0.
361	qMel(i,k,it)=0.
362	qTrdi(i,k,it)=0.
363	dtrcvMA(i,k,it)=0.
364	zmfd1a(i,k,it)=0.
365	zmfdam(i,k,it)=0.
366	zmfphi2(i,k,it)=0.
367	END DO
368	END DO
369	END DO
370
371	DO k = 1, klev
372	DO i = 1, klon
373	delp(i,k) = paprs(i,k)-paprs(i,k+1)
374	END DO
375	END DO
376
377	IF (debutphy) THEN
378	!!jyg
379	!$OMP BARRIER
380	ecrit_tra=86400. ! frequence de stokage en dur
381	! obsolete car remplace par des ecritures dans phys_output_write
382	!RomP >>>
383	!
384	!Config Key = convscav
385	!Config Desc = Convective scavenging switch: 0=off, 1=on.
386	!Config Def = .false.
387	!Config Help =
388	!
389	!$OMP MASTER
390	convscav_omp=.false.
391	call getin('convscav', convscav_omp)
392	iflag_vdf_trac_omp=1
393	call getin('iflag_vdf_trac', iflag_vdf_trac_omp)
394	iflag_con_trac_omp=1
395	call getin('iflag_con_trac', iflag_con_trac_omp)
396	iflag_the_trac_omp=1
397	call getin('iflag_the_trac', iflag_the_trac_omp)
398	!$OMP END MASTER
399	!$OMP BARRIER
400	convscav=convscav_omp
401	iflag_vdf_trac=iflag_vdf_trac_omp
402	iflag_con_trac=iflag_con_trac_omp
403	iflag_the_trac=iflag_the_trac_omp
404	write(lunout,*) 'phytrac passage dans routine conv avec lessivage', convscav
405	!
406	!Config Key = iflag_lscav
407	!Config Desc = Large scale scavenging parametrization: 0=none, 1=old(Genthon92),
408	! 2=1+PHeinrich, 3=Reddy_Boucher2004, 4=3+RPilon.
409	!Config Def = 1
410	!Config Help =
411	!
412	!$OMP MASTER
413	iflag_lscav_omp=1
414	call getin('iflag_lscav', iflag_lscav_omp)
415	ccntrAA_omp=1
416	ccntrENV_omp=1.
417	coefcoli_omp=0.001
418	call getin('ccntrAA', ccntrAA_omp)
419	call getin('ccntrENV', ccntrENV_omp)
420	call getin('coefcoli', coefcoli_omp)
421	!$OMP END MASTER
422	!$OMP BARRIER
423	iflag_lscav=iflag_lscav_omp
424	ccntrAA_in=ccntrAA_omp
425	ccntrENV_in=ccntrENV_omp
426	coefcoli_in=coefcoli_omp
427	!
428	SELECT CASE(iflag_lscav)
429	CASE(0)
430	WRITE(lunout,*) 'Large scale scavenging: none'
431	CASE(1)
432	WRITE(lunout,*) 'Large scale scavenging: C. Genthon, Tellus(1992), 44B, 371-389'
433	CASE(2)
434	WRITE(lunout,*) 'Large scale scavenging: C. Genthon, modified P. Heinrich'
435	CASE(3)
436	WRITE(lunout,*) 'Large scale scavenging: M. Shekkar Reddy and O. Boucher, JGR(2004), 109, D14202'
437	CASE(4)
438	WRITE(lunout,*) 'Large scale scavenging: Reddy and Boucher, modified R. Pilon'
439	END SELECT
440	!RomP <<<
441	WRITE(,) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
442	ALLOCATE( source(klon,nbtr), stat=ierr)
443	IF (ierr /= 0) CALL abort_physic('phytrac', 'pb in allocation 1',1)
444
445	ALLOCATE( aerosol(nbtr), stat=ierr)
446	IF (ierr /= 0) CALL abort_physic('phytrac', 'pb in allocation 2',1)
447
448
449	! Initialize module for specific tracers
450	SELECT CASE(type_trac)
451	CASE('lmdz')
452	CALL traclmdz_init(pctsrf,xlat,xlon,ftsol,tr_seri,t_seri,pplay,sh,pdtphys,aerosol,lessivage)
453	CASE('inca')
454	source(:,:)=0.
455	CALL tracinca_init(aerosol,lessivage)
456	CASE('repr')
457	source(:,:)=0.
458	END SELECT
459
460	!
461	!--initialising coefficients for scavenging in the case of NP
462	!
463	ALLOCATE(flag_cvltr(nbtr))
464	IF (iflag_con.EQ.3) THEN
465	!
466	ALLOCATE(ccntrAA(nbtr))
467	ALLOCATE(ccntrENV(nbtr))
468	ALLOCATE(coefcoli(nbtr))
469	!
470	DO it=1, nbtr
471	SELECT CASE(type_trac)
472	CASE('lmdz')
473	IF (convscav.and.aerosol(it)) THEN
474	flag_cvltr(it)=.true.
475	ccntrAA(it) =ccntrAA_in !--a modifier par JYG a lire depuis fichier
476	ccntrENV(it)=ccntrENV_in
477	coefcoli(it)=coefcoli_in
478	ELSE
479	flag_cvltr(it)=.false.
480	ENDIF
481
482	CASE('inca')
483	! IF ((it.EQ.id_Rn222) .OR. ((it.GE.id_SO2) .AND. (it.LE.id_NH3)) ) THEN
484	! !--gas-phase species
485	! flag_cvltr(it)=.false.
486	!
487	! ELSEIF ( (it.GE.id_CIDUSTM) .AND. (it.LE.id_AIN) ) THEN
488	! !--insoluble aerosol species
489	! flag_cvltr(it)=.true.
490	! ccntrAA(it)=0.7
491	! ccntrENV(it)=0.7
492	! coefcoli(it)=0.001
493	! ELSEIF ( (it.EQ.id_Pb210) .OR. ((it.GE.id_CSSSM) .AND. (it.LE.id_SSN))) THEN
494	! !--soluble aerosol species
495	! flag_cvltr(it)=.true.
496	! ccntrAA(it)=0.9
497	! ccntrENV(it)=0.9
498	! coefcoli(it)=0.001
499	! ELSE
500	! WRITE(lunout,*) 'pb it=', it
501	! CALL abort_physic('phytrac','pb it scavenging',1)
502	! ENDIF
503	!--test OB
504	!--for now we do not scavenge in cvltr
505	flag_cvltr(it)=.false.
506	END SELECT
507	ENDDO
508	!
509	ELSE ! iflag_con .ne. 3
510	flag_cvltr(:) = .false.
511	ENDIF
512	!
513	! Initialize diagnostic output
514	! ----------------------------
515	#ifdef CPP_IOIPSL
516	! INCLUDE "ini_histrac.h"
517	#endif
518	!
519	! print out all tracer flags
520	!
521	WRITE(lunout,*) 'print out all tracer flags'
522	WRITE(lunout,*) 'type_trac =', type_trac
523	WRITE(lunout,*) 'config_inca =', config_inca
524	WRITE(lunout,*) 'iflag_con_trac =', iflag_con_trac
525	WRITE(lunout,*) 'iflag_con =', iflag_con
526	WRITE(lunout,*) 'convscav =', convscav
527	WRITE(lunout,*) 'iflag_lscav =', iflag_lscav
528	WRITE(lunout,*) 'aerosol =', aerosol
529	WRITE(lunout,*) 'iflag_the_trac =', iflag_the_trac
530	WRITE(lunout,*) 'iflag_thermals =', iflag_thermals
531	WRITE(lunout,*) 'iflag_vdf_trac =', iflag_vdf_trac
532	WRITE(lunout,*) 'pbl_flg =', pbl_flg
533	WRITE(lunout,*) 'lessivage =', lessivage
534	write(lunout,*) 'flag_cvltr = ', flag_cvltr
535
536	IF (lessivage.AND.config_inca.EQ.'inca') THEN
537	CALL abort_physic('phytrac', 'lessivage=T config_inca=inca impossible',1)
538	STOP
539	ENDIF
540	!
541	END IF ! of IF (debutphy)
542	!############################################ END INITIALIZATION #######
543
544	DO k=1,klev
545	DO i=1,klon
546	zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
547	END DO
548	END DO
549	!
550	IF (id_be .GT. 0) THEN
551	DO k=1,klev
552	DO i=1,klon
553	sourceBE(i,k)=srcbe(i,k) !RomP -> pour sortie histrac
554	END DO
555	END DO
556	ENDIF
557
558	!===============================================================================
559	! -- Do specific treatment according to chemestry model or local LMDZ tracers
560	!
561	!===============================================================================
562	SELECT CASE(type_trac)
563	CASE('lmdz')
564	! -- Traitement des traceurs avec traclmdz
565	CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
566	cdragh, coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,iflag_vdf_trac>=0,sh, &
567	rh, pphi, ustar, wstar, ale_bl, ale_wake, u10m, v10m, &
568	tr_seri, source, d_tr_cl,d_tr_dec, zmasse) !RomP
569
570	CASE('inca')
571	! -- CHIMIE INCA config_inca = aero or chem --
572
573	CALL tracinca(&
574	nstep, julien, gmtime, lafin, &
575	pdtphys, t_seri, paprs, pplay, &
576	pmfu, upwd, ftsol, pctsrf, pphis, &
577	pphi, albsol, sh, rh, &
578	cldfra, rneb, diafra, cldliq, &
579	itop_con, ibas_con, pmflxr, pmflxs, &
580	prfl, psfl, aerosol_couple, flxmass_w, &
581	tau_aero, piz_aero, cg_aero, ccm, &
582	rfname, &
583	tr_seri, source)
584
585	CASE('repr')
586	! -- CHIMIE REPROBUS --
587
588	CALL tracreprobus(pdtphys, gmtime, debutphy, julien, &
589	presnivs, xlat, xlon, pphis, pphi, &
590	t_seri, pplay, paprs, sh , &
591	tr_seri)
592
593	END SELECT
594	!======================================================================
595	! -- Calcul de l'effet de la convection --
596	!======================================================================
597
598	IF (iflag_con_trac==1) THEN
599
600	DO it=1, nbtr
601	IF ( conv_flg(it) == 0 ) CYCLE
602	IF (iflag_con.LT.2) THEN
603	!--pas de transport convectif
604
605	d_tr_cv(:,:,it)=0.
606	ELSE IF (iflag_con.EQ.2) THEN
607	!--ancien transport convectif de Tiedtke
608
609	CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
610	pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
611	ELSE
612	!--nouveau transport convectif de Emanuel
613
614	IF (flag_cvltr(it)) THEN
615	!--nouveau transport convectif de Emanuel avec lessivage convectif
616	!
617	!
618	ccntrAA_3d(:,:) =ccntrAA(it)
619	ccntrENV_3d(:,:)=ccntrENV(it)
620	coefcoli_3d(:,:)=coefcoli(it)
621
622	!--beware this interface is a bit weird because it is called for each tracer
623	!--with the full array tr_seri even if only item it is processed
624
625	print*,'CV SCAV ',it,ccntrAA(it),ccntrENV(it)
626
627	CALL cvltr_scav(pdtphys, da, phi,phi2,d1a,dam, mp,ep, &
628	sigd,sij,wght_cvfd,clw,elij,epmlmMm,eplaMm, &
629	pmflxr,pmflxs,evap,t_seri,wdtrainA,wdtrainM, &
630	paprs,it,tr_seri,upwd,dnwd,itop_con,ibas_con, &
631	ccntrAA_3d,ccntrENV_3d,coefcoli_3d, &
632	d_tr_cv,d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qDi,qPr,&
633	qPa,qMel,qTrdi,dtrcvMA,Mint, &
634	zmfd1a,zmfphi2,zmfdam)
635
636
637	ELSE !---flag_cvltr(it).EQ.FALSE
638	!--nouveau transport convectif de Emanuel mais pas de lessivage convectif
639
640	!--beware this interface is a bit weird because it is called for each tracer
641	!--with the full array tr_seri even if only item it is processed
642	!
643	CALL cvltr_noscav(it,pdtphys, da, phi,mp,wght_cvfd,paprs,pplay, & !jyg
644	tr_seri,upwd,dnwd,d_tr_cv) !jyg
645
646	ENDIF
647
648	ENDIF !--iflag
649
650	!--on ajoute les tendances
651
652	DO k = 1, klev
653	DO i = 1, klon
654	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
655	END DO
656	END DO
657
658	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
659
660	END DO ! nbtr
661
662	END IF ! convection
663
664	!======================================================================
665	! -- Calcul de l'effet des thermiques --
666	!======================================================================
667
668	DO it=1,nbtr
669	DO k=1,klev
670	DO i=1,klon
671	d_tr_th(i,k,it)=0.
672	tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
673	tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
674	END DO
675	END DO
676	END DO
677
678	IF (iflag_thermals.GT.0.AND.iflag_the_trac>0) THEN
679
680	DO it=1, nbtr
681
682	CALL thermcell_dq(klon,klev,1,pdtphys,fm_therm,entr_therm, &
683	zmasse,tr_seri(1:klon,1:klev,it), &
684	d_tr_th(1:klon,1:klev,it),ztra_th,0 )
685
686	DO k=1,klev
687	DO i=1,klon
688	d_tr_th(i,k,it)=pdtphys*d_tr_th(i,k,it)
689	tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr_th(i,k,it),0.)
690	END DO
691	END DO
692
693	END DO ! it
694
695	END IF ! Thermiques
696
697	!======================================================================
698	! -- Calcul de l'effet de la couche limite --
699	!======================================================================
700
701	IF (iflag_vdf_trac==1) THEN
702
703	! Injection during BL mixing
704	!
705	DO it=1, nbtr
706	!
707	IF( pbl_flg(it) /= 0 ) THEN
708	!
709	CALL cltrac(pdtphys, coefh,t_seri, &
710	tr_seri(:,:,it), source(:,it), &
711	paprs, pplay, delp, &
712	d_tr_cl(:,:,it),d_tr_dry(:,it),flux_tr_dry(:,it))
713	!
714	tr_seri(:,:,it)=tr_seri(:,:,it)+d_tr_cl(:,:,it)
715	!
716	END IF
717	!
718	END DO
719	!
720	ELSE IF (iflag_vdf_trac==0) THEN
721	!
722	! Injection of source in the first model layer
723	!
724	DO it=1,nbtr
725	d_tr_cl(:,1,it)=source(:,it)RG/delp(:,1)pdtphys
726	tr_seri(:,1,it)=tr_seri(:,1,it)+d_tr_cl(:,1,it)
727	ENDDO
728	d_tr_cl(:,2:klev,1:nbtr)=0.
729	!
730	ELSE IF (iflag_vdf_trac==-1) THEN
731	!
732	! Nothing happens
733	!
734	d_tr_cl=0.
735	!
736	ELSE
737	!
738	CALL abort_physic('iflag_vdf_trac', 'cas non prevu',1)
739	!
740	END IF ! couche limite
741
742	!======================================================================
743	! Calcul de l'effet de la precipitation grande echelle
744	! POUR INCA le lessivage est fait directement dans INCA
745	!======================================================================
746
747	IF (lessivage) THEN
748
749	ql_incloud_ref = 10.e-4
750	ql_incloud_ref = 5.e-4
751
752
753	! calcul du contenu en eau liquide au sein du nuage
754	ql_incl = ql_incloud_ref
755	! choix du lessivage
756	!
757	IF (iflag_lscav .EQ. 3 .OR. iflag_lscav .EQ. 4) THEN
758	! ******** Olivier Boucher version (3) possibly with modified ql_incl (4)
759	!
760	DO it = 1, nbtr
761
762	IF (aerosol(it)) THEN
763	! incloud scavenging and removal by large scale rain ! orig : ql_incl was replaced by 0.5e-3 kg/kg
764	! the value 0.5e-3 kg/kg is from Giorgi and Chameides (1986), JGR
765	! Liu (2001) proposed to use 1.5e-3 kg/kg
766
767	!jyg<
768	!! CALL lsc_scav(pdtphys,it,iflag_lscav,ql_incl,prfl,psfl,rneb,beta_fisrt, &
769	CALL lsc_scav(pdtphys,it,iflag_lscav,aerosol,ql_incl,prfl,psfl,rneb,beta_fisrt, &
770	!>jyg
771	beta_v1,pplay,paprs,t_seri,tr_seri,d_tr_insc, &
772	d_tr_bcscav,d_tr_evapls,qPrls)
773
774	!large scale scavenging tendency
775	DO k = 1, klev
776	DO i = 1, klon
777	d_tr_ls(i,k,it)=d_tr_insc(i,k,it)+d_tr_bcscav(i,k,it)+d_tr_evapls(i,k,it)
778	tr_seri(i,k,it)=tr_seri(i,k,it)+d_tr_ls(i,k,it)
779	ENDDO
780	ENDDO
781	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'lsc scav it = '//solsym(it))
782	ENDIF
783
784	END DO !tr
785
786	ELSE IF (iflag_lscav .EQ. 2) THEN ! frac_impa, frac_nucl
787	! ********* modified old version
788
789	d_tr_lessi_nucl(:,:,:) = 0.
790	d_tr_lessi_impa(:,:,:) = 0.
791	flestottr(:,:,:) = 0.
792	! Tendance des aerosols nuclees et impactes
793	DO it = 1, nbtr
794	IF (aerosol(it)) THEN
795	his_dh(:)=0.
796	DO k = 1, klev
797	DO i = 1, klon
798	!PhH
799	zrho(i,k)=pplay(i,k)/t_seri(i,k)/RD
800	zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
801	!
802	END DO
803	END DO
804
805	DO k=klev-1, 1, -1
806	DO i=1, klon
807	! d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.)
808	dx=d_tr_ls(i,k,it)
809	his_dh(i)=his_dh(i)-dxzrho(i,k)zdz(i,k)/pdtphys ! kg/m2/s
810	evaplsc = prfl(i,k) - prfl(i,k+1) + psfl(i,k) - psfl(i,k+1)
811	! Evaporation Partielle -> Liberation Partielle 0.5*evap
812	IF ( evaplsc .LT.0..and.abs(prfl(i,k+1)+psfl(i,k+1)).gt.1.e-10) THEN
813	evaplsc = (-evaplsc)/(prfl(i,k+1)+psfl(i,k+1))
814	! evaplsc est donc positif, his_dh(i) est positif
815	!--------------
816	d_tr_evapls(i,k,it)=0.5evaplsc(d_tr_lessi_nucl(i,k+1,it) &
817	+d_tr_lessi_impa(i,k+1,it))
818	!------------- d_tr_evapls(i,k,it)=-0.5evaplsc(d_tr_lsc(i,k+1,it))
819	beta=0.5*evaplsc
820	if ((prfl(i,k)+psfl(i,k)).lt.1.e-10) THEN
821	beta=1.0*evaplsc
822	endif
823	dx=betahis_dh(i)/zrho(i,k)/zdz(i,k)pdtphys
824	his_dh(i)=(1.-beta)*his_dh(i) ! tracer from
825	d_tr_evapls(i,k,it)=dx
826	ENDIF
827	d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.) &
828	+d_tr_evapls(i,k,it)
829
830	!--------------
831	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
832	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
833	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
834	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
835	!
836	! Flux lessivage total
837	flestottr(i,k,it) = flestottr(i,k,it) - &
838	( d_tr_lessi_nucl(i,k,it) + &
839	d_tr_lessi_impa(i,k,it) ) * &
840	( paprs(i,k)-paprs(i,k+1) ) / &
841	(RG * pdtphys)
842	!! Mise a jour des traceurs due a l'impaction,nucleation
843	! tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
844	!! calcul de la tendance liee au lessivage stratiforme
845	! d_tr_ls(i,k,it)=tr_seri(i,k,it)*&
846	! (1.-1./(frac_impa(i,k)*frac_nucl(i,k)))
847	!--------------
848	END DO
849	END DO
850	END IF
851	END DO
852	! ********* end modified old version
853
854	ELSE IF (iflag_lscav .EQ. 1) THEN ! frac_impa, frac_nucl
855	! ********* old version
856
857	d_tr_lessi_nucl(:,:,:) = 0.
858	d_tr_lessi_impa(:,:,:) = 0.
859	flestottr(:,:,:) = 0.
860	!=========================
861	! LESSIVAGE LARGE SCALE :
862	!=========================
863
864	! Tendance des aerosols nuclees et impactes
865	! -----------------------------------------
866	DO it = 1, nbtr
867	IF (aerosol(it)) THEN
868	DO k = 1, klev
869	DO i = 1, klon
870	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
871	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
872	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
873	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
874
875	!
876	! Flux lessivage total
877	! ------------------------------------------------------------
878	flestottr(i,k,it) = flestottr(i,k,it) - &
879	( d_tr_lessi_nucl(i,k,it) + &
880	d_tr_lessi_impa(i,k,it) ) * &
881	( paprs(i,k)-paprs(i,k+1) ) / &
882	(RG * pdtphys)
883	!
884	! Mise a jour des traceurs due a l'impaction,nucleation
885	! ----------------------------------------------------------------------
886	tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
887	END DO
888	END DO
889	END IF
890	END DO
891
892	! ********* end old version
893	ENDIF ! iflag_lscav . EQ. 1, 2, 3 or 4
894	!
895	END IF ! lessivage
896
897	!=============================================================
898	! Ecriture des sorties
899	!=============================================================
900	#ifdef CPP_IOIPSL
901	! INCLUDE "write_histrac.h"
902	#endif
903
904	END SUBROUTINE phytrac
905
906	END MODULE

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

Download in other formats:

Original Format