Context Navigation

phytrac_mod.F90 @ 5456

Last change on this file since 5456 was 4795, checked in by nfevrier, 12 months ago
First save of N. Février's modifications
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 Property svn:keywords set to `Id`
File size: 50.5 KB

Rev	Line
[1877]	1	!$Id: phytrac_mod.F90 4795 2024-02-01 15:53:24Z aborella $
[1813]	2	MODULE phytrac_mod
	3	!=================================================================================
	4	! Interface between the LMDZ physical package and tracer computation.
	5	! Chemistry modules (INCA, Reprobus or the more specific traclmdz routine)
	6	! are called from phytrac.
	7	!
	8	!======================================================================
	9	! Auteur(s) FH
	10	! Objet: Moniteur general des tendances traceurs
	11	!
	12	! iflag_vdf_trac : Options for activating transport by vertical diffusion :
	13	! 1. notmal
	14	! 0. emission is injected in the first layer only, without diffusion
	15	! -1 no emission & no diffusion
	16	! Modification 2013/07/22 : transformed into a module to pass tendencies to
	17	! physics outputs. Additional keys for controling activation of sub processes.
	18	! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
	19	! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
	20	!=================================================================================
	21
	22	!
	23	! Tracer tendencies, for outputs
	24	!-------------------------------
	25	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_cl ! Td couche limite/traceur
	26	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_dec !RomP
	27	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_cv ! Td convection/traceur
	28	! RomP >>>
	29	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_insc
	30	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_bcscav
	31	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_evapls
	32	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_ls
	33	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_trsp
	34	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_sscav
	35	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_sat
	36	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_uscav
	37	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qPr,qDi ! concentration tra dans pluie,air descente insaturee
	38	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qPa,qMel
	39	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: qTrdi,dtrcvMA ! conc traceur descente air insaturee et td convective MA
	40	! RomP <<<
	41	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_th ! Td thermique
	42	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_lessi_impa ! Td du lessivage par impaction
	43	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: d_tr_lessi_nucl ! Td du lessivage par nucleation
	44	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: qPrls !jyg: concentration tra dans pluie LS a la surf.
	45	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: d_tr_dry ! Td depot sec/traceur (1st layer),ALLOCATABLE,SAVE jyg
	46	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: flux_tr_dry ! depot sec/traceur (surface),ALLOCATABLE,SAVE jyg
[4795]	47	!Nicolas Février
	48	INTEGER,DIMENSION(12),SAVE :: list_aero_wetdep
	49	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: beta_st,beta_cv
	50	!$OMP THREADPRIVATE(list_aero_wetdep,beta_st,beta_cv)
[1813]	51
[4795]	52
[1813]	53	!$OMP THREADPRIVATE(qPa,qMel,qTrdi,dtrcvMA,d_tr_th,d_tr_lessi_impa,d_tr_lessi_nucl)
	54	!$OMP THREADPRIVATE(d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qPr,qDi)
	55	!$OMP THREADPRIVATE(d_tr_insc,d_tr_bcscav,d_tr_evapls,d_tr_ls,qPrls)
[3043]	56	!$OMP THREADPRIVATE(d_tr_cl,d_tr_dry,flux_tr_dry,d_tr_dec,d_tr_cv)
[1813]	57
	58
	59	CONTAINS
	60
[3465]	61	SUBROUTINE phytrac_init()
	62	USE dimphy
[3649]	63	USE infotrac_phy, ONLY: nbtr, type_trac
	64	USE tracco2i_mod, ONLY: tracco2i_init
[3465]	65	IMPLICIT NONE
	66
	67	ALLOCATE(d_tr_cl(klon,klev,nbtr),d_tr_dry(klon,nbtr))
	68	ALLOCATE(flux_tr_dry(klon,nbtr),d_tr_dec(klon,klev,nbtr),d_tr_cv(klon,klev,nbtr))
	69	ALLOCATE(d_tr_insc(klon,klev,nbtr),d_tr_bcscav(klon,klev,nbtr))
	70	ALLOCATE(d_tr_evapls(klon,klev,nbtr),d_tr_ls(klon,klev,nbtr))
	71	ALLOCATE(qPrls(klon,nbtr),d_tr_trsp(klon,klev,nbtr))
	72	ALLOCATE(d_tr_sscav(klon,klev,nbtr),d_tr_sat(klon,klev,nbtr))
	73	ALLOCATE(d_tr_uscav(klon,klev,nbtr),qPr(klon,klev,nbtr),qDi(klon,klev,nbtr))
	74	ALLOCATE(qPa(klon,klev,nbtr),qMel(klon,klev,nbtr))
	75	ALLOCATE(qTrdi(klon,klev,nbtr),dtrcvMA(klon,klev,nbtr))
	76	ALLOCATE(d_tr_th(klon,klev,nbtr))
	77	ALLOCATE(d_tr_lessi_impa(klon,klev,nbtr),d_tr_lessi_nucl(klon,klev,nbtr))
[4795]	78	! Nicolas Février
	79	ALLOCATE(beta_st(klon,klev,12),beta_cv(klon,klev,12))
[3465]	80
[3649]	81
	82
	83	!===============================================================================
	84	! -- Do specific treatment according to chemestry model or local LMDZ tracers
	85	!
	86	!===============================================================================
	87	SELECT CASE(type_trac)
	88	CASE('co2i')
	89	! -- CO2 interactif --
	90	CALL tracco2i_init()
	91	END SELECT
	92
	93
[3465]	94	END SUBROUTINE phytrac_init
	95
[2146]	96	SUBROUTINE phytrac( &
	97	nstep, julien, gmtime, debutphy, &
	98	lafin, pdtphys, u, v, t_seri, &
	99	paprs, pplay, pmfu, pmfd, &
	100	pen_u, pde_u, pen_d, pde_d, &
	101	cdragh, coefh, fm_therm, entr_therm, &
	102	yu1, yv1, ftsol, pctsrf, &
	103	ustar, u10m, v10m, &
	104	wstar, ale_bl, ale_wake, &
	105	xlat, xlon, &
	106	frac_impa,frac_nucl,beta_fisrt,beta_v1, &
	107	presnivs, pphis, pphi, albsol, &
[2784]	108	sh, ch, rh, cldfra, rneb, &
[2146]	109	diafra, cldliq, itop_con, ibas_con, &
	110	pmflxr, pmflxs, prfl, psfl, &
	111	da, phi, mp, upwd, &
	112	phi2, d1a, dam, sij, wght_cvfd, & ! RomP +RL
	113	wdtrainA, wdtrainM, sigd, clw, elij, & ! RomP
	114	evap, ep, epmlmMm, eplaMm, & ! RomP
	115	dnwd, aerosol_couple, flxmass_w, &
	116	tau_aero, piz_aero, cg_aero, ccm, &
	117	rfname, &
	118	d_tr_dyn, & ! RomP
[3418]	119	tr_seri, init_source)
[2146]	120	!
	121	!======================================================================
	122	! Auteur(s) FH
	123	! Objet: Moniteur general des tendances traceurs
	124	! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
	125	! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
	126	!======================================================================
[1813]	127
[2146]	128	USE ioipsl
	129	USE phys_cal_mod, only : hour
	130	USE dimphy
[2320]	131	USE infotrac_phy, ONLY: nbtr, type_trac, conv_flg, solsym, pbl_flg
[2146]	132	USE mod_grid_phy_lmdz
	133	USE mod_phys_lmdz_para
	134	USE iophy
	135	USE traclmdz_mod
	136	USE tracinca_mod
	137	USE tracreprobus_mod
	138	USE indice_sol_mod
	139	USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root
[2311]	140	USE print_control_mod, ONLY: lunout
[2394]	141	USE aero_mod, ONLY : naero_grp
[1813]	142
[3361]	143	USE tracco2i_mod
	144
[2690]	145	#ifdef CPP_StratAer
	146	USE traccoag_mod
[2752]	147	USE phys_local_var_mod, ONLY: mdw
[3100]	148	USE phys_local_var_mod, ONLY: budg_dep_dry_ocs, budg_dep_wet_ocs
	149	USE phys_local_var_mod, ONLY: budg_dep_dry_so2, budg_dep_wet_so2
[2752]	150	USE phys_local_var_mod, ONLY: budg_dep_dry_h2so4, budg_dep_wet_h2so4
	151	USE phys_local_var_mod, ONLY: budg_dep_dry_part, budg_dep_wet_part
[3100]	152	USE infotrac, ONLY: nbtr_sulgas, id_OCS_strat, id_SO2_strat, id_H2SO4_strat
[2690]	153	USE aerophys
	154	#endif
	155
[2146]	156	IMPLICIT NONE
[1813]	157
[2146]	158	INCLUDE "YOMCST.h"
	159	INCLUDE "clesphys.h"
	160	INCLUDE "thermcell.h"
	161	!==========================================================================
	162	! -- ARGUMENT DESCRIPTION --
	163	!==========================================================================
[1813]	164
[2146]	165	! Input arguments
	166	!----------------
	167	!Configuration grille,temps:
	168	INTEGER,INTENT(IN) :: nstep ! Appel physique
	169	INTEGER,INTENT(IN) :: julien ! Jour julien
	170	REAL,INTENT(IN) :: gmtime ! Heure courante
	171	REAL,INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
	172	LOGICAL,INTENT(IN) :: debutphy ! le flag de l'initialisation de la physique
	173	LOGICAL,INTENT(IN) :: lafin ! le flag de la fin de la physique
[1813]	174
[2146]	175	REAL,DIMENSION(klon),INTENT(IN) :: xlat ! latitudes pour chaque point
	176	REAL,DIMENSION(klon),INTENT(IN) :: xlon ! longitudes pour chaque point
	177	!
	178	!Physique:
	179	!--------
	180	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
	181	REAL,DIMENSION(klon,klev),INTENT(IN) :: u ! variable not used
	182	REAL,DIMENSION(klon,klev),INTENT(IN) :: v ! variable not used
	183	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
	184	REAL,DIMENSION(klon,klev),INTENT(IN) :: rh ! humidite relative
[2784]	185	REAL,DIMENSION(klon,klev),INTENT(IN) :: ch ! eau liquide
[2146]	186	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
	187	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
	188	REAL,DIMENSION(klon,klev),INTENT(IN) :: pphi ! geopotentiel
	189	REAL,DIMENSION(klon),INTENT(IN) :: pphis
	190	REAL,DIMENSION(klev),INTENT(IN) :: presnivs
	191	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldliq ! eau liquide nuageuse
	192	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldfra ! fraction nuageuse (tous les nuages)
	193	REAL,DIMENSION(klon,klev),INTENT(IN) :: diafra ! fraction nuageuse (convection ou stratus artificiels)
	194	REAL,DIMENSION(klon,klev),INTENT(IN) :: rneb ! fraction nuageuse (grande echelle)
	195	!
	196	REAL :: ql_incl ! contenu en eau liquide nuageuse dans le nuage ! ql_incl=oliq/rneb
	197	REAL,DIMENSION(klon,klev),INTENT(IN) :: beta_fisrt ! taux de conversion de l'eau cond (de fisrtilp)
	198	REAL,DIMENSION(klon,klev),INTENT(out) :: beta_v1 ! -- (originale version)
[1813]	199
[2146]	200	!
	201	INTEGER,DIMENSION(klon),INTENT(IN) :: itop_con
	202	INTEGER,DIMENSION(klon),INTENT(IN) :: ibas_con
	203	REAL,DIMENSION(klon),INTENT(IN) :: albsol ! albedo surface
	204	!
	205	!Dynamique
	206	!--------
	207	REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: d_tr_dyn
	208	!
	209	!Convection:
	210	!----------
	211	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu ! flux de masse dans le panache montant
	212	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd ! flux de masse dans le panache descendant
	213	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
	214	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
	215	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
	216	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
[1813]	217
[2146]	218	!...Tiedke
	219	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
	220	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
[1813]	221
[2146]	222	LOGICAL,INTENT(IN) :: aerosol_couple
	223	REAL,DIMENSION(klon,klev),INTENT(IN) :: flxmass_w
[2394]	224	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: tau_aero
	225	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: piz_aero
	226	REAL,DIMENSION(klon,klev,naero_grp,2),INTENT(IN) :: cg_aero
	227	CHARACTER(len=4),DIMENSION(naero_grp),INTENT(IN) :: rfname
[2146]	228	REAL,DIMENSION(klon,klev,2),INTENT(IN) :: ccm
	229	!... K.Emanuel
	230	REAL,DIMENSION(klon,klev),INTENT(IN) :: da
	231	REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
	232	! RomP >>>
	233	REAL,DIMENSION(klon,klev),INTENT(IN) :: d1a,dam
	234	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2
	235	!
	236	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainA
	237	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainM
	238	REAL,DIMENSION(klon),INTENT(IN) :: sigd
	239	! ---- RomP flux entraine, detraine et precipitant kerry Emanuel
	240	REAL,DIMENSION(klon,klev),INTENT(IN) :: evap
	241	REAL,DIMENSION(klon,klev),INTENT(IN) :: ep
	242	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij
	243	REAL,DIMENSION(klon,klev),INTENT(IN) :: wght_cvfd !RL
	244	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij
	245	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm
	246	REAL,DIMENSION(klon,klev),INTENT(IN) :: eplaMm
	247	REAL,DIMENSION(klon,klev),INTENT(IN) :: clw
	248	! RomP <<<
[1813]	249
[2146]	250	!
	251	REAL,DIMENSION(klon,klev),INTENT(IN) :: mp
	252	REAL,DIMENSION(klon,klev),INTENT(IN) :: upwd ! saturated updraft mass flux
	253	REAL,DIMENSION(klon,klev),INTENT(IN) :: dnwd ! saturated downdraft mass flux
	254	!
	255	!Thermiques:
	256	!----------
	257	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: fm_therm
	258	REAL,DIMENSION(klon,klev),INTENT(IN) :: entr_therm
	259	!
	260	!Couche limite:
	261	!--------------
	262	!
	263	REAL,DIMENSION(:),INTENT(IN) :: cdragh ! (klon) coeff drag pour T et Q
	264	REAL,DIMENSION(:,:),INTENT(IN) :: coefh ! (klon,klev) coeff melange CL (m**2/s)
	265	REAL,DIMENSION(:),INTENT(IN) :: ustar,u10m,v10m ! (klon) u* & vent a 10m (m/s)
	266	REAL,DIMENSION(:),INTENT(IN) :: wstar,ale_bl,ale_wake ! (klon) w* and Avail. Lifting Ener.
	267	REAL,DIMENSION(:),INTENT(IN) :: yu1 ! (klon) vents au premier niveau
	268	REAL,DIMENSION(:),INTENT(IN) :: yv1 ! (klon) vents au premier niveau
	269	!
	270	!Lessivage:
	271	!----------
	272	!
	273	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrAA
	274	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrENV
	275	REAL, DIMENSION(:), ALLOCATABLE, SAVE :: coefcoli
	276	LOGICAL, DIMENSION(:), ALLOCATABLE, SAVE :: flag_cvltr
	277	!$OMP THREADPRIVATE(ccntrAA,ccntrENV,coefcoli,flag_cvltr)
	278	REAL, DIMENSION(klon,klev) :: ccntrAA_3d
	279	REAL, DIMENSION(klon,klev) :: ccntrENV_3d
	280	REAL, DIMENSION(klon,klev) :: coefcoli_3d
	281	!
	282	! pour le ON-LINE
	283	!
	284	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
	285	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
[1813]	286
[2146]	287	! Arguments necessaires pour les sources et puits de traceur:
	288	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
	289	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
[1813]	290
[2690]	291	#ifdef CPP_StratAer
	292	REAL,DIMENSION(klon) :: v_dep_dry !dry deposition velocity of stratospheric sulfate in m/s
	293	#endif
[2146]	294	! Output argument
	295	!----------------
	296	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
	297	REAL,DIMENSION(klon,klev) :: sourceBE
[3418]	298	REAL,DIMENSION(klon,nbtr), INTENT(IN) :: init_source
[2690]	299
[2146]	300	!=======================================================================================
	301	! -- LOCAL VARIABLES --
	302	!=======================================================================================
	303
[4795]	304	INTEGER :: i, k, it, j ! j rajouté par Nicolas Février
[2146]	305	INTEGER :: nsplit
	306
	307	!Sources et Reservoirs de traceurs (ex:Radon):
	308	!--------------------------------------------
	309	!
	310	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source ! a voir lorsque le flux de surface est prescrit
[1813]	311	!$OMP THREADPRIVATE(source)
	312
[2146]	313	!
	314	!Entrees/Sorties: (cf ini_histrac.h et write_histrac.h)
	315	!---------------
	316	INTEGER :: iiq, ierr
	317	INTEGER :: nhori, nvert
	318	REAL :: zsto, zout, zjulian
	319	INTEGER,SAVE :: nid_tra ! pointe vers le fichier histrac.nc
[1813]	320	!$OMP THREADPRIVATE(nid_tra)
[2146]	321	REAL,DIMENSION(klon) :: zx_tmp_fi2d ! variable temporaire grille physique
	322	INTEGER :: itau_w ! pas de temps ecriture = nstep + itau_phy
	323	LOGICAL,PARAMETER :: ok_sync=.TRUE.
	324	!
	325	! Nature du traceur
	326	!------------------
	327	LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol ! aerosol(it) = true => aerosol => lessivage
[1813]	328	!$OMP THREADPRIVATE(aerosol) ! aerosol(it) = false => gaz
[2146]	329	REAL,DIMENSION(klon,klev) :: delp ! epaisseur de couche (Pa)
	330	!
	331	! Tendances de traceurs (Td) et flux de traceurs:
	332	!------------------------
	333	REAL,DIMENSION(klon,klev) :: d_tr ! Td dans l'atmosphere
	334	REAL,DIMENSION(klon,klev) :: Mint
	335	REAL,DIMENSION(klon,klev,nbtr) :: zmfd1a
	336	REAL,DIMENSION(klon,klev,nbtr) :: zmfdam
	337	REAL,DIMENSION(klon,klev,nbtr) :: zmfphi2
	338	! Physique
	339	!----------
	340	REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
	341	REAL,DIMENSION(klon,klev) :: zmasse ! densité atmosphérique Kg/m2
	342	REAL,DIMENSION(klon,klev) :: ztra_th
	343	!PhH
	344	REAL,DIMENSION(klon,klev) :: zrho
	345	REAL,DIMENSION(klon,klev) :: zdz
	346	REAL :: evaplsc,dx,beta ! variable pour lessivage Genthon
	347	REAL,DIMENSION(klon) :: his_dh ! ---
	348	! in-cloud scav variables
	349	REAL :: ql_incloud_ref ! ref value of in-cloud condensed water content
	350
	351	!Controles:
	352	!---------
	353	INTEGER,SAVE :: iflag_vdf_trac,iflag_con_trac,iflag_the_trac
	354	INTEGER,SAVE :: iflag_con_trac_omp, iflag_vdf_trac_omp,iflag_the_trac_omp
[1820]	355	!$OMP THREADPRIVATE(iflag_vdf_trac,iflag_con_trac,iflag_the_trac)
[1813]	356
[2146]	357	LOGICAL,SAVE :: lessivage
[1813]	358	!$OMP THREADPRIVATE(lessivage)
	359
[2146]	360	!RomP >>>
	361	INTEGER,SAVE :: iflag_lscav_omp,iflag_lscav
[2210]	362	REAL, SAVE :: ccntrAA_in,ccntrAA_omp
	363	REAL, SAVE :: ccntrENV_in,ccntrENV_omp
	364	REAL, SAVE :: coefcoli_in,coefcoli_omp
	365
[2146]	366	LOGICAL,SAVE :: convscav_omp,convscav
[1822]	367	!$OMP THREADPRIVATE(iflag_lscav)
[2210]	368	!$OMP THREADPRIVATE(ccntrAA_in,ccntrENV_in,coefcoli_in)
[1822]	369	!$OMP THREADPRIVATE(convscav)
[4795]	370
	371
	372	! Nicolas Février
	373
	374	INTEGER,SAVE :: id_ASAPp1a,id_ASAPp2a,id_ASARp1a,id_ASARp2a,id_ASBCM,id_AIBCM,id_CIDUSTM,id_AIPOMM,id_ASPOMM,id_ASSSM,id_CSSSM,id_SSSSM
	375	!$OMP THREADPRIVATE(id_ASAPp1a,id_ASAPp2a,id_ASARp1a,id_ASARp2a,id_ASBCM,id_AIBCM,id_CIDUSTM,id_AIPOMM,id_ASPOMM,id_ASSSM,id_CSSSM,id_SSSSM)
	376
[2146]	377	!RomP <<<
	378	!######################################################################
	379	! -- INITIALIZATION --
	380	!######################################################################
[1813]	381
[4795]	382	! Nicolas Février
	383	! On force iflag_lscav à 4
	384	!iflag_lscav=4
	385
[2146]	386	DO k=1,klev
	387	DO i=1,klon
	388	sourceBE(i,k)=0.
	389	Mint(i,k)=0.
	390	zrho(i,k)=0.
	391	zdz(i,k)=0.
	392	END DO
	393	END DO
[1813]	394
[4795]	395	! Nicolas Février
	396	DO j=1, 12
	397	DO k=1,klev
	398	DO i=1,klon
	399	beta_st(i,k,j)=0.
	400	beta_cv(i,k,j)=0.
	401	END DO
	402	END DO
	403	END DO
	404
	405
[2146]	406	DO it=1, nbtr
	407	DO k=1,klev
	408	DO i=1,klon
	409	d_tr_insc(i,k,it)=0.
	410	d_tr_bcscav(i,k,it)=0.
	411	d_tr_evapls(i,k,it)=0.
	412	d_tr_ls(i,k,it)=0.
	413	d_tr_cv(i,k,it)=0.
	414	d_tr_cl(i,k,it)=0.
	415	d_tr_trsp(i,k,it)=0.
	416	d_tr_sscav(i,k,it)=0.
	417	d_tr_sat(i,k,it)=0.
	418	d_tr_uscav(i,k,it)=0.
	419	d_tr_lessi_impa(i,k,it)=0.
	420	d_tr_lessi_nucl(i,k,it)=0.
	421	qDi(i,k,it)=0.
	422	qPr(i,k,it)=0.
	423	qPa(i,k,it)=0.
	424	qMel(i,k,it)=0.
	425	qTrdi(i,k,it)=0.
	426	dtrcvMA(i,k,it)=0.
	427	zmfd1a(i,k,it)=0.
	428	zmfdam(i,k,it)=0.
	429	zmfphi2(i,k,it)=0.
	430	END DO
	431	END DO
[1813]	432	END DO
[2146]	433
[2710]	434	DO it=1, nbtr
	435	DO i=1,klon
	436	d_tr_dry(i,it)=0.
	437	flux_tr_dry(i,it)=0.
	438	END DO
	439	END DO
	440
[2146]	441	DO k = 1, klev
	442	DO i = 1, klon
	443	delp(i,k) = paprs(i,k)-paprs(i,k+1)
	444	END DO
	445	END DO
	446
	447	IF (debutphy) THEN
	448	!!jyg
[1813]	449	!$OMP BARRIER
[2146]	450	ecrit_tra=86400. ! frequence de stokage en dur
	451	! obsolete car remplace par des ecritures dans phys_output_write
	452	!RomP >>>
	453	!
	454	!Config Key = convscav
	455	!Config Desc = Convective scavenging switch: 0=off, 1=on.
[3450]	456	!Config Def = .FALSE.
[2146]	457	!Config Help =
	458	!
[1813]	459	!$OMP MASTER
[3450]	460	convscav_omp=.FALSE.
[2146]	461	call getin('convscav', convscav_omp)
	462	iflag_vdf_trac_omp=1
	463	call getin('iflag_vdf_trac', iflag_vdf_trac_omp)
	464	iflag_con_trac_omp=1
	465	call getin('iflag_con_trac', iflag_con_trac_omp)
	466	iflag_the_trac_omp=1
	467	call getin('iflag_the_trac', iflag_the_trac_omp)
[1813]	468	!$OMP END MASTER
	469	!$OMP BARRIER
[2146]	470	convscav=convscav_omp
	471	iflag_vdf_trac=iflag_vdf_trac_omp
	472	iflag_con_trac=iflag_con_trac_omp
	473	iflag_the_trac=iflag_the_trac_omp
	474	write(lunout,*) 'phytrac passage dans routine conv avec lessivage', convscav
	475	!
	476	!Config Key = iflag_lscav
	477	!Config Desc = Large scale scavenging parametrization: 0=none, 1=old(Genthon92),
	478	! 2=1+PHeinrich, 3=Reddy_Boucher2004, 4=3+RPilon.
	479	!Config Def = 1
	480	!Config Help =
	481	!
[1813]	482	!$OMP MASTER
[4795]	483	!iflag_lscav_omp=4 !Nicolas Février
[2146]	484	call getin('iflag_lscav', iflag_lscav_omp)
[2210]	485	ccntrAA_omp=1
	486	ccntrENV_omp=1.
	487	coefcoli_omp=0.001
	488	call getin('ccntrAA', ccntrAA_omp)
	489	call getin('ccntrENV', ccntrENV_omp)
	490	call getin('coefcoli', coefcoli_omp)
[1813]	491	!$OMP END MASTER
	492	!$OMP BARRIER
[2146]	493	iflag_lscav=iflag_lscav_omp
[2210]	494	ccntrAA_in=ccntrAA_omp
	495	ccntrENV_in=ccntrENV_omp
	496	coefcoli_in=coefcoli_omp
[2146]	497	!
	498	SELECT CASE(iflag_lscav)
	499	CASE(0)
	500	WRITE(lunout,*) 'Large scale scavenging: none'
	501	CASE(1)
	502	WRITE(lunout,*) 'Large scale scavenging: C. Genthon, Tellus(1992), 44B, 371-389'
	503	CASE(2)
	504	WRITE(lunout,*) 'Large scale scavenging: C. Genthon, modified P. Heinrich'
	505	CASE(3)
	506	WRITE(lunout,*) 'Large scale scavenging: M. Shekkar Reddy and O. Boucher, JGR(2004), 109, D14202'
	507	CASE(4)
	508	WRITE(lunout,*) 'Large scale scavenging: Reddy and Boucher, modified R. Pilon'
	509	END SELECT
	510	!RomP <<<
	511	WRITE(,) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
	512	ALLOCATE( source(klon,nbtr), stat=ierr)
[2311]	513	IF (ierr /= 0) CALL abort_physic('phytrac', 'pb in allocation 1',1)
[1813]	514
[2146]	515	ALLOCATE( aerosol(nbtr), stat=ierr)
[2311]	516	IF (ierr /= 0) CALL abort_physic('phytrac', 'pb in allocation 2',1)
[2146]	517
	518
	519	! Initialize module for specific tracers
	520	SELECT CASE(type_trac)
	521	CASE('lmdz')
	522	CALL traclmdz_init(pctsrf,xlat,xlon,ftsol,tr_seri,t_seri,pplay,sh,pdtphys,aerosol,lessivage)
	523	CASE('inca')
[3418]	524	source(:,:)=init_source(:,:)
[2146]	525	CALL tracinca_init(aerosol,lessivage)
[4795]	526	! Nicolas Février
	527	! On met le lessivage à TRUE tout en désactivant le lessivage de INCA
	528	! correspondant pour ne pas faire le lessivage deux fois.
	529	! Exemple : commenter CALL WETDEP dans aerosolmain.F90
	530	! lessivage=.TRUE.
[2146]	531	CASE('repr')
	532	source(:,:)=0.
[3361]	533	CASE('co2i')
	534	source(:,:)=0.
[3450]	535	lessivage = .FALSE.
	536	aerosol(:) = .FALSE.
	537	pbl_flg(:) = 1
	538	iflag_the_trac= 1
	539	iflag_vdf_trac= 1
	540	iflag_con_trac= 1
[2690]	541	#ifdef CPP_StratAer
	542	CASE('coag')
	543	source(:,:)=0.
	544	DO it= 1, nbtr_sulgas
	545	aerosol(it)=.FALSE.
	546	IF (it==id_H2SO4_strat) aerosol(it)=.TRUE.
	547	ENDDO
	548	DO it= nbtr_sulgas+1, nbtr
	549	aerosol(it)=.TRUE.
	550	ENDDO
	551	#endif
[2146]	552	END SELECT
	553
	554	!
	555	!--initialising coefficients for scavenging in the case of NP
	556	!
	557	ALLOCATE(flag_cvltr(nbtr))
	558	IF (iflag_con.EQ.3) THEN
	559	!
	560	ALLOCATE(ccntrAA(nbtr))
	561	ALLOCATE(ccntrENV(nbtr))
	562	ALLOCATE(coefcoli(nbtr))
	563	!
	564	DO it=1, nbtr
	565	SELECT CASE(type_trac)
	566	CASE('lmdz')
	567	IF (convscav.and.aerosol(it)) THEN
[3450]	568	flag_cvltr(it)=.TRUE.
[2210]	569	ccntrAA(it) =ccntrAA_in !--a modifier par JYG a lire depuis fichier
	570	ccntrENV(it)=ccntrENV_in
	571	coefcoli(it)=coefcoli_in
[2146]	572	ELSE
[3450]	573	flag_cvltr(it)=.FALSE.
[2146]	574	ENDIF
	575
[3125]	576	CASE('repr')
[3450]	577	flag_cvltr(it)=.FALSE.
[3125]	578
[2146]	579	CASE('inca')
	580	! IF ((it.EQ.id_Rn222) .OR. ((it.GE.id_SO2) .AND. (it.LE.id_NH3)) ) THEN
	581	! !--gas-phase species
[3450]	582	! flag_cvltr(it)=.FALSE.
[1813]	583	!
[2146]	584	! ELSEIF ( (it.GE.id_CIDUSTM) .AND. (it.LE.id_AIN) ) THEN
	585	! !--insoluble aerosol species
[3450]	586	! flag_cvltr(it)=.TRUE.
[2146]	587	! ccntrAA(it)=0.7
	588	! ccntrENV(it)=0.7
	589	! coefcoli(it)=0.001
	590	! ELSEIF ( (it.EQ.id_Pb210) .OR. ((it.GE.id_CSSSM) .AND. (it.LE.id_SSN))) THEN
	591	! !--soluble aerosol species
[3450]	592	! flag_cvltr(it)=.TRUE.
[2146]	593	! ccntrAA(it)=0.9
	594	! ccntrENV(it)=0.9
	595	! coefcoli(it)=0.001
	596	! ELSE
	597	! WRITE(lunout,*) 'pb it=', it
[2311]	598	! CALL abort_physic('phytrac','pb it scavenging',1)
[2146]	599	! ENDIF
	600	!--test OB
	601	!--for now we do not scavenge in cvltr
[3450]	602	flag_cvltr(it)=.FALSE.
[2690]	603
[3361]	604	CASE('co2i')
	605	!--co2 tracers are not scavenged
[3450]	606	flag_cvltr(it)=.FALSE.
[3361]	607
[2690]	608	#ifdef CPP_StratAer
	609	CASE('coag')
	610	IF (convscav.and.aerosol(it)) THEN
[3450]	611	flag_cvltr(it)=.TRUE.
[2690]	612	ccntrAA(it) =ccntrAA_in
	613	ccntrENV(it)=ccntrENV_in
	614	coefcoli(it)=coefcoli_in
	615	ELSE
[3450]	616	flag_cvltr(it)=.FALSE.
[2690]	617	ENDIF
	618	#endif
	619
[2146]	620	END SELECT
	621	ENDDO
	622	!
	623	ELSE ! iflag_con .ne. 3
[3450]	624	flag_cvltr(:) = .FALSE.
[2146]	625	ENDIF
	626	!
	627	! Initialize diagnostic output
	628	! ----------------------------
[1813]	629	#ifdef CPP_IOIPSL
[2146]	630	! INCLUDE "ini_histrac.h"
[1813]	631	#endif
[2146]	632	!
	633	! print out all tracer flags
	634	!
	635	WRITE(lunout,*) 'print out all tracer flags'
	636	WRITE(lunout,*) 'type_trac =', type_trac
	637	WRITE(lunout,*) 'config_inca =', config_inca
	638	WRITE(lunout,*) 'iflag_con_trac =', iflag_con_trac
	639	WRITE(lunout,*) 'iflag_con =', iflag_con
	640	WRITE(lunout,*) 'convscav =', convscav
	641	WRITE(lunout,*) 'iflag_lscav =', iflag_lscav
	642	WRITE(lunout,*) 'aerosol =', aerosol
	643	WRITE(lunout,*) 'iflag_the_trac =', iflag_the_trac
	644	WRITE(lunout,*) 'iflag_thermals =', iflag_thermals
	645	WRITE(lunout,*) 'iflag_vdf_trac =', iflag_vdf_trac
	646	WRITE(lunout,*) 'pbl_flg =', pbl_flg
	647	WRITE(lunout,*) 'lessivage =', lessivage
	648	write(lunout,*) 'flag_cvltr = ', flag_cvltr
[1813]	649
[4795]	650	! IF (lessivage .AND. type_trac .EQ. 'inca') THEN
	651	! CALL abort_physic('phytrac', 'lessivage=T config_inca=inca impossible',1)
	652	! ! STOP
	653	! ENDIF
[2146]	654	!
[3450]	655	ENDIF ! of IF (debutphy)
[2146]	656	!############################################ END INITIALIZATION #######
[1813]	657
[2146]	658	DO k=1,klev
	659	DO i=1,klon
	660	zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
	661	END DO
	662	END DO
	663	!
	664	IF (id_be .GT. 0) THEN
	665	DO k=1,klev
	666	DO i=1,klon
	667	sourceBE(i,k)=srcbe(i,k) !RomP -> pour sortie histrac
	668	END DO
	669	END DO
	670	ENDIF
[1813]	671
[2146]	672	!===============================================================================
	673	! -- Do specific treatment according to chemestry model or local LMDZ tracers
	674	!
	675	!===============================================================================
	676	SELECT CASE(type_trac)
	677	CASE('lmdz')
	678	! -- Traitement des traceurs avec traclmdz
	679	CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
	680	cdragh, coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,iflag_vdf_trac>=0,sh, &
	681	rh, pphi, ustar, wstar, ale_bl, ale_wake, u10m, v10m, &
[2180]	682	tr_seri, source, d_tr_cl,d_tr_dec, zmasse) !RomP
[1813]	683
[2146]	684	CASE('inca')
	685	! -- CHIMIE INCA config_inca = aero or chem --
[2637]	686	! Appel fait en fin de phytrac pour avoir les emissions modifiees par
	687	! la couche limite et la convection avant le calcul de la chimie
[2690]	688
[2146]	689	CASE('repr')
	690	! -- CHIMIE REPROBUS --
	691	CALL tracreprobus(pdtphys, gmtime, debutphy, julien, &
	692	presnivs, xlat, xlon, pphis, pphi, &
	693	t_seri, pplay, paprs, sh , &
[2180]	694	tr_seri)
[1813]	695
[3361]	696	CASE('co2i')
	697	! -- CO2 interactif --
	698	! -- source is updated with FF and BB emissions
[3450]	699	! -- and net fluxes from ocean and orchidee
[3361]	700	! -- sign convention : positive into the atmosphere
[3450]	701
[3361]	702	CALL tracco2i(pdtphys, debutphy, &
	703	xlat, xlon, pphis, pphi, &
	704	t_seri, pplay, paprs, tr_seri, source)
	705
[2690]	706	#ifdef CPP_StratAer
	707	CASE('coag')
	708	! --STRATOSPHERIC AER IN THE STRAT --
	709	CALL traccoag(pdtphys, gmtime, debutphy, julien, &
	710	presnivs, xlat, xlon, pphis, pphi, &
	711	t_seri, pplay, paprs, sh, rh , &
	712	tr_seri)
	713	#endif
	714
[2146]	715	END SELECT
[4795]	716
	717	! Nicolas Février
	718
	719	!IF (aerosol_couple) THEN !--we get aerosols from tr_seri array from INCA
	720	!
	721	!--copy fields from INCA tr_seri
	722	!--convert to ug m-3 unit for consistency with offline fields
	723	!
	724	DO it=1,nbtr
	725	SELECT CASE(trim(solsym(it)))
	726	CASE ("ASAPp1a")
	727	id_ASAPp1a = it
	728	aerosol(it)=.TRUE.
	729	CASE ("ASAPp2a")
	730	id_ASAPp2a = it
	731	aerosol(it)=.TRUE.
	732	CASE ("ASARp1a")
	733	id_ASARp1a = it
	734	aerosol(it)=.TRUE.
	735	CASE ("ASARp2a")
	736	id_ASARp2a = it
	737	aerosol(it)=.TRUE.
	738	CASE ("ASBCM")
	739	id_ASBCM = it
	740	aerosol(it)=.TRUE.
	741	CASE ("AIBCM")
	742	id_AIBCM = it
	743	aerosol(it)=.TRUE.
	744	CASE ("CIDUSTM")
	745	id_CIDUSTM = it
	746	aerosol(it)=.TRUE.
	747	CASE ("AIPOMM")
	748	id_AIPOMM = it
	749	aerosol(it)=.TRUE.
	750	CASE ("ASPOMM")
	751	id_ASPOMM = it
	752	aerosol(it)=.TRUE.
	753	CASE ("ASSSM")
	754	id_ASSSM = it
	755	aerosol(it)=.TRUE.
	756	CASE ("CSSSM")
	757	id_CSSSM = it
	758	aerosol(it)=.TRUE.
	759	CASE ("SSSSM")
	760	id_SSSSM = it
	761	aerosol(it)=.TRUE.
	762	END SELECT
	763	ENDDO
	764	!END IF
	765
	766	list_aero_wetdep=(/ id_ASAPp1a,id_ASAPp2a,id_ASARp1a,id_ASARp2a,id_ASBCM,id_AIBCM,id_CIDUSTM,id_AIPOMM,id_ASPOMM,id_ASSSM,id_CSSSM,id_SSSSM /)
	767
	768	!WRITE(lunout,*) 'id_AIBCM=',id_AIBCM
	769	!WRITE(lunout,*) 'id_ASBCM=',id_ASBCM
	770	!WRITE(lunout,*) 'shape_lmdz tr_seri=',shape(tr_seri)
	771
[2146]	772	!======================================================================
	773	! -- Calcul de l'effet de la convection --
	774	!======================================================================
[1813]	775
[4795]	776	! Nicolas Février
	777
	778	!WRITE(lunout,*) 'CHECK1 TR_SERI AIBCM',SUM(tr_seri(:,:,id_AIBCM))
	779	!WRITE(lunout,*) 'CHECK1 TR_SERI ASBCM',SUM(tr_seri(:,:,id_ASBCM))
	780
[2146]	781	IF (iflag_con_trac==1) THEN
[1813]	782
[2146]	783	DO it=1, nbtr
	784	IF ( conv_flg(it) == 0 ) CYCLE
	785	IF (iflag_con.LT.2) THEN
	786	!--pas de transport convectif
[2690]	787	d_tr_cv(:,:,it)=0.
[1813]	788
[2146]	789	ELSE IF (iflag_con.EQ.2) THEN
	790	!--ancien transport convectif de Tiedtke
[1813]	791
[2146]	792	CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
	793	pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
	794	ELSE
	795	!--nouveau transport convectif de Emanuel
[4795]	796
	797	! Nicolas Février
[1813]	798
[4795]	799	IF ( ANY( list_aero_wetdep==it ) ) THEN
	800	flag_cvltr(it)=.TRUE.
	801	ccntrAA(it) =ccntrAA_in
	802	ccntrENV(it)=ccntrENV_in
	803	coefcoli(it)=coefcoli_in
	804	END IF
	805
[2146]	806	IF (flag_cvltr(it)) THEN
	807	!--nouveau transport convectif de Emanuel avec lessivage convectif
	808	!
	809	!
	810	ccntrAA_3d(:,:) =ccntrAA(it)
	811	ccntrENV_3d(:,:)=ccntrENV(it)
	812	coefcoli_3d(:,:)=coefcoli(it)
[1813]	813
[2146]	814	!--beware this interface is a bit weird because it is called for each tracer
	815	!--with the full array tr_seri even if only item it is processed
[1813]	816
[2210]	817	print*,'CV SCAV ',it,ccntrAA(it),ccntrENV(it)
	818
[2146]	819	CALL cvltr_scav(pdtphys, da, phi,phi2,d1a,dam, mp,ep, &
	820	sigd,sij,wght_cvfd,clw,elij,epmlmMm,eplaMm, &
	821	pmflxr,pmflxs,evap,t_seri,wdtrainA,wdtrainM, &
	822	paprs,it,tr_seri,upwd,dnwd,itop_con,ibas_con, &
	823	ccntrAA_3d,ccntrENV_3d,coefcoli_3d, &
	824	d_tr_cv,d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qDi,qPr,&
	825	qPa,qMel,qTrdi,dtrcvMA,Mint, &
	826	zmfd1a,zmfphi2,zmfdam)
[1813]	827
[2146]	828
	829	ELSE !---flag_cvltr(it).EQ.FALSE
	830	!--nouveau transport convectif de Emanuel mais pas de lessivage convectif
	831
	832	!--beware this interface is a bit weird because it is called for each tracer
	833	!--with the full array tr_seri even if only item it is processed
	834	!
	835	CALL cvltr_noscav(it,pdtphys, da, phi,mp,wght_cvfd,paprs,pplay, & !jyg
	836	tr_seri,upwd,dnwd,d_tr_cv) !jyg
	837
	838	ENDIF
	839
	840	ENDIF !--iflag
	841
	842	!--on ajoute les tendances
	843
	844	DO k = 1, klev
	845	DO i = 1, klon
	846	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
	847	END DO
	848	END DO
	849
[4795]	850	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it),it) ! Nicolas Février it rajouté en argument
[2146]	851
[4795]	852	!WRITE(lunout,*) 'CHECK2 TR_SERI AIBCM',SUM(tr_seri(:,:,id_AIBCM))
	853	!WRITE(lunout,*) 'CHECK2 TR_SERI ASBCM',SUM(tr_seri(:,:,id_ASBCM))
	854
[2146]	855	END DO ! nbtr
	856
[2690]	857	#ifdef CPP_StratAer
	858	IF (type_trac=='coag') THEN
	859	! initialize wet deposition flux of sulfur
[3100]	860	budg_dep_wet_ocs(:)=0.0
	861	budg_dep_wet_so2(:)=0.0
[2752]	862	budg_dep_wet_h2so4(:)=0.0
	863	budg_dep_wet_part(:)=0.0
[2690]	864	! compute wet deposition flux of sulfur (sum over gases and particles)
	865	! and convert to kg(S)/m2/s
	866	DO i = 1, klon
	867	DO k = 1, klev
	868	DO it = 1, nbtr
	869	!do not include SO2 because most of it comes trom the troposphere
[3100]	870	IF (it==id_OCS_strat) THEN
	871	budg_dep_wet_ocs(i)=budg_dep_wet_ocs(i)+d_tr_cv(i,k,it)*(mSatom/mOCSmol) &
	872	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	873	ELSEIF (it==id_SO2_strat) THEN
	874	budg_dep_wet_so2(i)=budg_dep_wet_so2(i)+d_tr_cv(i,k,it)*(mSatom/mSO2mol) &
	875	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	876	ELSEIF (it==id_H2SO4_strat) THEN
[2752]	877	budg_dep_wet_h2so4(i)=budg_dep_wet_h2so4(i)+d_tr_cv(i,k,it)*(mSatom/mH2SO4mol) &
[2690]	878	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	879	ELSEIF (it.GT.nbtr_sulgas) THEN
[2752]	880	budg_dep_wet_part(i)=budg_dep_wet_part(i)+d_tr_cv(i,k,it)*(mSatom/mH2SO4mol) &
[2690]	881	& dens_aer_dry4./3.RPI(mdw(it-nbtr_sulgas)/2.)**3 &
	882	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	883	ENDIF
	884	ENDDO
	885	ENDDO
	886	ENDDO
	887	ENDIF
	888	#endif
	889
[3450]	890	ENDIF ! convection
[2146]	891
	892	!======================================================================
	893	! -- Calcul de l'effet des thermiques --
	894	!======================================================================
	895
[1813]	896	DO it=1,nbtr
[2146]	897	DO k=1,klev
	898	DO i=1,klon
	899	d_tr_th(i,k,it)=0.
	900	tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
[3099]	901	! the next safeguard causes some problem for stratospheric aerosol tracers (particle number)
[3109]	902	! and there is little justification for it so it is commented out (4 December 2017) by OB
	903	! if reinstated please keep the ifndef CPP_StratAer
	904	!#ifndef CPP_StratAer
	905	! tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
	906	!#endif
[2146]	907	END DO
	908	END DO
	909	END DO
[1813]	910
[2146]	911	IF (iflag_thermals.GT.0.AND.iflag_the_trac>0) THEN
[1813]	912
[2146]	913	DO it=1, nbtr
[1813]	914
[2146]	915	CALL thermcell_dq(klon,klev,1,pdtphys,fm_therm,entr_therm, &
	916	zmasse,tr_seri(1:klon,1:klev,it), &
	917	d_tr_th(1:klon,1:klev,it),ztra_th,0 )
[1813]	918
[2146]	919	DO k=1,klev
	920	DO i=1,klon
	921	d_tr_th(i,k,it)=pdtphys*d_tr_th(i,k,it)
	922	tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr_th(i,k,it),0.)
	923	END DO
	924	END DO
[1813]	925
[2146]	926	END DO ! it
[1813]	927
[3450]	928	ENDIF ! Thermiques
[1813]	929
[2146]	930	!======================================================================
	931	! -- Calcul de l'effet de la couche limite --
	932	!======================================================================
	933
	934	IF (iflag_vdf_trac==1) THEN
	935
	936	! Injection during BL mixing
	937	!
[2690]	938	#ifdef CPP_StratAer
	939	IF (type_trac=='coag') THEN
	940
	941	! initialize dry deposition flux of sulfur
[3100]	942	budg_dep_dry_ocs(:)=0.0
	943	budg_dep_dry_so2(:)=0.0
[2752]	944	budg_dep_dry_h2so4(:)=0.0
	945	budg_dep_dry_part(:)=0.0
[2690]	946
	947	! compute dry deposition velocity as function of surface type (numbers
	948	! from IPSL note 23, 2002)
	949	v_dep_dry(:) = pctsrf(:,is_ter) * 2.5e-3 &
	950	& + pctsrf(:,is_oce) * 0.5e-3 &
	951	& + pctsrf(:,is_lic) * 2.5e-3 &
	952	& + pctsrf(:,is_sic) * 2.5e-3
	953
	954	! compute surface dry deposition flux
	955	zrho(:,1)=pplay(:,1)/t_seri(:,1)/RD
	956
	957	DO it=1, nbtr
	958	source(:,it) = - v_dep_dry(:) * tr_seri(:,1,it) * zrho(:,1)
	959	ENDDO
	960
	961	ENDIF
	962	#endif
	963
[2146]	964	DO it=1, nbtr
[4795]	965	!WRITE(lunout,*) 'BOUM'
[2146]	966	IF( pbl_flg(it) /= 0 ) THEN
[4795]	967	!WRITE(lunout,*) 'BADABOUM'
[2146]	968	!
	969	CALL cltrac(pdtphys, coefh,t_seri, &
	970	tr_seri(:,:,it), source(:,it), &
	971	paprs, pplay, delp, &
	972	d_tr_cl(:,:,it),d_tr_dry(:,it),flux_tr_dry(:,it))
	973	!
	974	tr_seri(:,:,it)=tr_seri(:,:,it)+d_tr_cl(:,:,it)
	975	!
[2690]	976	#ifdef CPP_StratAer
	977	IF (type_trac=='coag') THEN
	978	! compute dry deposition flux of sulfur (sum over gases and particles)
[3100]	979	IF (it==id_OCS_strat) THEN
	980	budg_dep_dry_ocs(:)=budg_dep_dry_ocs(:)-source(:,it)*(mSatom/mOCSmol)
	981	ELSEIF (it==id_SO2_strat) THEN
	982	budg_dep_dry_so2(:)=budg_dep_dry_so2(:)-source(:,it)*(mSatom/mSO2mol)
	983	ELSEIF (it==id_H2SO4_strat) THEN
[2752]	984	budg_dep_dry_h2so4(:)=budg_dep_dry_h2so4(:)-source(:,it)*(mSatom/mH2SO4mol)
[2690]	985	ELSEIF (it.GT.nbtr_sulgas) THEN
[2752]	986	budg_dep_dry_part(:)=budg_dep_dry_part(:)-source(:,it)(mSatom/mH2SO4mol)dens_aer_dry &
[2690]	987	& 4./3.RPI(mdw(it-nbtr_sulgas)/2.)*3
	988	ENDIF
	989	ENDIF
	990	#endif
	991	!
	992	ENDIF
[2146]	993	!
[2690]	994	ENDDO
[2146]	995	!
	996	ELSE IF (iflag_vdf_trac==0) THEN
[4795]	997
	998	!WRITE(lunout,*) 'TRALALA'
[2146]	999	!
	1000	! Injection of source in the first model layer
	1001	!
	1002	DO it=1,nbtr
	1003	d_tr_cl(:,1,it)=source(:,it)RG/delp(:,1)pdtphys
	1004	tr_seri(:,1,it)=tr_seri(:,1,it)+d_tr_cl(:,1,it)
	1005	ENDDO
	1006	d_tr_cl(:,2:klev,1:nbtr)=0.
	1007	!
	1008	ELSE IF (iflag_vdf_trac==-1) THEN
[4795]	1009
	1010	!WRITE(lunout,*) 'BOFBOFBOF'
[2146]	1011	!
	1012	! Nothing happens
	1013	d_tr_cl=0.
	1014	!
	1015	ELSE
	1016	!
[2311]	1017	CALL abort_physic('iflag_vdf_trac', 'cas non prevu',1)
[2146]	1018	!
[3450]	1019	ENDIF ! couche limite
[2146]	1020
	1021	!======================================================================
	1022	! Calcul de l'effet de la precipitation grande echelle
	1023	! POUR INCA le lessivage est fait directement dans INCA
	1024	!======================================================================
[4795]	1025
	1026	! Nicolas Février
[2146]	1027
[4795]	1028	! id_ASAPp1a,id_ASAPp2a,id_ASARp1a,id_ASARp2a,id_ASBCM,id_AIBCM,id_CIDUSTM,id_AIPOMM,id_ASPOMM,id_ASSSM,id_CSSSM,id_SSSSM
	1029
	1030	! IF (aerosol_couple) THEN !--we get aerosols from tr_seri array from INCA
	1031	! !
	1032	! !--copy fields from INCA tr_seri
	1033	! !--convert to ug m-3 unit for consistency with offline fields
	1034	! !
	1035	! DO it=1,nbtr
	1036	! SELECT CASE(trim(solsym(it)))
	1037	! CASE ("ASAPp1a")
	1038	! id_ASAPp1a = it
	1039	! aerosol(it)=.TRUE.
	1040	! CASE ("ASAPp2a")
	1041	! id_ASAPp2a = it
	1042	! aerosol(it)=.TRUE.
	1043	! CASE ("ASARp1a")
	1044	! id_ASARp1a = it
	1045	! aerosol(it)=.TRUE.
	1046	! CASE ("ASARp2a")
	1047	! id_ASARp2a = it
	1048	! aerosol(it)=.TRUE.
	1049	! CASE ("ASBCM")
	1050	! id_ASBCM = it
	1051	! aerosol(it)=.TRUE.
	1052	! CASE ("AIBCM")
	1053	! id_AIBCM = it
	1054	! aerosol(it)=.TRUE.
	1055	! CASE ("CIDUSTM")
	1056	! id_CIDUSTM = it
	1057	! aerosol(it)=.TRUE.
	1058	! CASE ("AIPOMM")
	1059	! id_AIPOMM = it
	1060	! aerosol(it)=.TRUE.
	1061	! CASE ("ASPOMM")
	1062	! id_ASPOMM = it
	1063	! aerosol(it)=.TRUE.
	1064	! CASE ("ASSSM")
	1065	! id_ASSSM = it
	1066	! aerosol(it)=.TRUE.
	1067	! CASE ("CSSSM")
	1068	! id_CSSSM = it
	1069	! aerosol(it)=.TRUE.
	1070	! CASE ("SSSSM")
	1071	! id_SSSSM = it
	1072	! aerosol(it)=.TRUE.
	1073	! END SELECT
	1074	! ENDDO
	1075	! END IF
	1076
	1077	! list_aero_wetdep=(/ id_ASAPp1a,id_ASAPp2a,id_ASARp1a,id_ASARp2a,id_ASBCM,id_AIBCM,id_CIDUSTM,id_AIPOMM,id_ASPOMM,id_ASSSM,id_CSSSM,id_SSSSM /)
	1078
	1079	lessivage = .TRUE.
	1080
	1081	!WRITE(lunout,*) 'CHECK3 TR_SERI AIBCM',SUM(tr_seri(:,:,id_AIBCM))
	1082	!WRITE(lunout,*) 'CHECK3 TR_SERI ASBCM',SUM(tr_seri(:,:,id_ASBCM))
	1083
	1084
[2146]	1085	IF (lessivage) THEN
	1086
	1087	ql_incloud_ref = 10.e-4
	1088	ql_incloud_ref = 5.e-4
	1089
	1090
	1091	! calcul du contenu en eau liquide au sein du nuage
	1092	ql_incl = ql_incloud_ref
	1093	! choix du lessivage
	1094	!
[4795]	1095
	1096	!WRITE(lunout,*) 'Nico iflag_lscav=',iflag_lscav
	1097
[2146]	1098	IF (iflag_lscav .EQ. 3 .OR. iflag_lscav .EQ. 4) THEN
	1099	! ******** Olivier Boucher version (3) possibly with modified ql_incl (4)
	1100	!
[4795]	1101	DO j = 1,size(list_aero_wetdep)
	1102	it=list_aero_wetdep(j)
[2210]	1103
	1104	IF (aerosol(it)) THEN
[2146]	1105	! incloud scavenging and removal by large scale rain ! orig : ql_incl was replaced by 0.5e-3 kg/kg
	1106	! the value 0.5e-3 kg/kg is from Giorgi and Chameides (1986), JGR
	1107	! Liu (2001) proposed to use 1.5e-3 kg/kg
	1108
[4795]	1109	!WRITE(lunout,*) 'On_lessive_id = ',it
	1110
[2284]	1111	CALL lsc_scav(pdtphys,it,iflag_lscav,aerosol,ql_incl,prfl,psfl,rneb,beta_fisrt, &
[1813]	1112	beta_v1,pplay,paprs,t_seri,tr_seri,d_tr_insc, &
[4795]	1113	d_tr_bcscav,d_tr_evapls,qPrls,beta_st,j)
[1813]	1114
[2146]	1115	!large scale scavenging tendency
[4795]	1116	! Nicolas Février : commenté pour ne pas faire le lessivage ici. Il est à la place fait dans SETHET.
	1117	!DO k = 1, klev
	1118	! DO i = 1, klon
	1119	! d_tr_ls(i,k,it)=d_tr_insc(i,k,it)+d_tr_bcscav(i,k,it)+d_tr_evapls(i,k,it)
	1120	! tr_seri(i,k,it)=tr_seri(i,k,it)+d_tr_ls(i,k,it)
	1121	! ENDDO
	1122	!ENDDO
	1123	!CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'lsc scav it = '//solsym(it),it) ! Nicolas Février it rajouté en argument
[2210]	1124	ENDIF
	1125
[2146]	1126	END DO !tr
[1813]	1127
[4795]	1128	! Debug beta_st
	1129	print *,"BETA PHYTRAC",beta_st(5000,35,7)
	1130
	1131
[2690]	1132	#ifdef CPP_StratAer
	1133	IF (type_trac=='coag') THEN
	1134	! compute wet deposition flux of sulfur (sum over gases and
	1135	! particles) and convert to kg(S)/m2/s
	1136	! adding contribution of d_tr_ls to d_tr_cv (above)
	1137	DO i = 1, klon
	1138	DO k = 1, klev
	1139	DO it = 1, nbtr
[3100]	1140	IF (it==id_OCS_strat) THEN
	1141	budg_dep_wet_ocs(i)=budg_dep_wet_ocs(i)+d_tr_ls(i,k,it)*(mSatom/mOCSmol) &
	1142	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	1143	ELSEIF (it==id_SO2_strat) THEN
	1144	budg_dep_wet_so2(i)=budg_dep_wet_so2(i)+d_tr_ls(i,k,it)*(mSatom/mSO2mol) &
	1145	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	1146	ELSEIF (it==id_H2SO4_strat) THEN
[2752]	1147	budg_dep_wet_h2so4(i)=budg_dep_wet_h2so4(i)+d_tr_ls(i,k,it)*(mSatom/mH2SO4mol) &
[2690]	1148	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	1149	ELSEIF (it.GT.nbtr_sulgas) THEN
[2752]	1150	budg_dep_wet_part(i)=budg_dep_wet_part(i)+d_tr_ls(i,k,it)*(mSatom/mH2SO4mol) &
[2690]	1151	& dens_aer_dry4./3.RPI(mdw(it-nbtr_sulgas)/2.)**3 &
	1152	& *(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
	1153	ENDIF
	1154	ENDDO
	1155	ENDDO
	1156	ENDDO
	1157	ENDIF
	1158	#endif
	1159
[2146]	1160	ELSE IF (iflag_lscav .EQ. 2) THEN ! frac_impa, frac_nucl
	1161	! ********* modified old version
[1813]	1162
[2146]	1163	d_tr_lessi_nucl(:,:,:) = 0.
	1164	d_tr_lessi_impa(:,:,:) = 0.
	1165	flestottr(:,:,:) = 0.
	1166	! Tendance des aerosols nuclees et impactes
	1167	DO it = 1, nbtr
	1168	IF (aerosol(it)) THEN
	1169	his_dh(:)=0.
	1170	DO k = 1, klev
	1171	DO i = 1, klon
	1172	!PhH
	1173	zrho(i,k)=pplay(i,k)/t_seri(i,k)/RD
	1174	zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
	1175	!
[3450]	1176	ENDDO
	1177	ENDDO
[1813]	1178
[2146]	1179	DO k=klev-1, 1, -1
	1180	DO i=1, klon
	1181	! d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.)
	1182	dx=d_tr_ls(i,k,it)
	1183	his_dh(i)=his_dh(i)-dxzrho(i,k)zdz(i,k)/pdtphys ! kg/m2/s
	1184	evaplsc = prfl(i,k) - prfl(i,k+1) + psfl(i,k) - psfl(i,k+1)
	1185	! Evaporation Partielle -> Liberation Partielle 0.5*evap
	1186	IF ( evaplsc .LT.0..and.abs(prfl(i,k+1)+psfl(i,k+1)).gt.1.e-10) THEN
	1187	evaplsc = (-evaplsc)/(prfl(i,k+1)+psfl(i,k+1))
	1188	! evaplsc est donc positif, his_dh(i) est positif
	1189	!--------------
	1190	d_tr_evapls(i,k,it)=0.5evaplsc(d_tr_lessi_nucl(i,k+1,it) &
	1191	+d_tr_lessi_impa(i,k+1,it))
	1192	!------------- d_tr_evapls(i,k,it)=-0.5evaplsc(d_tr_lsc(i,k+1,it))
	1193	beta=0.5*evaplsc
	1194	if ((prfl(i,k)+psfl(i,k)).lt.1.e-10) THEN
	1195	beta=1.0*evaplsc
	1196	endif
	1197	dx=betahis_dh(i)/zrho(i,k)/zdz(i,k)pdtphys
	1198	his_dh(i)=(1.-beta)*his_dh(i) ! tracer from
	1199	d_tr_evapls(i,k,it)=dx
	1200	ENDIF
	1201	d_tr_ls(i,k,it)=tr_seri(i,k,it)(frac_impa(i,k)frac_nucl(i,k)-1.) &
	1202	+d_tr_evapls(i,k,it)
[1813]	1203
[2146]	1204	!--------------
	1205	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
	1206	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
	1207	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
	1208	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
	1209	!
	1210	! Flux lessivage total
	1211	flestottr(i,k,it) = flestottr(i,k,it) - &
	1212	( d_tr_lessi_nucl(i,k,it) + &
	1213	d_tr_lessi_impa(i,k,it) ) * &
	1214	( paprs(i,k)-paprs(i,k+1) ) / &
	1215	(RG * pdtphys)
	1216	!! Mise a jour des traceurs due a l'impaction,nucleation
	1217	! tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
	1218	!! calcul de la tendance liee au lessivage stratiforme
	1219	! d_tr_ls(i,k,it)=tr_seri(i,k,it)*&
	1220	! (1.-1./(frac_impa(i,k)*frac_nucl(i,k)))
	1221	!--------------
[3450]	1222	ENDDO
	1223	ENDDO
	1224	ENDIF
	1225	ENDDO
[2146]	1226	! ********* end modified old version
[1813]	1227
[2146]	1228	ELSE IF (iflag_lscav .EQ. 1) THEN ! frac_impa, frac_nucl
	1229	! ********* old version
[1813]	1230
[2146]	1231	d_tr_lessi_nucl(:,:,:) = 0.
	1232	d_tr_lessi_impa(:,:,:) = 0.
	1233	flestottr(:,:,:) = 0.
	1234	!=========================
	1235	! LESSIVAGE LARGE SCALE :
	1236	!=========================
[1813]	1237
[2146]	1238	! Tendance des aerosols nuclees et impactes
	1239	! -----------------------------------------
	1240	DO it = 1, nbtr
	1241	IF (aerosol(it)) THEN
	1242	DO k = 1, klev
	1243	DO i = 1, klon
	1244	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
	1245	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
	1246	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
	1247	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
[1813]	1248
[2146]	1249	!
	1250	! Flux lessivage total
	1251	! ------------------------------------------------------------
	1252	flestottr(i,k,it) = flestottr(i,k,it) - &
	1253	( d_tr_lessi_nucl(i,k,it) + &
	1254	d_tr_lessi_impa(i,k,it) ) * &
	1255	( paprs(i,k)-paprs(i,k+1) ) / &
	1256	(RG * pdtphys)
	1257	!
	1258	! Mise a jour des traceurs due a l'impaction,nucleation
	1259	! ----------------------------------------------------------------------
	1260	tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
[3450]	1261	ENDDO
	1262	ENDDO
	1263	ENDIF
	1264	ENDDO
[1813]	1265
[2146]	1266	! ********* end old version
	1267	ENDIF ! iflag_lscav . EQ. 1, 2, 3 or 4
	1268	!
[4795]	1269	!write(lunout,*) 'THATSALLFOLKS'
[3450]	1270	ENDIF ! lessivage
[2146]	1271
[4795]	1272	!WRITE(lunout,*) 'CHECK4 TR_SERI AIBCM',SUM(tr_seri(:,:,id_AIBCM))
	1273	!WRITE(lunout,*) 'CHECK4 TR_SERI ASBCM',SUM(tr_seri(:,:,id_ASBCM))
[2637]	1274
[4795]	1275
	1276
[2637]	1277	! -- CHIMIE INCA config_inca = aero or chem --
	1278	IF (type_trac == 'inca') THEN
	1279
	1280	CALL tracinca(&
	1281	nstep, julien, gmtime, lafin, &
	1282	pdtphys, t_seri, paprs, pplay, &
	1283	pmfu, upwd, ftsol, pctsrf, pphis, &
[2784]	1284	pphi, albsol, sh, ch, rh, &
[2637]	1285	cldfra, rneb, diafra, cldliq, &
	1286	itop_con, ibas_con, pmflxr, pmflxs, &
	1287	prfl, psfl, aerosol_couple, flxmass_w, &
	1288	tau_aero, piz_aero, cg_aero, ccm, &
	1289	rfname, &
[4795]	1290	tr_seri, source,beta_st)
[2637]	1291
	1292
	1293	ENDIF
[4795]	1294
	1295	!WRITE(lunout,*) 'CHECK5 TR_SERI AIBCM',SUM(tr_seri(:,:,id_AIBCM))
	1296	!WRITE(lunout,*) 'CHECK5 TR_SERI ASBCM',SUM(tr_seri(:,:,id_ASBCM))
	1297
	1298
	1299	!write(lunout,*) 'MERCI INCA POUR LES SERVICES'
[2146]	1300	!=============================================================
	1301	! Ecriture des sorties
	1302	!=============================================================
[1813]	1303	#ifdef CPP_IOIPSL
[2146]	1304	! INCLUDE "write_histrac.h"
[1813]	1305	#endif
	1306
[2146]	1307	END SUBROUTINE phytrac
[1813]	1308
	1309	END MODULE

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

Download in other formats:

Original Format