Context Navigation

phytrac.F90 @ 1191

Last change on this file since 1191 was 1191, checked in by jghattas, 15 years ago

Reecriture de phytrac et les routines concernes (Anthony Jamelot)

les suffix change de F -> F90 (nflxtr.F90,cltracrn.F90,initrrnpb.F90,cvltr.F90,minmaxqfi.F90,cltrac.F90,phytrac.F90)

Traitement d'un nouveau traceur berelium (optionel, toujours pour des
tests)(Anthony Jamelot)

radiornpb.F change du nom pour radio_decay.F90 car il traite maintenant tout les traceurs radioactives
ajoute init_be.F90

Nouveau interface dans phytrac pour serparer les calculs et appels
specifique a INCA avec les traitements des traceurs specifiques au LMDZ
(JG)

ajoute tracinca_mod.F90 pour les appeles a INCA
ajoute traclmdz_mod.F90 pour les calculs des traceurs specifiques a LMDZ
enleve fichier restartrac et ajoute la variable trs dans restartphy.nc

La convergence numerique a etait rompue uniquement pour les traceurs
LMDZ RN et PB.

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 14.9 KB

Line
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)
20	!
21	!======================================================================
22	! Auteur(s) FH
23	! Objet: Moniteur general des tendances traceurs
24	!======================================================================
25
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
35
36
37	IMPLICIT NONE
38
39	INCLUDE "YOMCST.h"
40	INCLUDE "dimensions.h"
41	INCLUDE "indicesol.h"
42	INCLUDE "clesphys.h"
43	INCLUDE "temps.h"
44	INCLUDE "paramet.h"
45	INCLUDE "control.h"
46	INCLUDE "thermcell.h"
47	!==========================================================================
48	! -- ARGUMENT DESCRIPTION --
49	!==========================================================================
50
51	! Input arguments
52	!----------------
53	!Configuration grille,temps:
54	INTEGER,INTENT(IN) :: nstep ! Appel physique
55	INTEGER,INTENT(IN) :: julien ! Jour julien
56	REAL,INTENT(IN) :: gmtime
57	REAL,INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
58	LOGICAL,INTENT(IN) :: debutphy ! le flag de l'initialisation de la physique
59	LOGICAL,INTENT(IN) :: lafin ! le flag de la fin de la physique
60
61	REAL,DIMENSION(klon),INTENT(IN) :: xlat ! latitudes pour chaque point
62	REAL,DIMENSION(klon),INTENT(IN) :: xlon ! longitudes pour chaque point
63	!
64	!Physique:
65	!--------
66	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
67	REAL,DIMENSION(klon,klev),INTENT(IN) :: u !
68	REAL,DIMENSION(klon,klev),INTENT(IN) :: v !
69	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
70	REAL,DIMENSION(klon,klev),INTENT(IN) :: rh ! humidite relative
71	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
72	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
73	REAL,DIMENSION(klon,klev),INTENT(IN) :: pphi ! geopotentiel
74	REAL,DIMENSION(klon),INTENT(IN) :: pphis
75	REAL,DIMENSION(klev),INTENT(IN) :: presnivs
76	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldliq ! eau liquide nuageuse
77	REAL,DIMENSION(klon,klev),INTENT(IN) :: cldfra ! fraction nuageuse (tous les nuages)
78	REAL,DIMENSION(klon,klev),INTENT(IN) :: diafra ! fraction nuageuse (convection ou stratus artificiels)
79	REAL,DIMENSION(klon,klev),INTENT(IN) :: rneb ! fraction nuageuse (grande echelle)
80	INTEGER,DIMENSION(klon),INTENT(IN) :: itop_con
81	INTEGER,DIMENSION(klon),INTENT(IN) :: ibas_con
82	REAL,DIMENSION(klon),INTENT(IN) :: albsol ! albedo surface
83	!
84	!Convection:
85	!----------
86	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu ! flux de masse dans le panache montant
87	REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd ! flux de masse dans le panache descendant
88	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
89	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
90	REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
91	REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
92
93	!...Tiedke
94	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
95	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
96
97	LOGICAL,INTENT(IN) :: aerosol_couple
98	REAL,DIMENSION(klon,klev),INTENT(IN) :: flxmass_w
99	REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: tau_aero
100	REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: piz_aero
101	REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: cg_aero
102	CHARACTER(len=4),DIMENSION(9),INTENT(IN) :: rfname
103	REAL,DIMENSION(klon,klev,2),INTENT(IN) :: ccm
104	!... K.Emanuel
105	REAL,DIMENSION(klon,klev),INTENT(IN) :: da
106	REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
107	REAL,DIMENSION(klon,klev),INTENT(IN) :: mp
108	REAL,DIMENSION(klon,klev),INTENT(IN) :: upwd ! saturated updraft mass flux
109	REAL,DIMENSION(klon,klev),INTENT(IN) :: dnwd ! saturated downdraft mass flux
110	!
111	!Thermiques:
112	!----------
113	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: fm_therm
114	REAL,DIMENSION(klon,klev),INTENT(IN) :: entr_therm
115	!
116	!Couche limite:
117	!--------------
118	!
119	REAL,DIMENSION(klon,klev),INTENT(IN) :: cdragh ! coeff drag pour T et Q
120	REAL,DIMENSION(klon,klev),INTENT(IN) :: coefh ! coeff melange CL (m**2/s)
121	REAL,DIMENSION(klon),INTENT(IN) :: yu1 ! vents au premier niveau
122	REAL,DIMENSION(klon),INTENT(IN) :: yv1 ! vents au premier niveau
123	!
124	!Lessivage:
125	!----------
126	!
127	! pour le ON-LINE
128	!
129	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
130	REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
131
132	! Arguments necessaires pour les sources et puits de traceur:
133	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol ! Temperature du sol (surf)(Kelvin)
134	REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
135
136
137	! Output argument
138	!----------------
139	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
140
141	!=======================================================================================
142	! -- LOCAL VARIABLES --
143	!=======================================================================================
144
145	INTEGER :: i, k, it
146	INTEGER :: nsplit
147
148	!Sources et Reservoirs de traceurs (ex:Radon):
149	!--------------------------------------------
150	!
151	REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source ! a voir lorsque le flux de surface est prescrit
152	!$OMP THREADPRIVATE(source)
153
154	!
155	!Entrees/Sorties: (cf ini_histrac.h et write_histrac.h)
156	!---------------
157	INTEGER :: iiq, ierr
158	INTEGER :: nhori, nvert
159	REAL :: zsto, zout, zjulian
160	INTEGER,SAVE :: nid_tra ! pointe vers le fichier histrac.nc
161	!$OMP THREADPRIVATE(nid_tra)
162	REAL,DIMENSION(klon) :: zx_tmp_fi2d ! variable temporaire grille physique
163	INTEGER :: itau_w ! pas de temps ecriture = nstep + itau_phy
164	LOGICAL,PARAMETER :: ok_sync=.TRUE.
165
166	!
167	! Nature du traceur
168	!------------------
169	LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol ! aerosol(it) = true => aerosol => lessivage
170	!$OMP THREADPRIVATE(aerosol) ! aerosol(it) = false => gaz
171	REAL,DIMENSION(klon,klev) :: delp ! epaisseur de couche (Pa)
172	!
173	! Tendances de traceurs (Td):
174	!------------------------
175	!
176	REAL,DIMENSION(klon,klev) :: d_tr ! Td dans l'atmosphere
177	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_cl ! Td couche limite/traceur
178	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_cv ! Td convection/traceur
179	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_th ! Td thermique
180	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_lessi_impa ! Td du lessivage par impaction
181	REAL,DIMENSION(klon,klev,nbtr) :: d_tr_lessi_nucl ! Td du lessivage par nucleation
182	!
183	! Physique
184	!----------
185	REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
186	REAL,DIMENSION(klon,klev) :: zmasse ! densité atmosphérique Kg/m2
187	REAL,DIMENSION(klon,klev) :: ztra_th
188
189	!Controles:
190	!---------
191	LOGICAL,SAVE :: couchelimite=.TRUE.
192	LOGICAL,SAVE :: convection=.TRUE.
193	LOGICAL,SAVE :: lessivage
194	!$OMP THREADPRIVATE(couchelimite,convection,lessivage)
195
196	CHARACTER(len=8),DIMENSION(nbtr) :: solsym
197
198
199	!######################################################################
200	! -- INITIALIZATION --
201	!######################################################################
202	IF (debutphy) THEN
203	WRITE(,) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
204	ALLOCATE( source(klon,nbtr), stat=ierr)
205	IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 1',1)
206
207	ALLOCATE( aerosol(nbtr), stat=ierr)
208	IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 2',1)
209
210
211	! Initialize module for specific tracers
212	SELECT CASE(type_trac)
213	CASE('lmdz')
214	CALL traclmdz_init(pctsrf, ftsol, aerosol, lessivage)
215	CASE('inca')
216	source(:,:)=0.
217	CALL tracinca_init(aerosol,lessivage)
218	END SELECT
219	!
220	! Initialize diagnostic output
221	! ----------------------------
222	INCLUDE "ini_histrac.h"
223
224	END IF
225	!############################################ END INITIALIZATION #######
226
227	!===============================================================================
228	! -- Do specific treatment according to chemestry model or local LMDZ tracers
229	!
230	!===============================================================================
231	SELECT CASE(type_trac)
232	CASE('lmdz')
233	! -- Traitement des traceurs avec traclmdz
234
235	CALL traclmdz(&
236	pdtphys, t_seri, &
237	paprs, pplay, cdragh, coefh, &
238	yu1, yv1, ftsol, pctsrf, &
239	xlat, couchelimite, &
240	tr_seri, source, solsym, d_tr_cl)
241
242	CASE('inca')
243	! -- CHIMIE INCA config_inca = aero or chem --
244
245	CALL tracinca(&
246	nstep, julien, gmtime, lafin, &
247	pdtphys, t_seri, paprs, pplay, &
248	pmfu, ftsol, pctsrf, pphis, &
249	pphi, albsol, sh, rh, &
250	cldfra, rneb, diafra, cldliq, &
251	itop_con, ibas_con, pmflxr, pmflxs, &
252	prfl, psfl, aerosol_couple, flxmass_w, &
253	tau_aero, piz_aero, cg_aero, ccm, &
254	rfname, &
255	tr_seri, source, solsym)
256	END SELECT
257
258	!======================================================================
259	! -- Calcul de l'effet de la convection --
260	!======================================================================
261	IF (convection) THEN
262	DO it=1, nbtr
263	IF ( conv_flg(it) == 0 ) CYCLE
264
265	IF (iflag_con.LT.2) THEN
266	d_tr_cv(:,:,:)=0.
267	ELSE IF (iflag_con.EQ.2) THEN
268	!..Tiedke
269	CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
270	pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
271	ELSE
272	!..K.Emanuel
273	CALL cvltr(pdtphys, da, phi, mp, paprs,pplay, tr_seri(:,:,it),&
274	upwd,dnwd,d_tr_cv(:,:,it))
275	END IF
276
277	DO k = 1, klev
278	DO i = 1, klon
279	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
280	END DO
281	END DO
282
283	CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
284
285	END DO ! nbtr
286	END IF ! convection
287
288	!======================================================================
289	! -- Calcul de l'effet des thermiques --
290	!======================================================================
291
292	DO k=1,klev
293	DO i=1,klon
294	zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
295	END DO
296	END DO
297
298	DO it=1,nbtr
299	DO k=1,klev
300	DO i=1,klon
301	d_tr_th(i,k,it)=0.
302	tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
303	tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
304	END DO
305	END DO
306	END DO
307
308	IF (iflag_thermals.GT.0) THEN
309	nsplit=10
310	DO it=1, nbtr
311	DO isplit=1,nsplit
312
313	CALL dqthermcell(klon,klev,pdtphys/nsplit, &
314	fm_therm,entr_therm,zmasse, &
315	tr_seri(1:klon,1:klev,it),d_tr,ztra_th)
316
317	DO k=1,klev
318	DO i=1,klon
319	d_tr(i,k)=pdtphys*d_tr(i,k)/nsplit
320	d_tr_th(i,k,it)=d_tr_th(i,k,it)+d_tr(i,k)
321	tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr(i,k),0.)
322	END DO
323	END DO
324	END DO ! nsplit
325	END DO ! it
326	END IF ! Thermiques
327
328	!======================================================================
329	! -- Calcul de l'effet de la couche limite --
330	!======================================================================
331
332	IF (couchelimite) THEN
333
334	DO k = 1, klev
335	DO i = 1, klon
336	delp(i,k) = paprs(i,k)-paprs(i,k+1)
337	END DO
338	END DO
339
340	DO it=1, nbtr
341
342	IF( pbl_flg(it) /= 0 ) THEN
343
344	CALL cltrac(pdtphys, coefh,t_seri, &
345	tr_seri(:,:,it), source(:,it), &
346	paprs, pplay, delp, &
347	d_tr_cl(:,:,it))
348
349	DO k = 1, klev
350	DO i = 1, klon
351	tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cl(i,k,it)
352	END DO
353	END DO
354	END IF
355
356	END DO
357
358	END IF ! couche limite
359
360
361	!======================================================================
362	! Calcul de l'effet de la precipitation
363	!======================================================================
364
365	IF (lessivage) THEN
366
367	d_tr_lessi_nucl(:,:,:) = 0.
368	d_tr_lessi_impa(:,:,:) = 0.
369	flestottr(:,:,:) = 0.
370	!=========================
371	! LESSIVAGE LARGE SCALE :
372	!=========================
373
374	! Tendance des aerosols nuclees et impactes
375	! -----------------------------------------
376	DO it = 1, nbtr
377	IF (aerosol(it)) THEN
378	DO k = 1, klev
379	DO i = 1, klon
380	d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) + &
381	( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
382	d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) + &
383	( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
384
385	!
386	! Flux lessivage total
387	! ------------------------------------------------------------
388	flestottr(i,k,it) = flestottr(i,k,it) - &
389	( d_tr_lessi_nucl(i,k,it) + &
390	d_tr_lessi_impa(i,k,it) ) * &
391	( paprs(i,k)-paprs(i,k+1) ) / &
392	(RG * pdtphys)
393	!
394	! Mise a jour des traceurs due a l'impaction,nucleation
395	! ----------------------------------------------------------------------
396	tr_seri(i,k,it)=tr_seri(i,k,it)frac_impa(i,k)frac_nucl(i,k)
397	END DO
398	END DO
399	END IF
400	END DO
401
402	END IF ! lessivage
403
404	!=============================================================
405	! Ecriture des sorties
406	!=============================================================
407	INCLUDE "write_histrac.h"
408
409
410	END SUBROUTINE phytrac

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

Download in other formats:

Original Format