Context Navigation

phytrac.F90 @ 5003

Last change on this file since 5003 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

Rev	Line
[1191]	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)
[524]	20	!
[1191]	21	!======================================================================
	22	! Auteur(s) FH
	23	! Objet: Moniteur general des tendances traceurs
	24	!======================================================================
[524]	25
[1191]	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
[1299]	35	USE control_mod
[524]	36
	37
[1299]	38
[1191]	39	IMPLICIT NONE
[959]	40
[1191]	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	!==========================================================================
[541]	51
[1191]	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
[619]	93
[1191]	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]
[524]	97
[1191]	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
[524]	132
[1191]	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)
[524]	136
	137
[1191]	138	! Output argument
	139	!----------------
	140	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
[524]	141
[1191]	142	!=======================================================================================
	143	! -- LOCAL VARIABLES --
	144	!=======================================================================================
[766]	145
[1191]	146	INTEGER :: i, k, it
	147	INTEGER :: nsplit
[524]	148
[1191]	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)
[524]	154
[1191]	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.
[524]	166
[1191]	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)
[766]	196
[1191]	197	CHARACTER(len=8),DIMENSION(nbtr) :: solsym
[524]	198
	199
[1191]	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
[524]	211
[1191]	212	! Initialize module for specific tracers
	213	SELECT CASE(type_trac)
	214	CASE('lmdz')
[1376]	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)
[1191]	217	CASE('inca')
	218	source(:,:)=0.
	219	CALL tracinca_init(aerosol,lessivage)
	220	END SELECT
	221	!
	222	! Initialize diagnostic output
	223	! ----------------------------
[1192]	224	#ifdef CPP_IOIPSL
[1191]	225	INCLUDE "ini_histrac.h"
[1192]	226	#endif
[1191]	227	END IF
	228	!############################################ END INITIALIZATION #######
[524]	229
[1379]	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
[1191]	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
[1379]	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)
[1191]	247
	248	CASE('inca')
	249	! -- CHIMIE INCA config_inca = aero or chem --
[776]	250
[1191]	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
[524]	263
[1191]	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
[959]	282
[1376]	283	!IM ajout traceurs RR
	284	! print*,'phytrac it,nseuil=',it,nseuil
	285	IF (it.lt.nseuil) THEN
[1191]	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
[1376]	291	END IF !(it.lt.nseuil) then
[959]	292
[1191]	293	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
	294
	295	END DO ! nbtr
	296	END IF ! convection
[959]	297
[524]	298	!======================================================================
[1191]	299	! -- Calcul de l'effet des thermiques --
[524]	300	!======================================================================
	301
[1191]	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)
[524]	308	END DO
[1191]	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
[524]	316
[1191]	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)
[524]	320
[1191]	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
[682]	331
[1191]	332	!======================================================================
	333	! -- Calcul de l'effet de la couche limite --
	334	!======================================================================
[682]	335
[1191]	336	IF (couchelimite) THEN
[524]	337
[1191]	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
[524]	343
[1191]	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
[959]	359
[1191]	360	END DO
	361
	362	END IF ! couche limite
[959]	363
[619]	364
[1191]	365	!======================================================================
	366	! Calcul de l'effet de la precipitation
	367	!======================================================================
[959]	368
[1191]	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	!=========================
[524]	377
[1191]	378	! Tendance des aerosols nuclees et impactes
	379	! -----------------------------------------
	380	DO it = 1, nbtr
	381	IF (aerosol(it)) THEN
	382	DO k = 1, klev
[524]	383	DO i = 1, klon
[1191]	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)
[524]	388
[1191]	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
[524]	407
[1191]	408	!=============================================================
	409	! Ecriture des sorties
	410	!=============================================================
[1192]	411	#ifdef CPP_IOIPSL
[1191]	412	INCLUDE "write_histrac.h"
[1192]	413	#endif
[524]	414
[1191]	415	END SUBROUTINE phytrac

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

Download in other formats:

Original Format