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

Contents of /trunk/Sources/phylmd/Radlwsw/sw2s.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: 15074 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 sw2s_m
2
3 IMPLICIT NONE
4
5 contains
6
7 SUBROUTINE sw2s(knu, flag_aer, tauae, pizae, cgae, paki, palbd, palbp, &
8 pcg, pcld, pclear, pdsig, pomega, poz, prmu, psec, ptau, pud, &
9 pwv, pqs, pfdown, pfup)
10 USE dimens_m
11 USE dimphy
12 USE raddim
13 USE radepsi
14 use swclr_m, only: swclr
15 use swr_m, only: swr
16
17 ! ------------------------------------------------------------------
18 ! PURPOSE.
19 ! --------
20
21 ! THIS ROUTINE COMPUTES THE SHORTWAVE RADIATION FLUXES IN THE
22 ! SECOND SPECTRAL INTERVAL FOLLOWING FOUQUART AND BONNEL (1980).
23
24 ! METHOD.
25 ! -------
26
27 ! 1. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING TO
28 ! CONTINUUM SCATTERING
29 ! 2. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING FOR
30 ! A GREY MOLECULAR ABSORPTION
31 ! 3. LAPLACE TRANSFORM ON THE PREVIOUS TO GET EFFECTIVE AMOUNTS
32 ! OF ABSORBERS
33 ! 4. APPLY H2O AND U.M.G. TRANSMISSION FUNCTIONS
34 ! 5. MULTIPLY BY OZONE TRANSMISSION FUNCTION
35
36 ! REFERENCE.
37 ! ----------
38
39 ! SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
40 ! DOCUMENTATION, AND FOUQUART AND BONNEL (1980)
41
42 ! AUTHOR.
43 ! -------
44 ! JEAN-JACQUES MORCRETTE *ECMWF*
45
46 ! MODIFICATIONS.
47 ! --------------
48 ! ORIGINAL : 89-07-14
49 ! 94-11-15 J.-J. MORCRETTE DIRECT/DIFFUSE ALBEDO
50 ! ------------------------------------------------------------------
51 ! * ARGUMENTS:
52
53 INTEGER knu
54 ! -OB
55 DOUBLE PRECISION flag_aer
56 DOUBLE PRECISION tauae(kdlon, kflev, 2)
57 DOUBLE PRECISION pizae(kdlon, kflev, 2)
58 DOUBLE PRECISION cgae(kdlon, kflev, 2)
59 DOUBLE PRECISION paki(kdlon, 2)
60 DOUBLE PRECISION palbd(kdlon, 2)
61 DOUBLE PRECISION palbp(kdlon, 2)
62 DOUBLE PRECISION pcg(kdlon, 2, kflev)
63 DOUBLE PRECISION pcld(kdlon, kflev)
64 DOUBLE PRECISION pclear(kdlon)
65 DOUBLE PRECISION pdsig(kdlon, kflev)
66 DOUBLE PRECISION pomega(kdlon, 2, kflev)
67 DOUBLE PRECISION poz(kdlon, kflev)
68 DOUBLE PRECISION pqs(kdlon, kflev)
69 DOUBLE PRECISION prmu(kdlon)
70 DOUBLE PRECISION psec(kdlon)
71 DOUBLE PRECISION ptau(kdlon, 2, kflev)
72 DOUBLE PRECISION pud(kdlon, 5, kflev+1)
73 DOUBLE PRECISION pwv(kdlon, kflev)
74
75 DOUBLE PRECISION pfdown(kdlon, kflev+1)
76 DOUBLE PRECISION pfup(kdlon, kflev+1)
77
78 ! * LOCAL VARIABLES:
79
80 INTEGER iind2(2), iind3(3)
81 DOUBLE PRECISION zcgaz(kdlon, kflev)
82 DOUBLE PRECISION zfd(kdlon, kflev+1)
83 DOUBLE PRECISION zfu(kdlon, kflev+1)
84 DOUBLE PRECISION zg(kdlon)
85 DOUBLE PRECISION zgg(kdlon)
86 DOUBLE PRECISION zpizaz(kdlon, kflev)
87 DOUBLE PRECISION zrayl(kdlon)
88 DOUBLE PRECISION zray1(kdlon, kflev+1)
89 DOUBLE PRECISION zray2(kdlon, kflev+1)
90 DOUBLE PRECISION zref(kdlon)
91 DOUBLE PRECISION zrefz(kdlon, 2, kflev+1)
92 DOUBLE PRECISION zre1(kdlon)
93 DOUBLE PRECISION zre2(kdlon)
94 DOUBLE PRECISION zrj(kdlon, 6, kflev+1)
95 DOUBLE PRECISION zrj0(kdlon, 6, kflev+1)
96 DOUBLE PRECISION zrk(kdlon, 6, kflev+1)
97 DOUBLE PRECISION zrk0(kdlon, 6, kflev+1)
98 DOUBLE PRECISION zrl(kdlon, 8)
99 DOUBLE PRECISION zrmue(kdlon, kflev+1)
100 DOUBLE PRECISION zrmu0(kdlon, kflev+1)
101 DOUBLE PRECISION zrmuz(kdlon)
102 DOUBLE PRECISION zrneb(kdlon)
103 DOUBLE PRECISION zr1(kdlon)
104 DOUBLE PRECISION zr2(kdlon, 2)
105 DOUBLE PRECISION zr3(kdlon, 3)
106 DOUBLE PRECISION zr4(kdlon)
107 DOUBLE PRECISION zr21(kdlon)
108 DOUBLE PRECISION zr22(kdlon)
109 DOUBLE PRECISION zs(kdlon)
110 DOUBLE PRECISION ztauaz(kdlon, kflev)
111 DOUBLE PRECISION zto1(kdlon)
112 DOUBLE PRECISION ztr(kdlon, 2, kflev+1)
113 DOUBLE PRECISION ztra1(kdlon, kflev+1)
114 DOUBLE PRECISION ztra2(kdlon, kflev+1)
115 DOUBLE PRECISION ztr1(kdlon)
116 DOUBLE PRECISION ztr2(kdlon)
117 DOUBLE PRECISION zw(kdlon)
118 DOUBLE PRECISION zw1(kdlon)
119 DOUBLE PRECISION zw2(kdlon, 2)
120 DOUBLE PRECISION zw3(kdlon, 3)
121 DOUBLE PRECISION zw4(kdlon)
122 DOUBLE PRECISION zw5(kdlon)
123
124 INTEGER jl, jk, k, jaj, ikm1, ikl, jn, jabs, jkm1
125 INTEGER jref, jkl, jklp1, jajp, jkki, jkkp4, jn2j, iabs
126 DOUBLE PRECISION zrmum1, zwh2o, zcneb, zaa, zbb, zrki, zre11
127
128 ! * Prescribed Data:
129
130 DOUBLE PRECISION rsun(2)
131 SAVE rsun
132 DOUBLE PRECISION rray(2, 6)
133 SAVE rray
134 DATA rsun(1)/0.441676d0/
135 DATA rsun(2)/0.558324d0/
136 DATA (rray(1,k), k=1, 6)/.428937d-01, .890743d+00, -.288555d+01, &
137 .522744d+01, -.469173d+01, .161645d+01/
138 DATA (rray(2,k), k=1, 6)/.697200d-02, .173297d-01, -.850903d-01, &
139 .248261d+00, -.302031d+00, .129662d+00/
140
141 ! ------------------------------------------------------------------
142
143 ! * 1. SECOND SPECTRAL INTERVAL (0.68-4.00 MICRON)
144 ! -------------------------------------------
145
146
147
148 ! * 1.1 OPTICAL THICKNESS FOR RAYLEIGH SCATTERING
149 ! -----------------------------------------
150
151
152 DO jl = 1, kdlon
153 zrmum1 = 1. - prmu(jl)
154 zrayl(jl) = rray(knu, 1) + zrmum1*(rray(knu,2)+zrmum1*(rray(knu, &
155 3)+zrmum1*(rray(knu,4)+zrmum1*(rray(knu,5)+zrmum1*rray(knu,6)))))
156 END DO
157
158
159 ! ------------------------------------------------------------------
160
161 ! * 2. CONTINUUM SCATTERING CALCULATIONS
162 ! ---------------------------------
163
164
165 ! * 2.1 CLEAR-SKY FRACTION OF THE COLUMN
166 ! --------------------------------
167
168
169 CALL swclr(knu, flag_aer, tauae, pizae, cgae, palbp, pdsig, zrayl, &
170 psec, zcgaz, zpizaz, zray1, zray2, zrefz, zrj0, zrk0, zrmu0, ztauaz, &
171 ztra1, ztra2)
172
173
174 ! * 2.2 CLOUDY FRACTION OF THE COLUMN
175 ! -----------------------------
176
177
178 CALL swr(knu, palbd, pcg, pcld, pomega, psec, ptau, zcgaz, &
179 zpizaz, zray1, zray2, zrefz, zrj, zrk, zrmue, ztauaz, ztra1, ztra2)
180
181
182 ! ------------------------------------------------------------------
183
184 ! * 3. SCATTERING CALCULATIONS WITH GREY MOLECULAR ABSORPTION
185 ! ------------------------------------------------------
186
187
188 jn = 2
189
190 DO jabs = 1, 2
191
192
193 ! * 3.1 SURFACE CONDITIONS
194 ! ------------------
195
196
197 DO jl = 1, kdlon
198 zrefz(jl, 2, 1) = palbd(jl, knu)
199 zrefz(jl, 1, 1) = palbd(jl, knu)
200 END DO
201
202
203 ! * 3.2 INTRODUCING CLOUD EFFECTS
204 ! -------------------------
205
206
207 DO jk = 2, kflev + 1
208 jkm1 = jk - 1
209 ikl = kflev + 1 - jkm1
210 DO jl = 1, kdlon
211 zrneb(jl) = pcld(jl, jkm1)
212 IF (jabs==1 .AND. zrneb(jl)>2.*zeelog) THEN
213 zwh2o = max(pwv(jl,jkm1), zeelog)
214 zcneb = max(zeelog, min(zrneb(jl),1.-zeelog))
215 zbb = pud(jl, jabs, jkm1)*pqs(jl, jkm1)/zwh2o
216 zaa = max((pud(jl,jabs,jkm1)-zcneb*zbb)/(1.-zcneb), zeelog)
217 ELSE
218 zaa = pud(jl, jabs, jkm1)
219 zbb = zaa
220 END IF
221 zrki = paki(jl, jabs)
222 zs(jl) = exp(-zrki*zaa*1.66)
223 zg(jl) = exp(-zrki*zaa/zrmue(jl,jk))
224 ztr1(jl) = 0.
225 zre1(jl) = 0.
226 ztr2(jl) = 0.
227 zre2(jl) = 0.
228
229 zw(jl) = pomega(jl, knu, jkm1)
230 zto1(jl) = ptau(jl, knu, jkm1)/zw(jl) + ztauaz(jl, jkm1)/zpizaz(jl, &
231 jkm1) + zbb*zrki
232
233 zr21(jl) = ptau(jl, knu, jkm1) + ztauaz(jl, jkm1)
234 zr22(jl) = ptau(jl, knu, jkm1)/zr21(jl)
235 zgg(jl) = zr22(jl)*pcg(jl, knu, jkm1) + (1.-zr22(jl))*zcgaz(jl, jkm1)
236 zw(jl) = zr21(jl)/zto1(jl)
237 zref(jl) = zrefz(jl, 1, jkm1)
238 zrmuz(jl) = zrmue(jl, jk)
239 END DO
240
241 CALL swde(zgg, zref, zrmuz, zto1, zw, zre1, zre2, ztr1, ztr2)
242
243 DO jl = 1, kdlon
244
245 zrefz(jl, 2, jk) = (1.-zrneb(jl))*(zray1(jl,jkm1)+zrefz(jl,2,jkm1)* &
246 ztra1(jl,jkm1)*ztra2(jl,jkm1))*zg(jl)*zs(jl) + zrneb(jl)*zre1(jl)
247
248 ztr(jl, 2, jkm1) = zrneb(jl)*ztr1(jl) + (ztra1(jl,jkm1))*zg(jl)*(1.- &
249 zrneb(jl))
250
251 zrefz(jl, 1, jk) = (1.-zrneb(jl))*(zray1(jl,jkm1)+zrefz(jl,1,jkm1)* &
252 ztra1(jl,jkm1)*ztra2(jl,jkm1)/(1.-zray2(jl,jkm1)*zrefz(jl,1, &
253 jkm1)))*zg(jl)*zs(jl) + zrneb(jl)*zre2(jl)
254
255 ztr(jl, 1, jkm1) = zrneb(jl)*ztr2(jl) + (ztra1(jl,jkm1)/(1.-zray2(jl, &
256 jkm1)*zrefz(jl,1,jkm1)))*zg(jl)*(1.-zrneb(jl))
257
258 END DO
259 END DO
260
261 ! * 3.3 REFLECT./TRANSMISSIVITY BETWEEN SURFACE AND LEVEL
262 ! -------------------------------------------------
263
264
265 DO jref = 1, 2
266
267 jn = jn + 1
268
269 DO jl = 1, kdlon
270 zrj(jl, jn, kflev+1) = 1.
271 zrk(jl, jn, kflev+1) = zrefz(jl, jref, kflev+1)
272 END DO
273
274 DO jk = 1, kflev
275 jkl = kflev + 1 - jk
276 jklp1 = jkl + 1
277 DO jl = 1, kdlon
278 zre11 = zrj(jl, jn, jklp1)*ztr(jl, jref, jkl)
279 zrj(jl, jn, jkl) = zre11
280 zrk(jl, jn, jkl) = zre11*zrefz(jl, jref, jkl)
281 END DO
282 END DO
283 END DO
284 END DO
285
286
287 ! ------------------------------------------------------------------
288
289 ! * 4. INVERT GREY AND CONTINUUM FLUXES
290 ! --------------------------------
291
292
293
294 ! * 4.1 UPWARD (ZRK) AND DOWNWARD (ZRJ) PSEUDO-FLUXES
295 ! ---------------------------------------------
296
297
298 DO jk = 1, kflev + 1
299 DO jaj = 1, 5, 2
300 jajp = jaj + 1
301 DO jl = 1, kdlon
302 zrj(jl, jaj, jk) = zrj(jl, jaj, jk) - zrj(jl, jajp, jk)
303 zrk(jl, jaj, jk) = zrk(jl, jaj, jk) - zrk(jl, jajp, jk)
304 zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
305 zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
306 END DO
307 END DO
308 END DO
309
310 DO jk = 1, kflev + 1
311 DO jaj = 2, 6, 2
312 DO jl = 1, kdlon
313 zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
314 zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
315 END DO
316 END DO
317 END DO
318
319 ! * 4.2 EFFECTIVE ABSORBER AMOUNTS BY INVERSE LAPLACE
320 ! ---------------------------------------------
321
322
323 DO jk = 1, kflev + 1
324 jkki = 1
325 DO jaj = 1, 2
326 iind2(1) = jaj
327 iind2(2) = jaj
328 DO jn = 1, 2
329 jn2j = jn + 2*jaj
330 jkkp4 = jkki + 4
331
332 ! * 4.2.1 EFFECTIVE ABSORBER AMOUNTS
333 ! --------------------------
334
335
336 DO jl = 1, kdlon
337 zw2(jl, 1) = log(zrj(jl,jn,jk)/zrj(jl,jn2j,jk))/paki(jl, jaj)
338 zw2(jl, 2) = log(zrk(jl,jn,jk)/zrk(jl,jn2j,jk))/paki(jl, jaj)
339 END DO
340
341 ! * 4.2.2 TRANSMISSION FUNCTION
342 ! ---------------------
343
344
345 CALL swtt1(knu, 2, iind2, zw2, zr2)
346
347 DO jl = 1, kdlon
348 zrl(jl, jkki) = zr2(jl, 1)
349 zrl(jl, jkkp4) = zr2(jl, 2)
350 END DO
351
352 jkki = jkki + 1
353 END DO
354 END DO
355
356 ! * 4.3 UPWARD AND DOWNWARD FLUXES WITH H2O AND UMG ABSORPTION
357 ! ------------------------------------------------------
358
359
360 DO jl = 1, kdlon
361 pfdown(jl, jk) = zrj(jl, 1, jk)*zrl(jl, 1)*zrl(jl, 3) + &
362 zrj(jl, 2, jk)*zrl(jl, 2)*zrl(jl, 4)
363 pfup(jl, jk) = zrk(jl, 1, jk)*zrl(jl, 5)*zrl(jl, 7) + &
364 zrk(jl, 2, jk)*zrl(jl, 6)*zrl(jl, 8)
365 END DO
366 END DO
367
368
369 ! ------------------------------------------------------------------
370
371 ! * 5. MOLECULAR ABSORPTION ON CLEAR-SKY FLUXES
372 ! ----------------------------------------
373
374
375
376 ! * 5.1 DOWNWARD FLUXES
377 ! ---------------
378
379
380 jaj = 2
381 iind3(1) = 1
382 iind3(2) = 2
383 iind3(3) = 3
384
385 DO jl = 1, kdlon
386 zw3(jl, 1) = 0.
387 zw3(jl, 2) = 0.
388 zw3(jl, 3) = 0.
389 zw4(jl) = 0.
390 zw5(jl) = 0.
391 zr4(jl) = 1.
392 zfd(jl, kflev+1) = zrj0(jl, jaj, kflev+1)
393 END DO
394 DO jk = 1, kflev
395 ikl = kflev + 1 - jk
396 DO jl = 1, kdlon
397 zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikl)/zrmu0(jl, ikl)
398 zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikl)/zrmu0(jl, ikl)
399 zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikl)/zrmu0(jl, ikl)
400 zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmu0(jl, ikl)
401 zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmu0(jl, ikl)
402 END DO
403
404 CALL swtt1(knu, 3, iind3, zw3, zr3)
405
406 DO jl = 1, kdlon
407 ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
408 zfd(jl, ikl) = zr3(jl, 1)*zr3(jl, 2)*zr3(jl, 3)*zr4(jl)* &
409 zrj0(jl, jaj, ikl)
410 END DO
411 END DO
412
413
414 ! * 5.2 UPWARD FLUXES
415 ! -------------
416
417
418 DO jl = 1, kdlon
419 zfu(jl, 1) = zfd(jl, 1)*palbp(jl, knu)
420 END DO
421
422 DO jk = 2, kflev + 1
423 ikm1 = jk - 1
424 DO jl = 1, kdlon
425 zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikm1)*1.66
426 zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikm1)*1.66
427 zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikm1)*1.66
428 zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
429 zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
430 END DO
431
432 CALL swtt1(knu, 3, iind3, zw3, zr3)
433
434 DO jl = 1, kdlon
435 ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
436 zfu(jl, jk) = zr3(jl, 1)*zr3(jl, 2)*zr3(jl, 3)*zr4(jl)* &
437 zrk0(jl, jaj, jk)
438 END DO
439 END DO
440
441
442 ! ------------------------------------------------------------------
443
444 ! * 6. INTRODUCTION OF OZONE AND H2O CONTINUUM ABSORPTION
445 ! --------------------------------------------------
446
447 iabs = 3
448
449 ! * 6.1 DOWNWARD FLUXES
450 ! ---------------
451
452 DO jl = 1, kdlon
453 zw1(jl) = 0.
454 zw4(jl) = 0.
455 zw5(jl) = 0.
456 zr1(jl) = 0.
457 pfdown(jl, kflev+1) = ((1.-pclear(jl))*pfdown(jl,kflev+1)+pclear(jl)*zfd( &
458 jl,kflev+1))*rsun(knu)
459 END DO
460
461 DO jk = 1, kflev
462 ikl = kflev + 1 - jk
463 DO jl = 1, kdlon
464 zw1(jl) = zw1(jl) + poz(jl, ikl)/zrmue(jl, ikl)
465 zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmue(jl, ikl)
466 zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmue(jl, ikl)
467 ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
468 END DO
469
470 CALL swtt(knu, iabs, zw1, zr1)
471
472 DO jl = 1, kdlon
473 pfdown(jl, ikl) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfdown(jl,ikl)+ &
474 pclear(jl)*zfd(jl,ikl))*rsun(knu)
475 END DO
476 END DO
477
478
479 ! * 6.2 UPWARD FLUXES
480 ! -------------
481
482 DO jl = 1, kdlon
483 pfup(jl, 1) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfup(jl,1)+pclear(jl)*zfu( &
484 jl,1))*rsun(knu)
485 END DO
486
487 DO jk = 2, kflev + 1
488 ikm1 = jk - 1
489 DO jl = 1, kdlon
490 zw1(jl) = zw1(jl) + poz(jl, ikm1)*1.66
491 zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
492 zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
493 ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
494 END DO
495
496 CALL swtt(knu, iabs, zw1, zr1)
497
498 DO jl = 1, kdlon
499 pfup(jl, jk) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfup(jl,jk)+pclear(jl)* &
500 zfu(jl,jk))*rsun(knu)
501 END DO
502 END DO
503
504 END SUBROUTINE sw2s
505
506 end module sw2s_m

  ViewVC Help
Powered by ViewVC 1.1.21