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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 52 - (hide annotations)
Fri Sep 23 12:28:01 2011 UTC (12 years, 7 months ago) by guez
Original Path: trunk/libf/phylmd/Radlwsw/radlwsw.f
File size: 14839 byte(s)
Split "conflx.f" into single-procedure files in directory "Conflx".

Split "cv_routines.f" into single-procedure files in directory
"CV_routines". Made module "cvparam" from included file
"cvparam.h". No included file other than "netcdf.inc" left in LMDZE.

1 guez 3 !
2     ! $Header: /home/cvsroot/LMDZ4/libf/phylmd/radlwsw.F,v 1.4 2005/06/06 13:16:33 fairhead Exp $
3     !
4     SUBROUTINE radlwsw(dist, rmu0, fract,
5     . paprs, pplay,tsol,albedo, alblw, t,q,wo,
6     . cldfra, cldemi, cldtaupd,
7     . heat,heat0,cool,cool0,radsol,albpla,
8     . topsw,toplw,solsw,sollw,
9     . sollwdown,
10     . topsw0,toplw0,solsw0,sollw0,
11     . lwdn0, lwdn, lwup0, lwup,
12     . swdn0, swdn, swup0, swup,
13     . ok_ade, ok_aie,
14     . tau_ae, piz_ae, cg_ae,
15     . topswad, solswad,
16     . cldtaupi, topswai, solswai)
17     c
18     use dimphy
19     use clesphys
20 guez 38 use SUPHEC_M
21 guez 3 use raddim, only: kflev, kdlon
22 guez 38 use yoethf_m
23 guez 3 IMPLICIT none
24     c======================================================================
25     c Auteur(s): Z.X. Li (LMD/CNRS) date: 19960719
26     c Objet: interface entre le modele et les rayonnements
27     c Arguments:
28     c dist-----input-R- distance astronomique terre-soleil
29     c rmu0-----input-R- cosinus de l'angle zenithal
30     c fract----input-R- duree d'ensoleillement normalisee
31     c co2_ppm--input-R- concentration du gaz carbonique (en ppm)
32     c solaire--input-R- constante solaire (W/m**2)
33     c paprs----input-R- pression a inter-couche (Pa)
34     c pplay----input-R- pression au milieu de couche (Pa)
35     c tsol-----input-R- temperature du sol (en K)
36     c albedo---input-R- albedo du sol (entre 0 et 1)
37     c t--------input-R- temperature (K)
38     c q--------input-R- vapeur d'eau (en kg/kg)
39     c wo-------input-R- contenu en ozone (en kg/kg) correction MPL 100505
40     c cldfra---input-R- fraction nuageuse (entre 0 et 1)
41     c cldtaupd---input-R- epaisseur optique des nuages dans le visible (present-day value)
42     c cldemi---input-R- emissivite des nuages dans l'IR (entre 0 et 1)
43     c ok_ade---input-L- apply the Aerosol Direct Effect or not?
44     c ok_aie---input-L- apply the Aerosol Indirect Effect or not?
45     c tau_ae, piz_ae, cg_ae-input-R- aerosol optical properties (calculated in aeropt.F)
46     c cldtaupi-input-R- epaisseur optique des nuages dans le visible
47     c calculated for pre-industrial (pi) aerosol concentrations, i.e. with smaller
48     c droplet concentration, thus larger droplets, thus generally cdltaupi cldtaupd
49     c it is needed for the diagnostics of the aerosol indirect radiative forcing
50     c
51     c heat-----output-R- echauffement atmospherique (visible) (K/jour)
52     c cool-----output-R- refroidissement dans l'IR (K/jour)
53     c radsol---output-R- bilan radiatif net au sol (W/m**2) (+ vers le bas)
54     c albpla---output-R- albedo planetaire (entre 0 et 1)
55     c topsw----output-R- flux solaire net au sommet de l'atm.
56     c toplw----output-R- ray. IR montant au sommet de l'atmosphere
57     c solsw----output-R- flux solaire net a la surface
58     c sollw----output-R- ray. IR montant a la surface
59     c solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
60     c topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
61     c solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
62     c topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
63     c
64     c ATTENTION: swai and swad have to be interpreted in the following manner:
65     c ---------
66     c ok_ade=F & ok_aie=F -both are zero
67     c ok_ade=T & ok_aie=F -aerosol direct forcing is F_{AD} = topsw-topswad
68     c indirect is zero
69     c ok_ade=F & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
70     c direct is zero
71     c ok_ade=T & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
72     c aerosol direct forcing is F_{AD} = topswai-topswad
73     c
74    
75     c======================================================================
76     c
77     real rmu0(klon), fract(klon), dist
78     cIM real co2_ppm
79     cIM real solaire
80     c
81     real, intent(in):: paprs(klon,klev+1)
82 guez 10 real, intent(in):: pplay(klon,klev)
83 guez 3 real albedo(klon), alblw(klon), tsol(klon)
84 guez 52 real, intent(in):: t(klon,klev)
85     real q(klon,klev)
86 guez 3 real, intent(in):: wo(klon,klev)
87     real cldfra(klon,klev), cldemi(klon,klev), cldtaupd(klon,klev)
88     real heat(klon,klev), cool(klon,klev)
89     real heat0(klon,klev), cool0(klon,klev)
90     real radsol(klon), topsw(klon), toplw(klon)
91     real solsw(klon), sollw(klon), albpla(klon)
92     real topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
93     real sollwdown(klon)
94     cIM output 3D
95     REAL*8 ZFSUP(KDLON,KFLEV+1)
96     REAL*8 ZFSDN(KDLON,KFLEV+1)
97     REAL*8 ZFSUP0(KDLON,KFLEV+1)
98     REAL*8 ZFSDN0(KDLON,KFLEV+1)
99     c
100     REAL*8 ZFLUP(KDLON,KFLEV+1)
101     REAL*8 ZFLDN(KDLON,KFLEV+1)
102     REAL*8 ZFLUP0(KDLON,KFLEV+1)
103     REAL*8 ZFLDN0(KDLON,KFLEV+1)
104     c
105     REAL*8 zx_alpha1, zx_alpha2
106     c
107     c
108     INTEGER k, kk, i, j, iof, nb_gr
109     EXTERNAL lw, sw
110     c
111     cIM ctes ds clesphys.h REAL*8 RCO2, RCH4, RN2O, RCFC11, RCFC12
112     REAL*8 PSCT
113     c
114     REAL*8 PALBD(kdlon,2), PALBP(kdlon,2)
115     REAL*8 PEMIS(kdlon), PDT0(kdlon), PVIEW(kdlon)
116     REAL*8 PPSOL(kdlon), PDP(kdlon,klev)
117     REAL*8 PTL(kdlon,kflev+1), PPMB(kdlon,kflev+1)
118     REAL*8 PTAVE(kdlon,kflev)
119     REAL*8 PWV(kdlon,kflev), PQS(kdlon,kflev), POZON(kdlon,kflev)
120     REAL*8 PAER(kdlon,kflev,5)
121     REAL*8 PCLDLD(kdlon,kflev)
122     REAL*8 PCLDLU(kdlon,kflev)
123     REAL*8 PCLDSW(kdlon,kflev)
124     REAL*8 PTAU(kdlon,2,kflev)
125     REAL*8 POMEGA(kdlon,2,kflev)
126     REAL*8 PCG(kdlon,2,kflev)
127     c
128     REAL*8 zfract(kdlon), zrmu0(kdlon), zdist
129     c
130     REAL*8 zheat(kdlon,kflev), zcool(kdlon,kflev)
131     REAL*8 zheat0(kdlon,kflev), zcool0(kdlon,kflev)
132     REAL*8 ztopsw(kdlon), ztoplw(kdlon)
133     REAL*8 zsolsw(kdlon), zsollw(kdlon), zalbpla(kdlon)
134     cIM
135     REAL*8 zsollwdown(kdlon)
136     c
137     REAL*8 ztopsw0(kdlon), ztoplw0(kdlon)
138     REAL*8 zsolsw0(kdlon), zsollw0(kdlon)
139     REAL*8 zznormcp
140     cIM output 3D : SWup, SWdn, LWup, LWdn
141     REAL swdn(klon,kflev+1),swdn0(klon,kflev+1)
142     REAL swup(klon,kflev+1),swup0(klon,kflev+1)
143     REAL lwdn(klon,kflev+1),lwdn0(klon,kflev+1)
144     REAL lwup(klon,kflev+1),lwup0(klon,kflev+1)
145     c-OB
146     cjq the following quantities are needed for the aerosol radiative forcings
147    
148     real topswad(klon), solswad(klon) ! output: aerosol direct forcing at TOA and surface
149     real topswai(klon), solswai(klon) ! output: aerosol indirect forcing atTOA and surface
150     real tau_ae(klon,klev,2), piz_ae(klon,klev,2), cg_ae(klon,klev,2) ! aerosol optical properties (see aeropt.F)
151     real cldtaupi(klon,klev) ! cloud optical thickness for pre-industrial aerosol concentrations
152     ! (i.e., with a smaller droplet concentrationand thus larger droplet radii)
153     logical ok_ade, ok_aie ! switches whether to use aerosol direct (indirect) effects or not
154     real*8 tauae(kdlon,kflev,2) ! aer opt properties
155     real*8 pizae(kdlon,kflev,2)
156     real*8 cgae(kdlon,kflev,2)
157     REAL*8 PTAUA(kdlon,2,kflev) ! present-day value of cloud opt thickness (PTAU is pre-industrial value), local use
158     REAL*8 POMEGAA(kdlon,2,kflev) ! dito for single scatt albedo
159     REAL*8 ztopswad(kdlon), zsolswad(kdlon) ! Aerosol direct forcing at TOAand surface
160     REAL*8 ztopswai(kdlon), zsolswai(kdlon) ! dito, indirect
161     cjq-end
162     !rv
163     tauae(:,:,:)=0.
164     pizae(:,:,:)=0.
165     cgae(:,:,:)=0.
166     !rv
167    
168     c
169     c-------------------------------------------
170     nb_gr = klon / kdlon
171     IF (nb_gr*kdlon .NE. klon) THEN
172     PRINT*, "kdlon mauvais:", klon, kdlon, nb_gr
173     stop 1
174     ENDIF
175     IF (kflev .NE. klev) THEN
176     PRINT*, "kflev differe de klev, kflev, klev"
177     stop 1
178     ENDIF
179     c-------------------------------------------
180     DO k = 1, klev
181     DO i = 1, klon
182     heat(i,k)=0.
183     cool(i,k)=0.
184     heat0(i,k)=0.
185     cool0(i,k)=0.
186     ENDDO
187     ENDDO
188     c
189     zdist = dist
190     c
191     cIM anciennes valeurs
192     c RCO2 = co2_ppm * 1.0e-06 * 44.011/28.97
193     c
194     cIM : on met RCO2, RCH4, RN2O, RCFC11 et RCFC12 dans clesphys.h /lecture ds conf_phys.F90
195     c RCH4 = 1.65E-06* 16.043/28.97
196     c RN2O = 306.E-09* 44.013/28.97
197     c RCFC11 = 280.E-12* 137.3686/28.97
198     c RCFC12 = 484.E-12* 120.9140/28.97
199     cIM anciennes valeurs
200     c RCH4 = 1.72E-06* 16.043/28.97
201     c RN2O = 310.E-09* 44.013/28.97
202     c
203     c PRINT*,'IMradlwsw : solaire, co2= ', solaire, co2_ppm
204     PSCT = solaire/zdist/zdist
205     c
206     DO 99999 j = 1, nb_gr
207     iof = kdlon*(j-1)
208     c
209     DO i = 1, kdlon
210     zfract(i) = fract(iof+i)
211     zrmu0(i) = rmu0(iof+i)
212     PALBD(i,1) = albedo(iof+i)
213     ! PALBD(i,2) = albedo(iof+i)
214     PALBD(i,2) = alblw(iof+i)
215     PALBP(i,1) = albedo(iof+i)
216     ! PALBP(i,2) = albedo(iof+i)
217     PALBP(i,2) = alblw(iof+i)
218     cIM cf. JLD pour etre en accord avec ORCHIDEE il faut mettre PEMIS(i) = 0.96
219     PEMIS(i) = 1.0
220     PVIEW(i) = 1.66
221     PPSOL(i) = paprs(iof+i,1)
222     zx_alpha1 = (paprs(iof+i,1)-pplay(iof+i,2))
223     . / (pplay(iof+i,1)-pplay(iof+i,2))
224     zx_alpha2 = 1.0 - zx_alpha1
225     PTL(i,1) = t(iof+i,1) * zx_alpha1 + t(iof+i,2) * zx_alpha2
226     PTL(i,klev+1) = t(iof+i,klev)
227     PDT0(i) = tsol(iof+i) - PTL(i,1)
228     ENDDO
229     DO k = 2, kflev
230     DO i = 1, kdlon
231     PTL(i,k) = (t(iof+i,k)+t(iof+i,k-1))*0.5
232     ENDDO
233     ENDDO
234     DO k = 1, kflev
235     DO i = 1, kdlon
236     PDP(i,k) = paprs(iof+i,k)-paprs(iof+i,k+1)
237     PTAVE(i,k) = t(iof+i,k)
238     PWV(i,k) = MAX (q(iof+i,k), 1.0e-12)
239     PQS(i,k) = PWV(i,k)
240     c wo: cm.atm (epaisseur en cm dans la situation standard)
241     c POZON: kg/kg
242     IF (bug_ozone) then
243     POZON(i,k) = MAX(wo(iof+i,k),1.0e-12)*RG/46.6968
244     . /(paprs(iof+i,k)-paprs(iof+i,k+1))
245     . *(paprs(iof+i,1)/101325.0)
246     ELSE
247     c le calcul qui suit est maintenant fait dans ozonecm (MPL)
248     POZON(i,k) = wo(i,k)
249     ENDIF
250     PCLDLD(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
251     PCLDLU(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
252     PCLDSW(i,k) = cldfra(iof+i,k)
253     PTAU(i,1,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! 1e-12 serait instable
254     PTAU(i,2,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! pour 32-bit machines
255     POMEGA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i,1,k))
256     POMEGA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i,2,k))
257     PCG(i,1,k) = 0.865
258     PCG(i,2,k) = 0.910
259     c-OB
260     cjq Introduced for aerosol indirect forcings.
261     cjq The following values use the cloud optical thickness calculated from
262     cjq present-day aerosol concentrations whereas the quantities without the
263     cjq "A" at the end are for pre-industial (natural-only) aerosol concentrations
264     cjq
265     PTAUA(i,1,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! 1e-12 serait instable
266     PTAUA(i,2,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! pour 32-bit machines
267     POMEGAA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i,1,k))
268     POMEGAA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i,2,k))
269     cjq-end
270     ENDDO
271     ENDDO
272     c
273     DO k = 1, kflev+1
274     DO i = 1, kdlon
275     PPMB(i,k) = paprs(iof+i,k)/100.0
276     ENDDO
277     ENDDO
278     c
279     DO kk = 1, 5
280     DO k = 1, kflev
281     DO i = 1, kdlon
282     PAER(i,k,kk) = 1.0E-15
283     ENDDO
284     ENDDO
285     ENDDO
286     c-OB
287     DO k = 1, kflev
288     DO i = 1, kdlon
289     tauae(i,k,1)=tau_ae(iof+i,k,1)
290     pizae(i,k,1)=piz_ae(iof+i,k,1)
291     cgae(i,k,1) =cg_ae(iof+i,k,1)
292     tauae(i,k,2)=tau_ae(iof+i,k,2)
293     pizae(i,k,2)=piz_ae(iof+i,k,2)
294     cgae(i,k,2) =cg_ae(iof+i,k,2)
295     ENDDO
296     ENDDO
297     c
298     c======================================================================
299     cIM ctes ds clesphys.h CALL LW(RCO2,RCH4,RN2O,RCFC11,RCFC12,
300     CALL LW(
301     . PPMB, PDP,
302     . PPSOL,PDT0,PEMIS,
303     . PTL, PTAVE, PWV, POZON, PAER,
304     . PCLDLD,PCLDLU,
305     . PVIEW,
306     . zcool, zcool0,
307     . ztoplw,zsollw,ztoplw0,zsollw0,
308     . zsollwdown,
309     . ZFLUP, ZFLDN, ZFLUP0,ZFLDN0)
310     cIM ctes ds clesphys.h CALL SW(PSCT, RCO2, zrmu0, zfract,
311     CALL SW(PSCT, zrmu0, zfract,
312     S PPMB, PDP,
313     S PPSOL, PALBD, PALBP,
314     S PTAVE, PWV, PQS, POZON, PAER,
315     S PCLDSW, PTAU, POMEGA, PCG,
316     S zheat, zheat0,
317     S zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,
318     S ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,
319     S tauae, pizae, cgae, ! aerosol optical properties
320     s PTAUA, POMEGAA,
321     s ztopswad,zsolswad,ztopswai,zsolswai, ! diagnosed aerosol forcing
322     J ok_ade, ok_aie) ! apply aerosol effects or not?
323    
324     c======================================================================
325     DO i = 1, kdlon
326     radsol(iof+i) = zsolsw(i) + zsollw(i)
327     topsw(iof+i) = ztopsw(i)
328     toplw(iof+i) = ztoplw(i)
329     solsw(iof+i) = zsolsw(i)
330     sollw(iof+i) = zsollw(i)
331     sollwdown(iof+i) = zsollwdown(i)
332     cIM
333     DO k = 1, kflev+1
334     lwdn0 ( iof+i,k) = ZFLDN0 ( i,k)
335     lwdn ( iof+i,k) = ZFLDN ( i,k)
336     lwup0 ( iof+i,k) = ZFLUP0 ( i,k)
337     lwup ( iof+i,k) = ZFLUP ( i,k)
338     ENDDO
339     c
340     topsw0(iof+i) = ztopsw0(i)
341     toplw0(iof+i) = ztoplw0(i)
342     solsw0(iof+i) = zsolsw0(i)
343     sollw0(iof+i) = zsollw0(i)
344     albpla(iof+i) = zalbpla(i)
345     cIM
346     DO k = 1, kflev+1
347     swdn0 ( iof+i,k) = ZFSDN0 ( i,k)
348     swdn ( iof+i,k) = ZFSDN ( i,k)
349     swup0 ( iof+i,k) = ZFSUP0 ( i,k)
350     swup ( iof+i,k) = ZFSUP ( i,k)
351     ENDDO !k=1, kflev+1
352     ENDDO
353     cjq-transform the aerosol forcings, if they have
354     cjq to be calculated
355     IF (ok_ade) THEN
356     DO i = 1, kdlon
357     topswad(iof+i) = ztopswad(i)
358     solswad(iof+i) = zsolswad(i)
359     ENDDO
360     ELSE
361     DO i = 1, kdlon
362     topswad(iof+i) = 0.0
363     solswad(iof+i) = 0.0
364     ENDDO
365     ENDIF
366     IF (ok_aie) THEN
367     DO i = 1, kdlon
368     topswai(iof+i) = ztopswai(i)
369     solswai(iof+i) = zsolswai(i)
370     ENDDO
371     ELSE
372     DO i = 1, kdlon
373     topswai(iof+i) = 0.0
374     solswai(iof+i) = 0.0
375     ENDDO
376     ENDIF
377     cjq-end
378     DO k = 1, kflev
379     c DO i = 1, kdlon
380     c heat(iof+i,k) = zheat(i,k)
381     c cool(iof+i,k) = zcool(i,k)
382     c heat0(iof+i,k) = zheat0(i,k)
383     c cool0(iof+i,k) = zcool0(i,k)
384     c ENDDO
385     DO i = 1, kdlon
386     C scale factor to take into account the difference between
387     C dry air and watter vapour scpecific heat capacity
388     zznormcp=1.0+RVTMP2*PWV(i,k)
389     heat(iof+i,k) = zheat(i,k)/zznormcp
390     cool(iof+i,k) = zcool(i,k)/zznormcp
391     heat0(iof+i,k) = zheat0(i,k)/zznormcp
392     cool0(iof+i,k) = zcool0(i,k)/zznormcp
393     ENDDO
394     ENDDO
395     c
396     99999 CONTINUE
397     RETURN
398     END

  ViewVC Help
Powered by ViewVC 1.1.21