Context Navigation

traclmdz_mod.F90 @ 4190

Last change on this file since 4190 was 3605, checked in by lguez, 5 years ago
Merge revisions 3427:3600 of trunk into branch `Ocean_skin`
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.7 KB

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

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

Download in other formats:

Original Format