Context Navigation

radlwsw_aero.F90 @ 1159

Last change on this file since 1159 was 1159, checked in by jghattas, 16 years ago

Recuperation de la variable solaire avec include clesphys.h au lieu de liste des arguments
Recuperation de qq modifications qui etaient fait sur radlwsw.F mais pas repercute sur la nouvelle version radlwsw_aero.F90 :
- Ajoute iflag_rrtm et 3 arguments qui sont necessaire pour iflag_rrtm=1 + Cosmetique

Marie-Pierre + Josefine

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 16.4 KB

Line
1	SUBROUTINE radlwsw_aero( &
2	dist, rmu0, fract, &
3	paprs, pplay,tsol,alb1, alb2, &
4	t,q,wo,&
5	cldfra, cldemi, cldtaupd,&
6	ok_ade, ok_aie,&
7	tau_aero, piz_aero, cg_aero,&
8	cldtaupi, new_aod, &
9	qsat, flwc, fiwc, &
10	heat,heat0,cool,cool0,radsol,albpla,&
11	topsw,toplw,solsw,sollw,&
12	sollwdown,&
13	topsw0,toplw0,solsw0,sollw0,&
14	lwdn0, lwdn, lwup0, lwup,&
15	swdn0, swdn, swup0, swup,&
16	topswad_aero, solswad_aero,&
17	topswai_aero, solswai_aero, &
18	topswad0_aero, solswad0_aero,&
19	topsw_aero, topsw0_aero,&
20	solsw_aero, solsw0_aero)
21
22
23
24	USE DIMPHY
25
26	IMPLICIT NONE
27	!======================================================================
28	! Auteur(s): Z.X. Li (LMD/CNRS) date: 19960719
29	! Objet: interface entre le modele et les rayonnements
30	! Arguments:
31	! dist-----input-R- distance astronomique terre-soleil
32	! rmu0-----input-R- cosinus de l'angle zenithal
33	! fract----input-R- duree d'ensoleillement normalisee
34	! co2_ppm--input-R- concentration du gaz carbonique (en ppm)
35	! paprs----input-R- pression a inter-couche (Pa)
36	! pplay----input-R- pression au milieu de couche (Pa)
37	! tsol-----input-R- temperature du sol (en K)
38	! alb1-----input-R- albedo du sol(entre 0 et 1) dans l'interval visible
39	! alb2-----input-R- albedo du sol(entre 0 et 1) dans l'interval proche infra-rouge
40	! t--------input-R- temperature (K)
41	! q--------input-R- vapeur d'eau (en kg/kg)
42	! wo-------input-R- contenu en ozone (en kg/kg) correction MPL 100505
43	! cldfra---input-R- fraction nuageuse (entre 0 et 1)
44	! cldtaupd---input-R- epaisseur optique des nuages dans le visible (present-day value)
45	! cldemi---input-R- emissivite des nuages dans l'IR (entre 0 et 1)
46	! ok_ade---input-L- apply the Aerosol Direct Effect or not?
47	! ok_aie---input-L- apply the Aerosol Indirect Effect or not?
48	! tau_ae, piz_ae, cg_ae-input-R- aerosol optical properties (calculated in aeropt.F)
49	! cldtaupi-input-R- epaisseur optique des nuages dans le visible
50	! calculated for pre-industrial (pi) aerosol concentrations, i.e. with smaller
51	! droplet concentration, thus larger droplets, thus generally cdltaupi cldtaupd
52	! it is needed for the diagnostics of the aerosol indirect radiative forcing
53	!
54	! heat-----output-R- echauffement atmospherique (visible) (K/jour)
55	! cool-----output-R- refroidissement dans l'IR (K/jour)
56	! radsol---output-R- bilan radiatif net au sol (W/m**2) (+ vers le bas)
57	! albpla---output-R- albedo planetaire (entre 0 et 1)
58	! topsw----output-R- flux solaire net au sommet de l'atm.
59	! toplw----output-R- ray. IR montant au sommet de l'atmosphere
60	! solsw----output-R- flux solaire net a la surface
61	! sollw----output-R- ray. IR montant a la surface
62	! solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
63	! topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
64	! solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
65	! topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
66	!
67	! ATTENTION: swai and swad have to be interpreted in the following manner:
68	! ---------
69	! ok_ade=F & ok_aie=F -both are zero
70	! ok_ade=T & ok_aie=F -aerosol direct forcing is F_{AD} = topsw-topswad
71	! indirect is zero
72	! ok_ade=F & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
73	! direct is zero
74	! ok_ade=T & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
75	! aerosol direct forcing is F_{AD} = topswai-topswad
76	!
77
78	!======================================================================
79
80	! ====================================================================
81	! Adapte au modele de chimie INCA par Celine Deandreis & Anne Cozic -- 2009
82	! 1 = ZERO
83	! 2 = AER total
84	! 3 = NAT
85	! 4 = BC
86	! 5 = SO4
87	! 6 = POM
88	! 7 = DUST
89	! 8 = SS
90	! 9 = NO3
91	!
92	! ====================================================================
93	include "YOETHF.h"
94	include "YOMCST.h"
95	include "clesphys.h"
96	include "iniprint.h"
97
98	! Input arguments
99	REAL, INTENT(in) :: dist
100	REAL, INTENT(in) :: rmu0(KLON), fract(KLON)
101	REAL, INTENT(in) :: paprs(KLON,KLEV+1), pplay(KLON,KLEV)
102	REAL, INTENT(in) :: alb1(KLON), alb2(KLON), tsol(KLON)
103	REAL, INTENT(in) :: t(KLON,KLEV), q(KLON,KLEV), wo(KLON,KLEV)
104	LOGICAL, INTENT(in) :: ok_ade, ok_aie ! switches whether to use aerosol direct (indirect) effects or not
105	REAL, INTENT(in) :: cldfra(KLON,KLEV), cldemi(KLON,KLEV), cldtaupd(KLON,KLEV)
106	REAL, INTENT(in) :: tau_aero(KLON,KLEV,9,2) ! aerosol optical properties (see aeropt.F)
107	REAL, INTENT(in) :: piz_aero(KLON,KLEV,9,2) ! aerosol optical properties (see aeropt.F)
108	REAL, INTENT(in) :: cg_aero(KLON,KLEV,9,2) ! aerosol optical properties (see aeropt.F)
109	REAL, INTENT(in) :: cldtaupi(KLON,KLEV) ! cloud optical thickness for pre-industrial aerosol concentrations
110	LOGICAL, INTENT(in) :: new_aod ! flag pour retrouver les resultats exacts de l'AR4 dans le cas ou l'on ne travaille qu'avec les sulfates
111	REAL, INTENT(in) :: qsat(klon,klev) ! Variable pour iflag_rrtm=1
112	REAL, INTENT(in) :: flwc(klon,klev) ! Variable pour iflag_rrtm=1
113	REAL, INTENT(in) :: fiwc(klon,klev) ! Variable pour iflag_rrtm=1
114
115	! Output arguments
116	REAL, INTENT(out) :: heat(KLON,KLEV), cool(KLON,KLEV)
117	REAL, INTENT(out) :: heat0(KLON,KLEV), cool0(KLON,KLEV)
118	REAL, INTENT(out) :: radsol(KLON), topsw(KLON), toplw(KLON)
119	REAL, INTENT(out) :: solsw(KLON), sollw(KLON), albpla(KLON)
120	REAL, INTENT(out) :: topsw0(KLON), toplw0(KLON), solsw0(KLON), sollw0(KLON)
121	REAL, INTENT(out) :: sollwdown(KLON)
122	REAL, INTENT(out) :: swdn(KLON,kflev+1),swdn0(KLON,kflev+1)
123	REAL, INTENT(out) :: swup(KLON,kflev+1),swup0(KLON,kflev+1)
124	REAL, INTENT(out) :: lwdn(KLON,kflev+1),lwdn0(KLON,kflev+1)
125	REAL, INTENT(out) :: lwup(KLON,kflev+1),lwup0(KLON,kflev+1)
126	REAL, INTENT(out) :: topswad_aero(KLON), solswad_aero(KLON) ! output: aerosol direct forcing at TOA and surface
127	REAL, INTENT(out) :: topswai_aero(KLON), solswai_aero(KLON) ! output: aerosol indirect forcing atTOA and surface
128	REAL, DIMENSION(klon), INTENT(out) :: topswad0_aero
129	REAL, DIMENSION(klon), INTENT(out) :: solswad0_aero
130	REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw_aero
131	REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw0_aero
132	REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw_aero
133	REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw0_aero
134
135	! Local variables
136	REAL*8 ZFSUP(KDLON,KFLEV+1)
137	REAL*8 ZFSDN(KDLON,KFLEV+1)
138	REAL*8 ZFSUP0(KDLON,KFLEV+1)
139	REAL*8 ZFSDN0(KDLON,KFLEV+1)
140	REAL*8 ZFLUP(KDLON,KFLEV+1)
141	REAL*8 ZFLDN(KDLON,KFLEV+1)
142	REAL*8 ZFLUP0(KDLON,KFLEV+1)
143	REAL*8 ZFLDN0(KDLON,KFLEV+1)
144	REAL*8 zx_alpha1, zx_alpha2
145	INTEGER k, kk, i, j, iof, nb_gr
146	REAL*8 PSCT
147	REAL*8 PALBD(kdlon,2), PALBP(kdlon,2)
148	REAL*8 PEMIS(kdlon), PDT0(kdlon), PVIEW(kdlon)
149	REAL*8 PPSOL(kdlon), PDP(kdlon,KLEV)
150	REAL*8 PTL(kdlon,kflev+1), PPMB(kdlon,kflev+1)
151	REAL*8 PTAVE(kdlon,kflev)
152	REAL*8 PWV(kdlon,kflev), PQS(kdlon,kflev), POZON(kdlon,kflev)
153	REAL*8 PAER(kdlon,kflev,5)
154	REAL*8 PCLDLD(kdlon,kflev)
155	REAL*8 PCLDLU(kdlon,kflev)
156	REAL*8 PCLDSW(kdlon,kflev)
157	REAL*8 PTAU(kdlon,2,kflev)
158	REAL*8 POMEGA(kdlon,2,kflev)
159	REAL*8 PCG(kdlon,2,kflev)
160	REAL*8 zfract(kdlon), zrmu0(kdlon), zdist
161	REAL*8 zheat(kdlon,kflev), zcool(kdlon,kflev)
162	REAL*8 zheat0(kdlon,kflev), zcool0(kdlon,kflev)
163	REAL*8 ztopsw(kdlon), ztoplw(kdlon)
164	REAL*8 zsolsw(kdlon), zsollw(kdlon), zalbpla(kdlon)
165	REAL*8 zsollwdown(kdlon)
166	REAL*8 ztopsw0(kdlon), ztoplw0(kdlon)
167	REAL*8 zsolsw0(kdlon), zsollw0(kdlon)
168	REAL*8 zznormcp
169	REAL*8 tauaero(kdlon,kflev,9,2) ! aer opt properties
170	REAL*8 pizaero(kdlon,kflev,9,2)
171	REAL*8 cgaero(kdlon,kflev,9,2)
172	REAL*8 PTAUA(kdlon,2,kflev) ! present-day value of cloud opt thickness (PTAU is pre-industrial value), local use
173	REAL*8 POMEGAA(kdlon,2,kflev) ! dito for single scatt albedo
174	REAL*8 ztopswadaero(kdlon), zsolswadaero(kdlon) ! Aerosol direct forcing at TOAand surface
175	REAL*8 ztopswad0aero(kdlon), zsolswad0aero(kdlon) ! Aerosol direct forcing at TOAand surface
176	REAL*8 ztopswaiaero(kdlon), zsolswaiaero(kdlon) ! dito, indirect
177	REAL*8 ztopsw_aero(kdlon,9), ztopsw0_aero(kdlon,9)
178	REAL*8 zsolsw_aero(kdlon,9), zsolsw0_aero(kdlon,9)
179
180	! initialisation
181	tauaero(:,:,:,:)=0.
182	pizaero(:,:,:,:)=0.
183	cgaero(:,:,:,:)=0.
184
185	!
186	!-------------------------------------------
187	nb_gr = KLON / kdlon
188	IF (nb_gr*kdlon .NE. KLON) THEN
189	PRINT*, "kdlon mauvais:", KLON, kdlon, nb_gr
190	CALL abort
191	ENDIF
192	IF (kflev .NE. KLEV) THEN
193	PRINT*, "kflev differe de KLEV, kflev, KLEV"
194	CALL abort
195	ENDIF
196	!-------------------------------------------
197	DO k = 1, KLEV
198	DO i = 1, KLON
199	heat(i,k)=0.
200	cool(i,k)=0.
201	heat0(i,k)=0.
202	cool0(i,k)=0.
203	ENDDO
204	ENDDO
205	!
206	zdist = dist
207	!
208	PSCT = solaire/zdist/zdist
209	DO j = 1, nb_gr
210	iof = kdlon*(j-1)
211	DO i = 1, kdlon
212	zfract(i) = fract(iof+i)
213	zrmu0(i) = rmu0(iof+i)
214	PALBD(i,1) = alb1(iof+i)
215	PALBD(i,2) = alb2(iof+i)
216	PALBP(i,1) = alb1(iof+i)
217	PALBP(i,2) = alb2(iof+i)
218	PEMIS(i) = 1.0
219	PVIEW(i) = 1.66
220	PPSOL(i) = paprs(iof+i,1)
221	zx_alpha1 = (paprs(iof+i,1)-pplay(iof+i,2))/(pplay(iof+i,1)-pplay(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, kflev
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, kflev
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	! wo: cm.atm (epaisseur en cm dans la situation standard)
239	! POZON: kg/kg
240	POZON(i,k) = MAX(wo(iof+i,k),1.0e-12)*RG/46.6968 &
241	/(paprs(iof+i,k)-paprs(iof+i,k+1))&
242	*(paprs(iof+i,1)/101325.0)
243	PCLDLD(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
244	PCLDLU(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
245	PCLDSW(i,k) = cldfra(iof+i,k)
246	PTAU(i,1,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! 1e-12 serait instable
247	PTAU(i,2,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! pour 32-bit machines
248	POMEGA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i,1,k))
249	POMEGA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i,2,k))
250	PCG(i,1,k) = 0.865
251	PCG(i,2,k) = 0.910
252	!-
253	! Introduced for aerosol indirect forcings.
254	! The following values use the cloud optical thickness calculated from
255	! present-day aerosol concentrations whereas the quantities without the
256	! "A" at the end are for pre-industial (natural-only) aerosol concentrations
257	!
258	PTAUA(i,1,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! 1e-12 serait instable
259	PTAUA(i,2,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! pour 32-bit machines
260	POMEGAA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i,1,k))
261	POMEGAA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i,2,k))
262	ENDDO
263	ENDDO
264	!
265	DO k = 1, kflev+1
266	DO i = 1, kdlon
267	PPMB(i,k) = paprs(iof+i,k)/100.0
268	ENDDO
269	ENDDO
270	!
271	DO kk = 1, 5
272	DO k = 1, kflev
273	DO i = 1, kdlon
274	PAER(i,k,kk) = 1.0E-15
275	ENDDO
276	ENDDO
277	ENDDO
278	DO k = 1, kflev
279	DO i = 1, kdlon
280	tauaero(i,k,:,1)=tau_aero(iof+i,k,:,1)
281	pizaero(i,k,:,1)=piz_aero(iof+i,k,:,1)
282	cgaero(i,k,:,1) =cg_aero(iof+i,k,:,1)
283	tauaero(i,k,:,2)=tau_aero(iof+i,k,:,2)
284	pizaero(i,k,:,2)=piz_aero(iof+i,k,:,2)
285	cgaero(i,k,:,2) =cg_aero(iof+i,k,:,2)
286	ENDDO
287	ENDDO
288
289	!
290	!===== iflag_rrtm ================================================
291	!
292	IF (iflag_rrtm == 0) THEN
293	! Old radiation scheme, used for AR4 runs
294	CALL LW_LMDAR4(&
295	PPMB, PDP,&
296	PPSOL,PDT0,PEMIS,&
297	PTL, PTAVE, PWV, POZON, PAER,&
298	PCLDLD,PCLDLU,&
299	PVIEW,&
300	zcool, zcool0,&
301	ztoplw,zsollw,ztoplw0,zsollw0,&
302	zsollwdown,&
303	ZFLUP, ZFLDN, ZFLUP0,ZFLDN0)
304
305
306	IF (.NOT. new_aod) THEN
307	! use old version
308	CALL SW_LMDAR4(PSCT, zrmu0, zfract,&
309	PPMB, PDP, &
310	PPSOL, PALBD, PALBP,&
311	PTAVE, PWV, PQS, POZON, PAER,&
312	PCLDSW, PTAU, POMEGA, PCG,&
313	zheat, zheat0,&
314	zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,&
315	ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,&
316	tau_aero(:,:,5,:), piz_aero(:,:,5,:), cg_aero(:,:,5,:),&
317	PTAUA, POMEGAA,&
318	ztopswadaero,zsolswadaero,&
319	ztopswaiaero,zsolswaiaero,&
320	ok_ade, ok_aie)
321
322	ELSE ! new_aod=T
323	CALL SW_AEROAR4(PSCT, zrmu0, zfract,&
324	PPMB, PDP,&
325	PPSOL, PALBD, PALBP,&
326	PTAVE, PWV, PQS, POZON, PAER,&
327	PCLDSW, PTAU, POMEGA, PCG,&
328	zheat, zheat0,&
329	zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,&
330	ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,&
331	tauaero, pizaero, cgaero, &
332	PTAUA, POMEGAA,&
333	ztopswadaero,zsolswadaero,&
334	ztopswad0aero,zsolswad0aero,&
335	ztopswaiaero,zsolswaiaero, &
336	ztopsw_aero,ztopsw0_aero,&
337	zsolsw_aero,zsolsw0_aero,&
338	ok_ade, ok_aie)
339
340	ENDIF
341
342	ELSE
343	!===== iflag_rrtm=1, on passe dans SW via RECMWFL ===============
344	WRITE(lunout,*) "Option iflag_rrtm=T ne fonctionne pas encore !!!"
345	CALL abort_gcm('radlwsw_aero','iflag_rrtm=T not valid',1)
346
347	ENDIF ! iflag_rrtm
348	!======================================================================
349
350	DO i = 1, kdlon
351	radsol(iof+i) = zsolsw(i) + zsollw(i)
352	topsw(iof+i) = ztopsw(i)
353	toplw(iof+i) = ztoplw(i)
354	solsw(iof+i) = zsolsw(i)
355	sollw(iof+i) = zsollw(i)
356	sollwdown(iof+i) = zsollwdown(i)
357	DO k = 1, kflev+1
358	lwdn0 ( iof+i,k) = ZFLDN0 ( i,k)
359	lwdn ( iof+i,k) = ZFLDN ( i,k)
360	lwup0 ( iof+i,k) = ZFLUP0 ( i,k)
361	lwup ( iof+i,k) = ZFLUP ( i,k)
362	ENDDO
363	topsw0(iof+i) = ztopsw0(i)
364	toplw0(iof+i) = ztoplw0(i)
365	solsw0(iof+i) = zsolsw0(i)
366	sollw0(iof+i) = zsollw0(i)
367	albpla(iof+i) = zalbpla(i)
368
369	DO k = 1, kflev+1
370	swdn0 ( iof+i,k) = ZFSDN0 ( i,k)
371	swdn ( iof+i,k) = ZFSDN ( i,k)
372	swup0 ( iof+i,k) = ZFSUP0 ( i,k)
373	swup ( iof+i,k) = ZFSUP ( i,k)
374	ENDDO
375	ENDDO
376	!-transform the aerosol forcings, if they have
377	! to be calculated
378	IF (ok_ade) THEN
379	DO i = 1, kdlon
380	topswad_aero(iof+i) = ztopswadaero(i)
381	topswad0_aero(iof+i) = ztopswad0aero(i)
382	solswad_aero(iof+i) = zsolswadaero(i)
383	solswad0_aero(iof+i) = zsolswad0aero(i)
384	topsw_aero(iof+i,:) = ztopsw_aero(iof+i,:)
385	topsw0_aero(iof+i,:) = ztopsw0_aero(iof+i,:)
386	solsw_aero(iof+i,:) = zsolsw_aero(iof+i,:)
387	solsw0_aero(iof+i,:) = zsolsw0_aero(iof+i,:)
388
389	ENDDO
390	ELSE
391	DO i = 1, kdlon
392	topswad_aero(iof+i) = 0.0
393	solswad_aero(iof+i) = 0.0
394	topswad0_aero(iof+i) = 0.0
395	solswad0_aero(iof+i) = 0.0
396	topsw_aero(iof+i,:) = 0.
397	topsw0_aero(iof+i,:) =0.
398	solsw_aero(iof+i,:) = 0.
399	solsw0_aero(iof+i,:) = 0.
400	ENDDO
401	ENDIF
402	IF (ok_aie) THEN
403	DO i = 1, kdlon
404	topswai_aero(iof+i) = ztopswaiaero(i)
405	solswai_aero(iof+i) = zsolswaiaero(i)
406	ENDDO
407	ELSE
408	DO i = 1, kdlon
409	topswai_aero(iof+i) = 0.0
410	solswai_aero(iof+i) = 0.0
411	ENDDO
412	ENDIF
413	DO k = 1, kflev
414	DO i = 1, kdlon
415	! scale factor to take into account the difference between
416	! dry air and watter vapour scpecifi! heat capacity
417	zznormcp=1.0+RVTMP2*PWV(i,k)
418	heat(iof+i,k) = zheat(i,k)/zznormcp
419	cool(iof+i,k) = zcool(i,k)/zznormcp
420	heat0(iof+i,k) = zheat0(i,k)/zznormcp
421	cool0(iof+i,k) = zcool0(i,k)/zznormcp
422	ENDDO
423	ENDDO
424
425	ENDDO ! j = 1, nb_gr
426
427
428	END SUBROUTINE radlwsw_aero
429
430

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format