Context Navigation

callradite_mod.F @ 3468

Last change on this file since 3468 was 3468, checked in by emillour, 4 weeks ago
Mars PCM: Remove obsolete/depreciated lwrite flag (which would trigger some very specific extra text outputs), in code and in reference callphys.def files. EM
File size: 24.0 KB

Line
1	MODULE callradite_mod
2
3	IMPLICIT NONE
4
5	CONTAINS
6
7	SUBROUTINE callradite(icount,ngrid,nlayer,nq,zday,ls,pq,albedo,
8	$ emis,mu0,pplev,pplay,pt,tsurf,fract,dist_sol,igout,
9	$ dtlw,dtsw,fluxsurf_lw,fluxsurf_dn_sw,fluxsurf_up_sw,
10	$ fluxtop_lw,fluxtop_dn_sw,fluxtop_up_sw,
11	& tau_pref_scenario,tau_pref_gcm,
12	& tau,aerosol,dsodust,tauscaling,dust_rad_adjust,IRtoVIScoef,
13	$ taucloudtes,rdust,rice,nuice,riceco2,nuiceco2,co2ice,
14	$ rstormdust,rtopdust,totstormfract,clearatm,dsords,dsotop,
15	$ nohmons,clearsky,totcloudfrac)
16
17	use aeropacity_mod, only: aeropacity
18	use updatereffrad_mod, only: updatereffrad
19	use dimradmars_mod, only: ndomainsz, nflev, nsun, nir
20	use dimradmars_mod, only: naerkind, name_iaer,
21	& iaer_dust_conrath,iaer_dust_doubleq,
22	& iaer_dust_submicron, iaer_h2o_ice, iaer_co2_ice,
23	& iaer_stormdust_doubleq,iaer_topdust_doubleq,flux_1AU
24	use yomlw_h, only: gcp, nlaylte
25	use comcstfi_h, only: g,cpp
26	use time_phylmdz_mod, only: daysec
27	use lwmain_mod, only: lwmain
28	use swmain_mod, only: swmain
29	use dust_param_mod, only: doubleq, active, submicron
30	IMPLICIT NONE
31	c=======================================================================
32	c subject:
33	c --------
34	c Subroutine designed to call the main canonic
35	c radiative transfer subroutine "lwmain" et "swmain"
36	c to compute radiative heating and cooling rate and
37	c radiative fluxes to the surface.
38	c
39	c These calculations are only valid on the part of the atmosphere
40	c where Local Thermal Equilibrium (LTE) is verified. In practice
41	c the calculations are only performed for the first "nlaylte"
42	c parameters (nlaylte is calculated by subroutine "nlthermeq"
43	c and stored in module "yomlw_h").
44	c
45	c The purpose of this subroutine is to:
46	c 1) Make some initial calculation at first call
47	c 2) Split the calculation in several sub-grid
48	c ("sub-domain") to save memory and
49	c be able run on a workstation at high resolution
50	c The sub-grid size is defined in dimradmars_mod
51	c 3) Compute the 3D scattering parameters depending on the
52	c size distribution of the different tracers (added by JBM)
53	c 4) call "lwmain" and "swmain"
54	c
55	c
56	c authors:
57	c ------
58	c Francois Forget / Christophe Hourdin / J.-B. Madeleine (2009)
59	c
60	c
61	c 3D scattering scheme user's guide (J.-B. Madeleine)
62	c ---------------------------------
63	c
64	c This routine has been modified to take into account 3D, time
65	c dependent scattering properties of the aerosols.
66	c---- The look-up tables that contain the scattering parameters
67	c of a given tracer, for different sizes, are read by SUAER.F90.
68	c The names of the corresponding ASCII files have to be set in
69	c this subroutine (file_id variable), and files must be in the
70	c directory specified in datafile_mod. Please make sure that the
71	c ASCII files are correctly written, and that the range
72	c of particle sizes is consistent with what you would expect.
73	c---- SUAER.F90 is in charge of reading the ASCII files and averaging
74	c the scattering parameters in each GCM channel, using the three last
75	c equations of Forget et al. 1998 (GRL 25, No.7, p.1105-1108).
76	c---- These look-up tables, loaded during the firstcall, are then
77	c constantly used by the subroutine "aeroptproperties.F" to compute,
78	c online, the 3D scattering parameters, based on the size distribution
79	c (reffrad and nueffrad) of the different tracers, in each grid box.
80	c These 3D size distributions are loaded by the "updatereffrad.F"
81	c subroutine. A log-normal distribution is then assumed in
82	c "aeroptproperties.F", along with a Gauss-Legendre integration.
83	c---- The optical depth at the visible reference wavelength (set in
84	c SUAER.F90, after the file_id variable) is then computed by
85	c the subroutine "aeropacity.F", by using the size and spatial
86	c distribution of the corresponding tracer. This connection has to
87	c be implemented in "aeropacity.F" when adding a new tracer. To do so,
88	c one can use equation 2 of Forget et al. 1998 (Icarus 131, p.302-316).
89	c---- The resulting variables "aerosol", "QVISsQREF3d", "omegaVIS3d" and
90	c "gVIS3d" (same in the infrared) are finally used by lwmain.F and
91	c swmain.F to solve the radiative transfer equation.
92	c
93	c changes:
94	c -------
95	c
96	c > SRL 7/2000
97	c
98	c This version has been modified to only calculate radiative tendencies
99	c over layers 1..NFLEV (set in dimradmars_mod). Returns zero for higher
100	c layers, if any.
101	c In other routines, nlayer -> nflev.
102	c Routines affected: lwflux, lwi, lwmain, lwxb, lwxd, lwxn.
103	c
104	c > J.-B. Madeleine 10W12
105	c This version uses the variable's splitting, which can be usefull
106	c when performing very high resolution simulation like LES.
107	c
108	c ----------
109	c Here, solar band#1 is spectral interval between "long1vis" and "long2vis"
110	c set in dimradmars_mod
111	c Here, solar band#2 is spectral interval between "long2vis" and "long3vis"
112	c set in dimradmars_mod
113	c
114	c input:
115	c -----
116	c icount counter of call to subroutine physic by gcm
117	c ngrid number of gridpoint of horizontal grid
118	c nlayer Number of layer
119	c nq Number of tracer
120	c ls Solar longitude (Ls) , radian
121	c zday Date (time since Ls=0, in martian days)
122	c pq(ngrid,nlayer,nq) Advected fields
123	c
124	c albedo (ngrid,2) hemispheric surface albedo
125	c albedo (i,1) : mean albedo for solar band#1
126	c (see below)
127	c albedo (i,2) : mean albedo for solar band#2
128	c (see below)
129	c emis Thermal IR surface emissivity (no unit)
130	c mu0(ngrid) cos of solar zenith angle
131	c (=1 when sun at zenith)
132	c pplay(ngrid,nlayer) pressure (Pa) in the middle of each layer
133	c pplev(ngrid,nlayer+1) pressure (Pa) at boundaries of each layer
134	c pt(ngrid,nlayer) atmospheric temperature in each layer (K)
135	c tsurf(ngrid) surface temperature (K)
136	c fract(ngrid) day fraction of the time interval
137	c =1 during the full day ; =0 during the night
138	c declin latitude of subsolar point
139	c dist_sol sun-Mars distance (AU)
140	c igout coordinate of analysed point for debugging
141	c reffrad(ngrid,nlayer,naerkind) Aerosol effective radius
142	c nueffrad(ngrid,nlayer,naerkind) Aerosol effective variance
143
144	c=======================================================================
145	c
146	c Declarations :
147	c -------------
148	c
149	include "callkeys.h"
150
151	c-----------------------------------------------------------------------
152	c Input/Output
153	c ------------
154	INTEGER,INTENT(IN) :: icount
155	INTEGER,INTENT(IN) :: ngrid,nlayer,nq
156	INTEGER,INTENT(IN) :: igout
157
158	REAL,INTENT(IN) :: pq(ngrid,nlayer,nq)
159	REAL,INTENT(INOUT) :: tauscaling(ngrid) ! Conversion factor for
160	! qdust and Ndust
161	REAL,INTENT(INOUT) :: dust_rad_adjust(ngrid) ! Radiative adjustment
162	! factor for dust
163	REAL,INTENT(INOUT) :: IRtoVIScoef(ngrid) ! conversion coefficient to apply on
164	! scenario absorption IR (9.3um) CDOD
165	! = tau_pref_gcm_VIS / tau_pref_gcm_IR
166	REAL,INTENT(IN) :: albedo(ngrid,2),emis(ngrid)
167	REAL,INTENT(IN) :: ls,zday
168
169	REAL,INTENT(IN) :: pplev(ngrid,nlayer+1),pplay(ngrid,nlayer)
170	REAL,INTENT(IN) :: pt(ngrid,nlayer)
171	REAL,INTENT(IN) :: tsurf(ngrid)
172	REAL,INTENT(IN) :: dist_sol,mu0(ngrid),fract(ngrid)
173	REAL,INTENT(OUT) :: dtlw(ngrid,nlayer) ! longwave (IR) heating rate (K/s)
174	REAL,INTENT(OUT) :: dtsw(ngrid,nlayer) ! shortwave (Solar) heating rate (K/s)
175	REAL,INTENT(OUT) :: fluxsurf_lw(ngrid) ! total LW (thermal IR) downward flux
176	! (W.m-2) at the surface
177	REAL,INTENT(OUT) :: fluxtop_lw(ngrid) ! outgoing total LW (thermal IR)
178	! upward flux (W.m-2) at the top of the atm.
179	REAL,INTENT(OUT) :: fluxsurf_dn_sw(ngrid,2) ! surface downward SW flux for
180	! solar bands #1 and #2 (W.m-2)
181	REAL,INTENT(OUT) :: fluxsurf_up_sw(ngrid,2) ! surface upward SW flux for
182	! solar bands #1 and #2 (W.m-2)
183	REAL,INTENT(OUT) :: fluxtop_dn_sw(ngrid,2) ! incoming downward SW flux for
184	! solar bands #1 and #2 (W.m-2) at top of atm.
185	REAL,INTENT(OUT) :: fluxtop_up_sw(ngrid,2) ! outgoing upward SW flux for
186	! solar bands #1 and #2 (W.m-2) at top of atm.
187	REAL,INTENT(OUT) :: tau_pref_scenario(ngrid) ! prescribed dust column
188	! visible opacity at odpref from scenario
189	REAL,INTENT(OUT) :: tau_pref_gcm(ngrid) ! computed dust column
190	! visible opacity at odpref in the GCM
191	REAL,INTENT(OUT) :: tau(ngrid,naerkind) ! Column visible optical depth
192	! for each aerosol
193	REAL,INTENT(OUT) :: taucloudtes(ngrid)! Cloud opacity at infrared
194	! reference wavelength using
195	! Qabs instead of Qext
196	! (direct comparison with TES)
197	REAL,INTENT(OUT) :: aerosol(ngrid,nlayer,naerkind) ! aerosol extinction
198	! optical depth at reference wavelength "longrefvis",
199	! set in dimradmars_h, for each kind of aerosol
200	REAL,INTENT(INOUT) :: dsodust(ngrid,nlayer)
201	REAL,INTENT(OUT) :: rdust(ngrid,nlayer) ! Dust geometric mean radius (m)
202	REAL,INTENT(OUT) :: rice(ngrid,nlayer) ! Ice geometric mean radius (m)
203	REAL,INTENT(OUT) :: nuice(ngrid,nlayer) ! Estimated effective variance
204	double precision,INTENT(OUT) :: riceco2(ngrid,nlayer) ! CO2 ice mean radius(m)
205	REAL,INTENT(OUT) :: nuiceco2(ngrid,nlayer) ! Effective variance
206	REAL,INTENT(IN) :: co2ice(ngrid) ! co2 ice surface layer (kg.m-2)
207
208	c rocket dust storm
209	LOGICAL,INTENT(IN) :: clearatm ! true for background dust
210	REAL,INTENT(IN) :: totstormfract(ngrid) ! dust storm mesh fraction
211	REAL,INTENT(OUT) :: rstormdust(ngrid,nlayer) ! Storm dust geometric mean radius (m)
212	REAL,INTENT(OUT) :: dsords(ngrid,nlayer) ! density scaled opacity for rocket dust storm dust
213
214	c entrainment by mountain top dust flows
215	LOGICAL, INTENT(IN) :: nohmons ! true for background dust
216	REAL,INTENT(OUT) :: rtopdust(ngrid,nlayer) ! Topdust geometric mean radius (m)
217	REAL,INTENT(OUT) :: dsotop(ngrid,nlayer) ! density scaled opacity for topmons dust
218
219	c sub-grid scale water ice clouds
220	LOGICAL,INTENT(IN) :: clearsky
221	REAL,INTENT(IN) :: totcloudfrac(ngrid)
222
223	c
224	c Local variables :
225	c -----------------
226
227	INTEGER j,l,ig,n,ich
228	INTEGER aer_count,iaer
229	INTEGER jd,ig0,nd
230
231	real cste_mars ! solar constant on Mars (Wm-2)
232	REAL ptlev(ngrid,nlayer+1)
233
234	INTEGER :: ndomain
235
236	c Thermal IR net radiative budget (W m-2)
237
238	real znetrad(ndomainsz,nflev)
239
240	real zfluxd_sw(ndomainsz,nflev+1,2)
241	real zfluxu_sw(ndomainsz,nflev+1,2)
242
243	REAL zplev(ndomainsz,nflev+1)
244	REAL zztlev(ndomainsz,nflev+1)
245	REAL zplay(ndomainsz,nflev)
246	REAL zt(ndomainsz,nflev)
247	REAL zaerosol(ndomainsz,nflev,naerkind)
248	REAL zalbedo(ndomainsz,2)
249	REAL zdp(ndomainsz,nflev)
250	REAL zdt0(ndomainsz)
251
252	REAL zzdtlw(ndomainsz,nflev)
253	REAL zzdtsw(ndomainsz,nflev)
254	REAL zzflux(ndomainsz,6)
255	real zrmuz
256
257	REAL :: zQVISsQREF3d(ndomainsz,nflev,nsun,naerkind)
258	REAL :: zomegaVIS3d(ndomainsz,nflev,nsun,naerkind)
259	REAL :: zgVIS3d(ndomainsz,nflev,nsun,naerkind)
260
261	REAL :: zQIRsQREF3d(ndomainsz,nflev,nir,naerkind)
262	REAL :: zomegaIR3d(ndomainsz,nflev,nir,naerkind)
263	REAL :: zgIR3d(ndomainsz,nflev,nir,naerkind)
264
265	c Aerosol size distribution
266	REAL :: reffrad(ngrid,nlayer,naerkind)
267	REAL :: nueffrad(ngrid,nlayer,naerkind)
268	c Aerosol optical properties
269	REAL :: QVISsQREF3d(ngrid,nlayer,nsun,naerkind)
270	REAL :: omegaVIS3d(ngrid,nlayer,nsun,naerkind)
271	REAL :: gVIS3d(ngrid,nlayer,nsun,naerkind)
272
273	REAL :: QIRsQREF3d(ngrid,nlayer,nir,naerkind)
274	REAL :: omegaIR3d(ngrid,nlayer,nir,naerkind)
275	REAL :: gIR3d(ngrid,nlayer,nir,naerkind)
276
277	REAL :: QREFvis3d(ngrid,nlayer,naerkind)
278	! QREFvis3d : Extinction efficiency at the VISible reference wavelength
279	REAL :: QREFir3d(ngrid,nlayer,naerkind)
280	! QREFir3d : Extinction efficiency at the InfraRed reference wavelength
281
282	REAL :: omegaREFvis3d(ngrid,nlayer,naerkind)
283	REAL :: omegaREFir3d(ngrid,nlayer,naerkind)
284
285	c local saved variables
286	c ---------------------
287
288	real zco2 ! volume fraction of CO2 in Mars atmosphere
289	!$OMP THREADPRIVATE(zco2)
290	DATA zco2/0.95/
291	SAVE zco2
292
293	LOGICAL firstcall
294	!$OMP THREADPRIVATE(firstcall)
295	DATA firstcall/.true./
296	SAVE firstcall
297
298
299
300	c----------------------------------------------------------------------
301
302	c Initialisation
303	c --------------
304
305	! compute ndomain
306	! AS: moved out of firstcall to allow nesting+evoluting domain
307	! ------------------------------------------------------------
308	ndomain= (ngrid-1) / ndomainsz + 1
309
310	IF (firstcall) THEN
311
312	write(,) 'Splitting radiative calculations: ',
313	$ ' ngrid,ndomainsz,ndomain',
314	$ ngrid,ndomainsz,ndomain
315
316	c Assign a number to the different scatterers
317	c -------------------------------------------
318
319	iaer_dust_conrath=0
320	iaer_dust_doubleq=0
321	iaer_dust_submicron=0
322	iaer_h2o_ice=0
323	iaer_co2_ice=0
324	iaer_stormdust_doubleq=0
325	iaer_topdust_doubleq=0
326
327	aer_count=0
328	if (.NOT.active) then
329	do iaer=1,naerkind
330	if (name_iaer(iaer).eq."dust_conrath") then
331	iaer_dust_conrath = iaer
332	aer_count = aer_count + 1
333	endif
334	enddo
335	endif
336	if (doubleq.AND.active) then
337	do iaer=1,naerkind
338	if (name_iaer(iaer).eq."dust_doubleq") then
339	iaer_dust_doubleq = iaer
340	aer_count = aer_count + 1
341	endif
342	enddo
343	endif
344	if (submicron.AND.active) then
345	do iaer=1,naerkind
346	if (name_iaer(iaer).eq."dust_submicron") then
347	iaer_dust_submicron = iaer
348	aer_count = aer_count + 1
349	endif
350	enddo
351	endif
352	if (water.AND.activice) then
353	do iaer=1,naerkind
354	if (name_iaer(iaer).eq."h2o_ice") then
355	iaer_h2o_ice = iaer
356	aer_count = aer_count + 1
357	endif
358	enddo
359	endif
360	if (co2clouds.AND.activeco2ice) then
361	do iaer=1,naerkind
362	if (name_iaer(iaer).eq."co2_ice") then
363	iaer_co2_ice = iaer
364	aer_count = aer_count + 1
365	endif
366	enddo
367	endif
368	if (rdstorm.AND.active) then
369	do iaer=1,naerkind
370	if (name_iaer(iaer).eq."stormdust_doubleq") then
371	iaer_stormdust_doubleq = iaer
372	aer_count = aer_count + 1
373	endif
374	enddo
375	end if
376	if (topflows.AND.active) then
377	do iaer=1,naerkind
378	if (name_iaer(iaer).eq."topdust_doubleq") then
379	iaer_topdust_doubleq = iaer
380	aer_count = aer_count + 1
381	endif
382	enddo
383	end if
384
385	c Check that we identified all tracers:
386	if (aer_count.ne.naerkind) then
387	write(,) "callradite: found only ",aer_count," scatterers"
388	write(,) " expected ",naerkind
389	write(,) "please make sure that the number of"
390	write(,) "scatterers in scatterers.h, the names"
391	write(,) "in callradite.F, and the flags in"
392	write(,) "callphys.def are all consistent!"
393	do iaer=1,naerkind
394	write(,)' ',iaer,' ',trim(name_iaer(iaer))
395	enddo
396	call abort_physic("callradite","incoherent scatterers",1)
397	else
398	write(,) "callradite: found all scatterers, namely:"
399	do iaer=1,naerkind
400	write(,)' ',iaer,' ',trim(name_iaer(iaer))
401	enddo
402	endif
403	c -------------------------------------------
404
405	gcp = g/cpp
406
407
408	c Loading the optical properties in external look-up tables:
409
410	CALL SUAER
411
412	! CALL SULW ! this step is now done in ini_yomlw_h
413
414	if (ngrid .EQ. 1) then
415	if (ndomainsz .NE. 1) then
416	print*
417	print*,'ATTENTION !!!'
418	print*,'pour tourner en 1D, '
419	print*,'fixer ndomainsz=1 dans phymars/dimradmars_h'
420	print*
421	call exit(1)
422	endif
423	endif
424
425	firstcall=.false.
426	END IF
427
428	c Computing aerosol optical properties and opacity
429	c ------------------------------------------------
430	c Updating aerosol size distributions:
431	CALL updatereffrad(ngrid,nlayer,
432	& rdust,rstormdust,rtopdust,rice,nuice,
433	& reffrad,nueffrad, riceco2, nuiceco2,
434	& pq,tauscaling,tau,pplay, pt)
435	c Computing 3D scattering parameters:
436	gVIS3d(:,:,:,:) = 0.
437	CALL aeroptproperties(ngrid,nlayer,reffrad,nueffrad,
438	& QVISsQREF3d,omegaVIS3d,gVIS3d,
439	& QIRsQREF3d,omegaIR3d,gIR3d,
440	& QREFvis3d,QREFir3d,
441	& omegaREFvis3d,omegaREFir3d)
442	c Computing aerosol optical depth in each layer:
443	CALL aeropacity(ngrid,nlayer,nq,zday,pplay,pplev,ls,
444	& pq,pt,tauscaling,dust_rad_adjust,IRtoVIScoef,
445	& tau_pref_scenario,tau_pref_gcm,tau,taucloudtes,
446	& aerosol,dsodust,reffrad,
447	& QREFvis3d,QREFir3d,omegaREFir3d,
448	& totstormfract,clearatm,dsords,dsotop,
449	& nohmons,clearsky,totcloudfrac)
450
451	c Starting loop on sub-domain
452	c ----------------------------
453	zgVIS3d(:,:,:,:) = 0.
454	zfluxd_sw(:,:,:) = 0.
455	zfluxu_sw(:,:,:) = 0.
456	zQVISsQREF3d(:,:,:,:) = 0.
457	zomegaVIS3d(:,:,:,:) = 0.
458	DO jd=1,ndomain
459	ig0=(jd-1)*ndomainsz
460	if (jd.eq.ndomain) then
461	nd=ngrid-ig0
462	else
463	nd=ndomainsz
464	endif
465
466	c Spliting input variable in sub-domain input variables
467	c ---------------------------------------------------
468
469	do l=1,nlaylte
470	do ig = 1,nd
471	do iaer = 1, naerkind
472	do ich = 1, nsun
473	zQVISsQREF3d(ig,l,ich,iaer) =
474	& QVISsQREF3d(ig0+ig,l,ich,iaer)
475	zomegaVIS3d(ig,l,ich,iaer) =
476	& omegaVIS3d(ig0+ig,l,ich,iaer)
477	zgVIS3d(ig,l,ich,iaer) =
478	& gVIS3d(ig0+ig,l,ich,iaer)
479	enddo
480	do ich = 1, nir
481	zQIRsQREF3d(ig,l,ich,iaer) =
482	& QIRsQREF3d(ig0+ig,l,ich,iaer)
483	zomegaIR3d(ig,l,ich,iaer) =
484	& omegaIR3d(ig0+ig,l,ich,iaer)
485	zgIR3d(ig,l,ich,iaer) =
486	& gIR3d(ig0+ig,l,ich,iaer)
487	enddo
488	enddo
489	enddo
490	enddo
491	zplev(:,:) = 0.
492	do l=1,nlaylte+1
493	do ig = 1,nd
494	zplev(ig,l) = pplev(ig0+ig,l)
495	enddo
496	enddo
497	zdp(:,:) = 0.
498
499	do l=1,nlaylte
500	do ig = 1,nd
501	zplay(ig,l) = pplay(ig0+ig,l)
502	zt(ig,l) = pt(ig0+ig,l)
503	c Thickness of each layer (Pa) :
504	zdp(ig,l)= pplev(ig0+ig,l) - pplev(ig0+ig,l+1)
505	enddo
506	enddo
507	zaerosol(:,:,:) = 0.
508	do n=1,naerkind
509	do l=1,nlaylte
510	do ig=1,nd
511	zaerosol(ig,l,n) = aerosol(ig0+ig,l,n)
512	enddo
513	enddo
514	enddo
515	zalbedo(:,:) = 0.
516	do j=1,2
517	do ig = 1,nd
518	zalbedo(ig,j) = albedo(ig0+ig,j)
519	enddo
520	enddo
521
522	c Intermediate levels: (computing tlev)
523	c ---------------------------------------
524	c Extrapolation for the air temperature above the surface
525	DO ig=1,nd
526	zztlev(ig,1)=zt(ig,1)+
527	s (zplev(ig,1)-zplay(ig,1))*
528	s (zt(ig,1)-zt(ig,2))/(zplay(ig,1)-zplay(ig,2))
529
530	zdt0(ig) = tsurf(ig0+ig) - zztlev(ig,1)
531	ENDDO
532
533	DO l=2,nlaylte
534	DO ig=1, nd
535	zztlev(ig,l)=0.5*(zt(ig,l-1)+zt(ig,l))
536	ENDDO
537	ENDDO
538
539	DO ig=1, nd
540	zztlev(ig,nlaylte+1)=zt(ig,nlaylte)
541	ENDDO
542
543
544	c Longwave ("lw") radiative transfer (= thermal infrared)
545	c -------------------------------------------------------
546	call lwmain (ig0,icount,nd,nflev
547	. ,zdp,zdt0,emis(ig0+1),zplev,zztlev,zt
548	. ,zaerosol,zzdtlw
549	. ,fluxsurf_lw(ig0+1),fluxtop_lw(ig0+1)
550	. ,znetrad
551	& ,zQIRsQREF3d,zomegaIR3d,zgIR3d
552	& ,co2ice(ig0+1))
553
554	c Shortwave ("sw") radiative transfer (= solar radiation)
555	c -------------------------------------------------------
556	c Mars solar constant (W m-2)
557	c flux_1AU = 1370 W.m-2 is the solar constant at 1 AU.
558	cste_mars=flux_1AU/(dist_sol*dist_sol)
559	zzdtsw(:,:) = 0.
560	call swmain ( nd, nflev,
561	S cste_mars, zalbedo,
562	S mu0(ig0+1), zdp, zplev, zaerosol, fract(ig0+1),
563	S zzdtsw, zfluxd_sw, zfluxu_sw,
564	& zQVISsQREF3d,zomegaVIS3d,zgVIS3d)
565	c ------------------------------------------------------------
566	c Un-spliting output variable from sub-domain input variables
567	c ------------------------------------------------------------
568
569	do l=1,nlaylte
570	do ig = 1,nd
571	dtlw(ig0+ig,l) = zzdtlw(ig,l)
572	dtsw(ig0+ig,l) = zzdtsw(ig,l)
573	enddo
574	enddo
575
576	ptlev(:, :) = 0.
577	do l=1,nlaylte+1
578	do ig = 1,nd
579	ptlev(ig0+ig,l) = zztlev(ig,l)
580	enddo
581	enddo
582
583	! copy SW fluxes at surface and TOA
584	do ig = 1,nd
585	! surface downward SW flux
586	fluxsurf_dn_sw(ig0+ig,1) = zfluxd_sw(ig,1,1)
587	fluxsurf_dn_sw(ig0+ig,2) = zfluxd_sw(ig,1,2)
588	! surface upward SW flux
589	fluxsurf_up_sw(ig0+ig,1) = zfluxu_sw(ig,1,1)
590	fluxsurf_up_sw(ig0+ig,2) = zfluxu_sw(ig,1,2)
591	! downward SW flux at top of atmosphere
592	fluxtop_dn_sw(ig0+ig,1) = zfluxd_sw(ig,nlaylte+1,1)
593	fluxtop_dn_sw(ig0+ig,2) = zfluxd_sw(ig,nlaylte+1,2)
594	! upward SW flux at top of atmosphere
595	fluxtop_up_sw(ig0+ig,1) = zfluxu_sw(ig,nlaylte+1,1)
596	fluxtop_up_sw(ig0+ig,2) = zfluxu_sw(ig,nlaylte+1,2)
597	enddo
598
599	ENDDO ! (boucle jd=1, ndomain)
600
601	c Zero tendencies for any remaining layers between nlaylte and nlayer
602	if (nlayer.gt.nlaylte) then
603	do l = nlaylte+1, nlayer
604	do ig = 1, ngrid
605	dtlw(ig, l) = 0.
606	dtsw(ig, l) = 0.
607	enddo
608	enddo
609	endif
610	c Output for debugging
611	c --------------------------------
612	c Write all nlayer layers, even though only nlaylte layers may have
613	c non-zero tendencies.
614
615	c PRINT*,'Diagnotique for the radiation'
616	c PRINT*,'albedo, emissiv, mu0,fract,fluxsurf_lw,fluxsurf_sw'
617	c PRINT*,albedo(igout,1),emis(igout),mu0(igout),
618	c s fract(igout), fluxsurf_lw(igout),
619	c $ fluxsurf_dn_sw(igout,1)+fluxsurf_dn_sw(igout,2)
620	c PRINT*,'Tlay Tlev Play Plev dT/dt SW dT/dt LW (K/s)'
621	c PRINT*,'daysec',daysec
622	c DO l=1,nlayer
623	c PRINT*,pt(igout,l),ptlev(igout,l),
624	c s pplay(igout,l),pplev(igout,l),
625	c s dtsw(igout,l),dtlw(igout,l)
626	c ENDDO
627
628
629	END SUBROUTINE callradite
630
631	END MODULE callradite_mod

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

Download in other formats:

Original Format