/[lmdze]/trunk/phylmd/Radlwsw/radlwsw.f
ViewVC logotype

Contents of /trunk/phylmd/Radlwsw/radlwsw.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 82 - (show annotations)
Wed Mar 5 14:57:53 2014 UTC (10 years, 2 months ago) by guez
File size: 13917 byte(s)
Changed all ".f90" suffixes to ".f".
1 module radlwsw_m
2
3 IMPLICIT none
4
5 contains
6
7 SUBROUTINE radlwsw(dist, rmu0, fract, paprs, play, tsol, albedo, alblw, &
8 t, q, wo, cldfra, cldemi, cldtaupd, heat, heat0, cool, cool0, radsol, &
9 albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, toplw0, solsw0, &
10 sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, swup0, swup, ok_ade, &
11 ok_aie, tau_ae, piz_ae, cg_ae, topswad, solswad, cldtaupi, topswai, &
12 solswai)
13
14 ! From LMDZ4/libf/phylmd/radlwsw.F, version 1.4 2005/06/06 13:16:33
15 ! Author: Z. X. Li (LMD/CNRS)
16 ! Date: 1996/07/19
17
18 ! Objet : interface entre le modèle et les rayonnements solaire et
19 ! infrarouge
20
21 ! ATTENTION: swai and swad have to be interpreted in the following manner:
22
23 ! not ok_ade and not ok_aie
24 ! both are zero
25
26 ! ok_ade and not ok_aie
27 ! aerosol direct forcing is F_{AD} = topsw - topswad
28 ! indirect is zero
29
30 ! not ok_ade and ok_aie
31 ! aerosol indirect forcing is F_{AI} = topsw - topswai
32 ! direct is zero
33
34 ! ok_ade and ok_aie
35 ! aerosol indirect forcing is F_{AI} = topsw - topswai
36 ! aerosol direct forcing is F_{AD} = topswai - topswad
37
38 USE clesphys, ONLY: bug_ozone, solaire
39 USE dimphy, ONLY: klev, klon
40 use lw_m, only: lw
41 USE raddim, ONLY: kdlon
42 USE suphec_m, ONLY: rg
43 use sw_m, only: sw
44 USE yoethf_m, ONLY: rvtmp2
45
46 ! Arguments:
47
48 real dist, rmu0(klon), fract(klon)
49 ! dist-----input-R- distance astronomique terre-soleil
50 ! rmu0-----input-R- cosinus de l'angle zenithal
51 ! fract----input-R- duree d'ensoleillement normalisee
52
53 real, intent(in):: paprs(klon, klev+1)
54 ! paprs----input-R- pression a inter-couche (Pa)
55 real, intent(in):: play(klon, klev)
56 ! play----input-R- pression au milieu de couche (Pa)
57 real tsol(klon), albedo(klon), alblw(klon)
58 ! albedo---input-R- albedo du sol (entre 0 et 1)
59 ! tsol-----input-R- temperature du sol (en K)
60 real, intent(in):: t(klon, klev)
61 ! t--------input-R- temperature (K)
62 real q(klon, klev)
63 ! q--------input-R- vapeur d'eau (en kg/kg)
64 real, intent(in):: wo(klon, klev)
65 ! wo-------input-R- contenu en ozone (en kg/kg) correction MPL 100505
66 real cldfra(klon, klev), cldemi(klon, klev)
67 ! cldfra---input-R- fraction nuageuse (entre 0 et 1)
68 ! cldemi---input-R- emissivite des nuages dans l'IR (entre 0 et 1)
69
70 real cldtaupd(klon, klev)
71 ! input-R- epaisseur optique des nuages dans le visible (present-day value)
72
73 real, intent(out):: heat(klon, klev)
74 ! échauffement atmosphérique (visible) (K/jour)
75
76 real heat0(klon, klev)
77 real cool(klon, klev)
78 ! cool-----output-R- refroidissement dans l'IR (K/jour)
79 real cool0(klon, klev)
80 real radsol(klon)
81 ! radsol---output-R- bilan radiatif net au sol (W/m**2) (+ vers le bas)
82 real albpla(klon)
83 ! albpla---output-R- albedo planetaire (entre 0 et 1)
84 real topsw(klon)
85 ! topsw----output-R- flux solaire net au sommet de l'atm.
86
87 real, intent(out):: toplw(klon)
88 ! rayonnement infrarouge montant au sommet de l'atmosphère
89
90 real, intent(out):: solsw(klon) ! flux solaire net à la surface
91
92 real, intent(out):: sollw(klon)
93 ! rayonnement infrarouge montant à la surface
94
95 real, intent(out):: sollwdown(klon)
96 real topsw0(klon)
97 real, intent(out):: toplw0(klon)
98 real solsw0(klon), sollw0(klon)
99 !IM output 3D: SWup, SWdn, LWup, LWdn
100 REAL lwdn0(klon, klev+1), lwdn(klon, klev+1)
101 REAL lwup0(klon, klev+1), lwup(klon, klev+1)
102 REAL swdn0(klon, klev+1), swdn(klon, klev+1)
103 REAL swup0(klon, klev+1), swup(klon, klev+1)
104
105 logical ok_ade, ok_aie
106 ! switches whether to use aerosol direct (indirect) effects or not
107 ! ok_ade---input-L- apply the Aerosol Direct Effect or not?
108 ! ok_aie---input-L- apply the Aerosol Indirect Effect or not?
109
110 real tau_ae(klon, klev, 2), piz_ae(klon, klev, 2), cg_ae(klon, klev, 2)
111 ! input-R- aerosol optical properties (calculated in aeropt.F)
112
113 real topswad(klon), solswad(klon)
114 ! output: aerosol direct forcing at TOA and surface
115 ! topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
116 ! solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
117
118 real cldtaupi(klon, klev)
119 ! cloud optical thickness for pre-industrial aerosol concentrations
120 ! (i.e. with a smaller droplet concentration and thus larger droplet radii)
121 ! -input-R- epaisseur optique des nuages dans le visible
122 ! calculated for pre-industrial (pi) aerosol concentrations,
123 ! i.e. with smaller droplet concentration, thus larger droplets,
124 ! thus generally cdltaupi cldtaupd it is needed for the
125 ! diagnostics of the aerosol indirect radiative forcing
126
127 real topswai(klon), solswai(klon)
128 ! output: aerosol indirect forcing atTOA and surface
129 ! topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
130 ! solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
131
132 ! Local:
133
134 double precision tauae(kdlon, klev, 2) ! aer opt properties
135 double precision pizae(kdlon, klev, 2)
136 double precision cgae(kdlon, klev, 2)
137
138 !IM output 3D
139 DOUBLE PRECISION ZFSUP(KDLON, KLEV+1)
140 DOUBLE PRECISION ZFSDN(KDLON, KLEV+1)
141 DOUBLE PRECISION ZFSUP0(KDLON, KLEV+1)
142 DOUBLE PRECISION ZFSDN0(KDLON, KLEV+1)
143
144 DOUBLE PRECISION ZFLUP(KDLON, KLEV+1)
145 DOUBLE PRECISION ZFLDN(KDLON, KLEV+1)
146 DOUBLE PRECISION ZFLUP0(KDLON, KLEV+1)
147 DOUBLE PRECISION ZFLDN0(KDLON, KLEV+1)
148
149 DOUBLE PRECISION zx_alpha1, zx_alpha2
150 INTEGER k, kk, i, iof, nb_gr
151 DOUBLE PRECISION PSCT
152
153 DOUBLE PRECISION PALBD(kdlon, 2), PALBP(kdlon, 2)
154 DOUBLE PRECISION PEMIS(kdlon), PDT0(kdlon), PVIEW(kdlon)
155 DOUBLE PRECISION PPSOL(kdlon), PDP(kdlon, klev)
156 DOUBLE PRECISION PTL(kdlon, klev+1), PPMB(kdlon, klev+1)
157 DOUBLE PRECISION PTAVE(kdlon, klev)
158 DOUBLE PRECISION PWV(kdlon, klev), PQS(kdlon, klev), POZON(kdlon, klev)
159 DOUBLE PRECISION PAER(kdlon, klev, 5)
160 DOUBLE PRECISION PCLDLD(kdlon, klev)
161 DOUBLE PRECISION PCLDLU(kdlon, klev)
162 DOUBLE PRECISION PCLDSW(kdlon, klev)
163 DOUBLE PRECISION PTAU(kdlon, 2, klev)
164 DOUBLE PRECISION POMEGA(kdlon, 2, klev)
165 DOUBLE PRECISION PCG(kdlon, 2, klev)
166
167 DOUBLE PRECISION zfract(kdlon), zrmu0(kdlon), zdist
168
169 DOUBLE PRECISION zheat(kdlon, klev), zcool(kdlon, klev)
170 DOUBLE PRECISION zheat0(kdlon, klev), zcool0(kdlon, klev)
171 DOUBLE PRECISION ztopsw(kdlon), ztoplw(kdlon)
172 DOUBLE PRECISION zsolsw(kdlon), zsollw(kdlon), zalbpla(kdlon)
173 DOUBLE PRECISION zsollwdown(kdlon)
174
175 DOUBLE PRECISION ztopsw0(kdlon), ztoplw0(kdlon)
176 DOUBLE PRECISION zsolsw0(kdlon), zsollw0(kdlon)
177 DOUBLE PRECISION zznormcp
178
179 !jq the following quantities are needed for the aerosol radiative forcings
180
181 DOUBLE PRECISION PTAUA(kdlon, 2, klev)
182 ! present-day value of cloud opt thickness (PTAU is pre-industrial
183 ! value), local use
184
185 DOUBLE PRECISION POMEGAA(kdlon, 2, klev) ! dito for single scatt albedo
186
187 DOUBLE PRECISION ztopswad(kdlon), zsolswad(kdlon)
188 ! Aerosol direct forcing at TOAand surface
189
190 DOUBLE PRECISION ztopswai(kdlon), zsolswai(kdlon) ! dito, indirect
191
192 !----------------------------------------------------------------------
193
194 tauae = 0.
195 pizae = 0.
196 cgae = 0.
197
198 nb_gr = klon / kdlon
199 IF (nb_gr * kdlon /= klon) THEN
200 PRINT *, "kdlon mauvais :", klon, kdlon, nb_gr
201 stop 1
202 ENDIF
203
204 heat = 0.
205 cool = 0.
206 heat0 = 0.
207 cool0 = 0.
208 zdist = dist
209 PSCT = solaire / zdist / zdist
210
211 loop_iof: DO iof = 0, klon - kdlon, kdlon
212 DO i = 1, kdlon
213 zfract(i) = fract(iof+i)
214 zrmu0(i) = rmu0(iof+i)
215 PALBD(i, 1) = albedo(iof+i)
216 PALBD(i, 2) = alblw(iof+i)
217 PALBP(i, 1) = albedo(iof+i)
218 PALBP(i, 2) = alblw(iof+i)
219 ! cf. JLD pour etre en accord avec ORCHIDEE il faut mettre
220 ! PEMIS(i) = 0.96
221 PEMIS(i) = 1.0
222 PVIEW(i) = 1.66
223 PPSOL(i) = paprs(iof+i, 1)
224 zx_alpha1 = (paprs(iof+i, 1)-play(iof+i, 2)) &
225 / (play(iof+i, 1)-play(iof+i, 2))
226 zx_alpha2 = 1.0 - zx_alpha1
227 PTL(i, 1) = t(iof+i, 1) * zx_alpha1 + t(iof+i, 2) * zx_alpha2
228 PTL(i, klev+1) = t(iof+i, klev)
229 PDT0(i) = tsol(iof+i) - PTL(i, 1)
230 ENDDO
231 DO k = 2, klev
232 DO i = 1, kdlon
233 PTL(i, k) = (t(iof+i, k)+t(iof+i, k-1))*0.5
234 ENDDO
235 ENDDO
236 DO k = 1, klev
237 DO i = 1, kdlon
238 PDP(i, k) = paprs(iof+i, k)-paprs(iof+i, k+1)
239 PTAVE(i, k) = t(iof+i, k)
240 PWV(i, k) = MAX (q(iof+i, k), 1.0e-12)
241 PQS(i, k) = PWV(i, k)
242 ! wo: cm.atm (epaisseur en cm dans la situation standard)
243 ! POZON: kg/kg
244 IF (bug_ozone) then
245 POZON(i, k) = MAX(wo(iof+i, k), 1.0e-12)*RG/46.6968 &
246 /(paprs(iof+i, k)-paprs(iof+i, k+1)) &
247 *(paprs(iof+i, 1)/101325.0)
248 ELSE
249 ! le calcul qui suit est maintenant fait dans ozonecm (MPL)
250 POZON(i, k) = wo(i, k)
251 ENDIF
252 PCLDLD(i, k) = cldfra(iof+i, k)*cldemi(iof+i, k)
253 PCLDLU(i, k) = cldfra(iof+i, k)*cldemi(iof+i, k)
254 PCLDSW(i, k) = cldfra(iof+i, k)
255 PTAU(i, 1, k) = MAX(cldtaupi(iof+i, k), 1.0e-05)
256 ! (1e-12 serait instable)
257 PTAU(i, 2, k) = MAX(cldtaupi(iof+i, k), 1.0e-05)
258 ! (pour 32-bit machines)
259 POMEGA(i, 1, k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i, 1, k))
260 POMEGA(i, 2, k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i, 2, k))
261 PCG(i, 1, k) = 0.865
262 PCG(i, 2, k) = 0.910
263
264 ! Introduced for aerosol indirect forcings. The
265 ! following values use the cloud optical thickness
266 ! calculated from present-day aerosol concentrations
267 ! whereas the quantities without the "A" at the end are
268 ! for pre-industial (natural-only) aerosol concentrations
269 PTAUA(i, 1, k) = MAX(cldtaupd(iof+i, k), 1.0e-05)
270 ! (1e-12 serait instable)
271 PTAUA(i, 2, k) = MAX(cldtaupd(iof+i, k), 1.0e-05)
272 ! (pour 32-bit machines)
273 POMEGAA(i, 1, k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i, 1, k))
274 POMEGAA(i, 2, k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i, 2, k))
275 !jq-end
276 ENDDO
277 ENDDO
278
279 DO k = 1, klev+1
280 DO i = 1, kdlon
281 PPMB(i, k) = paprs(iof+i, k)/100.0
282 ENDDO
283 ENDDO
284
285 DO kk = 1, 5
286 DO k = 1, klev
287 DO i = 1, kdlon
288 PAER(i, k, kk) = 1.0E-15
289 ENDDO
290 ENDDO
291 ENDDO
292
293 DO k = 1, klev
294 DO i = 1, kdlon
295 tauae(i, k, 1) = tau_ae(iof+i, k, 1)
296 pizae(i, k, 1) = piz_ae(iof+i, k, 1)
297 cgae(i, k, 1) =cg_ae(iof+i, k, 1)
298 tauae(i, k, 2) = tau_ae(iof+i, k, 2)
299 pizae(i, k, 2) = piz_ae(iof+i, k, 2)
300 cgae(i, k, 2) =cg_ae(iof+i, k, 2)
301 ENDDO
302 ENDDO
303
304 CALL LW(PPMB, PDP, PPSOL, PDT0, PEMIS, PTL, PTAVE, PWV, POZON, PAER, &
305 PCLDLD, PCLDLU, PVIEW, zcool, zcool0, ztoplw, zsollw, ztoplw0, &
306 zsollw0, zsollwdown, ZFLUP, ZFLDN, ZFLUP0, ZFLDN0)
307 CALL SW(PSCT, zrmu0, zfract, PPMB, PDP, PPSOL, PALBD, PALBP, PTAVE, &
308 PWV, PQS, POZON, PAER, PCLDSW, PTAU, POMEGA, PCG, zheat, zheat0, &
309 zalbpla, ztopsw, zsolsw, ztopsw0, zsolsw0, ZFSUP, ZFSDN, ZFSUP0, &
310 ZFSDN0, tauae, pizae, cgae, PTAUA, POMEGAA, ztopswad, zsolswad, &
311 ztopswai, zsolswai, ok_ade, ok_aie)
312
313 DO i = 1, kdlon
314 radsol(iof+i) = zsolsw(i) + zsollw(i)
315 topsw(iof+i) = ztopsw(i)
316 toplw(iof+i) = ztoplw(i)
317 solsw(iof+i) = zsolsw(i)
318 sollw(iof+i) = zsollw(i)
319 sollwdown(iof+i) = zsollwdown(i)
320
321 DO k = 1, klev+1
322 lwdn0 ( iof+i, k) = ZFLDN0 ( i, k)
323 lwdn ( iof+i, k) = ZFLDN ( i, k)
324 lwup0 ( iof+i, k) = ZFLUP0 ( i, k)
325 lwup ( iof+i, k) = ZFLUP ( i, k)
326 ENDDO
327
328 topsw0(iof+i) = ztopsw0(i)
329 toplw0(iof+i) = ztoplw0(i)
330 solsw0(iof+i) = zsolsw0(i)
331 sollw0(iof+i) = zsollw0(i)
332 albpla(iof+i) = zalbpla(i)
333
334 DO k = 1, klev+1
335 swdn0 ( iof+i, k) = ZFSDN0 ( i, k)
336 swdn ( iof+i, k) = ZFSDN ( i, k)
337 swup0 ( iof+i, k) = ZFSUP0 ( i, k)
338 swup ( iof+i, k) = ZFSUP ( i, k)
339 ENDDO
340 ENDDO
341 ! transform the aerosol forcings, if they have to be calculated
342 IF (ok_ade) THEN
343 DO i = 1, kdlon
344 topswad(iof+i) = ztopswad(i)
345 solswad(iof+i) = zsolswad(i)
346 ENDDO
347 ELSE
348 DO i = 1, kdlon
349 topswad(iof+i) = 0.0
350 solswad(iof+i) = 0.0
351 ENDDO
352 ENDIF
353 IF (ok_aie) THEN
354 DO i = 1, kdlon
355 topswai(iof+i) = ztopswai(i)
356 solswai(iof+i) = zsolswai(i)
357 ENDDO
358 ELSE
359 DO i = 1, kdlon
360 topswai(iof+i) = 0.0
361 solswai(iof+i) = 0.0
362 ENDDO
363 ENDIF
364
365 DO k = 1, klev
366 DO i = 1, kdlon
367 ! scale factor to take into account the difference
368 ! between dry air and water vapour specific heat capacity
369 zznormcp = 1. + RVTMP2 * PWV(i, k)
370 heat(iof+i, k) = zheat(i, k) / zznormcp
371 cool(iof+i, k) = zcool(i, k)/zznormcp
372 heat0(iof+i, k) = zheat0(i, k)/zznormcp
373 cool0(iof+i, k) = zcool0(i, k)/zznormcp
374 ENDDO
375 ENDDO
376 end DO loop_iof
377
378 END SUBROUTINE radlwsw
379
380 end module radlwsw_m

  ViewVC Help
Powered by ViewVC 1.1.21