Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

radiation_homogeneous_lw.F90

Last change on this file was 5185, checked in by abarral, 2 months ago
Replace REPROBUS CPP KEY by logical using handmade wonky wrapper
File size: 13.1 KB

Line
1	! radiation_homogeneous_lw.F90 - Longwave homogeneous-column (no cloud fraction) solver
2
3	! (C) Copyright 2016- 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
15	! Modifications
16	! 2017-04-11 R. Hogan Receive emission/albedo rather than planck/emissivity
17	! 2017-04-22 R. Hogan Store surface fluxes at all g-points
18	! 2017-10-23 R. Hogan Renamed single-character variables
19	! 2019-01-14 R. Hogan Save spectral flux profile if required
20
21	module radiation_homogeneous_lw
22
23	public
24
25	contains
26
27	!---------------------------------------------------------------------
28	! Longwave homogeneous solver, in which clouds are assumed to fill
29	! the gridbox horizontally
30	subroutine solver_homogeneous_lw(nlev,istartcol,iendcol, &
31	& config, cloud, &
32	& od, ssa, g, od_cloud, ssa_cloud, g_cloud, planck_hl, &
33	& emission, albedo, &
34	& flux)
35
36	use parkind1, only : jprb
37	use yomhook, only : lhook, dr_hook
38
39	use radiation_config, only : config_type
40	use radiation_cloud, only : cloud_type
41	use radiation_flux, only : flux_type, indexed_sum_profile
42	use radiation_two_stream, only : calc_two_stream_gammas_lw, &
43	& calc_reflectance_transmittance_lw, &
44	& calc_no_scattering_transmittance_lw
45	use radiation_adding_ica_lw, only : adding_ica_lw, calc_fluxes_no_scattering_lw
46	use radiation_lw_derivatives, only : calc_lw_derivatives_ica
47
48	implicit none
49
50	! Inputs
51	integer, intent(in) :: nlev ! number of model levels
52	integer, intent(in) :: istartcol, iendcol ! range of columns to process
53	type(config_type), intent(in) :: config
54	type(cloud_type), intent(in) :: cloud
55
56	! Gas and aerosol optical depth, single-scattering albedo and
57	! asymmetry factor at each longwave g-point
58	real(jprb), intent(in), dimension(config%n_g_lw, nlev, istartcol:iendcol) :: &
59	& od
60	real(jprb), intent(in), dimension(config%n_g_lw_if_scattering, nlev, istartcol:iendcol) :: &
61	& ssa, g
62
63	! Cloud and precipitation optical depth, single-scattering albedo and
64	! asymmetry factor in each longwave band
65	real(jprb), intent(in), dimension(config%n_bands_lw,nlev,istartcol:iendcol) :: &
66	& od_cloud
67	real(jprb), intent(in), dimension(config%n_bands_lw_if_scattering, &
68	& nlev,istartcol:iendcol) :: ssa_cloud, g_cloud
69
70	! Planck function at each half-level and the surface
71	real(jprb), intent(in), dimension(config%n_g_lw,nlev+1,istartcol:iendcol) :: &
72	& planck_hl
73
74	! Emission (Planck*emissivity) and albedo (1-emissivity) at the
75	! surface at each longwave g-point
76	real(jprb), intent(in), dimension(config%n_g_lw, istartcol:iendcol) &
77	& :: emission, albedo
78
79	! Output
80	type(flux_type), intent(inout):: flux
81
82	! Local variables
83
84	! Diffuse reflectance and transmittance for each layer in clear
85	! and all skies
86	real(jprb), dimension(config%n_g_lw, nlev) :: reflectance, transmittance
87
88	! Emission by a layer into the upwelling or downwelling diffuse
89	! streams, in clear and all skies
90	real(jprb), dimension(config%n_g_lw, nlev) :: source_up, source_dn
91
92	! Fluxes per g point
93	real(jprb), dimension(config%n_g_lw, nlev+1) :: flux_up, flux_dn
94
95	! Combined gas+aerosol+cloud optical depth, single scattering
96	! albedo and asymmetry factor
97	real(jprb), dimension(config%n_g_lw) :: od_total, ssa_total, g_total
98
99	! Two-stream coefficients
100	real(jprb), dimension(config%n_g_lw) :: gamma1, gamma2
101
102	! Optical depth of cloud in g-point space
103	real(jprb), dimension(config%n_g_lw) :: od_cloud_g
104
105	! Is there any cloud in the profile?
106	logical :: is_cloudy_profile
107
108	! Number of g points
109	integer :: ng
110
111	! Loop indices for level and column
112	integer :: jlev, jcol
113
114	real(jprb) :: hook_handle
115
116	if (lhook) call dr_hook('radiation_homogeneous_lw:solver_homogeneous_lw',0,hook_handle)
117
118	ng = config%n_g_lw
119
120	! Loop through columns
121	DO jcol = istartcol,iendcol
122
123	! Is there any cloud in the profile?
124	is_cloudy_profile = .false.
125	DO jlev = 1,nlev
126	if (cloud%fraction(jcol,jlev) >= config%cloud_fraction_threshold) then
127	is_cloudy_profile = .true.
128	exit
129	end if
130	end do
131
132	! If clear-sky fluxes need to be computed then we first compute
133	! the reflectance and transmittance of all layers, neglecting
134	! clouds. If clear-sky fluxes are not required then we only do
135	! the clear-sky layers since these will be needed when we come
136	! to do the total-sky fluxes.
137	DO jlev = 1,nlev
138	if (config%do_clear .or. cloud%fraction(jcol,jlev) &
139	& < config%cloud_fraction_threshold) then
140	if (config%do_lw_aerosol_scattering) then
141	! Scattering case: first compute clear-sky reflectance,
142	! transmittance etc at each model level
143	ssa_total = ssa(:,jlev,jcol)
144	g_total = g(:,jlev,jcol)
145	call calc_two_stream_gammas_lw(ng, ssa_total, g_total, &
146	& gamma1, gamma2)
147	call calc_reflectance_transmittance_lw(ng, &
148	& od(:,jlev,jcol), gamma1, gamma2, &
149	& planck_hl(:,jlev,jcol), planck_hl(:,jlev+1,jcol), &
150	& reflectance(:,jlev), transmittance(:,jlev), &
151	& source_up(:,jlev), source_dn(:,jlev))
152	else
153	! Non-scattering case: use simpler functions for
154	! transmission and emission
155	call calc_no_scattering_transmittance_lw(ng, od(:,jlev,jcol), &
156	& planck_hl(:,jlev,jcol), planck_hl(:,jlev+1, jcol), &
157	& transmittance(:,jlev), source_up(:,jlev), source_dn(:,jlev))
158	! Ensure that clear-sky reflectance is zero since it may be
159	! used in cloudy-sky case
160	reflectance(:,jlev) = 0.0_jprb
161	end if
162
163	end if
164	end do
165
166	if (config%do_clear) then
167	if (config%do_lw_aerosol_scattering) then
168	! Then use adding method to compute fluxes
169	call adding_ica_lw(ng, nlev, &
170	& reflectance, transmittance, source_up, source_dn, &
171	& emission(:,jcol), albedo(:,jcol), &
172	& flux_up, flux_dn)
173	else
174	! Simpler down-then-up method to compute fluxes
175	call calc_fluxes_no_scattering_lw(ng, nlev, &
176	& transmittance, source_up, source_dn, &
177	& emission(:,jcol), albedo(:,jcol), &
178	& flux_up, flux_dn)
179
180	end if
181
182	! Sum over g-points to compute broadband fluxes
183	flux%lw_up_clear(jcol,:) = sum(flux_up,1)
184	flux%lw_dn_clear(jcol,:) = sum(flux_dn,1)
185	! Store surface spectral downwelling fluxes
186	flux%lw_dn_surf_clear_g(:,jcol) = flux_dn(:,nlev+1)
187
188	! Save the spectral fluxes if required
189	if (config%do_save_spectral_flux) then
190	call indexed_sum_profile(flux_up, config%i_spec_from_reordered_g_lw, &
191	& flux%lw_up_clear_band(:,jcol,:))
192	call indexed_sum_profile(flux_dn, config%i_spec_from_reordered_g_lw, &
193	& flux%lw_dn_clear_band(:,jcol,:))
194	end if
195
196	end if ! Do clear-sky calculations
197
198	! Now the total-sky calculation. If this is a clear profile and
199	! clear-sky fluxes have been calculated then we can simply copy
200	! over the clear-sky fluxes, otherwise we need to compute fluxes
201	! now.
202	if (is_cloudy_profile .or. .not. config%do_clear) then
203	DO jlev = 1,nlev
204	! Compute combined gas+aerosol+cloud optical properties;
205	! note that for clear layers, the reflectance and
206	! transmittance have already been calculated
207	if (cloud%fraction(jcol,jlev) >= config%cloud_fraction_threshold) then
208	od_cloud_g = od_cloud(config%i_band_from_reordered_g_lw,jlev,jcol)
209	od_total = od(:,jlev,jcol) + od_cloud_g
210	ssa_total = 0.0_jprb
211	g_total = 0.0_jprb
212
213	if (config%do_lw_cloud_scattering) then
214	! Scattering case: calculate reflectance and
215	! transmittance at each model level
216	if (config%do_lw_aerosol_scattering) then
217	where (od_total > 0.0_jprb)
218	ssa_total = (ssa(:,jlev,jcol)*od(:,jlev,jcol) &
219	& + ssa_cloud(config%i_band_from_reordered_g_lw,jlev,jcol) &
220	& * od_cloud_g) &
221	& / od_total
222	end where
223	where (ssa_total > 0.0_jprb .AND. od_total > 0.0_jprb)
224	g_total = (g(:,jlev,jcol)ssa(:,jlev,jcol)od(:,jlev,jcol) &
225	& + g_cloud(config%i_band_from_reordered_g_lw,jlev,jcol) &
226	& * ssa_cloud(config%i_band_from_reordered_g_lw,jlev,jcol) &
227	& * od_cloud_g) &
228	& / (ssa_total*od_total)
229	end where
230	else
231	where (od_total > 0.0_jprb)
232	ssa_total = ssa_cloud(config%i_band_from_reordered_g_lw,jlev,jcol) &
233	& * od_cloud_g / od_total
234	end where
235	where (ssa_total > 0.0_jprb .AND. od_total > 0.0_jprb)
236	g_total = g_cloud(config%i_band_from_reordered_g_lw,jlev,jcol) &
237	& * ssa_cloud(config%i_band_from_reordered_g_lw,jlev,jcol) &
238	& * od_cloud_g / (ssa_total*od_total)
239	end where
240	end if
241
242	! Compute cloudy-sky reflectance, transmittance etc at
243	! each model level
244	call calc_two_stream_gammas_lw(ng, ssa_total, g_total, &
245	& gamma1, gamma2)
246	call calc_reflectance_transmittance_lw(ng, &
247	& od_total, gamma1, gamma2, &
248	& planck_hl(:,jlev,jcol), planck_hl(:,jlev+1,jcol), &
249	& reflectance(:,jlev), transmittance(:,jlev), &
250	& source_up(:,jlev), source_dn(:,jlev))
251	else
252	! No-scattering case: use simpler functions for
253	! transmission and emission
254	call calc_no_scattering_transmittance_lw(ng, od_total, &
255	& planck_hl(:,jlev,jcol), planck_hl(:,jlev+1, jcol), &
256	& transmittance(:,jlev), source_up(:,jlev), source_dn(:,jlev))
257	end if
258	end if ! is cloudy layer
259	end do
260
261	if (config%do_lw_cloud_scattering) then
262	! Use adding method to compute fluxes for an overcast sky
263	call adding_ica_lw(ng, nlev, reflectance, transmittance, source_up, source_dn, &
264	& emission(:,jcol), albedo(:,jcol), &
265	& flux_up, flux_dn)
266	else
267	! Simpler down-then-up method to compute fluxes
268	call calc_fluxes_no_scattering_lw(ng, nlev, &
269	& transmittance, source_up, source_dn, emission(:,jcol), albedo(:,jcol), &
270	& flux_up, flux_dn)
271	end if
272
273	! Store overcast broadband fluxes
274	flux%lw_up(jcol,:) = sum(flux_up,1)
275	flux%lw_dn(jcol,:) = sum(flux_dn,1)
276	! Store surface spectral downwelling fluxes
277	flux%lw_dn_surf_g(:,jcol) = flux_dn(:,nlev+1)
278
279	! Save the spectral fluxes if required
280	if (config%do_save_spectral_flux) then
281	call indexed_sum_profile(flux_up, config%i_spec_from_reordered_g_lw, &
282	& flux%lw_up_band(:,jcol,:))
283	call indexed_sum_profile(flux_dn, config%i_spec_from_reordered_g_lw, &
284	& flux%lw_dn_band(:,jcol,:))
285	end if
286
287	else
288	! No cloud in profile and clear-sky fluxes already
289	! calculated: copy them over
290	flux%lw_up(jcol,:) = flux%lw_up_clear(jcol,:)
291	flux%lw_dn(jcol,:) = flux%lw_dn_clear(jcol,:)
292	flux%lw_dn_surf_g(:,jcol) = flux%lw_dn_surf_clear_g(:,jcol)
293
294	if (config%do_save_spectral_flux) then
295	flux%lw_up_band(:,jcol,:) = flux%lw_up_clear_band(:,jcol,:)
296	flux%lw_dn_band(:,jcol,:) = flux%lw_dn_clear_band(:,jcol,:)
297	end if
298
299	end if
300
301	! Compute the longwave derivatives needed by Hogan and Bozzo
302	! (2015) approximate radiation update scheme, using clear-sky
303	! transmittance if no clouds were present in the profile,
304	! all-sky transmittance otherwise
305	if (config%do_lw_derivatives) then
306	call calc_lw_derivatives_ica(ng, nlev, jcol, transmittance, flux_up(:,nlev+1), &
307	& flux%lw_derivatives)
308
309	end if
310
311	end do
312
313	if (lhook) call dr_hook('radiation_homogeneous_lw:solver_homogeneous_lw',1,hook_handle)
314
315	end subroutine solver_homogeneous_lw
316
317	end module radiation_homogeneous_lw

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

Download in other formats:

Original Format