/[lmdze]/trunk/libf/phylmd/Radlwsw/radlwsw.f90
ViewVC logotype

Contents of /trunk/libf/phylmd/Radlwsw/radlwsw.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 71 - (show annotations)
Mon Jul 8 18:12:18 2013 UTC (10 years, 9 months ago) by guez
File size: 13846 byte(s)
No reason to call inidissip in ce0l.

In inidissip, set random seed to 1 beacuse PGI compiler does not
accept all zeros.

dq was computed needlessly in caladvtrac. Arguments masse and dq of
calfis not used.

Replaced real*8 by double precision.

Pass arrays with inverted order of vertical levels to conflx instead
of creating local variables for this inside conflx.

1 module radlwsw_m
2
3 IMPLICIT none
4
5 contains
6
7 SUBROUTINE radlwsw(dist, rmu0, fract, paprs, pplay, 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 rmu0(klon), fract(klon), dist
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):: pplay(klon, klev)
56 ! pplay----input-R- pression au milieu de couche (Pa)
57 real albedo(klon), alblw(klon), tsol(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 cool(klon, klev)
77 ! cool-----output-R- refroidissement dans l'IR (K/jour)
78 real heat0(klon, klev), cool0(klon, klev)
79 real radsol(klon), topsw(klon)
80 ! radsol---output-R- bilan radiatif net au sol (W/m**2) (+ vers le bas)
81 ! topsw----output-R- flux solaire net au sommet de l'atm.
82
83 real, intent(out):: toplw(klon)
84 ! rayonnement infrarouge montant au sommet de l'atmosphère
85
86 real solsw(klon), sollw(klon), albpla(klon)
87 ! solsw----output-R- flux solaire net a la surface
88 ! sollw----output-R- ray. IR montant a la surface
89 ! albpla---output-R- albedo planetaire (entre 0 et 1)
90 real topsw0(klon), solsw0(klon), sollw0(klon)
91 real, intent(out):: toplw0(klon)
92 real sollwdown(klon)
93 !IM output 3D
94 DOUBLE PRECISION ZFSUP(KDLON, KLEV+1)
95 DOUBLE PRECISION ZFSDN(KDLON, KLEV+1)
96 DOUBLE PRECISION ZFSUP0(KDLON, KLEV+1)
97 DOUBLE PRECISION ZFSDN0(KDLON, KLEV+1)
98
99 DOUBLE PRECISION ZFLUP(KDLON, KLEV+1)
100 DOUBLE PRECISION ZFLDN(KDLON, KLEV+1)
101 DOUBLE PRECISION ZFLUP0(KDLON, KLEV+1)
102 DOUBLE PRECISION ZFLDN0(KDLON, KLEV+1)
103
104 DOUBLE PRECISION zx_alpha1, zx_alpha2
105 INTEGER k, kk, i, iof, nb_gr
106 DOUBLE PRECISION PSCT
107
108 DOUBLE PRECISION PALBD(kdlon, 2), PALBP(kdlon, 2)
109 DOUBLE PRECISION PEMIS(kdlon), PDT0(kdlon), PVIEW(kdlon)
110 DOUBLE PRECISION PPSOL(kdlon), PDP(kdlon, klev)
111 DOUBLE PRECISION PTL(kdlon, klev+1), PPMB(kdlon, klev+1)
112 DOUBLE PRECISION PTAVE(kdlon, klev)
113 DOUBLE PRECISION PWV(kdlon, klev), PQS(kdlon, klev), POZON(kdlon, klev)
114 DOUBLE PRECISION PAER(kdlon, klev, 5)
115 DOUBLE PRECISION PCLDLD(kdlon, klev)
116 DOUBLE PRECISION PCLDLU(kdlon, klev)
117 DOUBLE PRECISION PCLDSW(kdlon, klev)
118 DOUBLE PRECISION PTAU(kdlon, 2, klev)
119 DOUBLE PRECISION POMEGA(kdlon, 2, klev)
120 DOUBLE PRECISION PCG(kdlon, 2, klev)
121
122 DOUBLE PRECISION zfract(kdlon), zrmu0(kdlon), zdist
123
124 DOUBLE PRECISION zheat(kdlon, klev), zcool(kdlon, klev)
125 DOUBLE PRECISION zheat0(kdlon, klev), zcool0(kdlon, klev)
126 DOUBLE PRECISION ztopsw(kdlon), ztoplw(kdlon)
127 DOUBLE PRECISION zsolsw(kdlon), zsollw(kdlon), zalbpla(kdlon)
128 DOUBLE PRECISION zsollwdown(kdlon)
129
130 DOUBLE PRECISION ztopsw0(kdlon), ztoplw0(kdlon)
131 DOUBLE PRECISION zsolsw0(kdlon), zsollw0(kdlon)
132 DOUBLE PRECISION zznormcp
133 !IM output 3D: SWup, SWdn, LWup, LWdn
134 REAL swdn(klon, klev+1), swdn0(klon, klev+1)
135 REAL swup(klon, klev+1), swup0(klon, klev+1)
136 REAL lwdn(klon, klev+1), lwdn0(klon, klev+1)
137 REAL lwup(klon, klev+1), lwup0(klon, klev+1)
138
139 !jq the following quantities are needed for the aerosol radiative forcings
140
141 real topswad(klon), solswad(klon)
142 ! output: aerosol direct forcing at TOA and surface
143 ! topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
144 ! solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
145
146 real topswai(klon), solswai(klon)
147 ! output: aerosol indirect forcing atTOA and surface
148 ! topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
149 ! solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
150
151 real tau_ae(klon, klev, 2), piz_ae(klon, klev, 2), cg_ae(klon, klev, 2)
152 ! input-R- aerosol optical properties (calculated in aeropt.F)
153
154 real cldtaupi(klon, klev)
155 ! cloud optical thickness for pre-industrial aerosol concentrations
156 ! (i.e. with a smaller droplet concentration and thus larger droplet radii)
157 ! -input-R- epaisseur optique des nuages dans le visible
158 ! calculated for pre-industrial (pi) aerosol concentrations,
159 ! i.e. with smaller droplet concentration, thus larger droplets,
160 ! thus generally cdltaupi cldtaupd it is needed for the
161 ! diagnostics of the aerosol indirect radiative forcing
162
163 logical ok_ade, ok_aie
164 ! switches whether to use aerosol direct (indirect) effects or not
165 ! ok_ade---input-L- apply the Aerosol Direct Effect or not?
166 ! ok_aie---input-L- apply the Aerosol Indirect Effect or not?
167
168 double precision tauae(kdlon, klev, 2) ! aer opt properties
169 double precision pizae(kdlon, klev, 2)
170 double precision cgae(kdlon, klev, 2)
171
172 DOUBLE PRECISION PTAUA(kdlon, 2, klev)
173 ! present-day value of cloud opt thickness (PTAU is pre-industrial
174 ! value), local use
175
176 DOUBLE PRECISION POMEGAA(kdlon, 2, klev) ! dito for single scatt albedo
177
178 DOUBLE PRECISION ztopswad(kdlon), zsolswad(kdlon)
179 ! Aerosol direct forcing at TOAand surface
180
181 DOUBLE PRECISION ztopswai(kdlon), zsolswai(kdlon) ! dito, indirect
182
183 !----------------------------------------------------------------------
184
185 tauae = 0.
186 pizae = 0.
187 cgae = 0.
188
189 nb_gr = klon / kdlon
190 IF (nb_gr * kdlon /= klon) THEN
191 PRINT *, "kdlon mauvais :", klon, kdlon, nb_gr
192 stop 1
193 ENDIF
194
195 heat = 0.
196 cool = 0.
197 heat0 = 0.
198 cool0 = 0.
199 zdist = dist
200 PSCT = solaire / zdist / zdist
201
202 loop_iof: DO iof = 0, klon - kdlon, kdlon
203 DO i = 1, kdlon
204 zfract(i) = fract(iof+i)
205 zrmu0(i) = rmu0(iof+i)
206 PALBD(i, 1) = albedo(iof+i)
207 PALBD(i, 2) = alblw(iof+i)
208 PALBP(i, 1) = albedo(iof+i)
209 PALBP(i, 2) = alblw(iof+i)
210 ! cf. JLD pour etre en accord avec ORCHIDEE il faut mettre
211 ! PEMIS(i) = 0.96
212 PEMIS(i) = 1.0
213 PVIEW(i) = 1.66
214 PPSOL(i) = paprs(iof+i, 1)
215 zx_alpha1 = (paprs(iof+i, 1)-pplay(iof+i, 2)) &
216 / (pplay(iof+i, 1)-pplay(iof+i, 2))
217 zx_alpha2 = 1.0 - zx_alpha1
218 PTL(i, 1) = t(iof+i, 1) * zx_alpha1 + t(iof+i, 2) * zx_alpha2
219 PTL(i, klev+1) = t(iof+i, klev)
220 PDT0(i) = tsol(iof+i) - PTL(i, 1)
221 ENDDO
222 DO k = 2, klev
223 DO i = 1, kdlon
224 PTL(i, k) = (t(iof+i, k)+t(iof+i, k-1))*0.5
225 ENDDO
226 ENDDO
227 DO k = 1, klev
228 DO i = 1, kdlon
229 PDP(i, k) = paprs(iof+i, k)-paprs(iof+i, k+1)
230 PTAVE(i, k) = t(iof+i, k)
231 PWV(i, k) = MAX (q(iof+i, k), 1.0e-12)
232 PQS(i, k) = PWV(i, k)
233 ! wo: cm.atm (epaisseur en cm dans la situation standard)
234 ! POZON: kg/kg
235 IF (bug_ozone) then
236 POZON(i, k) = MAX(wo(iof+i, k), 1.0e-12)*RG/46.6968 &
237 /(paprs(iof+i, k)-paprs(iof+i, k+1)) &
238 *(paprs(iof+i, 1)/101325.0)
239 ELSE
240 ! le calcul qui suit est maintenant fait dans ozonecm (MPL)
241 POZON(i, k) = wo(i, k)
242 ENDIF
243 PCLDLD(i, k) = cldfra(iof+i, k)*cldemi(iof+i, k)
244 PCLDLU(i, k) = cldfra(iof+i, k)*cldemi(iof+i, k)
245 PCLDSW(i, k) = cldfra(iof+i, k)
246 PTAU(i, 1, k) = MAX(cldtaupi(iof+i, k), 1.0e-05)
247 ! (1e-12 serait instable)
248 PTAU(i, 2, k) = MAX(cldtaupi(iof+i, k), 1.0e-05)
249 ! (pour 32-bit machines)
250 POMEGA(i, 1, k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i, 1, k))
251 POMEGA(i, 2, k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i, 2, k))
252 PCG(i, 1, k) = 0.865
253 PCG(i, 2, k) = 0.910
254
255 ! Introduced for aerosol indirect forcings. The
256 ! following values use the cloud optical thickness
257 ! calculated from present-day aerosol concentrations
258 ! whereas the quantities without the "A" at the end are
259 ! for pre-industial (natural-only) aerosol concentrations
260 PTAUA(i, 1, k) = MAX(cldtaupd(iof+i, k), 1.0e-05)
261 ! (1e-12 serait instable)
262 PTAUA(i, 2, k) = MAX(cldtaupd(iof+i, k), 1.0e-05)
263 ! (pour 32-bit machines)
264 POMEGAA(i, 1, k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i, 1, k))
265 POMEGAA(i, 2, k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i, 2, k))
266 !jq-end
267 ENDDO
268 ENDDO
269
270 DO k = 1, klev+1
271 DO i = 1, kdlon
272 PPMB(i, k) = paprs(iof+i, k)/100.0
273 ENDDO
274 ENDDO
275
276 DO kk = 1, 5
277 DO k = 1, klev
278 DO i = 1, kdlon
279 PAER(i, k, kk) = 1.0E-15
280 ENDDO
281 ENDDO
282 ENDDO
283
284 DO k = 1, klev
285 DO i = 1, kdlon
286 tauae(i, k, 1) = tau_ae(iof+i, k, 1)
287 pizae(i, k, 1) = piz_ae(iof+i, k, 1)
288 cgae(i, k, 1) =cg_ae(iof+i, k, 1)
289 tauae(i, k, 2) = tau_ae(iof+i, k, 2)
290 pizae(i, k, 2) = piz_ae(iof+i, k, 2)
291 cgae(i, k, 2) =cg_ae(iof+i, k, 2)
292 ENDDO
293 ENDDO
294
295 CALL LW(PPMB, PDP, PPSOL, PDT0, PEMIS, PTL, PTAVE, PWV, POZON, PAER, &
296 PCLDLD, PCLDLU, PVIEW, zcool, zcool0, ztoplw, zsollw, ztoplw0, &
297 zsollw0, zsollwdown, ZFLUP, ZFLDN, ZFLUP0, ZFLDN0)
298 CALL SW(PSCT, zrmu0, zfract, PPMB, PDP, PPSOL, PALBD, PALBP, PTAVE, &
299 PWV, PQS, POZON, PAER, PCLDSW, PTAU, POMEGA, PCG, zheat, zheat0, &
300 zalbpla, ztopsw, zsolsw, ztopsw0, zsolsw0, ZFSUP, ZFSDN, ZFSUP0, &
301 ZFSDN0, tauae, pizae, cgae, PTAUA, POMEGAA, ztopswad, zsolswad, &
302 ztopswai, zsolswai, ok_ade, ok_aie)
303
304 DO i = 1, kdlon
305 radsol(iof+i) = zsolsw(i) + zsollw(i)
306 topsw(iof+i) = ztopsw(i)
307 toplw(iof+i) = ztoplw(i)
308 solsw(iof+i) = zsolsw(i)
309 sollw(iof+i) = zsollw(i)
310 sollwdown(iof+i) = zsollwdown(i)
311
312 DO k = 1, klev+1
313 lwdn0 ( iof+i, k) = ZFLDN0 ( i, k)
314 lwdn ( iof+i, k) = ZFLDN ( i, k)
315 lwup0 ( iof+i, k) = ZFLUP0 ( i, k)
316 lwup ( iof+i, k) = ZFLUP ( i, k)
317 ENDDO
318
319 topsw0(iof+i) = ztopsw0(i)
320 toplw0(iof+i) = ztoplw0(i)
321 solsw0(iof+i) = zsolsw0(i)
322 sollw0(iof+i) = zsollw0(i)
323 albpla(iof+i) = zalbpla(i)
324
325 DO k = 1, klev+1
326 swdn0 ( iof+i, k) = ZFSDN0 ( i, k)
327 swdn ( iof+i, k) = ZFSDN ( i, k)
328 swup0 ( iof+i, k) = ZFSUP0 ( i, k)
329 swup ( iof+i, k) = ZFSUP ( i, k)
330 ENDDO
331 ENDDO
332 ! transform the aerosol forcings, if they have to be calculated
333 IF (ok_ade) THEN
334 DO i = 1, kdlon
335 topswad(iof+i) = ztopswad(i)
336 solswad(iof+i) = zsolswad(i)
337 ENDDO
338 ELSE
339 DO i = 1, kdlon
340 topswad(iof+i) = 0.0
341 solswad(iof+i) = 0.0
342 ENDDO
343 ENDIF
344 IF (ok_aie) THEN
345 DO i = 1, kdlon
346 topswai(iof+i) = ztopswai(i)
347 solswai(iof+i) = zsolswai(i)
348 ENDDO
349 ELSE
350 DO i = 1, kdlon
351 topswai(iof+i) = 0.0
352 solswai(iof+i) = 0.0
353 ENDDO
354 ENDIF
355
356 DO k = 1, klev
357 DO i = 1, kdlon
358 ! scale factor to take into account the difference
359 ! between dry air and water vapour specific heat capacity
360 zznormcp = 1. + RVTMP2 * PWV(i, k)
361 heat(iof+i, k) = zheat(i, k) / zznormcp
362 cool(iof+i, k) = zcool(i, k)/zznormcp
363 heat0(iof+i, k) = zheat0(i, k)/zznormcp
364 cool0(iof+i, k) = zcool0(i, k)/zznormcp
365 ENDDO
366 ENDDO
367 end DO loop_iof
368
369 END SUBROUTINE radlwsw
370
371 end module radlwsw_m

  ViewVC Help
Powered by ViewVC 1.1.21