Context Navigation

traclmdz_mod.F90 @ 5136

Last change on this file since 5136 was 5117, checked in by abarral, 2 months ago
rename modules properly lmdz_* move some unused files to obsolete/ (lint) uppercase fortran keywords
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: 22.1 KB

Line
1	!$Id $
2
3	MODULE traclmdz_mod
4
5	! In this module all tracers specific to LMDZ are treated. This module is used
6	! only if running without any other chemestry model as INCA or REPROBUS.
7	USE lmdz_abort_physic, ONLY: abort_physic
8	IMPLICIT NONE
9
10	REAL, DIMENSION(:, :), ALLOCATABLE, SAVE :: masktr ! Masque reservoir de sol traceur
11	!$OMP THREADPRIVATE(masktr) ! Masque de l'echange avec la surface (1 = reservoir) ou (possible >= 1 )
12	REAL, DIMENSION(:, :), ALLOCATABLE, SAVE :: fshtr ! Flux surfacique dans le reservoir de sol
13	!$OMP THREADPRIVATE(fshtr)
14	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: hsoltr ! Epaisseur equivalente du reservoir de sol
15	!$OMP THREADPRIVATE(hsoltr)
16
17	!Radioelements:
18	!--------------
19
20	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: tautr ! Constante de decroissance radioactive
21	!$OMP THREADPRIVATE(tautr)
22	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: vdeptr ! Vitesse de depot sec dans la couche Brownienne
23	!$OMP THREADPRIVATE(vdeptr)
24	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: scavtr ! Coefficient de lessivage
25	!$OMP THREADPRIVATE(scavtr)
26	REAL, DIMENSION(:, :), ALLOCATABLE, SAVE :: srcbe ! Production du beryllium7 dans l atmosphere (U/s/kgA)
27	!$OMP THREADPRIVATE(srcbe)
28
29	LOGICAL, DIMENSION(:), ALLOCATABLE, SAVE :: radio ! radio(it) = true => decroisssance radioactive
30	!$OMP THREADPRIVATE(radio)
31
32	REAL, DIMENSION(:, :), ALLOCATABLE, SAVE :: trs ! Conc. radon ds le sol
33	!$OMP THREADPRIVATE(trs)
34
35	INTEGER, SAVE :: id_aga ! Identification number for tracer : Age of stratospheric air
36	!$OMP THREADPRIVATE(id_aga)
37	INTEGER, SAVE :: lev_1p5km ! Approximative vertical layer number at 1.5km above surface, used for calculation of the age of air. The result shouldn't be that sensible to the exactness of this value as long as it is in the lower troposphere.
38	!$OMP THREADPRIVATE(lev_1p5km)
39
40	INTEGER, SAVE :: id_rn, id_pb ! Identification number for tracer : radon (Rn222), lead (Pb210)
41	!$OMP THREADPRIVATE(id_rn, id_pb)
42
43	INTEGER, SAVE :: id_be ! Activation et position du traceur Be7 [ id_be=0 -> desactive ]
44	!$OMP THREADPRIVATE(id_be)
45
46	INTEGER, SAVE :: id_pcsat, id_pcocsat, id_pcq ! traceurs pseudo-vapeur CL qsat, qsat_oc, q
47	!$OMP THREADPRIVATE(id_pcsat, id_pcocsat, id_pcq)
48	INTEGER, SAVE :: id_pcs0, id_pcos0, id_pcq0 ! traceurs pseudo-vapeur CL qsat, qsat_oc, q
49	! ! qui ne sont pas transportes par la convection
50	!$OMP THREADPRIVATE(id_pcs0, id_pcos0, id_pcq0)
51
52	INTEGER, SAVE :: id_o3
53	!$OMP THREADPRIVATE(id_o3)
54	! index of ozone tracer with Cariolle parameterization
55	! 0 means no ozone tracer
56
57	LOGICAL, SAVE :: rnpb = .FALSE. ! Presence du couple Rn222, Pb210
58	!$OMP THREADPRIVATE(rnpb)
59
60
61	CONTAINS
62
63	SUBROUTINE traclmdz_from_restart(trs_in)
64	! This SUBROUTINE initialize the module saved variable trs with values from restart file (startphy.nc).
65	! This SUBROUTINE is called from phyetat0 after the field trs_in has been read.
66
67	USE dimphy
68	USE infotrac_phy, ONLY: nbtr
69
70	! Input argument
71	REAL, DIMENSION(klon, nbtr), INTENT(IN) :: trs_in
72
73	! Local variables
74	INTEGER :: ierr
75
76	! Allocate restart variables trs
77	ALLOCATE(trs(klon, nbtr), stat = ierr)
78	IF (ierr /= 0) CALL abort_physic('traclmdz_from_restart', 'pb in allocation 1', 1)
79
80	! Initialize trs with values read from restart file
81	trs(:, :) = trs_in(:, :)
82
83	END SUBROUTINE traclmdz_from_restart
84
85
86	SUBROUTINE traclmdz_init(pctsrf, xlat, xlon, ftsol, tr_seri, t_seri, pplay, sh, pdtphys, aerosol, lessivage)
87	! This SUBROUTINE allocates and initialize module variables and control variables.
88	! Initialization of the tracers should be done here only for those not found in the restart file.
89	USE dimphy
90	USE infotrac_phy, ONLY: nbtr, nqtot, tracers, pbl_flg, conv_flg
91	USE regr_pr_comb_coefoz_m, ONLY: alloc_coefoz
92	USE press_coefoz_m, ONLY: press_coefoz
93	USE lmdz_grid_phy
94	USE lmdz_phys_para
95	USE indice_sol_mod
96	USE lmdz_print_control, ONLY: lunout
97
98	! Input variables
99	REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: pctsrf ! Pourcentage de sol f(nature du sol)
100	REAL, DIMENSION(klon), INTENT(IN) :: xlat ! latitudes en degres pour chaque point
101	REAL, DIMENSION(klon), INTENT(IN) :: xlon ! longitudes en degres pour chaque point
102	REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
103	REAL, DIMENSION(klon, klev, nbtr), INTENT(INOUT) :: tr_seri! Concentration Traceur [U/KgA]
104	REAL, DIMENSION(klon, klev), INTENT(IN) :: t_seri ! Temperature
105	REAL, DIMENSION(klon, klev), INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
106	REAL, DIMENSION(klon, klev), INTENT(IN) :: sh ! humidite specifique
107	REAL, INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
108
109	! Output variables
110	LOGICAL, DIMENSION(nbtr), INTENT(OUT) :: aerosol
111	LOGICAL, INTENT(OUT) :: lessivage
112
113	! Local variables
114	INTEGER :: ierr, it, iq, i, k
115	REAL, DIMENSION(klon_glo, klev) :: varglo ! variable temporaire sur la grille global
116	REAL, DIMENSION(klev) :: mintmp, maxtmp
117	LOGICAL :: zero
118	! RomP >>> profil initial Be7
119	INTEGER ilesfil
120	parameter (ilesfil = 1)
121	INTEGER irr, kradio
122	real beryllium(klon, klev)
123	! profil initial Pb210
124	INTEGER ilesfil2
125	parameter (ilesfil2 = 1)
126	INTEGER irr2, kradio2
127	real plomb(klon, klev)
128	!! RomP <<<
129	! --------------------------------------------
130	! Allocation
131	! --------------------------------------------
132	ALLOCATE(scavtr(nbtr), stat = ierr)
133	IF (ierr /= 0) CALL abort_physic('traclmdz_init', 'pb in allocation 9', 1)
134	scavtr(:) = 1.
135
136	ALLOCATE(radio(nbtr), stat = ierr)
137	IF (ierr /= 0) CALL abort_physic('traclmdz_init', 'pb in allocation 11', 1)
138	radio(:) = .FALSE. ! Par defaut pas decroissance radioactive
139
140	ALLOCATE(masktr(klon, nbtr), stat = ierr)
141	IF (ierr /= 0) CALL abort_physic('traclmdz_init', 'pb in allocation 2', 1)
142
143	ALLOCATE(fshtr(klon, nbtr), stat = ierr)
144	IF (ierr /= 0) CALL abort_physic('traclmdz_init', 'pb in allocation 3', 1)
145
146	ALLOCATE(hsoltr(nbtr), stat = ierr)
147	IF (ierr /= 0) CALL abort_physic('traclmdz_init', 'pb in allocation 4', 1)
148
149	ALLOCATE(tautr(nbtr), stat = ierr)
150	IF (ierr /= 0) CALL abort_physic('traclmdz_init', 'pb in allocation 5', 1)
151	tautr(:) = 0.
152
153	ALLOCATE(vdeptr(nbtr), stat = ierr)
154	IF (ierr /= 0) CALL abort_physic('traclmdz_init', 'pb in allocation 6', 1)
155	vdeptr(:) = 0.
156
157	lessivage = .TRUE.
158	!!jyg(20130206) : le choix d activation du lessivage est fait dans phytrac avec iflag_lscav
159	!! CALL getin('lessivage',lessivage)
160	!! IF(lessivage) THEN
161	!! PRINT*,'lessivage lsc ON'
162	!! else
163	!! PRINT*,'lessivage lsc OFF'
164	!! endif
165	aerosol(:) = .FALSE. ! Tous les traceurs sont des gaz par defaut
166
167	! Recherche des traceurs connus : Be7, O3, CO2,...
168	! --------------------------------------------
169	id_rn = 0; id_pb = 0; id_aga = 0; id_be = 0; id_o3 = 0
170	id_pcsat = 0; id_pcocsat = 0; id_pcq = 0; id_pcs0 = 0; id_pcos0 = 0; id_pcq0 = 0
171	it = 0
172	DO iq = 1, nqtot
173	IF(.NOT.(tracers(iq)%isInPhysics)) CYCLE
174	it = it + 1
175	SELECT CASE(tracers(iq)%name)
176	CASE("RN"); id_rn = it ! radon
177	CASE("PB"); id_pb = it ! plomb
178	CASE("Aga", "AGA"); id_aga = it ! Age of stratospheric air
179	CASE("Be", "BE", "Be7", "BE7"); id_be = it ! Recherche du Beryllium 7
180	CASE("o3", "O3"); id_o3 = it ! Recherche de l'ozone
181	CASE("pcsat", "Pcsat"); id_pcsat = it
182	CASE("pcocsat", "Pcocsat"); id_pcocsat = it
183	CASE("pcq", "Pcq"); id_pcq = it
184	CASE("pcs0", "Pcs0"); id_pcs0 = it
185	CASE("pcos0", "Pcos0"); id_pcos0 = it
186	CASE("pcq0", "Pcq0"); id_pcq0 = it
187	CASE DEFAULT
188	WRITE(lunout, *) 'This is an unknown tracer in LMDZ : ', trim(tracers(iq)%name)
189	END SELECT
190
191	SELECT CASE(tracers(iq)%name)
192	CASE("PB") !--- RomP >>> profil initial de PB210
193	OPEN(ilesfil2, file = 'prof.pb210', status = 'old', iostat = irr2)
194	IF(irr2 == 0) THEN
195	READ(ilesfil2, *) kradio2
196	WRITE(lunout, *)'number of levels for pb210 profile ', kradio2
197	DO k = kradio2, 1, -1; READ (ilesfil2, *) plomb(:, k);
198	END DO
199	CLOSE(ilesfil2)
200	tr_seri(:, :, id_pb) = plomb(:, :)
201	ELSE
202	WRITE(lunout, *)'prof.pb210 does not exist: use restart values'
203	END IF
204	CASE("Aga", "AGA")
205	radio(id_aga) = .FALSE.
206	aerosol(id_aga) = .FALSE.
207	pbl_flg(id_aga) = 0
208	! Find the first model layer above 1.5km from the surface
209	IF (klev>=30) THEN
210	lev_1p5km = 6 !--- NB: This value is for klev=39
211	ELSE IF (klev>=10) THEN
212	lev_1p5km = 5 !--- NB: This value is for klev=19
213	ELSE
214	lev_1p5km = klev / 2
215	END IF
216	CASE("Be", "BE", "Be7", "BE7")
217	ALLOCATE(srcbe(klon, klev))
218	radio(id_be) = .TRUE.
219	aerosol(id_be) = .TRUE. !--- Le Be est un aerosol
220	CALL init_be(pctsrf, pplay, masktr(:, id_be), tautr(id_be), vdeptr(id_be), scavtr(id_be), srcbe)
221	WRITE(lunout, *) 'Initialisation srcBe: OK'
222	!--- RomP >>> profil initial de Be7
223	OPEN(ilesfil, file = 'prof.be7', status = 'old', iostat = irr)
224	IF(irr == 0) THEN
225	READ(ilesfil, *) kradio
226	WRITE(lunout, *)'number of levels for Be7 profile ', kradio
227	DO k = kradio, 1, -1; READ(ilesfil, *) beryllium(:, k);
228	END DO
229	CLOSE(ilesfil)
230	tr_seri(:, :, id_be) = beryllium(:, :)
231	ELSE
232	WRITE(lunout, *)'Prof. Be7 does not exist: use restart values'
233	END IF
234	CASE("o3", "O3") !--- Parametrisation par la chimie de Cariolle
235	CALL alloc_coefoz !--- Allocate ozone coefficients
236	CALL press_coefoz !--- Read input pressure levels
237	CASE("pcs0", "Pcs0", "pcos0", "Pcos0", "pcq0", "Pcq0")
238	conv_flg(it) = 0 !--- No transport by convection for this tracer
239	END SELECT
240	END DO
241
242	! Valeurs specifiques pour les traceurs Rn222 et Pb210
243	! ----------------------------------------------
244	IF (id_rn/=0 .AND. id_pb/=0) THEN
245	rnpb = .TRUE.
246	radio(id_rn) = .TRUE.
247	radio(id_pb) = .TRUE.
248	pbl_flg(id_rn) = 0 ! au lieu de clsol=true ! CL au sol calcule
249	pbl_flg(id_pb) = 0 ! au lieu de clsol=true
250	aerosol(id_rn) = .FALSE.
251	aerosol(id_pb) = .TRUE. ! le Pb est un aerosol
252
253	CALL initrrnpb (ftsol, pctsrf, masktr, fshtr, hsoltr, tautr, vdeptr, scavtr)
254	END IF
255
256	! Check if all tracers have restart values
257	! ----------------------------------------------
258	it = 0
259	DO iq = 1, nqtot
260	IF(.NOT.(tracers(iq)%isAdvected .AND. tracers(iq)%isInPhysics)) CYCLE
261	it = it + 1
262	! Test if tracer is zero everywhere.
263	! Done by master process MPI and master thread OpenMP
264	CALL gather(tr_seri(:, :, it), varglo)
265	IF (is_mpi_root .AND. is_omp_root) THEN
266	mintmp = MINVAL(varglo, dim = 1)
267	maxtmp = MAXVAL(varglo, dim = 1)
268	IF (MINVAL(mintmp, dim = 1)==0. .AND. MAXVAL(maxtmp, dim = 1)==0.) THEN
269	! Tracer is zero everywhere
270	zero = .TRUE.
271	ELSE
272	zero = .FALSE.
273	END IF
274	END IF
275
276	! Distribute variable at all processes
277	CALL bcast(zero)
278
279	! Initalize tracer that was not found in restart file.
280	IF (zero) THEN
281	! The tracer was not found in restart file or it was equal zero everywhere.
282	WRITE(lunout, *) "The tracer ", trim(tracers(iq)%name), " will be initialized"
283	IF (it==id_pcsat .OR. it==id_pcq .OR. &
284	it==id_pcs0 .OR. it==id_pcq0) THEN
285	tr_seri(:, :, it) = 100.
286	ELSE IF (it==id_pcocsat .OR. it==id_pcos0) THEN
287	DO i = 1, klon
288	IF (pctsrf (i, is_oce) == 0.) THEN
289	tr_seri(i, :, it) = 0.
290	ELSE
291	tr_seri(i, :, it) = 100.
292	END IF
293	END DO
294	ELSE
295	! No specific initialization exist for this tracer
296	tr_seri(:, :, it) = 0.
297	END IF
298	END IF
299	END DO
300
301	END SUBROUTINE traclmdz_init
302
303	SUBROUTINE traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
304	cdragh, coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon, couchelimite, sh, &
305	rh, pphi, ustar, wstar, ale_bl, ale_wake, zu10m, zv10m, &
306	tr_seri, source, d_tr_cl, d_tr_dec, zmasse) !RomP
307
308	USE dimphy
309	USE infotrac_phy, ONLY: nbtr, pbl_flg
310	USE lmdz_strings, ONLY: int2str
311	USE regr_pr_comb_coefoz_m, ONLY: regr_pr_comb_coefoz
312	USE o3_chem_m, ONLY: o3_chem
313	USE indice_sol_mod
314	USE lmdz_q_sat, ONLY: q_sat
315
316	INCLUDE "YOMCST.h"
317
318	!==========================================================================
319	! -- DESCRIPTION DES ARGUMENTS --
320	!==========================================================================
321
322	! Input arguments
323
324	!Configuration grille,temps:
325	INTEGER, INTENT(IN) :: nstep ! nombre d'appels de la physiq
326	INTEGER, INTENT(IN) :: julien ! Jour julien
327	REAL, INTENT(IN) :: gmtime
328	REAL, INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
329	REAL, DIMENSION(klon), INTENT(IN) :: xlat ! latitudes pour chaque point
330	REAL, INTENT(IN) :: xlon(:) ! dim(klon) longitude
331
332	!Physique:
333	!--------
334	REAL, DIMENSION(klon, klev), INTENT(IN) :: t_seri ! Temperature
335	REAL, DIMENSION(klon, klev + 1), INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
336	REAL, DIMENSION(klon, klev), INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
337	REAL, INTENT(IN) :: zmasse (:, :) ! dim(klon,klev) density of air, in kg/m2
338
339
340	!Couche limite:
341	!--------------
342
343	REAL, DIMENSION(klon), INTENT(IN) :: cdragh ! coeff drag pour T et Q
344	REAL, DIMENSION(klon, klev), INTENT(IN) :: coefh ! coeff melange CL (m**2/s)
345	REAL, DIMENSION(klon), INTENT(IN) :: yu1 ! vents au premier niveau
346	REAL, DIMENSION(klon), INTENT(IN) :: yv1 ! vents au premier niveau
347	LOGICAL, INTENT(IN) :: couchelimite
348	REAL, DIMENSION(klon, klev), INTENT(IN) :: sh ! humidite specifique
349	REAL, DIMENSION(klon, klev), INTENT(IN) :: rh ! Humidite relative
350	REAL, DIMENSION(klon, klev), INTENT(IN) :: pphi ! geopotentie
351	REAL, DIMENSION(klon), INTENT(IN) :: ustar ! ustar (m/s)
352	REAL, DIMENSION(klon), INTENT(IN) :: wstar, ale_bl, ale_wake ! wstar (m/s) and Avail. Lifti. Energ.
353	REAL, DIMENSION(klon), INTENT(IN) :: zu10m ! vent zonal 10m (m/s)
354	REAL, DIMENSION(klon), INTENT(IN) :: zv10m ! vent zonal 10m (m/s)
355
356	! Arguments necessaires pour les sources et puits de traceur:
357	REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
358	REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: pctsrf ! Pourcentage de sol f(nature du sol)
359
360	! InOutput argument
361	REAL, DIMENSION(klon, klev, nbtr), INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
362
363	! Output argument
364	REAL, DIMENSION(klon, nbtr), INTENT(OUT) :: source ! a voir lorsque le flux de surface est prescrit
365	REAL, DIMENSION(klon, klev, nbtr), INTENT(OUT) :: d_tr_cl ! Td couche limite/traceur
366
367	!=======================================================================================
368	! -- VARIABLES LOCALES TRACEURS --
369	!=======================================================================================
370
371	INTEGER :: i, k, it
372	INTEGER :: lmt_pas ! number of time steps of "physics" per day
373
374	REAL, DIMENSION(klon) :: d_trs ! Td dans le reservoir
375	REAL, DIMENSION(klon, klev) :: qsat ! pression de la vapeur a saturation
376	REAL, DIMENSION(klon, klev, nbtr) :: d_tr_dec ! Td radioactive
377	REAL :: zrho ! Masse Volumique de l'air KgA/m3
378	REAL :: amn, amx
379
380	!=================================================================
381	! Ajout de la production en Be7 (Beryllium) srcbe U/s/kgA
382	!=================================================================
383
384	IF (id_be /= 0) THEN
385	DO k = 1, klev
386	DO i = 1, klon
387	tr_seri(i, k, id_be) = tr_seri(i, k, id_be) + srcbe(i, k) * pdtphys
388	END DO
389	END DO
390	WRITE(, ) 'Ajout srcBe dans tr_seri: OK'
391	END IF
392
393
394	!=================================================================
395	! Update pseudo-vapor tracers
396	!=================================================================
397
398	CALL q_sat(klon * klev, t_seri, pplay, qsat)
399
400	IF (id_pcsat /= 0) THEN
401	DO k = 1, klev
402	DO i = 1, klon
403	IF (pplay(i, k)>=85000.) THEN
404	tr_seri(i, k, id_pcsat) = qsat(i, k)
405	ELSE
406	tr_seri(i, k, id_pcsat) = MIN (qsat(i, k), tr_seri(i, k, id_pcsat))
407	END IF
408	END DO
409	END DO
410	END IF
411
412	IF (id_pcocsat /= 0) THEN
413	DO k = 1, klev
414	DO i = 1, klon
415	IF (pplay(i, k)>=85000.) THEN
416	IF (pctsrf (i, is_oce) > 0.) THEN
417	tr_seri(i, k, id_pcocsat) = qsat(i, k)
418	ELSE
419	tr_seri(i, k, id_pcocsat) = 0.
420	END IF
421	ELSE
422	tr_seri(i, k, id_pcocsat) = MIN (qsat(i, k), tr_seri(i, k, id_pcocsat))
423	END IF
424	END DO
425	END DO
426	END IF
427
428	IF (id_pcq /= 0) THEN
429	DO k = 1, klev
430	DO i = 1, klon
431	IF (pplay(i, k)>=85000.) THEN
432	tr_seri(i, k, id_pcq) = sh(i, k)
433	ELSE
434	tr_seri(i, k, id_pcq) = MIN (qsat(i, k), tr_seri(i, k, id_pcq))
435	END IF
436	END DO
437	END DO
438	END IF
439
440	IF (id_pcs0 /= 0) THEN
441	DO k = 1, klev
442	DO i = 1, klon
443	IF (pplay(i, k)>=85000.) THEN
444	tr_seri(i, k, id_pcs0) = qsat(i, k)
445	ELSE
446	tr_seri(i, k, id_pcs0) = MIN (qsat(i, k), tr_seri(i, k, id_pcs0))
447	END IF
448	END DO
449	END DO
450	END IF
451
452	IF (id_pcos0 /= 0) THEN
453	DO k = 1, klev
454	DO i = 1, klon
455	IF (pplay(i, k)>=85000.) THEN
456	IF (pctsrf (i, is_oce) > 0.) THEN
457	tr_seri(i, k, id_pcos0) = qsat(i, k)
458	ELSE
459	tr_seri(i, k, id_pcos0) = 0.
460	END IF
461	ELSE
462	tr_seri(i, k, id_pcos0) = MIN (qsat(i, k), tr_seri(i, k, id_pcos0))
463	END IF
464	END DO
465	END DO
466	END IF
467
468	IF (id_pcq0 /= 0) THEN
469	DO k = 1, klev
470	DO i = 1, klon
471	IF (pplay(i, k)>=85000.) THEN
472	tr_seri(i, k, id_pcq0) = sh(i, k)
473	ELSE
474	tr_seri(i, k, id_pcq0) = MIN (qsat(i, k), tr_seri(i, k, id_pcq0))
475	END IF
476	END DO
477	END DO
478	END IF
479
480	!=================================================================
481	! Update tracer : Age of stratospheric air
482	!=================================================================
483	IF (id_aga/=0) THEN
484
485	! Bottom layers
486	DO k = 1, lev_1p5km
487	tr_seri(:, k, id_aga) = 0.0
488	END DO
489
490	! Layers above 1.5km
491	DO k = lev_1p5km + 1, klev - 1
492	tr_seri(:, k, id_aga) = tr_seri(:, k, id_aga) + pdtphys
493	END DO
494
495	! Top layer
496	tr_seri(:, klev, id_aga) = tr_seri(:, klev - 1, id_aga)
497
498	END IF
499
500	!======================================================================
501	! -- Calcul de l'effet de la couche limite --
502	!======================================================================
503
504	IF (couchelimite) THEN
505	source(:, :) = 0.0
506
507	IF (id_be /=0) THEN
508	DO i = 1, klon
509	zrho = pplay(i, 1) / t_seri(i, 1) / RD
510	source(i, id_be) = - vdeptr(id_be) * tr_seri(i, 1, id_be) * zrho
511	END DO
512	END IF
513
514	END IF
515
516	DO it = 1, nbtr
517	IF (couchelimite .AND. pbl_flg(it) == 0 .AND. (it==id_rn .OR. it==id_pb)) THEN
518	! couche limite avec quantite dans le sol calculee
519	CALL cltracrn(it, pdtphys, yu1, yv1, &
520	cdragh, coefh, t_seri, ftsol, pctsrf, &
521	tr_seri(:, :, it), trs(:, it), &
522	paprs, pplay, zmasse * rg, &
523	masktr(:, it), fshtr(:, it), hsoltr(it), &
524	tautr(it), vdeptr(it), &
525	xlat, d_tr_cl(:, :, it), d_trs)
526
527	DO k = 1, klev
528	DO i = 1, klon
529	tr_seri(i, k, it) = tr_seri(i, k, it) + d_tr_cl(i, k, it)
530	END DO
531	END DO
532
533	! Traceur dans le reservoir sol
534	DO i = 1, klon
535	trs(i, it) = trs(i, it) + d_trs(i)
536	END DO
537	END IF
538	END DO
539
540
541	!======================================================================
542	! Calcul de l'effet du puits radioactif
543	!======================================================================
544	CALL radio_decay (radio, rnpb, pdtphys, tautr, tr_seri, d_tr_dec)
545
546	DO it = 1, nbtr
547	IF(radio(it)) THEN
548	DO k = 1, klev
549	DO i = 1, klon
550	tr_seri(i, k, it) = tr_seri(i, k, it) + d_tr_dec(i, k, it)
551	END DO
552	END DO
553	CALL minmaxqfi(tr_seri(:, :, it), 0., 1.e33, 'puits rn it=' // TRIM(int2str(it)))
554	END IF
555	END DO
556
557	!======================================================================
558	! Parameterization of ozone chemistry
559	!======================================================================
560
561	IF (id_o3 /= 0) THEN
562	lmt_pas = NINT(86400. / pdtphys)
563	IF (MOD(nstep - 1, lmt_pas) == 0) THEN
564	! Once per day, update the coefficients for ozone chemistry:
565	CALL regr_pr_comb_coefoz(julien, xlat, paprs, pplay)
566	END IF
567	CALL o3_chem(julien, gmtime, t_seri, zmasse, pdtphys, xlat, &
568	xlon, tr_seri(:, :, id_o3))
569	END IF
570
571	END SUBROUTINE traclmdz
572
573
574	SUBROUTINE traclmdz_to_restart(trs_out)
575	! This SUBROUTINE is called from phyredem.F where the module
576	! variable trs is written to restart file (restartphy.nc)
577	USE dimphy
578	USE infotrac_phy, ONLY: nbtr
579
580	REAL, DIMENSION(klon, nbtr), INTENT(OUT) :: trs_out
581	INTEGER :: ierr
582
583	IF (ALLOCATED(trs)) THEN
584	trs_out(:, :) = trs(:, :)
585	ELSE
586	! No previous allocate of trs. This is the case for create_etat0_limit.
587	trs_out(:, :) = 0.0
588	END IF
589
590	END SUBROUTINE traclmdz_to_restart
591
592
593	END MODULE traclmdz_mod

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

Download in other formats:

Original Format