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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 178 - (show annotations)
Fri Mar 11 18:47:26 2016 UTC (8 years, 2 months ago) by guez
File size: 13538 byte(s)
Moved variables date0, deltat, datasz_max, ncvar_ids, point, buff_pos,
buffer, regular from module histcom_var to modules where they are
defined.

Removed procedure ioipslmpp, useless for a sequential program.

Added argument datasz_max to histwrite_real (to avoid circular
dependency with histwrite).

Removed useless variables and computations everywhere.

Changed real litteral constants from default kind to double precision
in lwb, lwu, lwvn, sw1s, swtt, swtt1, swu.

Removed unused arguments: paer of sw, sw1s, sw2s, swclr; pcldsw of
sw1s, sw2s; pdsig, prayl of swr; co2_ppm of clmain, clqh; tsol of
transp_lay; nsrf of screenp; kcrit and kknu of gwstress; pstd of
orosetup.

Added output of relative humidity.

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

  ViewVC Help
Powered by ViewVC 1.1.21