Context Navigation

fisrtilp.F @ 1338

Last change on this file since 1338 was 1299, checked in by Laurent Fairhead, 15 years ago

Nettoyage general pour se rapprocher des normes et éviter des erreurs a la
compilation:

tous les FLOAT() sont remplacés par des REAL()
tous les STOP dans phylmd sont remplacés par des appels à abort_gcm
le common control défini dans le fichier control.h est remplacé par le module control_mod pour éviter des messages sur l'alignement des variables dans les déclarations
des $Header$ remplacés par des $Id$ pour svn

Quelques remplacements à faire ont pu m'échapper

General cleanup of the code to try and adhere to norms and to prevent some
compilation errors:

all FLOAT() instructions have been replaced by REAL() instructions
all STOP instructions in phylmd have been replaced by calls to abort_gcm
the common block control defined in the control.h file has been replaced by the control_mod to prevent compilation warnings on the alignement of declared variables
$Header$ replaced by $Id$ for svn

Some changes which should have been made might have escaped me

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

File size: 17.0 KB

Line
1	!
2	! $Id: fisrtilp.F 1299 2010-01-20 14:27:21Z idelkadi $
3	!
4	c
5	SUBROUTINE fisrtilp(dtime,paprs,pplay,t,q,ptconv,ratqs,
6	s d_t, d_q, d_ql, rneb, radliq, rain, snow,
7	s pfrac_impa, pfrac_nucl, pfrac_1nucl,
8	s frac_impa, frac_nucl,
9	s prfl, psfl, rhcl)
10
11	c
12	USE dimphy
13	IMPLICIT none
14	c======================================================================
15	c Auteur(s): Z.X. Li (LMD/CNRS)
16	c Date: le 20 mars 1995
17	c Objet: condensation et precipitation stratiforme.
18	c schema de nuage
19	c======================================================================
20	c======================================================================
21	cym#include "dimensions.h"
22	cym#include "dimphy.h"
23	#include "YOMCST.h"
24	#include "tracstoke.h"
25	#include "fisrtilp.h"
26	c
27	c Arguments:
28	c
29	REAL dtime ! intervalle du temps (s)
30	REAL paprs(klon,klev+1) ! pression a inter-couche
31	REAL pplay(klon,klev) ! pression au milieu de couche
32	REAL t(klon,klev) ! temperature (K)
33	REAL q(klon,klev) ! humidite specifique (kg/kg)
34	REAL d_t(klon,klev) ! incrementation de la temperature (K)
35	REAL d_q(klon,klev) ! incrementation de la vapeur d'eau
36	REAL d_ql(klon,klev) ! incrementation de l'eau liquide
37	REAL rneb(klon,klev) ! fraction nuageuse
38	REAL radliq(klon,klev) ! eau liquide utilisee dans rayonnements
39	REAL rhcl(klon,klev) ! humidite relative en ciel clair
40	REAL rain(klon) ! pluies (mm/s)
41	REAL snow(klon) ! neige (mm/s)
42	REAL prfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
43	REAL psfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
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	real zct ,zcl
57	cAA
58	c
59	c Options du programme:
60	c
61	REAL seuil_neb ! un nuage existe vraiment au-dela
62	PARAMETER (seuil_neb=0.001)
63
64	INTEGER ninter ! sous-intervals pour la precipitation
65	INTEGER ncoreczq
66	PARAMETER (ninter=5)
67	LOGICAL evap_prec ! evaporation de la pluie
68	PARAMETER (evap_prec=.TRUE.)
69	REAL ratqs(klon,klev) ! determine la largeur de distribution de vapeur
70	logical ptconv(klon,klev) ! determine la largeur de distribution de vapeur
71
72	real zpdf_sig(klon),zpdf_k(klon),zpdf_delta(klon)
73	real Zpdf_a(klon),zpdf_b(klon),zpdf_e1(klon),zpdf_e2(klon)
74	real erf
75	c
76	LOGICAL cpartiel ! condensation partielle
77	PARAMETER (cpartiel=.TRUE.)
78	REAL t_coup
79	PARAMETER (t_coup=234.0)
80	c
81	c Variables locales:
82	c
83	INTEGER i, k, n, kk
84	REAL zqs(klon), zdqs(klon), zdelta, zcor, zcvm5
85	REAL zrfl(klon), zrfln(klon), zqev, zqevt
86	REAL zoliq(klon), zcond(klon), zq(klon), zqn(klon), zdelq
87	REAL ztglace, zt(klon)
88	INTEGER nexpo ! exponentiel pour glace/eau
89	REAL zdz(klon),zrho(klon),ztot , zrhol(klon)
90	REAL zchau ,zfroi ,zfice(klon),zneb(klon)
91	c
92	LOGICAL appel1er
93	SAVE appel1er
94	c$OMP THREADPRIVATE(appel1er)
95	c
96	c---------------------------------------------------------------
97	c
98	cAA Variables traceurs:
99	cAA Provisoire !!! Parametres alpha du lessivage
100	cAA A priori on a 4 scavenging # possibles
101	c
102	REAL a_tr_sca(4)
103	save a_tr_sca
104	c$OMP THREADPRIVATE(a_tr_sca)
105	c
106	c Variables intermediaires
107	c
108	REAL zalpha_tr
109	REAL zfrac_lessi
110	REAL zprec_cond(klon)
111	cAA
112	REAL zmair, zcpair, zcpeau
113	C Pour la conversion eau-neige
114	REAL zlh_solid(klon), zm_solid
115	cIM
116	cym INTEGER klevm1
117	c---------------------------------------------------------------
118	c
119	c Fonctions en ligne:
120	c
121	REAL fallvs,fallvc ! vitesse de chute pour crystaux de glace
122	REAL zzz
123	#include "YOETHF.h"
124	#include "FCTTRE.h"
125	fallvc (zzz) = 3.29/2.0 * ((zzz)*0.16) ffallv_con
126	fallvs (zzz) = 3.29/2.0 * ((zzz)*0.16) ffallv_lsc
127	c
128	DATA appel1er /.TRUE./
129	cym
130	zdelq=0.0
131
132	IF (appel1er) THEN
133	c
134	PRINT*, 'fisrtilp, ninter:', ninter
135	PRINT*, 'fisrtilp, evap_prec:', evap_prec
136	PRINT*, 'fisrtilp, cpartiel:', cpartiel
137	IF (ABS(dtime/REAL(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	c 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	!cdir collapse
153	DO k = 1, klev
154	DO i = 1, klon
155	pfrac_nucl(i,k)=1.
156	pfrac_1nucl(i,k)=1.
157	pfrac_impa(i,k)=1.
158	ENDDO
159	ENDDO
160
161	ENDIF ! test sur appel1er
162	c
163	cMAf Initialisation a 0 de zoliq
164	c DO i = 1, klon
165	c zoliq(i)=0.
166	c ENDDO
167	c Determiner les nuages froids par leur temperature
168	c nexpo regle la raideur de la transition eau liquide / eau glace.
169	c
170	ztglace = RTT - 15.0
171	nexpo = 6
172	ccc nexpo = 1
173	c
174	c Initialiser les sorties:
175	c
176	!cdir collapse
177	DO k = 1, klev+1
178	DO i = 1, klon
179	prfl(i,k) = 0.0
180	psfl(i,k) = 0.0
181	ENDDO
182	ENDDO
183
184	!cdir collapse
185	DO k = 1, klev
186	DO i = 1, klon
187	d_t(i,k) = 0.0
188	d_q(i,k) = 0.0
189	d_ql(i,k) = 0.0
190	rneb(i,k) = 0.0
191	radliq(i,k) = 0.0
192	frac_nucl(i,k) = 1.
193	frac_impa(i,k) = 1.
194	ENDDO
195	ENDDO
196	DO i = 1, klon
197	rain(i) = 0.0
198	snow(i) = 0.0
199	zoliq(i)=0.
200	c ENDDO
201	c
202	c Initialiser le flux de precipitation a zero
203	c
204	c DO i = 1, klon
205	zrfl(i) = 0.0
206	zneb(i) = seuil_neb
207	ENDDO
208	c
209	c
210	cAA Pour plus de securite
211
212	zalpha_tr = 0.
213	zfrac_lessi = 0.
214
215	cAA----------------------------------------------------------
216	c
217	ncoreczq=0
218	c Boucle verticale (du haut vers le bas)
219	c
220	cIM : klevm1
221	cym klevm1=klev-1
222	DO 9999 k = klev, 1, -1
223	c
224	cAA----------------------------------------------------------
225	c
226	DO i = 1, klon
227	zt(i)=t(i,k)
228	zq(i)=q(i,k)
229	ENDDO
230	c
231	c Calculer la varition de temp. de l'air du a la chaleur sensible
232	C transporter par la pluie.
233	C Il resterait a rajouter cet effet de la chaleur sensible sur les
234	C flux de surface, du a la diff. de temp. entre le 1er niveau et la
235	C surface.
236	C
237	IF(k.LE.klevm1) THEN
238	DO i = 1, klon
239	cIM
240	zmair=(paprs(i,k)-paprs(i,k+1))/RG
241	zcpair=RCPD(1.0+RVTMP2zq(i))
242	zcpeau=RCPD*RVTMP2
243	zt(i) = ( (t(i,k+1)+d_t(i,k+1))zrfl(i)dtime*zcpeau
244	$ + zmairzcpairzt(i) )
245	$ / (zmairzcpair + zrfl(i)dtime*zcpeau)
246	C C WRITE (6,*) 'cppluie ', zt(i)-(t(i,k+1)+d_t(i,k+1))
247	ENDDO
248	ENDIF
249	c
250	c
251	c Calculer l'evaporation de la precipitation
252	c
253
254
255	IF (evap_prec) THEN
256	DO i = 1, klon
257	IF (zrfl(i) .GT.0.) THEN
258	IF (thermcep) THEN
259	zdelta=MAX(0.,SIGN(1.,RTT-zt(i)))
260	zqs(i)= R2ES*FOEEW(zt(i),zdelta)/pplay(i,k)
261	zqs(i)=MIN(0.5,zqs(i))
262	zcor=1./(1.-RETV*zqs(i))
263	zqs(i)=zqs(i)*zcor
264	ELSE
265	IF (zt(i) .LT. t_coup) THEN
266	zqs(i) = qsats(zt(i)) / pplay(i,k)
267	ELSE
268	zqs(i) = qsatl(zt(i)) / pplay(i,k)
269	ENDIF
270	ENDIF
271	zqev = MAX (0.0, (zqs(i)-zq(i))*zneb(i) )
272	zqevt = coef_eva * (1.0-zq(i)/zqs(i)) * SQRT(zrfl(i))
273	. * (paprs(i,k)-paprs(i,k+1))/pplay(i,k)zt(i)RD/RG
274	zqevt = MAX(0.0,MIN(zqevt,zrfl(i)))
275	. * RG*dtime/(paprs(i,k)-paprs(i,k+1))
276	zqev = MIN (zqev, zqevt)
277	zrfln(i) = zrfl(i) - zqev*(paprs(i,k)-paprs(i,k+1))
278	. /RG/dtime
279
280	c pour la glace, on r�vapore toute la pr�ip dans la couche du dessous
281	c la glace venant de la couche du dessus est simplement dans la couche
282	c du dessous.
283
284	IF (zt(i) .LT. t_coup.and.reevap_ice) zrfln(i)=0.
285
286	zq(i) = zq(i) - (zrfln(i)-zrfl(i))
287	. * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime
288	zt(i) = zt(i) + (zrfln(i)-zrfl(i))
289	. * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime
290	. * RLVTT/RCPD/(1.0+RVTMP2*zq(i))
291	zrfl(i) = zrfln(i)
292	ENDIF
293	ENDDO
294	ENDIF
295	c
296	c Calculer Qs et L/Cp*dQs/dT:
297	c
298	IF (thermcep) THEN
299	DO i = 1, klon
300	zdelta = MAX(0.,SIGN(1.,RTT-zt(i)))
301	zcvm5 = R5LESRLVTT(1.-zdelta) + R5IESRLSTTzdelta
302	zcvm5 = zcvm5 /RCPD/(1.0+RVTMP2*zq(i))
303	zqs(i) = R2ES*FOEEW(zt(i),zdelta)/pplay(i,k)
304	zqs(i) = MIN(0.5,zqs(i))
305	zcor = 1./(1.-RETV*zqs(i))
306	zqs(i) = zqs(i)*zcor
307	zdqs(i) = FOEDE(zt(i),zdelta,zcvm5,zqs(i),zcor)
308	ENDDO
309	ELSE
310	DO i = 1, klon
311	IF (zt(i).LT.t_coup) THEN
312	zqs(i) = qsats(zt(i))/pplay(i,k)
313	zdqs(i) = dqsats(zt(i),zqs(i))
314	ELSE
315	zqs(i) = qsatl(zt(i))/pplay(i,k)
316	zdqs(i) = dqsatl(zt(i),zqs(i))
317	ENDIF
318	ENDDO
319	ENDIF
320	c
321	c Determiner la condensation partielle et calculer la quantite
322	c de l'eau condensee:
323	c
324	IF (cpartiel) THEN
325
326	c print*,'Dans partiel k=',k
327	c
328	c Calcul de l'eau condensee et de la fraction nuageuse et de l'eau
329	c nuageuse a partir des PDF de Sandrine Bony.
330	c rneb : fraction nuageuse
331	c zqn : eau totale dans le nuage
332	c zcond : eau condensee moyenne dans la maille.
333	c on prend en compte le r�hauffement qui diminue la partie condensee
334	c
335	c Version avec les raqts
336
337	if (iflag_pdf.eq.0) then
338
339	do i=1,klon
340	zdelq = min(ratqs(i,k),0.99) * zq(i)
341	rneb(i,k) = (zq(i)+zdelq-zqs(i)) / (2.0*zdelq)
342	zqn(i) = (zq(i)+zdelq+zqs(i))/2.0
343	enddo
344
345	else
346	c
347	c Version avec les nouvelles PDFs.
348	do i=1,klon
349	if(zq(i).lt.1.e-15) then
350	ncoreczq=ncoreczq+1
351	zq(i)=1.e-15
352	endif
353	enddo
354	do i=1,klon
355	zpdf_sig(i)=ratqs(i,k)*zq(i)
356	zpdf_k(i)=-sqrt(log(1.+(zpdf_sig(i)/zq(i))**2))
357	zpdf_delta(i)=log(zq(i)/zqs(i))
358	zpdf_a(i)=zpdf_delta(i)/(zpdf_k(i)*sqrt(2.))
359	zpdf_b(i)=zpdf_k(i)/(2.*sqrt(2.))
360	zpdf_e1(i)=zpdf_a(i)-zpdf_b(i)
361	zpdf_e1(i)=sign(min(abs(zpdf_e1(i)),5.),zpdf_e1(i))
362	zpdf_e1(i)=1.-erf(zpdf_e1(i))
363	zpdf_e2(i)=zpdf_a(i)+zpdf_b(i)
364	zpdf_e2(i)=sign(min(abs(zpdf_e2(i)),5.),zpdf_e2(i))
365	zpdf_e2(i)=1.-erf(zpdf_e2(i))
366	if (zpdf_e1(i).lt.1.e-10) then
367	rneb(i,k)=0.
368	zqn(i)=zqs(i)
369	else
370	rneb(i,k)=0.5*zpdf_e1(i)
371	zqn(i)=zq(i)*zpdf_e2(i)/zpdf_e1(i)
372	endif
373
374	enddo
375
376	endif ! iflag_pdf
377
378	DO i=1,klon
379	IF (rneb(i,k) .LE. 0.0) THEN
380	zqn(i) = 0.0
381	rneb(i,k) = 0.0
382	zcond(i) = 0.0
383	rhcl(i,k)=zq(i)/zqs(i)
384	ELSE IF (rneb(i,k) .GE. 1.0) THEN
385	zqn(i) = zq(i)
386	rneb(i,k) = 1.0
387	zcond(i) = MAX(0.0,zqn(i)-zqs(i))
388	rhcl(i,k)=1.0
389	ELSE
390	zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)
391	rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i)
392	ENDIF
393	ENDDO
394	! do i=1,klon
395	! IF (rneb(i,k) .LE. 0.0) zqn(i) = 0.0
396	! IF (rneb(i,k) .GE. 1.0) zqn(i) = zq(i)
397	! rneb(i,k) = MAX(0.0,MIN(1.0,rneb(i,k)))
398	!c zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)/(1.+zdqs(i))
399	!c On ne divise pas par 1+zdqs pour forcer a avoir l'eau predite par
400	!c la convection.
401	!c ATTENTION !!! Il va falloir verifier tout ca.
402	! zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)
403	!c print*,'ZDQS ',zdqs(i)
404	!c--Olivier
405	! rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i)
406	! IF (rneb(i,k) .LE. 0.0) rhcl(i,k)=zq(i)/zqs(i)
407	! IF (rneb(i,k) .GE. 1.0) rhcl(i,k)=1.0
408	!c--fin
409	! ENDDO
410	ELSE
411	DO i = 1, klon
412	IF (zq(i).GT.zqs(i)) THEN
413	rneb(i,k) = 1.0
414	ELSE
415	rneb(i,k) = 0.0
416	ENDIF
417	zcond(i) = MAX(0.0,zq(i)-zqs(i))/(1.+zdqs(i))
418	ENDDO
419	ENDIF
420	c
421	DO i = 1, klon
422	zq(i) = zq(i) - zcond(i)
423	c zt(i) = zt(i) + zcond(i) * RLVTT/RCPD
424	zt(i) = zt(i) + zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*zq(i))
425	ENDDO
426	c
427	c Partager l'eau condensee en precipitation et eau liquide nuageuse
428	c
429	DO i = 1, klon
430	IF (rneb(i,k).GT.0.0) THEN
431	zoliq(i) = zcond(i)
432	zrho(i) = pplay(i,k) / zt(i) / RD
433	zdz(i) = (paprs(i,k)-paprs(i,k+1)) / (zrho(i)*RG)
434	zfice(i) = 1.0 - (zt(i)-ztglace) / (273.13-ztglace)
435	zfice(i) = MIN(MAX(zfice(i),0.0),1.0)
436	zfice(i) = zfice(i)**nexpo
437	zneb(i) = MAX(rneb(i,k), seuil_neb)
438	radliq(i,k) = zoliq(i)/REAL(ninter+1)
439	ENDIF
440	ENDDO
441	c
442	DO n = 1, ninter
443	DO i = 1, klon
444	IF (rneb(i,k).GT.0.0) THEN
445	zrhol(i) = zrho(i) * zoliq(i) / zneb(i)
446
447	IF (zneb(i).EQ.seuil_neb) THEN
448	ztot = 0.0
449	ELSE
450	c quantite d'eau a eliminer: zchau
451	c meme chose pour la glace: zfroi
452	if (ptconv(i,k)) then
453	zcl =cld_lc_con
454	zct =1./cld_tau_con
455	zfroi = dtime/REAL(ninter)/zdz(i)*zoliq(i)
456	. fallvc(zrhol(i)) zfice(i)
457	else
458	zcl =cld_lc_lsc
459	zct =1./cld_tau_lsc
460	zfroi = dtime/REAL(ninter)/zdz(i)*zoliq(i)
461	. fallvs(zrhol(i)) zfice(i)
462	endif
463	zchau = zct dtime/REAL(ninter) zoliq(i)
464	. (1.0-EXP(-(zoliq(i)/zneb(i)/zcl )2)) (1.-zfice(i))
465	ztot = zchau + zfroi
466	ztot = MAX(ztot ,0.0)
467	ENDIF
468	ztot = MIN(ztot,zoliq(i))
469	zoliq(i) = MAX(zoliq(i)-ztot , 0.0)
470	radliq(i,k) = radliq(i,k) + zoliq(i)/REAL(ninter+1)
471	ENDIF
472	ENDDO
473	ENDDO
474	c
475	DO i = 1, klon
476	IF (rneb(i,k).GT.0.0) THEN
477	d_ql(i,k) = zoliq(i)
478	zrfl(i) = zrfl(i)+ MAX(zcond(i)-zoliq(i),0.0)
479	. * (paprs(i,k)-paprs(i,k+1))/(RG*dtime)
480	ENDIF
481	IF (zt(i).LT.RTT) THEN
482	psfl(i,k)=zrfl(i)
483	ELSE
484	prfl(i,k)=zrfl(i)
485	ENDIF
486	ENDDO
487	c
488	c Calculer les tendances de q et de t:
489	c
490	DO i = 1, klon
491	d_q(i,k) = zq(i) - q(i,k)
492	d_t(i,k) = zt(i) - t(i,k)
493	ENDDO
494	c
495	cAA--------------- Calcul du lessivage stratiforme -------------
496
497	DO i = 1,klon
498	c
499	zprec_cond(i) = MAX(zcond(i)-zoliq(i),0.0)
500	. * (paprs(i,k)-paprs(i,k+1))/RG
501	IF (rneb(i,k).GT.0.0.and.zprec_cond(i).gt.0.) THEN
502	cAA lessivage nucleation LMD5 dans la couche elle-meme
503	if (t(i,k) .GE. ztglace) THEN
504	zalpha_tr = a_tr_sca(3)
505	else
506	zalpha_tr = a_tr_sca(4)
507	endif
508	zfrac_lessi = 1. - EXP(zalpha_tr*zprec_cond(i)/zneb(i))
509	pfrac_nucl(i,k)=pfrac_nucl(i,k)(1.-zneb(i)zfrac_lessi)
510	frac_nucl(i,k)= 1.-zneb(i)*zfrac_lessi
511	c
512	c nucleation avec un facteur -1 au lieu de -0.5
513	zfrac_lessi = 1. - EXP(-zprec_cond(i)/zneb(i))
514	pfrac_1nucl(i,k)=pfrac_1nucl(i,k)(1.-zneb(i)zfrac_lessi)
515	ENDIF
516	c
517	ENDDO ! boucle sur i
518	c
519	cAA Lessivage par impaction dans les couches en-dessous
520	DO kk = k-1, 1, -1
521	DO i = 1, klon
522	IF (rneb(i,k).GT.0.0.and.zprec_cond(i).gt.0.) THEN
523	if (t(i,kk) .GE. ztglace) THEN
524	zalpha_tr = a_tr_sca(1)
525	else
526	zalpha_tr = a_tr_sca(2)
527	endif
528	zfrac_lessi = 1. - EXP(zalpha_tr*zprec_cond(i)/zneb(i))
529	pfrac_impa(i,kk)=pfrac_impa(i,kk)(1.-zneb(i)zfrac_lessi)
530	frac_impa(i,kk)= 1.-zneb(i)*zfrac_lessi
531	ENDIF
532	ENDDO
533	ENDDO
534	c
535	cAA----------------------------------------------------------
536	c FIN DE BOUCLE SUR K
537	9999 CONTINUE
538	c
539	cAA-----------------------------------------------------------
540	c
541	c Pluie ou neige au sol selon la temperature de la 1ere couche
542	c
543	DO i = 1, klon
544	IF ((t(i,1)+d_t(i,1)) .LT. RTT) THEN
545	snow(i) = zrfl(i)
546	zlh_solid(i) = RLSTT-RLVTT
547	ELSE
548	rain(i) = zrfl(i)
549	zlh_solid(i) = 0.
550	ENDIF
551	ENDDO
552	C
553	C For energy conservation : when snow is present, the solification
554	c latent heat is considered.
555	DO k = 1, klev
556	DO i = 1, klon
557	zcpair=RCPD(1.0+RVTMP2(q(i,k)+d_q(i,k)))
558	zmair=(paprs(i,k)-paprs(i,k+1))/RG
559	zm_solid = (prfl(i,k)-prfl(i,k+1)+psfl(i,k)-psfl(i,k+1))*dtime
560	d_t(i,k) = d_t(i,k) + zlh_solid(i) zm_solid / (zcpairzmair)
561	END DO
562	END DO
563	c
564
565	if (ncoreczq>0) then
566	print*,'WARNING : ZQ dans fisrtilp ',ncoreczq,' val < 1.e-15.'
567	endif
568	RETURN
569	END

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

Download in other formats:

Original Format