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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 145 - (hide annotations)
Tue Jun 16 15:23:29 2015 UTC (8 years, 11 months ago) by guez
File size: 13571 byte(s)
Renamed bibio to misc.

In procedure fxhyp, use the fact that xf is an odd function of xtild.

In procedure invert_zoom_x, replace linear search in xf by
bisection. Also, use result from previous loop iteration as initial
guess. Variable "it" cannot be equal to 2 * nmax after search.

Unused arguments: hm of cv3_feed; ph, qnk, tv,tvp of cv3_mixing; ppsol
of lw; rconst, temp of vdif_kcay; rconst, plev, temp, ustar, l_mix of
yamada.

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

  ViewVC Help
Powered by ViewVC 1.1.21