Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

radiation_ecckd_interface.F90

Last change on this file was 4853, checked in by idelkadi, 8 months ago

Ecrad update in LMDZ / Implementation of Ecrad double call in LMDZ

version 1.6.1 (commit 210d7911380f53a788c3cad73b3cf9b4e022ef87)
interface routines between lmdz and ecrad grouped in a new "lmdz" directory
double call method redesigned so as to go through the Ecrad initialization and configuration part only once for the entire simulation
clean-up in the read.F routine: delete unitules arguments
modification of default gas model in namelist (default: ECCKD)

File size: 13.2 KB

Line
1	! radiation_ecckd_interface.F90 - Interface to ecCKD gas optics model
2	!
3	! (C) Copyright 2020- ECMWF.
4	!
5	! This software is licensed under the terms of the Apache Licence Version 2.0
6	! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
7	!
8	! In applying this licence, ECMWF does not waive the privileges and immunities
9	! granted to it by virtue of its status as an intergovernmental organisation
10	! nor does it submit to any jurisdiction.
11	!
12	! Author: Robin Hogan
13	! Email: r.j.hogan@ecmwf.int
14	! License: see the COPYING file for details
15	!
16
17	module radiation_ecckd_interface
18
19	implicit none
20
21	public :: setup_gas_optics, set_gas_units, gas_optics !, planck_function
22
23	contains
24
25	!---------------------------------------------------------------------
26	! Setup the ecCKD generalized gas optics model
27	subroutine setup_gas_optics(config)
28
29	use parkind1, only : jprb
30	use radiation_config
31	use yomhook, only : lhook, dr_hook, jphook
32
33	type(config_type), intent(inout), target :: config
34
35	integer :: jj
36
37	real(jphook) :: hook_handle
38
39	if (lhook) call dr_hook('radiation_ecckd_interface:setup_gas_optics',0,hook_handle)
40
41	if (config%do_sw .and. config%i_gas_model_sw == IGasModelECCKD) then
42
43	! Read shortwave ecCKD gas optics NetCDF file
44	call config%gas_optics_sw%read(trim(config%gas_optics_sw_file_name), &
45	& config%iverbosesetup)
46
47	! Copy over relevant properties
48	config%n_g_sw = config%gas_optics_sw%ng
49
50	if (config%do_cloud_aerosol_per_sw_g_point) then
51	! Bands and g points are the same
52	config%n_bands_sw = config%n_g_sw
53	else
54	! Bands are groups of g points and span a continuous region of
55	! wavenumber space
56	config%n_bands_sw = config%gas_optics_sw%spectral_def%nband
57	end if
58
59	allocate(config%i_band_from_g_sw (config%n_g_sw))
60	allocate(config%i_band_from_reordered_g_sw(config%n_g_sw))
61	allocate(config%i_g_from_reordered_g_sw (config%n_g_sw))
62
63	if (config%do_cloud_aerosol_per_sw_g_point) then
64	config%i_band_from_g_sw = [ (jj, jj = 1,config%n_g_sw) ]
65	config%i_band_from_reordered_g_sw = [ (jj, jj = 1,config%n_g_sw) ]
66	else
67	config%i_band_from_g_sw &
68	& = config%gas_optics_sw%spectral_def%i_band_number
69	config%i_band_from_reordered_g_sw &
70	& = config%gas_optics_sw%spectral_def%i_band_number
71	end if
72	config%i_g_from_reordered_g_sw = [ (jj, jj = 1,config%n_g_sw) ]
73
74	if (config%iverbosesetup >= 2) then
75	call config%gas_optics_sw%print()
76	end if
77
78	! Override solar spectral irradiance, if filename provided
79	if (config%use_spectral_solar_cycle) then
80	call config%gas_optics_sw%read_spectral_solar_cycle(config%ssi_file_name, &
81	& config%iverbosesetup, config%use_updated_solar_spectrum)
82	end if
83
84	end if
85
86	if (config%do_lw .and. config%i_gas_model_lw == IGasModelECCKD) then
87
88	! Read longwave ecCKD gas optics NetCDF file
89	call config%gas_optics_lw%read(trim(config%gas_optics_lw_file_name), &
90	& config%iverbosesetup)
91
92	! Copy over relevant properties
93	config%n_g_lw = config%gas_optics_lw%ng
94
95	if (config%do_cloud_aerosol_per_lw_g_point) then
96	! Bands and g points are the same
97	config%n_bands_lw = config%n_g_lw
98	else
99	! Bands are groups of g points and span a continuous region of
100	! wavenumber space
101	config%n_bands_lw = config%gas_optics_lw%spectral_def%nband
102	end if
103
104	allocate(config%i_band_from_g_lw (config%n_g_lw))
105	allocate(config%i_band_from_reordered_g_lw(config%n_g_lw))
106	allocate(config%i_g_from_reordered_g_lw (config%n_g_lw))
107
108	if (config%do_cloud_aerosol_per_lw_g_point) then
109	config%i_band_from_g_lw = [ (jj, jj = 1,config%n_g_lw) ]
110	config%i_band_from_reordered_g_lw = [ (jj, jj = 1,config%n_g_lw) ]
111	else
112	config%i_band_from_g_lw &
113	& = config%gas_optics_lw%spectral_def%i_band_number
114	config%i_band_from_reordered_g_lw &
115	& = config%gas_optics_lw%spectral_def%i_band_number
116	end if
117	config%i_g_from_reordered_g_lw = [ (jj, jj = 1,config%n_g_lw) ]
118
119	if (config%iverbosesetup >= 2) then
120	call config%gas_optics_lw%print()
121	end if
122
123	end if
124
125	! The i_spec_* variables are used solely for storing spectral
126	! data, and this can either be by band or by g-point
127	if (config%do_save_spectral_flux) then
128	if (config%do_save_gpoint_flux) then
129	config%n_spec_sw = config%n_g_sw
130	config%n_spec_lw = config%n_g_lw
131	config%i_spec_from_reordered_g_sw => config%i_g_from_reordered_g_sw
132	config%i_spec_from_reordered_g_lw => config%i_g_from_reordered_g_lw
133	else
134	config%n_spec_sw = config%n_bands_sw
135	config%n_spec_lw = config%n_bands_lw
136	config%i_spec_from_reordered_g_sw => config%i_band_from_reordered_g_sw
137	config%i_spec_from_reordered_g_lw => config%i_band_from_reordered_g_lw
138	end if
139	else
140	config%n_spec_sw = 0
141	config%n_spec_lw = 0
142	nullify(config%i_spec_from_reordered_g_sw)
143	nullify(config%i_spec_from_reordered_g_lw)
144	end if
145
146	if (lhook) call dr_hook('radiation_ecckd_interface:setup_gas_optics',1,hook_handle)
147
148	end subroutine setup_gas_optics
149
150
151	!---------------------------------------------------------------------
152	! Scale gas mixing ratios according to required units
153	subroutine set_gas_units(gas)
154
155	use radiation_gas, only : gas_type, IVolumeMixingRatio
156	use yomhook, only : lhook, dr_hook, jphook
157
158	type(gas_type), intent(inout) :: gas
159
160	real(jphook) :: hook_handle
161
162	if (lhook) call dr_hook('radiation_ecckd_interface:set_gas_units',0,hook_handle)
163
164	call gas%set_units(IVolumeMixingRatio)
165
166	if (lhook) call dr_hook('radiation_ecckd_interface:set_gas_units',1,hook_handle)
167
168	end subroutine set_gas_units
169
170
171	!---------------------------------------------------------------------
172	! Compute gas optical depths, shortwave scattering, Planck function
173	! and incoming shortwave radiation at top-of-atmosphere
174	subroutine gas_optics(ncol,nlev,istartcol,iendcol, &
175	& config, single_level, thermodynamics, gas, &
176	& od_lw, od_sw, ssa_sw, lw_albedo, planck_hl, lw_emission, &
177	& incoming_sw)
178
179	use parkind1, only : jprb
180	use yomhook, only : lhook, dr_hook, jphook
181
182	use radiation_config, only : config_type, IGasModelECCKD
183	use radiation_thermodynamics, only : thermodynamics_type
184	use radiation_single_level, only : single_level_type
185	use radiation_gas_constants, only : NMaxGases
186	use radiation_gas
187
188	integer, intent(in) :: ncol ! number of columns
189	integer, intent(in) :: nlev ! number of levels
190	integer, intent(in) :: istartcol, iendcol ! range of cols to process
191	type(config_type), intent(in) :: config
192	type(single_level_type), intent(in) :: single_level
193	type(thermodynamics_type),intent(in) :: thermodynamics
194	type(gas_type), intent(in) :: gas
195
196	! Longwave albedo of the surface
197	real(jprb), dimension(config%n_g_lw,istartcol:iendcol), &
198	& intent(in), optional :: lw_albedo
199
200	! Gaseous layer optical depth in longwave and shortwave, and
201	! shortwave single scattering albedo (i.e. fraction of extinction
202	! due to Rayleigh scattering) at each g-point
203	real(jprb), dimension(config%n_g_lw,nlev,istartcol:iendcol), intent(out) :: &
204	& od_lw
205	real(jprb), dimension(config%n_g_sw,nlev,istartcol:iendcol), intent(out) :: &
206	& od_sw, ssa_sw
207
208	! The Planck function (emitted flux from a black body) at half
209	! levels at each longwave g-point
210	real(jprb), dimension(config%n_g_lw,nlev+1,istartcol:iendcol), &
211	& intent(out), optional :: planck_hl
212	! Planck function for the surface (W m-2)
213	real(jprb), dimension(config%n_g_lw,istartcol:iendcol), &
214	& intent(out), optional :: lw_emission
215
216	! The incoming shortwave flux into a plane perpendicular to the
217	! incoming radiation at top-of-atmosphere in each of the shortwave
218	! g-points
219	real(jprb), dimension(config%n_g_sw,istartcol:iendcol), &
220	& intent(out), optional :: incoming_sw
221
222	! Temperature at full levels (K)
223	real(jprb) :: temperature_fl(istartcol:iendcol,nlev)
224
225	real(jprb) :: concentration_scaling(NMaxGases)
226
227	logical :: is_volume_mixing_ratio
228
229	integer :: jcol, jlev, jg
230
231	real(jphook) :: hook_handle
232
233	if (lhook) call dr_hook('radiation_ecckd_interface:gas_optics',0,hook_handle)
234
235	!temperature_fl(istartcol:iendcol,:) &
236	! & = 0.5_jprb * (thermodynamics%temperature_hl(istartcol:iendcol,1:nlev) &
237	! & +thermodynamics%temperature_hl(istartcol:iendcol,2:nlev+1))
238
239	temperature_fl(istartcol:iendcol,:) &
240	& = (thermodynamics%temperature_hl(istartcol:iendcol,1:nlev) &
241	& *thermodynamics%pressure_hl(istartcol:iendcol,1:nlev) &
242	& +thermodynamics%temperature_hl(istartcol:iendcol,2:nlev+1) &
243	& *thermodynamics%pressure_hl(istartcol:iendcol,2:nlev+1)) &
244	& / (thermodynamics%pressure_hl(istartcol:iendcol,1:nlev) &
245	& +thermodynamics%pressure_hl(istartcol:iendcol,2:nlev+1))
246
247	! Check that the gas concentrations are stored in volume mixing
248	! ratio with no scaling; if not, return a vector of scalings
249	call gas%assert_units(IVolumeMixingRatio, scale_factor=1.0_jprb, &
250	& istatus=is_volume_mixing_ratio)
251	if (.not. is_volume_mixing_ratio) then
252	call gas%get_scaling(IVolumeMixingRatio, concentration_scaling)
253	else
254	concentration_scaling = 1.0_jprb
255	end if
256
257	if (config%do_sw .and. config%i_gas_model_sw == IGasModelECCKD) then
258
259	if (is_volume_mixing_ratio) then
260	call config%gas_optics_sw%calc_optical_depth(ncol,nlev,istartcol,iendcol, &
261	& NMaxGases, thermodynamics%pressure_hl, &
262	& temperature_fl, gas%mixing_ratio, &
263	& od_sw, rayleigh_od_fl=ssa_sw)
264	else
265	call config%gas_optics_sw%calc_optical_depth(ncol,nlev,istartcol,iendcol, &
266	& NMaxGases, thermodynamics%pressure_hl, &
267	& temperature_fl, gas%mixing_ratio, &
268	& od_sw, rayleigh_od_fl=ssa_sw, concentration_scaling=concentration_scaling)
269	end if
270
271	! At this point od_sw = absorption optical depth and ssa_sw =
272	! rayleigh optical depth: convert to total optical depth and
273	! single-scattering albedo
274	do jcol = istartcol,iendcol
275	do jlev = 1, nlev
276	do jg = 1, config%n_g_sw
277	od_sw(jg,jlev,jcol) = od_sw(jg,jlev,jcol) + ssa_sw(jg,jlev,jcol)
278	ssa_sw(jg,jlev,jcol) = ssa_sw(jg,jlev,jcol) / od_sw(jg,jlev,jcol)
279	end do
280	end do
281	end do
282
283	if (present(incoming_sw)) then
284	if (single_level%spectral_solar_cycle_multiplier == 0.0_jprb) then
285	call config%gas_optics_sw%calc_incoming_sw(single_level%solar_irradiance, &
286	& incoming_sw)
287	else
288	call config%gas_optics_sw%calc_incoming_sw(single_level%solar_irradiance, &
289	& incoming_sw, single_level%spectral_solar_cycle_multiplier)
290	end if
291	end if
292
293	end if
294
295	if (config%do_lw .and. config%i_gas_model_lw == IGasModelECCKD) then
296
297	if (is_volume_mixing_ratio) then
298	call config%gas_optics_lw%calc_optical_depth(ncol,nlev,istartcol,iendcol, &
299	& NMaxGases, thermodynamics%pressure_hl, &
300	& temperature_fl, gas%mixing_ratio, &
301	& od_lw)
302	else
303	call config%gas_optics_lw%calc_optical_depth(ncol,nlev,istartcol,iendcol, &
304	& NMaxGases, thermodynamics%pressure_hl, &
305	& temperature_fl, gas%mixing_ratio, &
306	& od_lw, concentration_scaling=concentration_scaling)
307	end if
308
309	! Calculate the Planck function for each g point
310	do jcol = istartcol,iendcol
311	call config%gas_optics_lw%calc_planck_function(nlev+1, &
312	& thermodynamics%temperature_hl(jcol,:), planck_hl(:,:,jcol))
313	end do
314	call config%gas_optics_lw%calc_planck_function(iendcol+1-istartcol, &
315	& single_level%skin_temperature(istartcol:iendcol), &
316	& lw_emission(:,:))
317	!NEC$ forced_collapse
318	lw_emission = lw_emission * (1.0_jprb - lw_albedo)
319
320	end if
321
322	if (lhook) call dr_hook('radiation_ecckd_interface:gas_optics',1,hook_handle)
323
324	end subroutine gas_optics
325
326	! !---------------------------------------------------------------------
327	! ! Externally facing function for computing the Planck function
328	! ! without reference to any gas profile; typically this would be used
329	! ! for computing the emission by a surface.
330	! subroutine planck_function(config, temperature, planck_surf)
331
332	! use parkind1, only : jprb
333	! use radiation_config, only : config_type
334
335	! type(config_type), intent(in) :: config
336	! real(jprb), intent(in) :: temperature
337
338	! ! Planck function of the surface (W m-2)
339	! real(jprb), dimension(config%n_g_lw), intent(out) :: planck_surf
340
341	! end subroutine planck_function
342
343	end module radiation_ecckd_interface

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

Download in other formats:

Original Format