Context Navigation

phytrac.F90 @ 1742

Last change on this file since 1742 was 1742, checked in by idelkadi, 11 years ago

1- Inclusion des developpements de la these de Romain Pilon sur le
lessivage des aerosols :

a/ par les pluies convectives (modifs cv30_routines et cv3_routines pour

sortir les champs nécessaires au calcul off-line ; modif cvltr)

b/ par les pluies stratiformes (modifs phytrac et introduction

lsc_scav).

2- Choix entre plusieurs schemas pour les pluies stratiformes, commande
par iflag_lscav.

3- Quelques corrections dans la convection "Nouvelle Physique" pour
assurer la conservation des traceurs (cv3p1_mixing et cva_driver) (travail
de Robin Locatelli).

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 25.5 KB

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

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

Download in other formats:

Original Format