/[lmdze]/trunk/libf/phylmd/newmicro.f
ViewVC logotype

Contents of /trunk/libf/phylmd/newmicro.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 38 - (show annotations)
Thu Jan 6 17:52:19 2011 UTC (13 years, 3 months ago) by guez
File size: 11829 byte(s)
Extracted ASCII art from "inigeom" into a separate text file in the
documentation.

"test_disvert" now creates a separate file for layer thicknesses.

Moved variables from module "yomcst" to module "suphec_m" because this
is where those variables are defined. Kept in "yomcst" only parameters
of Earth orbit. Gave the attribute "parameter" to some variables of
module "suphec_m".

Variables of module "yoethf" were defined in procedure "suphec". Moved
these definitions to a new procedure "yoethf" in module "yoethf_m".

1 !
2 ! $Header: /home/cvsroot/LMDZ4/libf/phylmd/newmicro.F,v 1.2 2004/06/03 09:22:43 lmdzadmin Exp $
3 !
4 SUBROUTINE newmicro (paprs, pplay,ok_newmicro,
5 . t, pqlwp, pclc, pcltau, pclemi,
6 . pch, pcl, pcm, pct, pctlwp,
7 s xflwp, xfiwp, xflwc, xfiwc,
8 e ok_aie,
9 e sulfate, sulfate_pi,
10 e bl95_b0, bl95_b1,
11 s cldtaupi, re, fl)
12 use dimens_m
13 use dimphy
14 use SUPHEC_M
15 use nuagecom
16 IMPLICIT none
17 c======================================================================
18 c Auteur(s): Z.X. Li (LMD/CNRS) date: 19930910
19 c Objet: Calculer epaisseur optique et emmissivite des nuages
20 c======================================================================
21 c Arguments:
22 c t-------input-R-temperature
23 c pqlwp---input-R-eau liquide nuageuse dans l'atmosphere (kg/kg)
24 c pclc----input-R-couverture nuageuse pour le rayonnement (0 a 1)
25 c
26 c ok_aie--input-L-apply aerosol indirect effect or not
27 c sulfate-input-R-sulfate aerosol mass concentration [um/m^3]
28 c sulfate_pi-input-R-dito, pre-industrial value
29 c bl95_b0-input-R-a parameter, may be varied for tests (s-sea, l-land)
30 c bl95_b1-input-R-a parameter, may be varied for tests ( -"- )
31 c
32 c cldtaupi-output-R-pre-industrial value of cloud optical thickness,
33 c needed for the diagnostics of the aerosol indirect
34 c radiative forcing (see radlwsw)
35 c re------output-R-Cloud droplet effective radius multiplied by fl [um]
36 c fl------output-R-Denominator to re, introduced to avoid problems in
37 c the averaging of the output. fl is the fraction of liquid
38 c water clouds within a grid cell
39 c pcltau--output-R-epaisseur optique des nuages
40 c pclemi--output-R-emissivite des nuages (0 a 1)
41 c======================================================================
42 C
43 c
44 REAL, intent(in):: paprs(klon,klev+1)
45 real, intent(in):: pplay(klon,klev)
46 REAL t(klon,klev)
47 c
48 REAL pclc(klon,klev)
49 REAL pqlwp(klon,klev)
50 REAL pcltau(klon,klev), pclemi(klon,klev)
51 c
52 REAL pct(klon), pctlwp(klon), pch(klon), pcl(klon), pcm(klon)
53 c
54 LOGICAL lo
55 c
56 REAL cetahb, cetamb
57 PARAMETER (cetahb = 0.45, cetamb = 0.80)
58 C
59 INTEGER i, k
60 cIM: 091003 REAL zflwp, zradef, zfice, zmsac
61 REAL zflwp(klon), zradef, zfice, zmsac
62 cIM: 091003 rajout
63 REAL xflwp(klon), xfiwp(klon)
64 REAL xflwc(klon,klev), xfiwc(klon,klev)
65 c
66 REAL radius, rad_chaud
67 cc PARAMETER (rad_chau1=13.0, rad_chau2=9.0, rad_froid=35.0)
68 ccc PARAMETER (rad_chaud=15.0, rad_froid=35.0)
69 c sintex initial PARAMETER (rad_chaud=10.0, rad_froid=30.0)
70 REAL coef, coef_froi, coef_chau
71 PARAMETER (coef_chau=0.13, coef_froi=0.09)
72 REAL seuil_neb, t_glace
73 PARAMETER (seuil_neb=0.001, t_glace=273.0-15.0)
74 INTEGER nexpo ! exponentiel pour glace/eau
75 PARAMETER (nexpo=6)
76 ccc PARAMETER (nexpo=1)
77
78 c -- sb:
79 logical ok_newmicro
80 c parameter (ok_newmicro=.FALSE.)
81 cIM: 091003 real rel, tc, rei, zfiwp
82 real rel, tc, rei, zfiwp(klon)
83 real k_liq, k_ice0, k_ice, DF
84 parameter (k_liq=0.0903, k_ice0=0.005) ! units=m2/g
85 parameter (DF=1.66) ! diffusivity factor
86 c sb --
87 cjq for the aerosol indirect effect
88 cjq introduced by Johannes Quaas (quaas@lmd.jussieu.fr), 27/11/2003
89 cjq
90 LOGICAL ok_aie ! Apply AIE or not?
91 LOGICAL ok_a1lwpdep ! a1 LWP dependent?
92
93 REAL sulfate(klon, klev) ! sulfate aerosol mass concentration [ug m-3]
94 REAL cdnc(klon, klev) ! cloud droplet number concentration [m-3]
95 REAL re(klon, klev) ! cloud droplet effective radius [um]
96 REAL sulfate_pi(klon, klev) ! sulfate aerosol mass concentration [ug m-3] (pre-industrial value)
97 REAL cdnc_pi(klon, klev) ! cloud droplet number concentration [m-3] (pi value)
98 REAL re_pi(klon, klev) ! cloud droplet effective radius [um] (pi value)
99
100 REAL fl(klon, klev) ! xliq * rneb (denominator to re; fraction of liquid water clouds within the grid cell)
101
102 REAL bl95_b0, bl95_b1 ! Parameter in B&L 95-Formula
103
104 REAL cldtaupi(klon, klev) ! pre-industrial cloud opt thickness for diag
105 cjq-end
106 c
107 c Calculer l'epaisseur optique et l'emmissivite des nuages
108 c
109 cIM inversion des DO
110 DO i = 1, klon
111 xflwp(i)=0.
112 xfiwp(i)=0.
113 DO k = 1, klev
114 c
115 xflwc(i,k)=0.
116 xfiwc(i,k)=0.
117 c
118 rad_chaud = rad_chau1
119 IF (k.LE.3) rad_chaud = rad_chau2
120 pclc(i,k) = MAX(pclc(i,k), seuil_neb)
121 zflwp(i) = 1000.*pqlwp(i,k)/RG/pclc(i,k)
122 . *(paprs(i,k)-paprs(i,k+1))
123 zfice = 1.0 - (t(i,k)-t_glace) / (273.13-t_glace)
124 zfice = MIN(MAX(zfice,0.0),1.0)
125 zfice = zfice**nexpo
126 radius = rad_chaud * (1.-zfice) + rad_froid * zfice
127 coef = coef_chau * (1.-zfice) + coef_froi * zfice
128 pcltau(i,k) = 3.0/2.0 * zflwp(i) / radius
129 pclemi(i,k) = 1.0 - EXP( - coef * zflwp(i))
130
131 if (ok_newmicro) then
132
133 c -- liquid/ice cloud water paths:
134
135 zfice = 1.0 - (t(i,k)-t_glace) / (273.13-t_glace)
136 zfice = MIN(MAX(zfice,0.0),1.0)
137
138 zflwp(i) = 1000.*(1.-zfice)*pqlwp(i,k)/pclc(i,k)
139 : *(paprs(i,k)-paprs(i,k+1))/RG
140 zfiwp(i) = 1000.*zfice*pqlwp(i,k)/pclc(i,k)
141 : *(paprs(i,k)-paprs(i,k+1))/RG
142
143 xflwp(i) = xflwp(i)+ (1.-zfice)*pqlwp(i,k)
144 : *(paprs(i,k)-paprs(i,k+1))/RG
145 xfiwp(i) = xfiwp(i)+ zfice*pqlwp(i,k)
146 : *(paprs(i,k)-paprs(i,k+1))/RG
147
148 cIM Total Liquid/Ice water content
149 xflwc(i,k) = xflwc(i,k)+(1.-zfice)*pqlwp(i,k)
150 xfiwc(i,k) = xfiwc(i,k)+zfice*pqlwp(i,k)
151 cIM In-Cloud Liquid/Ice water content
152 c xflwc(i,k) = xflwc(i,k)+(1.-zfice)*pqlwp(i,k)/pclc(i,k)
153 c xfiwc(i,k) = xfiwc(i,k)+zfice*pqlwp(i,k)/pclc(i,k)
154
155 c -- effective cloud droplet radius (microns):
156
157 c for liquid water clouds:
158 IF (ok_aie) THEN
159 ! Formula "D" of Boucher and Lohmann, Tellus, 1995
160 !
161 cdnc(i,k) = 10.**(bl95_b0+bl95_b1*
162 . log(MAX(sulfate(i,k),1.e-4))/log(10.))*1.e6 !-m-3
163 ! Cloud droplet number concentration (CDNC) is restricted
164 ! to be within [20, 1000 cm^3]
165 !
166 cdnc(i,k)=MIN(1000.e6,MAX(20.e6,cdnc(i,k)))
167 !
168 !
169 cdnc_pi(i,k) = 10.**(bl95_b0+bl95_b1*
170 . log(MAX(sulfate_pi(i,k),1.e-4))/log(10.))*1.e6 !-m-3
171 cdnc_pi(i,k)=MIN(1000.e6,MAX(20.e6,cdnc_pi(i,k)))
172 !
173 !
174 ! air density: pplay(i,k) / (RD * zT(i,k))
175 ! factor 1.1: derive effective radius from volume-mean radius
176 ! factor 1000 is the water density
177 ! _chaud means that this is the CDR for liquid water clouds
178 !
179 rad_chaud =
180 . 1.1 * ( (pqlwp(i,k) * pplay(i,k) / (RD * T(i,k)) )
181 . / (4./3. * RPI * 1000. * cdnc(i,k)) )**(1./3.)
182 !
183 ! Convert to um. CDR shall be at least 3 um.
184 !
185 c rad_chaud = MAX(rad_chaud*1.e6, 3.)
186 rad_chaud = MAX(rad_chaud*1.e6, 5.)
187
188 ! Pre-industrial cloud opt thickness
189 !
190 ! "radius" is calculated as rad_chaud above (plus the
191 ! ice cloud contribution) but using cdnc_pi instead of
192 ! cdnc.
193 radius =
194 . 1.1 * ( (pqlwp(i,k) * pplay(i,k) / (RD * T(i,k)) )
195 . / (4./3. * RPI * 1000. * cdnc_pi(i,k)) )**(1./3.)
196 radius = MAX(radius*1.e6, 5.)
197
198 tc = t(i,k)-273.15
199 rei = 0.71*tc + 61.29
200 if (tc.le.-81.4) rei = 3.5
201 if (zflwp(i).eq.0.) radius = 1.
202 if (zfiwp(i).eq.0. .or. rei.le.0.) rei = 1.
203 cldtaupi(i,k) = 3.0/2.0 * zflwp(i) / radius
204 . + zfiwp(i) * (3.448e-03 + 2.431/rei)
205 ENDIF ! ok_aie
206 ! For output diagnostics
207 !
208 ! Cloud droplet effective radius [um]
209 !
210 ! we multiply here with f * xl (fraction of liquid water
211 ! clouds in the grid cell) to avoid problems in the
212 ! averaging of the output.
213 ! In the output of IOIPSL, derive the real cloud droplet
214 ! effective radius as re/fl
215 !
216 fl(i,k) = pclc(i,k)*(1.-zfice)
217 re(i,k) = rad_chaud*fl(i,k)
218
219 c-jq end
220
221 rel = rad_chaud
222 c for ice clouds: as a function of the ambiant temperature
223 c [formula used by Iacobellis and Somerville (2000), with an
224 c asymptotical value of 3.5 microns at T<-81.4 C added to be
225 c consistent with observations of Heymsfield et al. 1986]:
226 tc = t(i,k)-273.15
227 rei = 0.71*tc + 61.29
228 if (tc.le.-81.4) rei = 3.5
229
230 c -- cloud optical thickness :
231
232 c [for liquid clouds, traditional formula,
233 c for ice clouds, Ebert & Curry (1992)]
234
235 if (zflwp(i).eq.0.) rel = 1.
236 if (zfiwp(i).eq.0. .or. rei.le.0.) rei = 1.
237 pcltau(i,k) = 3.0/2.0 * ( zflwp(i)/rel )
238 . + zfiwp(i) * (3.448e-03 + 2.431/rei)
239
240 c -- cloud infrared emissivity:
241
242 c [the broadband infrared absorption coefficient is parameterized
243 c as a function of the effective cld droplet radius]
244
245 c Ebert and Curry (1992) formula as used by Kiehl & Zender (1995):
246 k_ice = k_ice0 + 1.0/rei
247
248 pclemi(i,k) = 1.0
249 . - EXP( - coef_chau*zflwp(i) - DF*k_ice*zfiwp(i) )
250
251 endif ! ok_newmicro
252
253 lo = (pclc(i,k) .LE. seuil_neb)
254 IF (lo) pclc(i,k) = 0.0
255 IF (lo) pcltau(i,k) = 0.0
256 IF (lo) pclemi(i,k) = 0.0
257
258 IF (lo) cldtaupi(i,k) = 0.0
259 IF (.NOT.ok_aie) cldtaupi(i,k)=pcltau(i,k)
260 ENDDO
261 ENDDO
262 ccc DO k = 1, klev
263 ccc DO i = 1, klon
264 ccc t(i,k) = t(i,k)
265 ccc pclc(i,k) = MAX( 1.e-5 , pclc(i,k) )
266 ccc lo = pclc(i,k) .GT. (2.*1.e-5)
267 ccc zflwp = pqlwp(i,k)*1000.*(paprs(i,k)-paprs(i,k+1))
268 ccc . /(rg*pclc(i,k))
269 ccc zradef = 10.0 + (1.-sigs(k))*45.0
270 ccc pcltau(i,k) = 1.5 * zflwp / zradef
271 ccc zfice=1.0-MIN(MAX((t(i,k)-263.)/(273.-263.),0.0),1.0)
272 ccc zmsac = 0.13*(1.0-zfice) + 0.08*zfice
273 ccc pclemi(i,k) = 1.-EXP(-zmsac*zflwp)
274 ccc if (.NOT.lo) pclc(i,k) = 0.0
275 ccc if (.NOT.lo) pcltau(i,k) = 0.0
276 ccc if (.NOT.lo) pclemi(i,k) = 0.0
277 ccc ENDDO
278 ccc ENDDO
279 cccccc print*, 'pas de nuage dans le rayonnement'
280 cccccc DO k = 1, klev
281 cccccc DO i = 1, klon
282 cccccc pclc(i,k) = 0.0
283 cccccc pcltau(i,k) = 0.0
284 cccccc pclemi(i,k) = 0.0
285 cccccc ENDDO
286 cccccc ENDDO
287 C
288 C COMPUTE CLOUD LIQUID PATH AND TOTAL CLOUDINESS
289 C
290 DO i = 1, klon
291 pct(i)=1.0
292 pch(i)=1.0
293 pcm(i) = 1.0
294 pcl(i) = 1.0
295 pctlwp(i) = 0.0
296 ENDDO
297 C
298 DO k = klev, 1, -1
299 DO i = 1, klon
300 pctlwp(i) = pctlwp(i)
301 . + pqlwp(i,k)*(paprs(i,k)-paprs(i,k+1))/RG
302 pct(i) = pct(i)*(1.0-pclc(i,k))
303 if (pplay(i,k).LE.cetahb*paprs(i,1))
304 . pch(i) = pch(i)*(1.0-pclc(i,k))
305 if (pplay(i,k).GT.cetahb*paprs(i,1) .AND.
306 . pplay(i,k).LE.cetamb*paprs(i,1))
307 . pcm(i) = pcm(i)*(1.0-pclc(i,k))
308 if (pplay(i,k).GT.cetamb*paprs(i,1))
309 . pcl(i) = pcl(i)*(1.0-pclc(i,k))
310 ENDDO
311 ENDDO
312 C
313 DO i = 1, klon
314 pct(i)=1.-pct(i)
315 pch(i)=1.-pch(i)
316 pcm(i)=1.-pcm(i)
317 pcl(i)=1.-pcl(i)
318 ENDDO
319 C
320 RETURN
321 END

  ViewVC Help
Powered by ViewVC 1.1.21