Context Navigation

phytrac.F90 @ 4138

Last change on this file since 4138 was 1379, checked in by lguez, 15 years ago

Added optional ozone tracer with chemistry parameterized by Daniel
Cariolle. This tracer is passive: it has no influence on the rest of
the simulation.

Added variable "zmasse" in file "histrac.nc". Corrected long name of
variable "pplay" in "histrac.nc". Changed name of variable "t" to "T"
in "histrac.nc", corrected long name and unit.

In "phytrac", moved definition of "zmasse" toward the beginning of the
procedure, so that "zmasse" can be given as argument to
"traclmdz". Also added arguments "julien", "gmtime" and "xlon" to
"traclmdz". The four additional arguments are required for the ozone
tracer.

In module "traclmdz_mod", factorized declaration "implicit none" that
was in each procedure. There is now an equivalent single declaration
at the module level.

In procedure "traclmdz", removed variable "delp". Use "zmasse * rg"
instead since we now have "zmasse" as an argument.

Tests. Compilations on Brodie only, with optimization options "-debug"
and "-dev", parallelization options "none", "mpi", "omp" and
"mpi_omp", this revision and revision 1373. Run cases with and without
ozone tracer, 1 and 2 MPI processes, 1 and 2 OpenMP
threads. Comparisons of all cases are ok, except for strange
variations in variables "d_tr_cl_RN" and "d_tr_cl_PB" of file
"histrac.nc", variables "RN" and "PB" of "restart.nc", variables
"trs_RN" and "trs_PB" of "restartphy.nc". Relative variations of these
variables between cases are of order 1e-7 or less, after one day of
simulation.

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

File size: 15.2 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	xlat, frac_impa,frac_nucl,xlon, &
11	presnivs, pphis, pphi, albsol, &
12	sh, rh, cldfra, rneb, &
13	diafra, cldliq, itop_con, ibas_con, &
14	pmflxr, pmflxs, prfl, psfl, &
15	da, phi, mp, upwd, &
16	dnwd, aerosol_couple, flxmass_w, &
17	tau_aero, piz_aero, cg_aero, ccm, &
18	rfname, &
19	tr_seri)
20	!
21	!======================================================================
22	! Auteur(s) FH
23	! Objet: Moniteur general des tendances traceurs
24	!======================================================================
25
26	USE ioipsl
27	USE dimphy
28	USE infotrac
29	USE mod_grid_phy_lmdz
30	USE mod_phys_lmdz_para
31	USE comgeomphy
32	USE iophy
33	USE traclmdz_mod
34	USE tracinca_mod
35	USE control_mod
36
37
38
39	IMPLICIT NONE
40
41	INCLUDE "YOMCST.h"
42	INCLUDE "dimensions.h"
43	INCLUDE "indicesol.h"
44	INCLUDE "clesphys.h"
45	INCLUDE "temps.h"
46	INCLUDE "paramet.h"
47	INCLUDE "thermcell.h"
48	!==========================================================================
49	! -- ARGUMENT DESCRIPTION --
50	!==========================================================================
51
52	! Input arguments
53	!----------------
54	!Configuration grille,temps:
55	INTEGER,INTENT(IN) :: nstep ! Appel physique
56	INTEGER,INTENT(IN) :: julien ! Jour julien
57	REAL,INTENT(IN) :: gmtime
58	REAL,INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
59	LOGICAL,INTENT(IN) :: debutphy ! le flag de l'initialisation de la physique
60	LOGICAL,INTENT(IN) :: lafin ! le flag de la fin de la physique
61
62	REAL,DIMENSION(klon),INTENT(IN) :: xlat ! latitudes pour chaque point
63	REAL,DIMENSION(klon),INTENT(IN) :: xlon ! longitudes pour chaque point
64	!
65	!Physique:
66	!--------
67	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
68	REAL,DIMENSION(klon,klev),INTENT(IN) :: u !
69	REAL,DIMENSION(klon,klev),INTENT(IN) :: v !
70	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
71	REAL,DIMENSION(klon,klev),INTENT(IN) :: rh ! humidite relative
72	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
73	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
74	REAL,DIMENSION(klon,klev),INTENT(IN) :: pphi ! geopotentiel
75	REAL,DIMENSION(klon),INTENT(IN) :: pphis
76	REAL,DIMENSION(klev),INTENT(IN) :: presnivs
77	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldliq ! eau liquide nuageuse
78	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldfra ! fraction nuageuse (tous les nuages)
79	REAL,DIMENSION(klon,klev),INTENT(IN) :: diafra ! fraction nuageuse (convection ou stratus artificiels)
80	REAL,DIMENSION(klon,klev),INTENT(IN) :: rneb ! fraction nuageuse (grande echelle)
81	INTEGER,DIMENSION(klon),INTENT(IN) :: itop_con
82	INTEGER,DIMENSION(klon),INTENT(IN) :: ibas_con
83	REAL,DIMENSION(klon),INTENT(IN) :: albsol ! albedo surface
84	!
85	!Convection:
86	!----------
87	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu ! flux de masse dans le panache montant
88	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd ! flux de masse dans le panache descendant
89	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
90	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
91	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
92	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
93
94	!...Tiedke
95	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
96	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
97
98	LOGICAL,INTENT(IN) :: aerosol_couple
99	REAL,DIMENSION(klon,klev),INTENT(IN) :: flxmass_w
100	REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: tau_aero
101	REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: piz_aero
102	REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: cg_aero
103	CHARACTER(len=4),DIMENSION(9),INTENT(IN) :: rfname
104	REAL,DIMENSION(klon,klev,2),INTENT(IN) :: ccm
105	!... K.Emanuel
106	REAL,DIMENSION(klon,klev),INTENT(IN) :: da
107	REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
108	REAL,DIMENSION(klon,klev),INTENT(IN) :: mp
109	REAL,DIMENSION(klon,klev),INTENT(IN) :: upwd ! saturated updraft mass flux
110	REAL,DIMENSION(klon,klev),INTENT(IN) :: dnwd ! saturated downdraft mass flux
111	!
112	!Thermiques:
113	!----------
114	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: fm_therm
115	REAL,DIMENSION(klon,klev),INTENT(IN) :: entr_therm
116	!
117	!Couche limite:
118	!--------------
119	!
120	REAL,DIMENSION(klon,klev),INTENT(IN) :: cdragh ! coeff drag pour T et Q
121	REAL,DIMENSION(klon,klev),INTENT(IN) :: coefh ! coeff melange CL (m**2/s)
122	REAL,DIMENSION(klon),INTENT(IN) :: yu1 ! vents au premier niveau
123	REAL,DIMENSION(klon),INTENT(IN) :: yv1 ! vents au premier niveau
124	!
125	!Lessivage:
126	!----------
127	!
128	! pour le ON-LINE
129	!
130	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
131	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
132
133	! Arguments necessaires pour les sources et puits de traceur:
134	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
135	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
136
137
138	! Output argument
139	!----------------
140	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
141
142	!=======================================================================================
143	! -- LOCAL VARIABLES --
144	!=======================================================================================
145
146	INTEGER :: i, k, it
147	INTEGER :: nsplit
148
149	!Sources et Reservoirs de traceurs (ex:Radon):
150	!--------------------------------------------
151	!
152	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source ! a voir lorsque le flux de surface est prescrit
153	!$OMP THREADPRIVATE(source)
154
155	!
156	!Entrees/Sorties: (cf ini_histrac.h et write_histrac.h)
157	!---------------
158	INTEGER :: iiq, ierr
159	INTEGER :: nhori, nvert
160	REAL :: zsto, zout, zjulian
161	INTEGER,SAVE :: nid_tra ! pointe vers le fichier histrac.nc
162	!$OMP THREADPRIVATE(nid_tra)
163	REAL,DIMENSION(klon) :: zx_tmp_fi2d ! variable temporaire grille physique
164	INTEGER :: itau_w ! pas de temps ecriture = nstep + itau_phy
165	LOGICAL,PARAMETER :: ok_sync=.TRUE.
166
167	!
168	! Nature du traceur
169	!------------------
170	LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol ! aerosol(it) = true => aerosol => lessivage
171	!$OMP THREADPRIVATE(aerosol) ! aerosol(it) = false => gaz
172	REAL,DIMENSION(klon,klev) :: delp ! epaisseur de couche (Pa)
173	!
174	! Tendances de traceurs (Td):
175	!------------------------
176	!
177	REAL,DIMENSION(klon,klev) :: d_tr ! Td dans l'atmosphere
178	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_cl ! Td couche limite/traceur
179	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_cv ! Td convection/traceur
180	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_th ! Td thermique
181	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_lessi_impa ! Td du lessivage par impaction
182	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_lessi_nucl ! Td du lessivage par nucleation
183	!
184	! Physique
185	!----------
186	REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
187	REAL,DIMENSION(klon,klev) :: zmasse ! densité atmosphérique Kg/m2
188	REAL,DIMENSION(klon,klev) :: ztra_th
189
190	!Controles:
191	!---------
192	LOGICAL,SAVE :: couchelimite=.TRUE.
193	LOGICAL,SAVE :: convection=.TRUE.
194	LOGICAL,SAVE :: lessivage
195	!$OMP THREADPRIVATE(couchelimite,convection,lessivage)
196
197	CHARACTER(len=8),DIMENSION(nbtr) :: solsym
198
199
200	!######################################################################
201	! -- INITIALIZATION --
202	!######################################################################
203	IF (debutphy) THEN
204	WRITE(,) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
205	ALLOCATE( source(klon,nbtr), stat=ierr)
206	IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 1',1)
207
208	ALLOCATE( aerosol(nbtr), stat=ierr)
209	IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 2',1)
210
211
212	! Initialize module for specific tracers
213	SELECT CASE(type_trac)
214	CASE('lmdz')
215	!IM ajout t_seri, pplay, sh CALL traclmdz_init(pctsrf, ftsol, tr_seri, aerosol, lessivage)
216	CALL traclmdz_init(pctsrf, ftsol, tr_seri, t_seri, pplay, sh, aerosol, lessivage)
217	CASE('inca')
218	source(:,:)=0.
219	CALL tracinca_init(aerosol,lessivage)
220	END SELECT
221	!
222	! Initialize diagnostic output
223	! ----------------------------
224	#ifdef CPP_IOIPSL
225	INCLUDE "ini_histrac.h"
226	#endif
227	END IF
228	!############################################ END INITIALIZATION #######
229
230	DO k=1,klev
231	DO i=1,klon
232	zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
233	END DO
234	END DO
235
236	!===============================================================================
237	! -- Do specific treatment according to chemestry model or local LMDZ tracers
238	!
239	!===============================================================================
240	SELECT CASE(type_trac)
241	CASE('lmdz')
242	! -- Traitement des traceurs avec traclmdz
243
244	CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
245	cdragh, coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon, couchelimite, &
246	sh, tr_seri, source, solsym, d_tr_cl, zmasse)
247
248	CASE('inca')
249	! -- CHIMIE INCA config_inca = aero or chem --
250
251	CALL tracinca(&
252	nstep, julien, gmtime, lafin, &
253	pdtphys, t_seri, paprs, pplay, &
254	pmfu, ftsol, pctsrf, pphis, &
255	pphi, albsol, sh, rh, &
256	cldfra, rneb, diafra, cldliq, &
257	itop_con, ibas_con, pmflxr, pmflxs, &
258	prfl, psfl, aerosol_couple, flxmass_w, &
259	tau_aero, piz_aero, cg_aero, ccm, &
260	rfname, &
261	tr_seri, source, solsym)
262	END SELECT
263
264	!======================================================================
265	! -- Calcul de l'effet de la convection --
266	!======================================================================
267	IF (convection) THEN
268	DO it=1, nbtr
269	IF ( conv_flg(it) == 0 ) CYCLE
270
271	IF (iflag_con.LT.2) THEN
272	d_tr_cv(:,:,:)=0.
273	ELSE IF (iflag_con.EQ.2) THEN
274	!..Tiedke
275	CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
276	pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
277	ELSE
278	!..K.Emanuel
279	CALL cvltr(pdtphys, da, phi, mp, paprs,pplay, tr_seri(:,:,it),&
280	upwd,dnwd,d_tr_cv(:,:,it))
281	END IF
282
283	!IM ajout traceurs RR
284	! print*,'phytrac it,nseuil=',it,nseuil
285	IF (it.lt.nseuil) THEN
286	DO k = 1, klev
287	DO i = 1, klon
288	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
289	END DO
290	END DO
291	END IF !(it.lt.nseuil) then
292
293	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
294
295	END DO ! nbtr
296	END IF ! convection
297
298	!======================================================================
299	! -- Calcul de l'effet des thermiques --
300	!======================================================================
301
302	DO it=1,nbtr
303	DO k=1,klev
304	DO i=1,klon
305	d_tr_th(i,k,it)=0.
306	tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
307	tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
308	END DO
309	END DO
310	END DO
311
312	IF (iflag_thermals.GT.0) THEN
313	nsplit=10
314	DO it=1, nbtr
315	DO isplit=1,nsplit
316
317	CALL dqthermcell(klon,klev,pdtphys/nsplit, &
318	fm_therm,entr_therm,zmasse, &
319	tr_seri(1:klon,1:klev,it),d_tr,ztra_th)
320
321	DO k=1,klev
322	DO i=1,klon
323	d_tr(i,k)=pdtphys*d_tr(i,k)/nsplit
324	d_tr_th(i,k,it)=d_tr_th(i,k,it)+d_tr(i,k)
325	tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr(i,k),0.)
326	END DO
327	END DO
328	END DO ! nsplit
329	END DO ! it
330	END IF ! Thermiques
331
332	!======================================================================
333	! -- Calcul de l'effet de la couche limite --
334	!======================================================================
335
336	IF (couchelimite) THEN
337
338	DO k = 1, klev
339	DO i = 1, klon
340	delp(i,k) = paprs(i,k)-paprs(i,k+1)
341	END DO
342	END DO
343
344	DO it=1, nbtr
345
346	IF( pbl_flg(it) /= 0 ) THEN
347
348	CALL cltrac(pdtphys, coefh,t_seri, &
349	tr_seri(:,:,it), source(:,it), &
350	paprs, pplay, delp, &
351	d_tr_cl(:,:,it))
352
353	DO k = 1, klev
354	DO i = 1, klon
355	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cl(i,k,it)
356	END DO
357	END DO
358	END IF
359
360	END DO
361
362	END IF ! couche limite
363
364
365	!======================================================================
366	! Calcul de l'effet de la precipitation
367	!======================================================================
368
369	IF (lessivage) THEN
370
371	d_tr_lessi_nucl(:,:,:) = 0.
372	d_tr_lessi_impa(:,:,:) = 0.
373	flestottr(:,:,:) = 0.
374	!=========================
375	! LESSIVAGE LARGE SCALE :
376	!=========================
377
378	! Tendance des aerosols nuclees et impactes
379	! -----------------------------------------
380	DO it = 1, nbtr
381	IF (aerosol(it)) THEN
382	DO k = 1, klev
383	DO i = 1, klon
384	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
385	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
386	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
387	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
388
389	!
390	! Flux lessivage total
391	! ------------------------------------------------------------
392	flestottr(i,k,it) = flestottr(i,k,it) - &
393	( d_tr_lessi_nucl(i,k,it) + &
394	d_tr_lessi_impa(i,k,it) ) * &
395	( paprs(i,k)-paprs(i,k+1) ) / &
396	(RG * pdtphys)
397	!
398	! Mise a jour des traceurs due a l'impaction,nucleation
399	! ----------------------------------------------------------------------
400	tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
401	END DO
402	END DO
403	END IF
404	END DO
405
406	END IF ! lessivage
407
408	!=============================================================
409	! Ecriture des sorties
410	!=============================================================
411	#ifdef CPP_IOIPSL
412	INCLUDE "write_histrac.h"
413	#endif
414
415	END SUBROUTINE phytrac

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

Download in other formats:

Original Format