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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 72 - (hide annotations)
Tue Jul 23 13:00:07 2013 UTC (10 years, 9 months ago) by guez
Original Path: trunk/libf/phylmd/Radlwsw/radlwsw.f90
File size: 13917 byte(s)
NaN to signalling NaN in gfortran_debug.mk.

Removed unused procedures in getincom and getincom2. In procedure
conf_interface, replaced call to getincom by new namelist. Moved
procedure conf_interface into module interface_surf.

Added variables sig1 and w01 to startphy.nc and restartphy.nc, for
procedure cv_driver. Renamed (ema_)?work1 and (ema_)?work2 to sig1 and
w01 in concvl and physiq.

Deleted unused arguments of clmain, clqh and intersurf_hq, among which
(y)?sollwdown. Following LMDZ, in physiq, read sollw instead of
sollwdown from startphy.nc, write sollw instead of sollwdown to
restartphy.nc.

In procedure sw, initialized zfs[ud][pn]a[di], for runs where ok_ade
and ok_aie are false. (Following LMDZ.)

Added dimension klev to startphy.nc and restartphy.nc, and deleted
dimension horizon_vertical. Made t_ancien and q_ancien two-dimensional
NetCDF variables. Bug fix: in phyetat0, define ratqs, clwcon and
rnebcon for vertical levels >=2.

Bug fix: set mfg, p[de]n_[ud] to 0. when iflag_con >= 3. (Following LMDZ.)

1 guez 53 module radlwsw_m
2 guez 3
3 guez 53 IMPLICIT none
4 guez 3
5 guez 53 contains
6    
7 guez 72 SUBROUTINE radlwsw(dist, rmu0, 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 71 USE clesphys, ONLY: bug_ozone, 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 72 real dist, rmu0(klon), fract(klon)
49 guez 53 ! 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 guez 69
53     real, intent(in):: paprs(klon, klev+1)
54 guez 53 ! paprs----input-R- pression a inter-couche (Pa)
55 guez 72 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 guez 69 ! albedo---input-R- albedo du sol (entre 0 et 1)
59 guez 53 ! tsol-----input-R- temperature du sol (en K)
60 guez 69 real, intent(in):: t(klon, klev)
61 guez 53 ! t--------input-R- temperature (K)
62 guez 69 real q(klon, klev)
63 guez 53 ! q--------input-R- vapeur d'eau (en kg/kg)
64 guez 69 real, intent(in):: wo(klon, klev)
65 guez 53 ! wo-------input-R- contenu en ozone (en kg/kg) correction MPL 100505
66 guez 69 real cldfra(klon, klev), cldemi(klon, klev)
67 guez 53 ! 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 guez 69 real cldtaupd(klon, klev)
71     ! input-R- epaisseur optique des nuages dans le visible (present-day value)
72 guez 53
73     real, intent(out):: heat(klon, klev)
74     ! échauffement atmosphérique (visible) (K/jour)
75    
76 guez 72 real heat0(klon, klev)
77 guez 53 real cool(klon, klev)
78 guez 69 ! cool-----output-R- refroidissement dans l'IR (K/jour)
79 guez 72 real cool0(klon, klev)
80     real radsol(klon)
81 guez 69 ! radsol---output-R- bilan radiatif net au sol (W/m**2) (+ vers le bas)
82 guez 72 real albpla(klon)
83     ! albpla---output-R- albedo planetaire (entre 0 et 1)
84     real topsw(klon)
85 guez 69 ! topsw----output-R- flux solaire net au sommet de l'atm.
86 guez 62
87     real, intent(out):: toplw(klon)
88     ! rayonnement infrarouge montant au sommet de l'atmosphère
89    
90 guez 72 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 guez 62 real, intent(out):: toplw0(klon)
98 guez 72 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 guez 53 !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 guez 62 INTEGER k, kk, i, iof, nb_gr
151 guez 53 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 guez 62 loop_iof: DO iof = 0, klon - kdlon, kdlon
212 guez 53 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 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     ! 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 guez 62 end DO loop_iof
377 guez 53
378     END SUBROUTINE radlwsw
379    
380     end module radlwsw_m

  ViewVC Help
Powered by ViewVC 1.1.21