source: LMDZ4/trunk/libf/phylmd/traclmdz_mod.F90 @ 1404

Last change on this file since 1404 was 1403, checked in by Laurent Fairhead, 14 years ago

Merged LMDZ4V5.0-dev branch changes r1292:r1399 to trunk.

Validation:
Validation consisted in compiling the HEAD revision of the trunk,
LMDZ4V5.0-dev branch and the merged sources and running different
configurations on local and SX8 machines comparing results.

Local machine: bench configuration, 32x24x11, gfortran

  • IPSLCM5A configuration (comparison between trunk and merged sources):
    • numerical convergence on dynamical fields over 3 days
    • start files are equivalent (except for RN and PB fields)
    • daily history files equivalent
  • MH07 configuration, new physics package (comparison between LMDZ4V5.0-dev branch and merged sources):
    • numerical convergence on dynamical fields over 3 days
    • start files are equivalent (except for RN and PB fields)
    • daily history files equivalent

SX8 machine (brodie), 96x95x39 on 4 processors:

  • IPSLCM5A configuration:
    • start files are equivalent (except for RN and PB fields)
    • monthly history files equivalent
  • MH07 configuration:
    • start files are equivalent (except for RN and PB fields)
    • monthly history files equivalent

Changes to the makegcm and create_make_gcm scripts to take into account
main programs in F90 files


Fusion de la branche LMDZ4V5.0-dev (r1292:r1399) au tronc principal

Validation:
La validation a consisté à compiler la HEAD de le trunk et de la banche
LMDZ4V5.0-dev et les sources fusionnées et de faire tourner le modéle selon
différentes configurations en local et sur SX8 et de comparer les résultats

En local: 32x24x11, config bench/gfortran

  • pour une config IPSLCM5A (comparaison tronc/fusion):
    • convergence numérique sur les champs dynamiques après 3 jours
    • restart et restartphy égaux (à part sur RN et Pb)
    • fichiers histoire égaux
  • pour une config nlle physique (MH07) (comparaison LMDZ4v5.0-dev/fusion):
    • convergence numérique sur les champs dynamiques après 3 jours
    • restart et restartphy égaux
    • fichiers histoire équivalents

Sur brodie, 96x95x39 sur 4 proc:

  • pour une config IPSLCM5A:
    • restart et restartphy égaux (à part sur RN et PB)
    • pas de différence dans les fichiers histmth.nc
  • pour une config MH07
    • restart et restartphy égaux (à part sur RN et PB)
    • pas de différence dans les fichiers histmth.nc

Changement sur makegcm et create_make-gcm pour pouvoir prendre en compte des
programmes principaux en *F90

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