Context Navigation

phytrac.F90 @ 1763

Last change on this file since 1763 was 1759, checked in by Ehouarn Millour, 12 years ago
IOIPSL routine getin is not threadsafe, so when running in OpenMP, it should be called by only one thread (and the result copied to other threads in the case of threadprivate variables). EM
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 25.8 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),save :: chtratimestep,chtratimestep_omp
206	!$OMP THREADPRIVATE(chtratimestep)
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_omp,iflag_lscav
267	LOGICAL,SAVE :: convscav_omp,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	!$OMP MASTER
312	chtratimestep_omp='DefFreq'
313	CALL getin('tra_time_step',chtratimestep_omp)
314	!$OMP END MASTER
315	!$OMP BARRIER
316	chtratimestep=chtratimestep_omp
317	IF (chtratimestep .NE. 'DefFreq') THEN
318	call convers_timesteps(chtratimestep,pdtphys,ecrit_tra)
319	ENDIF
320	!RomP >>>
321	!
322	!Config Key = convscav
323	!Config Desc = Convective scavenging switch: 0=off, 1=on.
324	!Config Def = .false.
325	!Config Help =
326	!
327	!$OMP MASTER
328	convscav_omp=.false.
329	call getin('convscav', convscav_omp)
330	!$OMP END MASTER
331	!$OMP BARRIER
332	convscav=convscav_omp
333	print*,'phytrac passage dans routine conv avec lessivage', convscav
334	!
335	!Config Key = iflag_lscav
336	!Config Desc = Large scale scavenging parametrization: 0=none, 1=old(Genthon92),
337	! 2=1+PHeinrich, 3=Reddy_Boucher2004, 4=3+RPilon.
338	!Config Def = 1
339	!Config Help =
340	!
341	!$OMP MASTER
342	iflag_lscav_omp=1
343	call getin('iflag_lscav', iflag_lscav_omp)
344	!$OMP END MASTER
345	!$OMP BARRIER
346	iflag_lscav=iflag_lscav_omp
347	!
348	SELECT CASE(iflag_lscav)
349	CASE(0)
350	PRINT*, 'Large scale scavenging: none'
351	CASE(1)
352	PRINT*, 'Large scale scavenging: C. Genthon, Tellus(1992), 44B, 371-389'
353	CASE(2)
354	PRINT*, 'Large scale scavenging: C. Genthon, modified P. Heinrich'
355	CASE(3)
356	PRINT*, 'Large scale scavenging: M. Shekkar Reddy and O. Boucher, JGR(2004), 109, D14202'
357	CASE(4)
358	PRINT*, 'Large scale scavenging: Reddy and Boucher, modified R. Pilon'
359	END SELECT
360	!RomP <<<
361	WRITE(,) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
362	ALLOCATE( source(klon,nbtr), stat=ierr)
363	IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 1',1)
364
365	ALLOCATE( aerosol(nbtr), stat=ierr)
366	IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 2',1)
367
368
369	! Initialize module for specific tracers
370	SELECT CASE(type_trac)
371	CASE('lmdz')
372	CALL traclmdz_init(pctsrf, xlat, xlon, ftsol, tr_seri, t_seri, pplay, sh, pdtphys, aerosol, lessivage)
373	CASE('inca')
374	source(:,:)=0.
375	CALL tracinca_init(aerosol,lessivage)
376	CASE('repr')
377	source(:,:)=0.
378	END SELECT
379	!
380	! Initialize diagnostic output
381	! ----------------------------
382	#ifdef CPP_IOIPSL
383	INCLUDE "ini_histrac.h"
384	#endif
385	END IF ! of IF (debutphy)
386	!############################################ END INITIALIZATION #######
387
388	DO k=1,klev
389	DO i=1,klon
390	zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
391	END DO
392	END DO
393	!
394	IF (id_be .GT. 0) THEN
395	DO k=1,klev
396	DO i=1,klon
397	sourceBE(i,k)=srcbe(i,k) !RomP -> pour sortie histrac
398	END DO
399	END DO
400	ENDIF
401
402	!===============================================================================
403	! -- Do specific treatment according to chemestry model or local LMDZ tracers
404	!
405	!===============================================================================
406	SELECT CASE(type_trac)
407	CASE('lmdz')
408	! -- Traitement des traceurs avec traclmdz
409	CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
410	cdragh, coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,couchelimite,sh, &
411	rh, pphi, ustar, u10m, v10m, &
412	!! tr_seri, source, solsym, d_tr_cl, zmasse) !RomP
413	tr_seri, source, solsym, d_tr_cl,d_tr_dec, zmasse) !RomP
414	CASE('inca')
415	! -- CHIMIE INCA config_inca = aero or chem --
416
417	CALL tracinca(&
418	nstep, julien, gmtime, lafin, &
419	pdtphys, t_seri, paprs, pplay, &
420	pmfu, ftsol, pctsrf, pphis, &
421	pphi, albsol, sh, rh, &
422	cldfra, rneb, diafra, cldliq, &
423	itop_con, ibas_con, pmflxr, pmflxs, &
424	prfl, psfl, aerosol_couple, flxmass_w, &
425	tau_aero, piz_aero, cg_aero, ccm, &
426	rfname, &
427	tr_seri, source, solsym)
428
429	CASE('repr')
430	! -- CHIMIE REPROBUS --
431
432	CALL tracreprobus(pdtphys, gmtime, debutphy, julien, &
433	presnivs, xlat, xlon, pphis, pphi, &
434	t_seri, pplay, paprs, sh , &
435	tr_seri, solsym)
436
437	END SELECT
438	!======================================================================
439	! -- Calcul de l'effet de la convection --
440	!======================================================================
441
442	IF (convection) THEN
443	DO it=1, nbtr
444	IF ( conv_flg(it) == 0 ) CYCLE
445	IF (iflag_con.LT.2) THEN
446	d_tr_cv(:,:,it)=0.
447	ELSE IF (iflag_con.EQ.2) THEN
448	!..Tiedke
449	CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
450	pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
451	! RomP >>>
452	ELSE
453	!..K.Emanuel !RomP modif arg
454	if (convscav.and.aerosol(it)) then ! lessivage convectif pour aerosol
455	!
456	CALL cvltr(pdtphys, da, phi,phi2,d1a,dam, mp,ep, &
457	sigd,sij,clw,elij,epmlmMm,eplaMm, &
458	pmflxr,pmflxs,evap,t_seri,wdtrainA,wdtrainM, &
459	paprs,it,tr_seri,upwd,dnwd,itop_con,ibas_con, &
460	d_tr_cv,d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qDi,qPr,&
461	qPa,qMel,qTrdi,dtrcvMA,Mint, &
462	zmfd1a,zmfphi2,zmfdam)
463	else !pas de lessivage convectif ou n'est pas un aerosol
464	CALL cvltrorig(it,pdtphys, da, phi,mp,paprs,pplay,tr_seri,&
465	upwd,dnwd,d_tr_cv)
466	endif
467	END IF
468	! RomP <<<
469
470	DO k = 1, klev
471	DO i = 1, klon
472	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
473	END DO
474	END DO
475
476	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
477
478	END DO ! nbtr
479	END IF ! convection
480
481	!======================================================================
482	! -- Calcul de l'effet des thermiques --
483	!======================================================================
484
485	DO it=1,nbtr
486	DO k=1,klev
487	DO i=1,klon
488	d_tr_th(i,k,it)=0.
489	tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
490	tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
491	END DO
492	END DO
493	END DO
494
495	IF (iflag_thermals.GT.0) THEN
496	nsplit=10
497	DO it=1, nbtr
498	DO isplit=1,nsplit
499
500	CALL dqthermcell(klon,klev,pdtphys/nsplit, &
501	fm_therm,entr_therm,zmasse, &
502	tr_seri(1:klon,1:klev,it),d_tr,ztra_th)
503
504	DO k=1,klev
505	DO i=1,klon
506	d_tr(i,k)=pdtphys*d_tr(i,k)/nsplit
507	d_tr_th(i,k,it)=d_tr_th(i,k,it)+d_tr(i,k)
508	tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr(i,k),0.)
509	END DO
510	END DO
511	END DO ! nsplit
512	END DO ! it
513	END IF ! Thermiques
514
515	!======================================================================
516	! -- Calcul de l'effet de la couche limite --
517	!======================================================================
518
519	IF (couchelimite) THEN
520
521	DO k = 1, klev
522	DO i = 1, klon
523	delp(i,k) = paprs(i,k)-paprs(i,k+1)
524	END DO
525	END DO
526
527	DO it=1, nbtr
528
529	IF( pbl_flg(it) /= 0 ) THEN
530
531	CALL cltrac(pdtphys, coefh,t_seri, &
532	tr_seri(:,:,it), source(:,it), &
533	paprs, pplay, delp, &
534	d_tr_cl(:,:,it),d_tr_dry(:,it),flux_tr_dry(:,it))
535
536	DO k = 1, klev
537	DO i = 1, klon
538	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cl(i,k,it)
539	END DO
540	END DO
541	END IF
542
543	END DO
544
545	END IF ! couche limite
546
547
548	!======================================================================
549	! Calcul de l'effet de la precipitation grande echelle
550	!======================================================================
551	IF (lessivage) THEN
552
553	ql_incloud_ref = 10.e-4
554	ql_incloud_ref = 5.e-4
555
556
557	! calcul du contenu en eau liquide au sein du nuage
558	ql_incl = ql_incloud_ref
559	! choix du lessivage
560	!
561	IF (iflag_lscav .EQ. 3 .OR. iflag_lscav .EQ. 4) THEN
562	! ******** Olivier Boucher version (3) possibly with modified ql_incl (4)
563	!
564	DO it = 1, nbtr
565	! incloud scavenging and removal by large scale rain ! orig : ql_incl was replaced by 0.5e-3 kg/kg
566	! the value 0.5e-3 kg/kg is from Giorgi and Chameides (1986), JGR
567	! Liu (2001) proposed to use 1.5e-3 kg/kg
568
569	CALL lsc_scav(pdtphys,it,iflag_lscav,ql_incl,prfl,psfl,rneb,beta_fisrt, &
570	beta_v1,pplay,paprs,t_seri,tr_seri,d_tr_insc, &
571	d_tr_bcscav,d_tr_evapls,qPrls)
572
573	!large scale scavenging tendency
574	DO k = 1, klev
575	DO i = 1, klon
576	d_tr_ls(i,k,it)=d_tr_insc(i,k,it)+d_tr_bcscav(i,k,it)+d_tr_evapls(i,k,it)
577	tr_seri(i,k,it)=tr_seri(i,k,it)+d_tr_ls(i,k,it)
578	ENDDO
579	ENDDO
580	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'lsc scav it = '//solsym(it))
581	END DO !tr
582
583	ELSE IF (iflag_lscav .EQ. 2) THEN ! frac_impa, frac_nucl
584	! ********* modified old version
585
586	d_tr_lessi_nucl(:,:,:) = 0.
587	d_tr_lessi_impa(:,:,:) = 0.
588	flestottr(:,:,:) = 0.
589	! Tendance des aerosols nuclees et impactes
590	DO it = 1, nbtr
591	IF (aerosol(it)) THEN
592	his_dh(:)=0.
593	DO k = 1, klev
594	DO i = 1, klon
595	!PhH
596	zrho(i,k)=pplay(i,k)/t_seri(i,k)/RD
597	zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
598	!
599	END DO
600	END DO
601
602	DO k=klev-1, 1, -1
603	DO i=1, klon
604	! d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.)
605	dx=d_tr_ls(i,k,it)
606	his_dh(i)=his_dh(i)-dxzrho(i,k)zdz(i,k)/pdtphys ! kg/m2/s
607	evaplsc = prfl(i,k) - prfl(i,k+1) + psfl(i,k) - psfl(i,k+1)
608	! Evaporation Partielle -> Liberation Partielle 0.5*evap
609	IF ( evaplsc .LT.0..and.abs(prfl(i,k+1)+psfl(i,k+1)).gt.1.e-10) THEN
610	evaplsc = (-evaplsc)/(prfl(i,k+1)+psfl(i,k+1))
611	! evaplsc est donc positif, his_dh(i) est positif
612	!--------------
613	d_tr_evapls(i,k,it)=0.5evaplsc(d_tr_lessi_nucl(i,k+1,it) &
614	+d_tr_lessi_impa(i,k+1,it))
615	!------------- d_tr_evapls(i,k,it)=-0.5evaplsc(d_tr_lsc(i,k+1,it))
616	beta=0.5*evaplsc
617	if ((prfl(i,k)+psfl(i,k)).lt.1.e-10) THEN
618	beta=1.0*evaplsc
619	endif
620	dx=betahis_dh(i)/zrho(i,k)/zdz(i,k)pdtphys
621	his_dh(i)=(1.-beta)*his_dh(i) ! tracer from
622	d_tr_evapls(i,k,it)=dx
623	ENDIF
624	d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.) &
625	+d_tr_evapls(i,k,it)
626
627	!--------------
628	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
629	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
630	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
631	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
632	!
633	! Flux lessivage total
634	flestottr(i,k,it) = flestottr(i,k,it) - &
635	( d_tr_lessi_nucl(i,k,it) + &
636	d_tr_lessi_impa(i,k,it) ) * &
637	( paprs(i,k)-paprs(i,k+1) ) / &
638	(RG * pdtphys)
639	!! Mise a jour des traceurs due a l'impaction,nucleation
640	! tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
641	!! calcul de la tendance liee au lessivage stratiforme
642	! d_tr_ls(i,k,it)=tr_seri(i,k,it)*&
643	! (1.-1./(frac_impa(i,k)*frac_nucl(i,k)))
644	!--------------
645	END DO
646	END DO
647	END IF
648	END DO
649	! ********* end modified old version
650
651	ELSE IF (iflag_lscav .EQ. 1) THEN ! frac_impa, frac_nucl
652	! ********* old version
653
654	d_tr_lessi_nucl(:,:,:) = 0.
655	d_tr_lessi_impa(:,:,:) = 0.
656	flestottr(:,:,:) = 0.
657	!=========================
658	! LESSIVAGE LARGE SCALE :
659	!=========================
660
661	! Tendance des aerosols nuclees et impactes
662	! -----------------------------------------
663	DO it = 1, nbtr
664	IF (aerosol(it)) THEN
665	DO k = 1, klev
666	DO i = 1, klon
667	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
668	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
669	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
670	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
671
672	!
673	! Flux lessivage total
674	! ------------------------------------------------------------
675	flestottr(i,k,it) = flestottr(i,k,it) - &
676	( d_tr_lessi_nucl(i,k,it) + &
677	d_tr_lessi_impa(i,k,it) ) * &
678	( paprs(i,k)-paprs(i,k+1) ) / &
679	(RG * pdtphys)
680	!
681	! Mise a jour des traceurs due a l'impaction,nucleation
682	! ----------------------------------------------------------------------
683	tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
684	END DO
685	END DO
686	END IF
687	END DO
688
689	! ********* end old version
690	ENDIF ! iflag_lscav . EQ. 1, 2, 3 or 4
691	!
692	END IF ! lessivage
693
694	!=============================================================
695	! Ecriture des sorties
696	!=============================================================
697	#ifdef CPP_IOIPSL
698	INCLUDE "write_histrac.h"
699	#endif
700
701	END SUBROUTINE phytrac

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

Download in other formats:

Original Format