Context Navigation

traclmdz_mod.F90 @ 1409

Last change on this file since 1409 was 1409, checked in by jghattas, 14 years ago

Added tracer "Age of stratospheric air", activated by putting Aga in tracer.def.
The logical rnpb is now controled during execution. rnpb is true only if both tracers RN and PB existe in tracer.def. RN and PB can now be removed from tracer.def
In tracer.def, the 2 water traceurs (H2Ov and H2Ol) must still be the first 2 tracers. The following tracers have no specific order(RN and PB can now change places). Still a minimum of 3 tracers are required.

File size: 21.3 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	!
8
9	IMPLICIT NONE
10
11	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: masktr ! Masque reservoir de sol traceur
12	!$OMP THREADPRIVATE(masktr) ! Masque de l'echange avec la surface (1 = reservoir) ou (possible >= 1 )
13	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: fshtr ! Flux surfacique dans le reservoir de sol
14	!$OMP THREADPRIVATE(fshtr)
15	REAL,DIMENSION(:),ALLOCATABLE,SAVE :: hsoltr ! Epaisseur equivalente du reservoir de sol
16	!$OMP THREADPRIVATE(hsoltr)
17	!
18	!Radioelements:
19	!--------------
20	!
21	REAL,DIMENSION(:),ALLOCATABLE,SAVE :: tautr ! Constante de decroissance radioactive
22	!$OMP THREADPRIVATE(tautr)
23	REAL,DIMENSION(:),ALLOCATABLE,SAVE :: vdeptr ! Vitesse de depot sec dans la couche Brownienne
24	!$OMP THREADPRIVATE(vdeptr)
25	REAL,DIMENSION(:),ALLOCATABLE,SAVE :: scavtr ! Coefficient de lessivage
26	!$OMP THREADPRIVATE(scavtr)
27	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: srcbe ! Production du beryllium7 dans l atmosphere (U/s/kgA)
28	!$OMP THREADPRIVATE(srcbe)
29
30	LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: radio ! radio(it) = true => decroisssance radioactive
31	!$OMP THREADPRIVATE(radio)
32
33	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: trs ! Conc. radon ds le sol
34	!$OMP THREADPRIVATE(trs)
35
36	INTEGER,SAVE :: id_aga ! Identification number for tracer : Age of stratospheric air
37	!$OMP THREADPRIVATE(id_aga)
38	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.
39	!$OMP THREADPRIVATE(lev_1p5km)
40
41	INTEGER,SAVE :: id_rn, id_pb ! Identification number for tracer : radon (Rn222), lead (Pb210)
42	!$OMP THREADPRIVATE(id_rn, id_pb)
43
44	INTEGER,SAVE :: id_be ! Activation et position du traceur Be7 [ id_be=0 -> desactive ]
45	!$OMP THREADPRIVATE(id_be)
46
47	!IM ajout traceurs RR
48	INTEGER,SAVE :: id_dry !traceur dry intrusions
49	!$OMP THREADPRIVATE(id_dry)
50	INTEGER,SAVE :: id_pcsat, id_pcocsat, id_pcq ! traceurs pseudo-vapeur CL qsat, qsat_oc, q
51	!$OMP THREADPRIVATE(id_pcsat, id_pcocsat, id_pcq)
52	INTEGER,SAVE :: id_pcs0, id_pcos0, id_pcq0 ! traceurs pseudo-vapeur CL qsat, qsat_oc, q
53	! ! qui ne sont pas transportes par la convection
54	!$OMP THREADPRIVATE(id_pcs0, id_pcos0, id_pcq0)
55
56	INTEGER, SAVE:: id_o3
57	!$OMP THREADPRIVATE(id_o3)
58	! index of ozone tracer with Cariolle parameterization
59	! 0 means no ozone tracer
60
61	LOGICAL,SAVE :: rnpb=.FALSE. ! Presence du couple Rn222, Pb210
62	!$OMP THREADPRIVATE(rnpb)
63
64
65	CONTAINS
66
67	SUBROUTINE traclmdz_from_restart(trs_in)
68	! This subroutine initialize the module saved variable trs with values from restart file (startphy.nc).
69	! This subroutine is called from phyetat0 after the field trs_in has been read.
70
71	USE dimphy
72	USE infotrac
73
74	! Input argument
75	REAL,DIMENSION(klon,nbtr), INTENT(IN) :: trs_in
76
77	! Local variables
78	INTEGER :: ierr
79
80	! Allocate restart variables trs
81	ALLOCATE( trs(klon,nbtr), stat=ierr)
82	IF (ierr /= 0) CALL abort_gcm('traclmdz_from_restart', 'pb in allocation 1',1)
83
84	! Initialize trs with values read from restart file
85	trs(:,:) = trs_in(:,:)
86
87	END SUBROUTINE traclmdz_from_restart
88
89
90	SUBROUTINE traclmdz_init(pctsrf, ftsol, tr_seri, t_seri, pplay, sh, aerosol, lessivage)
91	! This subroutine allocates and initialize module variables and control variables.
92	USE dimphy
93	USE infotrac
94	USE regr_pr_comb_coefoz_m, ONLY: alloc_coefoz
95	USE press_coefoz_m, ONLY: press_coefoz
96	USE carbon_cycle_mod, ONLY : carbon_cycle_init, carbon_cycle_tr, carbon_cycle_cpl
97
98	INCLUDE "indicesol.h"
99
100	! Input variables
101	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol f(nature du sol)
102	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
103	!IM traceurs RR REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: tr_seri! Concentration Traceur [U/KgA]
104	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri! Concentration Traceur [U/KgA]
105	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
106	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
107	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
108
109	! Output variables
110	LOGICAL,DIMENSION(nbtr), INTENT(OUT) :: aerosol
111	LOGICAL,INTENT(OUT) :: lessivage
112
113	! Local variables
114	INTEGER :: ierr, it, iiq, i, k
115	REAL,DIMENSION(klon,klev) :: qsat ! pression de la vapeur a saturation
116
117	! --------------------------------------------
118	! Allocation
119	! --------------------------------------------
120
121	ALLOCATE( scavtr(nbtr), stat=ierr)
122	IF (ierr /= 0) CALL abort_gcm('traclmdz_init', 'pb in allocation 9',1)
123	scavtr(:)=1.
124
125	ALLOCATE( radio(nbtr), stat=ierr)
126	IF (ierr /= 0) CALL abort_gcm('traclmdz_init', 'pb in allocation 11',1)
127	radio(:) = .false. ! Par defaut pas decroissance radioactive
128
129	ALLOCATE( masktr(klon,nbtr), stat=ierr)
130	IF (ierr /= 0) CALL abort_gcm('traclmdz_init', 'pb in allocation 2',1)
131
132	ALLOCATE( fshtr(klon,nbtr), stat=ierr)
133	IF (ierr /= 0) CALL abort_gcm('traclmdz_init', 'pb in allocation 3',1)
134
135	ALLOCATE( hsoltr(nbtr), stat=ierr)
136	IF (ierr /= 0) CALL abort_gcm('traclmdz_init', 'pb in allocation 4',1)
137
138	ALLOCATE( tautr(nbtr), stat=ierr)
139	IF (ierr /= 0) CALL abort_gcm('traclmdz_init', 'pb in allocation 5',1)
140	tautr(:) = 0.
141
142	ALLOCATE( vdeptr(nbtr), stat=ierr)
143	IF (ierr /= 0) CALL abort_gcm('traclmdz_init', 'pb in allocation 6',1)
144	vdeptr(:) = 0.
145
146
147	lessivage = .TRUE.
148	aerosol(:) = .FALSE. ! Tous les traceurs sont des gaz par defaut
149
150	!
151	! Recherche des traceurs connus : Be7, O3, CO2,...
152	! --------------------------------------------
153	id_rn=0; id_pb=0; id_aga=0; id_be=0; id_o3=0
154	DO it=1,nbtr
155	iiq=niadv(it+2)
156	IF ( tname(iiq) == "RN" ) THEN
157	id_rn=it ! radon
158	ELSE IF ( tname(iiq) == "PB") THEN
159	id_pb=it ! plomb
160	ELSE IF ( tname(iiq) == "Aga" .OR. tname(iiq) == "AGA" ) THEN
161	! Age of stratospheric air
162	id_aga=it
163	radio(id_aga) = .FALSE.
164	aerosol(id_aga) = .FALSE.
165	pbl_flg(id_aga) = 0
166
167	! Find the first model layer above 1.5km from the surface
168	IF (klev>=30) THEN
169	lev_1p5km=6 ! NB! This value is for klev=39
170	ELSE IF (klev>=10) THEN
171	lev_1p5km=5 ! NB! This value is for klev=19
172	ELSE
173	lev_1p5km=klev/2
174	END IF
175	ELSE IF ( tname(iiq) == "BE" .OR. tname(iiq) == "Be" .OR. &
176	tname(iiq) == "BE7" .OR. tname(iiq) == "Be7" ) THEN
177	! Recherche du Beryllium 7
178	id_be=it
179	ALLOCATE( srcbe(klon,klev) )
180	radio(id_be) = .TRUE.
181	aerosol(id_be) = .TRUE. ! le Be est un aerosol
182	CALL init_be(pctsrf,masktr(:,id_be),tautr(id_be),vdeptr(id_be),scavtr(id_be),srcbe)
183	WRITE(,) 'Initialisation srcBe: OK'
184	ELSE IF (tname(iiq)=="O3" .OR. tname(iiq)=="o3") THEN
185	! Recherche de l'ozone : parametrization de la chimie par Cariolle
186	id_o3=it
187	CALL alloc_coefoz ! allocate ozone coefficients
188	CALL press_coefoz ! read input pressure levels
189	END IF
190	END DO
191
192	id_dry=0
193	DO it=1,nbtr
194	iiq=niadv(it+2)
195	IF ( tname(iiq) == "dry" .OR. tname(iiq) == "Dry" ) THEN
196	id_dry=it
197	END IF
198	END DO
199
200	id_pcsat=0
201	DO it=1,nbtr
202	iiq=niadv(it+2)
203	IF ( tname(iiq) == "pcsat" .OR. tname(iiq) == "Pcsat" ) THEN
204	id_pcsat=it
205	END IF
206	END DO
207
208	id_pcocsat=0
209	DO it=1,nbtr
210	iiq=niadv(it+2)
211	IF ( tname(iiq) == "pcocsat" .OR. tname(iiq) == "Pcocsat" ) THEN
212	id_pcocsat=it
213	END IF
214	END DO
215
216	id_pcq=0
217	DO it=1,nbtr
218	iiq=niadv(it+2)
219	IF ( tname(iiq) == "pcq" .OR. tname(iiq) == "Pcq" ) THEN
220	id_pcq=it
221	END IF
222	END DO
223
224	id_pcs0=0
225	DO it=1,nbtr
226	iiq=niadv(it+2)
227	IF ( tname(iiq) == "pcs0" .OR. tname(iiq) == "Pcs0" ) THEN
228	id_pcs0=it
229	END IF
230	END DO
231
232	id_pcos0=0
233	DO it=1,nbtr
234	iiq=niadv(it+2)
235	IF ( tname(iiq) == "pcos0" .OR. tname(iiq) == "Pcos0" ) THEN
236	id_pcos0=it
237	END IF
238	END DO
239
240	id_pcq0=0
241	DO it=1,nbtr
242	iiq=niadv(it+2)
243	IF ( tname(iiq) == "pcq0" .OR. tname(iiq) == "Pcq0" ) THEN
244	id_pcq0=it
245	END IF
246	END DO
247
248	!
249	! Valeurs specifiques pour les traceurs Rn222 et Pb210
250	! ----------------------------------------------
251	IF ( id_rn/=0 .AND. id_pb/=0 ) THEN
252	rnpb = .TRUE.
253	radio(id_rn)= .TRUE.
254	radio(id_pb)= .TRUE.
255	pbl_flg(id_rn) = 0 ! au lieu de clsol=true ! CL au sol calcule
256	pbl_flg(id_pb) = 0 ! au lieu de clsol=true
257	aerosol(id_rn) = .FALSE.
258	aerosol(id_pb) = .TRUE. ! le Pb est un aerosol
259
260	CALL initrrnpb (ftsol,pctsrf,masktr,fshtr,hsoltr,tautr,vdeptr,scavtr)
261	END IF
262
263	!
264	! Initialisation de module carbon_cycle_mod
265	! ----------------------------------------------
266	IF (carbon_cycle_tr .OR. carbon_cycle_cpl) THEN
267	CALL carbon_cycle_init(tr_seri, aerosol, radio)
268	END IF
269
270	!IM initialisation traceurs pseudo-vapeurs
271	call q_sat(klon*klev,t_seri,pplay,qsat)
272	IF ( id_pcsat /= 0 ) THEN
273	DO k = 1, klev
274	DO i = 1, klon
275	IF ( pplay(i,k).GE.85000.) THEN
276	tr_seri(i,k,id_pcsat) = qsat(i,k)
277	ELSE
278	tr_seri(i,k,id_pcsat) = 100.
279	END IF
280	END DO
281	END DO
282	END IF
283
284	IF ( id_pcocsat /= 0 ) THEN
285	DO k = 1, klev
286	DO i = 1, klon
287	IF ( pplay(i,k).GE.85000.) THEN
288	IF ( pctsrf (i, is_oce) > 0. ) THEN
289	tr_seri(i,k,id_pcocsat) = qsat(i,k)
290	ELSE
291	tr_seri(i,k,id_pcocsat) = 100.
292	END IF
293	END IF
294	END DO
295	END DO
296	END IF
297
298	IF ( id_pcq /= 0 ) THEN
299	DO k = 1, klev
300	DO i = 1, klon
301	IF ( pplay(i,k).GE.85000.) THEN
302	tr_seri(i,k,id_pcq) = sh(i,k)
303	ELSE
304	tr_seri(i,k,id_pcq) = 100.
305	END IF
306	END DO
307	END DO
308	END IF
309
310	IF ( id_pcs0 /= 0 ) THEN
311	DO k = 1, klev
312	DO i = 1, klon
313	IF ( pplay(i,k).GE.85000.) THEN
314	tr_seri(i,k,id_pcs0) = qsat(i,k)
315	ELSE
316	tr_seri(i,k,id_pcs0) = 100.
317	END IF
318	END DO
319	END DO
320	END IF
321
322	IF ( id_pcos0 /= 0 ) THEN
323	DO k = 1, klev
324	DO i = 1, klon
325	IF ( pplay(i,k).GE.85000.) THEN
326	IF ( pctsrf (i, is_oce) > 0. ) THEN
327	tr_seri(i,k,id_pcos0) = qsat(i,k)
328	ELSE
329	tr_seri(i,k,id_pcos0) = 100.
330	END IF
331	END IF
332	END DO
333	END DO
334	END IF
335
336	IF ( id_pcq0 /= 0 ) THEN
337	DO k = 1, klev
338	DO i = 1, klon
339	IF ( pplay(i,k).GE.85000.) THEN
340	tr_seri(i,k,id_pcq0) = sh(i,k)
341	ELSE
342	tr_seri(i,k,id_pcq0) = 100.
343	END IF
344	END DO
345	END DO
346	END IF
347
348	END SUBROUTINE traclmdz_init
349
350	SUBROUTINE traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
351	cdragh, coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon, couchelimite, sh, &
352	tr_seri, source, solsym, d_tr_cl, zmasse)
353
354	USE dimphy
355	USE infotrac
356	USE regr_pr_comb_coefoz_m, ONLY: regr_pr_comb_coefoz
357	USE o3_chem_m, ONLY: o3_chem
358	USE carbon_cycle_mod, ONLY : carbon_cycle, carbon_cycle_tr, carbon_cycle_cpl
359	INCLUDE "YOMCST.h"
360	INCLUDE "indicesol.h"
361
362	!==========================================================================
363	! -- DESCRIPTION DES ARGUMENTS --
364	!==========================================================================
365
366	! Input arguments
367	!
368	!Configuration grille,temps:
369	INTEGER,INTENT(IN) :: nstep ! nombre d'appels de la physiq
370	INTEGER,INTENT(IN) :: julien ! Jour julien
371	REAL,INTENT(IN) :: gmtime
372	REAL,INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
373	REAL,DIMENSION(klon),INTENT(IN) :: xlat ! latitudes pour chaque point
374	REAL, INTENT(IN):: xlon(:) ! dim(klon) longitude
375
376	!
377	!Physique:
378	!--------
379	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
380	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
381	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
382	REAL,intent(in):: zmasse (:, :) ! dim(klon,klev) density of air, in kg/m2
383
384
385	!Couche limite:
386	!--------------
387	!
388	REAL,DIMENSION(klon,klev),INTENT(IN) :: cdragh ! coeff drag pour T et Q
389	REAL,DIMENSION(klon,klev),INTENT(IN) :: coefh ! coeff melange CL (m**2/s)
390	REAL,DIMENSION(klon),INTENT(IN) :: yu1 ! vents au premier niveau
391	REAL,DIMENSION(klon),INTENT(IN) :: yv1 ! vents au premier niveau
392	LOGICAL,INTENT(IN) :: couchelimite
393	!IM traceurs RR
394	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
395
396	! Arguments necessaires pour les sources et puits de traceur:
397	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
398	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol f(nature du sol)
399
400
401	! InOutput argument
402	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
403
404	! Output argument
405	CHARACTER(len=8),DIMENSION(nbtr), INTENT(OUT) :: solsym
406	REAL,DIMENSION(klon,nbtr), INTENT(OUT) :: source ! a voir lorsque le flux de surface est prescrit
407	REAL,DIMENSION(klon,klev,nbtr), INTENT(OUT) :: d_tr_cl ! Td couche limite/traceur
408
409	!=======================================================================================
410	! -- VARIABLES LOCALES TRACEURS --
411	!=======================================================================================
412
413	INTEGER :: i, k, it
414	INTEGER lmt_pas ! number of time steps of "physics" per day
415
416	REAL,DIMENSION(klon) :: d_trs ! Td dans le reservoir
417	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_dec ! Td radioactive
418	REAL :: zrho ! Masse Volumique de l'air KgA/m3
419
420	!IM traceurs RR
421	REAL,DIMENSION(klon,klev) :: qsat ! pression de la vapeur a saturation
422	REAL :: amn, amx
423	!
424	!=================================================================
425	! Ajout de la production en Be7 (Beryllium) srcbe U/s/kgA
426	!=================================================================
427	!
428	IF ( id_be /= 0 ) THEN
429	DO k = 1, klev
430	DO i = 1, klon
431	tr_seri(i,k,id_be) = tr_seri(i,k,id_be)+srcbe(i,k)*pdtphys
432	END DO
433	END DO
434	WRITE(,) 'Ajout srcBe dans tr_seri: OK'
435	END IF
436
437	!IM ajout traceurs RR
438	call q_sat(klon*klev,t_seri,pplay,qsat)
439
440	IF ( id_pcsat /= 0 ) THEN
441	DO k = 1, klev
442	DO i = 1, klon
443	IF ( pplay(i,k).GE.85000.) THEN
444	tr_seri(i,k,id_pcsat) = qsat(i,k)
445	END IF
446	END DO
447	END DO
448	END IF
449
450	IF ( id_pcocsat /= 0 ) THEN
451	DO k = 1, klev
452	DO i = 1, klon
453	IF ( pplay(i,k).GE.85000.) THEN
454	IF ( pctsrf (i, is_oce) > 0. ) THEN
455	tr_seri(i,k,id_pcocsat) = qsat(i,k)
456	END IF
457	END IF
458	END DO
459	END DO
460	END IF
461
462	IF ( id_pcq /= 0 ) THEN
463	DO k = 1, klev
464	DO i = 1, klon
465	IF ( pplay(i,k).GE.85000.) THEN
466	tr_seri(i,k,id_pcq) = sh(i,k)
467	END IF
468	END DO
469	END DO
470	END IF
471
472	IF ( id_pcs0 /= 0 ) THEN
473	DO k = 1, klev
474	DO i = 1, klon
475	IF ( pplay(i,k).GE.85000.) THEN
476	tr_seri(i,k,id_pcs0) = qsat(i,k)
477	END IF
478	END DO
479	END DO
480	END IF
481
482	IF ( id_pcos0 /= 0 ) THEN
483	DO k = 1, klev
484	DO i = 1, klon
485	IF ( pplay(i,k).GE.85000.) THEN
486	IF ( pctsrf (i, is_oce) > 0. ) THEN
487	tr_seri(i,k,id_pcos0) = qsat(i,k)
488	END IF
489	END IF
490	END DO
491	END DO
492	END IF
493
494	IF ( id_pcq0 /= 0 ) THEN
495	DO k = 1, klev
496	DO i = 1, klon
497	IF ( pplay(i,k).GE.85000.) THEN
498	tr_seri(i,k,id_pcq0) = sh(i,k)
499	END IF
500	END DO
501	END DO
502	END IF
503
504	DO it=1,nbtr
505	WRITE(solsym(it),'(i2)') it
506	END DO
507
508	!=================================================================
509	! Update tracer : Age of stratospheric air
510	!=================================================================
511	IF (id_aga/=0) THEN
512
513	! Bottom layers
514	DO k = 1, lev_1p5km
515	tr_seri(:,k,id_aga) = 0.0
516	END DO
517
518	! Layers above 1.5km
519	DO k = lev_1p5km+1,klev-1
520	tr_seri(:,k,id_aga) = tr_seri(:,k,id_aga) + pdtphys
521	END DO
522
523	! Top layer
524	tr_seri(:,klev,id_aga) = tr_seri(:,klev-1,id_aga)
525
526	END IF
527
528	!======================================================================
529	! -- Calcul de l'effet de la couche limite --
530	!======================================================================
531
532	IF (couchelimite) THEN
533	source(:,:) = 0.0
534
535	IF (id_be /=0) THEN
536	DO i=1, klon
537	zrho = pplay(i,1)/t_seri(i,1)/RD
538	source(i,id_be) = - vdeptr(id_be)tr_seri(i,1,id_be)zrho
539	END DO
540	END IF
541
542	END IF
543
544	DO it=1, nbtr
545	IF (couchelimite .AND. pbl_flg(it) == 0 .AND. (it==id_rn .OR. it==id_pb)) THEN
546	! couche limite avec quantite dans le sol calculee
547	CALL cltracrn(it, pdtphys, yu1, yv1, &
548	cdragh, coefh,t_seri,ftsol,pctsrf, &
549	tr_seri(:,:,it),trs(:,it), &
550	paprs, pplay, zmasse * rg, &
551	masktr(:,it),fshtr(:,it),hsoltr(it),&
552	tautr(it),vdeptr(it), &
553	xlat,d_tr_cl(:,:,it),d_trs)
554
555	DO k = 1, klev
556	DO i = 1, klon
557	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cl(i,k,it)
558	END DO
559	END DO
560
561	! Traceur dans le reservoir sol
562	DO i = 1, klon
563	trs(i,it) = trs(i,it) + d_trs(i)
564	END DO
565	END IF
566	END DO
567
568	!IM traceurs RR
569	IF ( id_pcsat /= 0 ) THEN
570	DO k = 1, klev
571	DO i = 1, klon
572	IF ( pplay(i,k).LT.85000.) THEN
573	tr_seri(i,k,id_pcsat) = MIN (qsat(i,k), tr_seri(i,k,id_pcsat))
574	END IF
575	END DO
576	END DO
577	END IF
578
579	IF ( id_pcocsat /= 0 ) THEN
580	DO k = 1, klev
581	DO i = 1, klon
582	IF ( pplay(i,k).LT.85000.) THEN
583	tr_seri(i,k,id_pcocsat) = MIN (qsat(i,k), tr_seri(i,k,id_pcocsat))
584	END IF
585	END DO
586	END DO
587	END IF
588
589	IF ( id_pcq /= 0 ) THEN
590	DO k = 1, klev
591	DO i = 1, klon
592	IF ( pplay(i,k).LT.85000.) THEN
593	tr_seri(i,k,id_pcq) = MIN (qsat(i,k), tr_seri(i,k,id_pcq))
594	END IF
595	END DO
596	END DO
597	END IF
598
599	IF ( id_pcs0 /= 0 ) THEN
600	DO k = 1, klev
601	DO i = 1, klon
602	IF ( pplay(i,k).LT.85000.) THEN
603	tr_seri(i,k,id_pcs0) = MIN (qsat(i,k), tr_seri(i,k,id_pcs0))
604	END IF
605	END DO
606	END DO
607	END IF
608
609	IF ( id_pcos0 /= 0 ) THEN
610	DO k = 1, klev
611	DO i = 1, klon
612	IF ( pplay(i,k).LT.85000.) THEN
613	tr_seri(i,k,id_pcos0) = MIN (qsat(i,k), tr_seri(i,k,id_pcos0))
614	END IF
615	END DO
616	END DO
617	END IF
618
619	IF ( id_pcq0 /= 0 ) THEN
620	DO k = 1, klev
621	DO i = 1, klon
622	IF ( pplay(i,k).LT.85000.) THEN
623	tr_seri(i,k,id_pcq0) = MIN (qsat(i,k), tr_seri(i,k,id_pcq0))
624	END IF
625	END DO
626	END DO
627	END IF
628	!======================================================================
629	! Calcul de l'effet du puits radioactif
630	!======================================================================
631	CALL radio_decay (radio,rnpb,pdtphys,tautr,tr_seri,d_tr_dec)
632
633	DO it=1,nbtr
634	IF(radio(it)) then
635	DO k = 1, klev
636	DO i = 1, klon
637	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_dec(i,k,it)
638	END DO
639	END DO
640	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'puits rn it='//solsym(it))
641	END IF
642	END DO
643
644	!======================================================================
645	! Parameterization of ozone chemistry
646	!======================================================================
647
648	IF (id_o3 /= 0) then
649	lmt_pas = NINT(86400./pdtphys)
650	IF (MOD(nstep - 1, lmt_pas) == 0) THEN
651	! Once per day, update the coefficients for ozone chemistry:
652	CALL regr_pr_comb_coefoz(julien, xlat, paprs, pplay)
653	END IF
654	CALL o3_chem(julien, gmtime, t_seri, zmasse, pdtphys, xlat, &
655	xlon, tr_seri(:, :, id_o3))
656	END IF
657
658	!======================================================================
659	! Calcul de cycle de carbon
660	!======================================================================
661	IF (carbon_cycle_tr .OR. carbon_cycle_cpl) THEN
662	CALL carbon_cycle(nstep, pdtphys, pctsrf, tr_seri)
663	END IF
664
665	END SUBROUTINE traclmdz
666
667
668	SUBROUTINE traclmdz_to_restart(trs_out)
669	! This subroutine is called from phyredem.F where the module
670	! variable trs is written to restart file (restartphy.nc)
671	USE dimphy
672	USE infotrac
673
674	REAL,DIMENSION(klon,nbtr), INTENT(OUT) :: trs_out
675	INTEGER :: ierr
676
677	IF ( ALLOCATED(trs) ) THEN
678	trs_out(:,:) = trs(:,:)
679	ELSE
680	! No previous allocate of trs. This is the case for create_etat0_limit.
681	trs_out(:,:) = 0.0
682	END IF
683
684	END SUBROUTINE traclmdz_to_restart
685
686
687	END MODULE traclmdz_mod

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

Download in other formats:

Original Format