/[lmdze]/trunk/libf/phylmd/fisrtilp.f90
ViewVC logotype

Annotation of /trunk/libf/phylmd/fisrtilp.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 12 - (hide annotations)
Mon Jul 21 16:05:07 2008 UTC (15 years, 10 months ago) by guez
Original Path: trunk/libf/phylmd/fisrtilp.f
File size: 16541 byte(s)
-- Minor modification of input/output:

Created procedure "read_logic". Variables of module "logic" are read
by "read_logic" instead of "conf_gcm". Variable "offline" of module
"conf_gcm" is read from namelist instead of "*.def".

Deleted arguments "dtime", "co2_ppm_etat0", "solaire_etat0",
"tabcntr0" and local variables "radpas", "tab_cntrl" of
"phyetat0". "phyetat0" does not read "controle" in "startphy.nc" any
longer. "phyetat0" now reads global attribute "itau_phy" from
"startphy.nc". "phyredem" does not create variable "controle" in
"startphy.nc" any longer. "phyredem" now writes global attribute
"itau_phy" of "startphy.nc". Deleted argument "tabcntr0" of
"printflag". Removed diagnostic messages written by "printflag" for
comparison of the variable "controle" of "startphy.nc" and the
variables read from "*.def" or namelist input.

-- Removing unwanted functionality:

Removed variable "lunout" from module "iniprint", replaced everywhere
by standard output.

Removed case "ocean == 'couple'" in "clmain", "interfsurf_hq" and
"physiq". Removed procedure "interfoce_cpl".

-- Should not change anything at run time:

Automated creation of graphs in documentation. More documentation on
input files.

Converted Fortran files to free format: "phyredem.f90", "printflag.f90".

Split module "clesphy" into "clesphys" and "clesphys2".

Removed variables "conser", "leapf", "forward", "apphys", "apdiss" and
"statcl" from module "logic". Added arguments "conser" to "advect",
"leapf" to "integrd". Added local variables "forward", "leapf",
"apphys", "conser", "apdiss" in "leapfrog".

Added intent attributes.

Deleted arguments "dtime" of "phyredem", "pdtime" of "flxdtdq", "sh"
of "phytrac", "dt" of "yamada".

Deleted local variables "dtime", "co2_ppm_etat0", "solaire_etat0",
"length", "tabcntr0" in "physiq". Replaced all references to "dtime"
by references to "pdtphys".

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

  ViewVC Help
Powered by ViewVC 1.1.21