Context Navigation

fisrtilp_tr.F @ 3133

Last change on this file since 3133 was 517, checked in by lmdzadmin, 21 years ago
Inclusion des modifications de O. Boucher et de J. Quaas pour le calcul des premiers effets directs et indirects dus aux aerosols LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 12.7 KB

Line
1	c $Header$
2	c
3	SUBROUTINE fisrtilp_tr(dtime,paprs,pplay,t,q,ratqs,
4	s d_t, d_q, d_ql, rneb, radliq, rain, snow,
5	s pfrac_impa, pfrac_nucl, pfrac_1nucl,
6	s frac_impa, frac_nucl,
7	s prfl, psfl,
8	s RHcl) ! relative humidity in clear sky (needed for aer optical properties; aeropt.F)
9
10	c
11	IMPLICIT none
12	c======================================================================
13	c Auteur(s): Z.X. Li (LMD/CNRS)
14	c Date: le 20 mars 1995
15	c Objet: condensation et precipitation stratiforme.
16	c schema de nuage
17	c======================================================================
18	c======================================================================
19	#include "dimensions.h"
20	#include "dimphy.h"
21	#include "YOMCST.h"
22	#include "tracstoke.h"
23	c
24	c Arguments:
25	c
26	REAL dtime ! intervalle du temps (s)
27	REAL paprs(klon,klev+1) ! pression a inter-couche
28	REAL pplay(klon,klev) ! pression au milieu de couche
29	REAL t(klon,klev) ! temperature (K)
30	REAL q(klon,klev) ! humidite specifique (kg/kg)
31	REAL d_t(klon,klev) ! incrementation de la temperature (K)
32	REAL d_q(klon,klev) ! incrementation de la vapeur d'eau
33	REAL d_ql(klon,klev) ! incrementation de l'eau liquide
34	REAL rneb(klon,klev) ! fraction nuageuse
35	REAL radliq(klon,klev) ! eau liquide utilisee dans rayonnements
36	REAL rain(klon) ! pluies (mm/s)
37	REAL snow(klon) ! neige (mm/s)
38	REAL prfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
39	REAL psfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
40
41	Cjq For aerosol opt properties needed (see aeropt.F)
42	REAL RHcl(klon,klev)
43
44	cAA
45	c Coeffients de fraction lessivee : pour OFF-LINE
46	c
47	REAL pfrac_nucl(klon,klev)
48	REAL pfrac_1nucl(klon,klev)
49	REAL pfrac_impa(klon,klev)
50	c
51	c Fraction d'aerosols lessivee par impaction et par nucleation
52	c POur ON-LINE
53	c
54	REAL frac_impa(klon,klev)
55	REAL frac_nucl(klon,klev)
56	cAA
57	c
58	c Options du programme:
59	c
60	REAL seuil_neb ! un nuage existe vraiment au-dela
61	PARAMETER (seuil_neb=0.001)
62	REAL ct ! inverse du temps pour qu'un nuage precipite
63	PARAMETER (ct=1./1800.)
64	REAL cl ! seuil de precipitation
65	PARAMETER (cl=2.6e-4)
66	ccc PARAMETER (cl=2.3e-4)
67	ccc PARAMETER (cl=2.0e-4)
68	INTEGER ninter ! sous-intervals pour la precipitation
69	PARAMETER (ninter=5)
70	LOGICAL evap_prec ! evaporation de la pluie
71	PARAMETER (evap_prec=.TRUE.)
72	REAL coef_eva
73	PARAMETER (coef_eva=2.0E-05)
74	LOGICAL calcrat ! calculer ratqs au lieu de fixer sa valeur
75	REAL ratqs(klon,klev) ! determine la largeur de distribution de vapeur
76	PARAMETER (calcrat=.TRUE.)
77	REAL zx_min, rat_max
78	PARAMETER (zx_min=1.0, rat_max=0.01)
79	REAL zx_max, rat_min
80	PARAMETER (zx_max=0.1, rat_min=0.3)
81	REAL zx
82	c
83	LOGICAL cpartiel ! condensation partielle
84	PARAMETER (cpartiel=.TRUE.)
85	REAL t_coup
86	PARAMETER (t_coup=234.0)
87	c
88	c Variables locales:
89	c
90	INTEGER i, k, n, kk
91	REAL zqs(klon), zdqs(klon), zdelta, zcor, zcvm5
92	REAL zrfl(klon), zrfln(klon), zqev, zqevt
93	REAL zoliq(klon), zcond(klon), zq(klon), zqn(klon), zdelq
94	REAL ztglace, zt(klon)
95	INTEGER nexpo ! exponentiel pour glace/eau
96	REAL zdz(klon),zrho(klon),ztot(klon), zrhol(klon)
97	REAL zchau(klon),zfroi(klon),zfice(klon),zneb(klon)
98	c
99	LOGICAL appel1er
100	SAVE appel1er
101	c
102	c---------------------------------------------------------------
103	c
104	cAA Variables traceurs:
105	cAA Provisoire !!! Parametres alpha du lessivage
106	cAA A priori on a 4 scavenging # possibles
107	c
108	REAL a_tr_sca(4)
109	save a_tr_sca
110	c
111	c Variables intermediaires
112	c
113	REAL zalpha_tr
114	REAL zfrac_lessi
115	REAL zprec_cond(klon)
116	cAA
117	c---------------------------------------------------------------
118	c
119	c Fonctions en ligne:
120	c
121	REAL fallv ! vitesse de chute pour crystaux de glace
122	REAL zzz
123	#include "YOETHF.h"
124	#include "FCTTRE.h"
125	fallv (zzz) = 3.29/2.0 * ((zzz)**0.16)
126	ccc fallv (zzz) = 3.29/3.0 * ((zzz)**0.16)
127	ccc fallv (zzz) = 3.29 * ((zzz)**0.16)
128	c
129	DATA appel1er /.TRUE./
130	c
131	IF (appel1er) THEN
132	c
133	PRINT*, 'fisrtilp, calcrat:', calcrat
134	PRINT*, 'fisrtilp, ninter:', ninter
135	PRINT*, 'fisrtilp, evap_prec:', evap_prec
136	PRINT*, 'fisrtilp, cpartiel:', cpartiel
137	IF (ABS(dtime/FLOAT(ninter)-360.0).GT.0.001) THEN
138	PRINT*, 'fisrtilp: Ce n est pas prevu, voir Z.X.Li', dtime
139	PRINT*, 'Je prefere un sous-intervalle de 6 minutes'
140	CALL abort
141	ENDIF
142	appel1er = .FALSE.
143	c
144	cAA initialiation provisoire
145	a_tr_sca(1) = -0.5
146	a_tr_sca(2) = -0.5
147	a_tr_sca(3) = -0.5
148	a_tr_sca(4) = -0.5
149	c
150	cAA Initialisation a 1 des coefs des fractions lessivees
151	c
152	DO k = 1, klev
153	DO i = 1, klon
154	pfrac_nucl(i,k)=1.
155	pfrac_1nucl(i,k)=1.
156	pfrac_impa(i,k)=1.
157	ENDDO
158	ENDDO
159
160	ENDIF ! test sur appel1er
161	c
162	cMAf Initialisation a 0 de zoliq
163	DO i = 1, klon
164	zoliq(i)=0.
165	ENDDO
166	c Determiner les nuages froids par leur temperature
167	c
168	ztglace = RTT - 15.0
169	nexpo = 6
170	ccc nexpo = 1
171	c
172	c Initialiser les sorties:
173	c
174	DO k = 1, klev+1
175	DO i = 1, klon
176	prfl(i,k) = 0.0
177	psfl(i,k) = 0.0
178	ENDDO
179	ENDDO
180
181	DO k = 1, klev
182	DO i = 1, klon
183	d_t(i,k) = 0.0
184	d_q(i,k) = 0.0
185	d_ql(i,k) = 0.0
186	rneb(i,k) = 0.0
187	radliq(i,k) = 0.0
188	frac_nucl(i,k) = 1.
189	frac_impa(i,k) = 1.
190	ENDDO
191	ENDDO
192	DO i = 1, klon
193	rain(i) = 0.0
194	snow(i) = 0.0
195	ENDDO
196	c
197	c Initialiser le flux de precipitation a zero
198	c
199	DO i = 1, klon
200	zrfl(i) = 0.0
201	zneb(i) = seuil_neb
202	ENDDO
203	c
204	c
205	cAA Pour plus de securite
206
207	zalpha_tr = 0.
208	zfrac_lessi = 0.
209
210	cAA----------------------------------------------------------
211	c
212	c Boucle verticale (du haut vers le bas)
213	c
214	DO 9999 k = klev, 1, -1
215	c
216	cAA----------------------------------------------------------
217	c
218	DO i = 1, klon
219	zt(i)=t(i,k)
220	zq(i)=q(i,k)
221	ENDDO
222	c
223	c Calculer l'evaporation de la precipitation
224	c
225	IF (evap_prec) THEN
226	DO i = 1, klon
227	IF (zrfl(i) .GT.0.) THEN
228	IF (thermcep) THEN
229	zdelta=MAX(0.,SIGN(1.,RTT-zt(i)))
230	zqs(i)= R2ES*FOEEW(zt(i),zdelta)/pplay(i,k)
231	zqs(i)=MIN(0.5,zqs(i))
232	zcor=1./(1.-RETV*zqs(i))
233	zqs(i)=zqs(i)*zcor
234	ELSE
235	IF (zt(i) .LT. t_coup) THEN
236	zqs(i) = qsats(zt(i)) / pplay(i,k)
237	ELSE
238	zqs(i) = qsatl(zt(i)) / pplay(i,k)
239	ENDIF
240	ENDIF
241	zqev = MAX (0.0, (zqs(i)-zq(i))*zneb(i) )
242	zqevt = coef_eva * (1.0-zq(i)/zqs(i)) * SQRT(zrfl(i))
243	. * (paprs(i,k)-paprs(i,k+1))/pplay(i,k)zt(i)RD/RG
244	zqevt = MAX(0.0,MIN(zqevt,zrfl(i)))
245	. * RG*dtime/(paprs(i,k)-paprs(i,k+1))
246	zqev = MIN (zqev, zqevt)
247	zrfln(i) = zrfl(i) - zqev*(paprs(i,k)-paprs(i,k+1))
248	. /RG/dtime
249	zq(i) = zq(i) - (zrfln(i)-zrfl(i))
250	. * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime
251	zt(i) = zt(i) + (zrfln(i)-zrfl(i))
252	. * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime
253	. * RLVTT/RCPD/(1.0+RVTMP2*zq(i))
254	zrfl(i) = zrfln(i)
255	ENDIF
256	ENDDO
257	ENDIF
258	c
259	c Calculer Qs et L/Cp*dQs/dT:
260	c
261	IF (thermcep) THEN
262	DO i = 1, klon
263	zdelta = MAX(0.,SIGN(1.,RTT-zt(i)))
264	zcvm5 = R5LESRLVTT(1.-zdelta) + R5IESRLSTTzdelta
265	zcvm5 = zcvm5 /RCPD/(1.0+RVTMP2*zq(i))
266	zqs(i) = R2ES*FOEEW(zt(i),zdelta)/pplay(i,k)
267	zqs(i) = MIN(0.5,zqs(i))
268	zcor = 1./(1.-RETV*zqs(i))
269	zqs(i) = zqs(i)*zcor
270	zdqs(i) = FOEDE(zt(i),zdelta,zcvm5,zqs(i),zcor)
271	ENDDO
272	ELSE
273	DO i = 1, klon
274	IF (zt(i).LT.t_coup) THEN
275	zqs(i) = qsats(zt(i))/pplay(i,k)
276	zdqs(i) = dqsats(zt(i),zqs(i))
277	ELSE
278	zqs(i) = qsatl(zt(i))/pplay(i,k)
279	zdqs(i) = dqsatl(zt(i),zqs(i))
280	ENDIF
281	ENDDO
282	ENDIF
283	c
284	c Determiner la condensation partielle et calculer la quantite
285	c de l'eau condensee:
286	c
287	IF (cpartiel) THEN
288	DO i = 1, klon
289	c
290	zdelq = ratqs(i,k) * zq(i)
291	rneb(i,k) = (zq(i)+zdelq-zqs(i)) / (2.0*zdelq)
292	zqn(i) = (zq(i)+zdelq+zqs(i))/2.0
293	IF (rneb(i,k) .LE. 0.0) zqn(i) = 0.0
294	IF (rneb(i,k) .GE. 1.0) zqn(i) = zq(i)
295	rneb(i,k) = MAX(0.0,MIN(1.0,rneb(i,k)))
296	zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)/(1.+zdqs(i))
297
298	c--Olivier
299	RHcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i)
300	IF (rneb(i,k) .LE. 0.0) RHcl(i,k)=zq(i)/zqs(i)
301	IF (rneb(i,k) .GE. 1.0) RHcl(i,k)=1.0
302	c--fin
303
304	ENDDO
305	ELSE
306	DO i = 1, klon
307	IF (zq(i).GT.zqs(i)) THEN
308	rneb(i,k) = 1.0
309	ELSE
310	rneb(i,k) = 0.0
311	ENDIF
312	zcond(i) = MAX(0.0,zq(i)-zqs(i))/(1.+zdqs(i))
313	ENDDO
314	ENDIF
315	c
316	DO i = 1, klon
317	zq(i) = zq(i) - zcond(i)
318	zt(i) = zt(i) + zcond(i) * RLVTT/RCPD
319	ENDDO
320	c
321	c Partager l'eau condensee en precipitation et eau liquide nuageuse
322	c
323	DO i = 1, klon
324	IF (rneb(i,k).GT.0.0) THEN
325	zoliq(i) = zcond(i)
326	zrho(i) = pplay(i,k) / zt(i) / RD
327	zdz(i) = (paprs(i,k)-paprs(i,k+1)) / (zrho(i)*RG)
328	zfice(i) = 1.0 - (zt(i)-ztglace) / (273.13-ztglace)
329	zfice(i) = MIN(MAX(zfice(i),0.0),1.0)
330	zfice(i) = zfice(i)**nexpo
331	zneb(i) = MAX(rneb(i,k), seuil_neb)
332	radliq(i,k) = zoliq(i)/FLOAT(ninter+1)
333	ENDIF
334	ENDDO
335	c
336	DO n = 1, ninter
337	DO i = 1, klon
338	IF (rneb(i,k).GT.0.0) THEN
339	zchau(i) = ctdtime/FLOAT(ninter) zoliq(i)
340	. * (1.0-EXP(-(zoliq(i)/zneb(i)/cl)*2)) (1.-zfice(i))
341	zrhol(i) = zrho(i) * zoliq(i) / zneb(i)
342	zfroi(i) = dtime/FLOAT(ninter)/zdz(i)*zoliq(i)
343	. fallv(zrhol(i)) zfice(i)
344	ztot(i) = zchau(i) + zfroi(i)
345	IF (zneb(i).EQ.seuil_neb) ztot(i) = 0.0
346	ztot(i) = MIN(MAX(ztot(i),0.0),zoliq(i))
347	zoliq(i) = MAX(zoliq(i)-ztot(i), 0.0)
348	radliq(i,k) = radliq(i,k) + zoliq(i)/FLOAT(ninter+1)
349	ENDIF
350	ENDDO
351	ENDDO
352	c
353	DO i = 1, klon
354	IF (rneb(i,k).GT.0.0) THEN
355	d_ql(i,k) = zoliq(i)
356	zrfl(i) = zrfl(i)+ MAX(zcond(i)-zoliq(i),0.0)
357	. * (paprs(i,k)-paprs(i,k+1))/(RG*dtime)
358	ENDIF
359	IF (zt(i).LT.RTT) THEN
360	psfl(i,k)=zrfl(i)
361	ELSE
362	prfl(i,k)=zrfl(i)
363	ENDIF
364	ENDDO
365	c
366	c Calculer les tendances de q et de t:
367	c
368	DO i = 1, klon
369	d_q(i,k) = zq(i) - q(i,k)
370	d_t(i,k) = zt(i) - t(i,k)
371	ENDDO
372	c
373	cAA--------------- Calcul du lessivage stratiforme -------------
374
375	DO i = 1,klon
376	c
377	zprec_cond(i) = MAX(zcond(i)-zoliq(i),0.0)
378	. * (paprs(i,k)-paprs(i,k+1))/RG
379	IF (rneb(i,k).GT.0.0.and.zprec_cond(i).gt.0.) THEN
380	cAA lessivage nucleation LMD5 dans la couche elle-meme
381	if (t(i,k) .GE. ztglace) THEN
382	zalpha_tr = a_tr_sca(3)
383	else
384	zalpha_tr = a_tr_sca(4)
385	endif
386	zfrac_lessi = 1. - EXP(zalpha_tr*zprec_cond(i)/zneb(i))
387	pfrac_nucl(i,k)=pfrac_nucl(i,k)(1.-zneb(i)zfrac_lessi)
388	frac_nucl(i,k)= 1.-zneb(i)*zfrac_lessi
389	c
390	c nucleation avec un facteur -1 au lieu de -0.5
391	zfrac_lessi = 1. - EXP(-zprec_cond(i)/zneb(i))
392	pfrac_1nucl(i,k)=pfrac_1nucl(i,k)(1.-zneb(i)zfrac_lessi)
393	ENDIF
394	c
395	ENDDO ! boucle sur i
396	c
397	cAA Lessivage par impaction dans les couches en-dessous
398	DO kk = k-1, 1, -1
399	DO i = 1, klon
400	IF (rneb(i,k).GT.0.0.and.zprec_cond(i).gt.0.) THEN
401	if (t(i,kk) .GE. ztglace) THEN
402	zalpha_tr = a_tr_sca(1)
403	else
404	zalpha_tr = a_tr_sca(2)
405	endif
406	zfrac_lessi = 1. - EXP(zalpha_tr*zprec_cond(i)/zneb(i))
407	pfrac_impa(i,kk)=pfrac_impa(i,kk)(1.-zneb(i)zfrac_lessi)
408	frac_impa(i,kk)= 1.-zneb(i)*zfrac_lessi
409	ENDIF
410	ENDDO
411	ENDDO
412	c
413	cAA----------------------------------------------------------
414	c FIN DE BOUCLE SUR K
415	9999 CONTINUE
416	c
417	cAA-----------------------------------------------------------
418	c
419	c Pluie ou neige au sol selon la temperature de la 1ere couche
420	c
421	DO i = 1, klon
422	IF ((t(i,1)+d_t(i,1)) .LT. RTT) THEN
423	snow(i) = zrfl(i)
424	ELSE
425	rain(i) = zrfl(i)
426	ENDIF
427	ENDDO
428	c
429	RETURN
430	END

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

Download in other formats:

Original Format