CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
rrtmgp_sw_gas_optics.F90
1
4
6 use machine, only: kind_phys
7 use mo_rte_kind, only: wl
8 use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp
9 use mo_gas_concentrations, only: ty_gas_concs
10 use radiation_tools, only: check_error_msg
11 use netcdf
12#ifdef MPI
13 use mpi_f08
14#endif
15
16 implicit none
17 real(kind_phys),parameter :: &
18 tsi_default = 1360.85767381726, &
19 mg_default = 0.1567652, &
20 sb_default = 902.7126
21
22 ! RRTMGP k-distribution LUTs.
23 type(ty_gas_optics_rrtmgp) :: sw_gas_props
24 integer :: &
25 ntempssw, npresssw, ngptssw, nabsorberssw, nextrabsorberssw, nminorabsorberssw, &
26 nmixingfracssw, nlayerssw, nbndssw, npairssw, nminor_absorber_intervals_lowersw,&
27 nminor_absorber_intervals_uppersw, ncontributors_lowersw, ncontributors_uppersw
28 integer, dimension(:), allocatable :: &
29 kminor_start_lowersw, & !< Starting index in the [1, nContributors] vector for a contributor
30
31 kminor_start_uppersw
33 integer, dimension(:,:), allocatable :: &
34 band2gptsw, & !< Beginning and ending gpoint for each band
35 minor_limits_gpt_lowersw, & !< Beginning and ending gpoint for each minor interval in lower atmosphere
36 minor_limits_gpt_uppersw
37 integer, dimension(:,:,:), allocatable :: &
38 key_speciessw
39 real(kind_phys) :: &
40 press_ref_tropsw, & !< Reference pressure separating the lower and upper atmosphere [Pa]
41 temp_ref_psw, & !< Standard spectroscopic reference pressure [Pa]
42 temp_ref_tsw, & !< Standard spectroscopic reference temperature [K]
43 tsi_defaultsw, & !<
44 mg_defaultsw, & !< Mean value of Mg2 index over the average solar cycle from the NRLSSI2 model of solar variability
45 sb_defaultsw
46 real(kind_phys), dimension(:), allocatable :: &
47 press_refsw, & !< Pressures for reference atmosphere; press_ref(# reference layers) [Pa]
48 temp_refsw, & !< Temperatures for reference atmosphere; temp_ref(# reference layers) [K]
49 solar_quietsw, & !< Spectrally-dependent quiet sun irradiance from the NRLSSI2 model of solar variability
50 solar_facularsw, & !< Spectrally-dependent facular term from the NRLSSI2 model of solar variability
51 solar_sunspotsw
52 real(kind_phys), dimension(:,:), allocatable :: &
53 band_limssw
54 real(kind_phys), dimension(:,:,:), allocatable :: &
55 vmr_refsw, & !< Volume mixing ratios for reference atmosphere
56 kminor_lowersw, & !< (transformed from [nTemp x nEta x nGpt x nAbsorbers] array to
57
58 kminor_uppersw, &
60 rayl_lowersw, &
61 rayl_uppersw
62 real(kind_phys), dimension(:,:,:,:), allocatable :: &
63 kmajorsw
64 character(len=32), dimension(:), allocatable :: &
65 gas_namessw, & !< Names of absorbing gases
66 gas_minorsw, & !< Name of absorbing minor gas
67 identifier_minorsw, & !< Unique string identifying minor gas
68 minor_gases_lowersw, & !< Names of minor absorbing gases in lower atmosphere
69 minor_gases_uppersw, & !< Names of minor absorbing gases in upper atmosphere
70 scaling_gas_lowersw, & !< Absorption also depends on the concentration of this gas
71 scaling_gas_uppersw
72 logical(wl), dimension(:), allocatable :: &
73 minor_scales_with_density_lowersw, & !< Density scaling is applied to minor absorption coefficients
74 minor_scales_with_density_uppersw, & !< Density scaling is applied to minor absorption coefficients
75 scale_by_complement_lowersw, & !< Absorption is scaled by concentration of scaling_gas (F) or its complement (T)
76 scale_by_complement_uppersw
77contains
78
87 subroutine rrtmgp_sw_gas_optics_init(rrtmgp_root_dir, rrtmgp_sw_file_gas, &
88 active_gases_array, mpicomm, mpirank, mpiroot, errmsg, errflg)
89
90 ! Inputs
91 character(len=128),intent(in) :: &
92 rrtmgp_root_dir, & !< RTE-RRTMGP root directory
93 rrtmgp_sw_file_gas
94 character(len=*), dimension(:), intent(in) :: &
95 active_gases_array
96 type(mpi_comm),intent(in) :: &
97 mpicomm
98 integer,intent(in) :: &
99 mpirank, & !< Current MPI rank
100 mpiroot
101
102 ! Outputs
103 character(len=*), intent(out) :: &
104 errmsg
105 integer, intent(out) :: &
106 errflg
107
108 ! Local variables
109 integer :: status, ncid, dimid, varID, mpierr, iChar
110 integer,dimension(:),allocatable :: temp1, temp2, temp3, temp4
111 character(len=264) :: sw_gas_props_file
112 type(ty_gas_concs) :: gas_concs ! RRTMGP DDT containing active trace gases
113
114 ! Initialize
115 errmsg = ''
116 errflg = 0
117
118 ! Filenames are set in the gfphysics_nml
119 sw_gas_props_file = trim(rrtmgp_root_dir)//trim(rrtmgp_sw_file_gas)
120
121 ! #######################################################################################
122 !
123 ! Read dimensions for k-distribution fields...
124 ! (ONLY master processor(0), if MPI enabled)
125 !
126 ! #######################################################################################
127#ifdef MPI
128 if (mpirank .eq. mpiroot) then
129#endif
130 write (*,*) 'Reading RRTMGP shortwave k-distribution metadata ... '
131
132 ! Open file
133 status = nf90_open(trim(sw_gas_props_file), nf90_nowrite, ncid)
134
135 ! Read dimensions for k-distribution fields
136 status = nf90_inq_dimid(ncid, 'temperature', dimid)
137 status = nf90_inquire_dimension(ncid, dimid, len=ntempssw)
138 status = nf90_inq_dimid(ncid, 'pressure', dimid)
139 status = nf90_inquire_dimension(ncid, dimid, len=npresssw)
140 status = nf90_inq_dimid(ncid, 'absorber', dimid)
141 status = nf90_inquire_dimension(ncid, dimid, len=nabsorberssw)
142 status = nf90_inq_dimid(ncid, 'minor_absorber',dimid)
143 status = nf90_inquire_dimension(ncid, dimid, len=nminorabsorberssw)
144 status = nf90_inq_dimid(ncid, 'absorber_ext', dimid)
145 status = nf90_inquire_dimension(ncid, dimid, len=nextrabsorberssw)
146 status = nf90_inq_dimid(ncid, 'mixing_fraction', dimid)
147 status = nf90_inquire_dimension(ncid, dimid, len=nmixingfracssw)
148 status = nf90_inq_dimid(ncid, 'atmos_layer', dimid)
149 status = nf90_inquire_dimension(ncid, dimid, len=nlayerssw)
150 status = nf90_inq_dimid(ncid, 'bnd', dimid)
151 status = nf90_inquire_dimension(ncid, dimid, len=nbndssw)
152 status = nf90_inq_dimid(ncid, 'gpt', dimid)
153 status = nf90_inquire_dimension(ncid, dimid, len=ngptssw)
154 status = nf90_inq_dimid(ncid, 'pair', dimid)
155 status = nf90_inquire_dimension(ncid, dimid, len=npairssw)
156 status = nf90_inq_dimid(ncid, 'contributors_lower',dimid)
157 status = nf90_inquire_dimension(ncid, dimid, len=ncontributors_lowersw)
158 status = nf90_inq_dimid(ncid, 'contributors_upper', dimid)
159 status = nf90_inquire_dimension(ncid, dimid, len=ncontributors_uppersw)
160 status = nf90_inq_dimid(ncid, 'minor_absorber_intervals_lower', dimid)
161 status = nf90_inquire_dimension(ncid, dimid, len=nminor_absorber_intervals_lowersw)
162 status = nf90_inq_dimid(ncid, 'minor_absorber_intervals_upper', dimid)
163 status = nf90_inquire_dimension(ncid, dimid, len=nminor_absorber_intervals_uppersw)
164
165#ifdef MPI
166 endif ! On master processor
167
168 ! Other processors waiting...
169 call mpi_barrier(mpicomm, mpierr)
170
171 ! #######################################################################################
172 !
173 ! Broadcast dimensions...
174 ! (ALL processors)
175 !
176 ! #######################################################################################
177 call mpi_bcast(nbndssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
178 call mpi_bcast(ngptssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
179 call mpi_bcast(nmixingfracssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
180 call mpi_bcast(ntempssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
181 call mpi_bcast(npresssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
182 call mpi_bcast(nabsorberssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
183 call mpi_bcast(nextrabsorberssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
184 call mpi_bcast(nminorabsorberssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
185 call mpi_bcast(nlayerssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
186 call mpi_bcast(npairssw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
187 call mpi_bcast(ncontributors_uppersw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
188 call mpi_bcast(ncontributors_lowersw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
189 call mpi_bcast(nminor_absorber_intervals_uppersw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
190 call mpi_bcast(nminor_absorber_intervals_lowersw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
191#endif
192
193 ! #######################################################################################
194 !
195 ! Allocate space for arrays...
196 ! (ALL processors)
197 !
198 ! #######################################################################################
199 if (.not. allocated(gas_namessw)) &
200 allocate(gas_namessw(nabsorberssw))
201 if (.not. allocated(scaling_gas_lowersw)) &
202 allocate(scaling_gas_lowersw(nminor_absorber_intervals_lowersw))
203 if (.not. allocated(scaling_gas_uppersw)) &
204 allocate(scaling_gas_uppersw(nminor_absorber_intervals_uppersw))
205 if (.not. allocated(gas_minorsw)) &
206 allocate(gas_minorsw(nminorabsorberssw))
207 if (.not. allocated(identifier_minorsw)) &
208 allocate(identifier_minorsw(nminorabsorberssw))
209 if (.not. allocated(minor_gases_lowersw)) &
210 allocate(minor_gases_lowersw(nminor_absorber_intervals_lowersw))
211 if (.not. allocated(minor_gases_uppersw)) &
212 allocate(minor_gases_uppersw(nminor_absorber_intervals_uppersw))
213 if (.not. allocated(minor_limits_gpt_lowersw)) &
214 allocate(minor_limits_gpt_lowersw(npairssw,nminor_absorber_intervals_lowersw))
215 if (.not. allocated(minor_limits_gpt_uppersw)) &
216 allocate(minor_limits_gpt_uppersw(npairssw,nminor_absorber_intervals_uppersw))
217 if (.not. allocated(band2gptsw)) &
218 allocate(band2gptsw(2,nbndssw))
219 if (.not. allocated(key_speciessw)) &
220 allocate(key_speciessw(2,nlayerssw,nbndssw))
221 if (.not. allocated(band_limssw)) &
222 allocate(band_limssw(2,nbndssw))
223 if (.not. allocated(press_refsw)) &
224 allocate(press_refsw(npresssw))
225 if (.not. allocated(temp_refsw)) &
226 allocate(temp_refsw(ntempssw))
227 if (.not. allocated(vmr_refsw)) &
228 allocate(vmr_refsw(nlayerssw, nextrabsorberssw, ntempssw))
229 if (.not. allocated(kminor_lowersw)) &
230 allocate(kminor_lowersw(ncontributors_lowersw, nmixingfracssw, ntempssw))
231 if (.not. allocated(kmajorsw)) &
232 allocate(kmajorsw(ngptssw, nmixingfracssw, npresssw+1, ntempssw))
233 if (.not. allocated(kminor_start_lowersw)) &
234 allocate(kminor_start_lowersw(nminor_absorber_intervals_lowersw))
235 if (.not. allocated(kminor_uppersw)) &
236 allocate(kminor_uppersw(ncontributors_uppersw, nmixingfracssw, ntempssw))
237 if (.not. allocated(kminor_start_uppersw)) &
238 allocate(kminor_start_uppersw(nminor_absorber_intervals_uppersw))
239 if (.not. allocated(minor_scales_with_density_lowersw)) &
240 allocate(minor_scales_with_density_lowersw(nminor_absorber_intervals_lowersw))
241 if (.not. allocated(minor_scales_with_density_uppersw)) &
242 allocate(minor_scales_with_density_uppersw(nminor_absorber_intervals_uppersw))
243 if (.not. allocated(scale_by_complement_lowersw)) &
244 allocate(scale_by_complement_lowersw(nminor_absorber_intervals_lowersw))
245 if (.not. allocated(scale_by_complement_uppersw)) &
246 allocate(scale_by_complement_uppersw(nminor_absorber_intervals_uppersw))
247 if (.not. allocated(rayl_uppersw)) &
248 allocate(rayl_uppersw(ngptssw, nmixingfracssw, ntempssw))
249 if (.not. allocated(rayl_lowersw)) &
250 allocate(rayl_lowersw(ngptssw, nmixingfracssw, ntempssw))
251 if (.not. allocated(solar_quietsw)) &
252 allocate(solar_quietsw(ngptssw))
253 if (.not. allocated(solar_facularsw)) &
254 allocate(solar_facularsw(ngptssw))
255 if (.not. allocated(solar_sunspotsw)) &
256 allocate(solar_sunspotsw(ngptssw))
257 if (.not. allocated(temp1)) &
258 allocate(temp1(nminor_absorber_intervals_lowersw))
259 if (.not. allocated(temp2)) &
260 allocate(temp2(nminor_absorber_intervals_uppersw))
261 if (.not. allocated(temp3)) &
262 allocate(temp3(nminor_absorber_intervals_lowersw))
263 if (.not. allocated(temp4)) &
264 allocate(temp4(nminor_absorber_intervals_uppersw))
265
266 ! #######################################################################################
267 !
268 ! Read in data ...
269 ! (ONLY master processor(0), if MPI enabled)
270 !
271 ! #######################################################################################
272#ifdef MPI
273 if (mpirank .eq. mpiroot) then
274#endif
275 write (*,*) 'Reading RRTMGP shortwave k-distribution data ... '
276 status = nf90_inq_varid(ncid, 'gas_names', varid)
277 status = nf90_get_var( ncid, varid, gas_namessw)
278 status = nf90_inq_varid(ncid, 'scaling_gas_lower', varid)
279 status = nf90_get_var( ncid, varid, scaling_gas_lowersw)
280 status = nf90_inq_varid(ncid, 'scaling_gas_upper', varid)
281 status = nf90_get_var( ncid, varid, scaling_gas_uppersw)
282 status = nf90_inq_varid(ncid, 'gas_minor', varid)
283 status = nf90_get_var( ncid, varid, gas_minorsw)
284 status = nf90_inq_varid(ncid, 'identifier_minor', varid)
285 status = nf90_get_var( ncid, varid, identifier_minorsw)
286 status = nf90_inq_varid(ncid, 'minor_gases_lower', varid)
287 status = nf90_get_var( ncid, varid, minor_gases_lowersw)
288 status = nf90_inq_varid(ncid, 'minor_gases_upper', varid)
289 status = nf90_get_var( ncid, varid, minor_gases_uppersw)
290 status = nf90_inq_varid(ncid, 'minor_limits_gpt_lower', varid)
291 status = nf90_get_var( ncid, varid, minor_limits_gpt_lowersw)
292 status = nf90_inq_varid(ncid, 'minor_limits_gpt_upper', varid)
293 status = nf90_get_var( ncid, varid, minor_limits_gpt_uppersw)
294 status = nf90_inq_varid(ncid, 'bnd_limits_gpt', varid)
295 status = nf90_get_var( ncid, varid, band2gptsw)
296 status = nf90_inq_varid(ncid, 'key_species', varid)
297 status = nf90_get_var( ncid, varid, key_speciessw)
298 status = nf90_inq_varid(ncid,'bnd_limits_wavenumber', varid)
299 status = nf90_get_var( ncid, varid, band_limssw)
300 status = nf90_inq_varid(ncid, 'press_ref', varid)
301 status = nf90_get_var( ncid, varid, press_refsw)
302 status = nf90_inq_varid(ncid, 'temp_ref', varid)
303 status = nf90_get_var( ncid, varid, temp_refsw)
304 status = nf90_inq_varid(ncid, 'absorption_coefficient_ref_P', varid)
305 status = nf90_get_var( ncid, varid, temp_ref_psw)
306 status = nf90_inq_varid(ncid, 'absorption_coefficient_ref_T', varid)
307 status = nf90_get_var( ncid, varid, temp_ref_tsw)
308 status = nf90_inq_varid(ncid, 'tsi_default', varid)
309 if (status .eq. 0) then
310 status = nf90_get_var( ncid, varid, tsi_defaultsw)
311 else
312 tsi_defaultsw = tsi_default
313 endif
314 status = nf90_inq_varid(ncid, 'mg_default', varid)
315 if (status .eq. 0) then
316 status = nf90_get_var( ncid, varid, mg_defaultsw)
317 else
318 mg_defaultsw = mg_default
319 endif
320 status = nf90_inq_varid(ncid, 'sb_default', varid)
321 if (status .eq. 0) then
322 status = nf90_get_var( ncid, varid, sb_defaultsw)
323 else
324 sb_defaultsw = sb_default
325 endif
326 status = nf90_inq_varid(ncid, 'press_ref_trop', varid)
327 status = nf90_get_var( ncid, varid, press_ref_tropsw)
328 status = nf90_inq_varid(ncid, 'kminor_lower', varid)
329 status = nf90_get_var( ncid, varid, kminor_lowersw)
330 status = nf90_inq_varid(ncid, 'kminor_upper', varid)
331 status = nf90_get_var( ncid, varid, kminor_uppersw)
332 status = nf90_inq_varid(ncid, 'vmr_ref', varid)
333 status = nf90_get_var( ncid, varid, vmr_refsw)
334 status = nf90_inq_varid(ncid, 'kmajor', varid)
335 status = nf90_get_var( ncid, varid, kmajorsw)
336 status = nf90_inq_varid(ncid, 'kminor_start_lower', varid)
337 status = nf90_get_var( ncid, varid, kminor_start_lowersw)
338 status = nf90_inq_varid(ncid, 'kminor_start_upper', varid)
339 status = nf90_get_var( ncid, varid, kminor_start_uppersw)
340 status = nf90_inq_varid(ncid, 'solar_source_quiet', varid)
341 status = nf90_get_var( ncid, varid, solar_quietsw)
342 status = nf90_inq_varid(ncid, 'solar_source_facular', varid)
343 status = nf90_get_var( ncid, varid, solar_facularsw)
344 status = nf90_inq_varid(ncid, 'solar_source_sunspot', varid)
345 status = nf90_get_var( ncid, varid, solar_sunspotsw)
346 status = nf90_inq_varid(ncid, 'rayl_lower', varid)
347 status = nf90_get_var( ncid, varid, rayl_lowersw)
348 status = nf90_inq_varid(ncid, 'rayl_upper', varid)
349 status = nf90_get_var( ncid, varid, rayl_uppersw)
350
351 ! Logical fields are read in as integers and then converted to logicals.
352 status = nf90_inq_varid(ncid,'minor_scales_with_density_lower', varid)
353 status = nf90_get_var( ncid, varid,temp1)
354 minor_scales_with_density_lowersw(:) = .false.
355 where(temp1 .eq. 1) minor_scales_with_density_lowersw(:) = .true.
356 status = nf90_inq_varid(ncid,'minor_scales_with_density_upper', varid)
357 status = nf90_get_var( ncid, varid,temp2)
358 minor_scales_with_density_uppersw(:) = .false.
359 where(temp2 .eq. 1) minor_scales_with_density_uppersw(:) = .true.
360 status = nf90_inq_varid(ncid,'scale_by_complement_lower', varid)
361 status = nf90_get_var( ncid, varid,temp3)
362 scale_by_complement_lowersw(:) = .false.
363 where(temp3 .eq. 1) scale_by_complement_lowersw(:) = .true.
364 status = nf90_inq_varid(ncid,'scale_by_complement_upper', varid)
365 status = nf90_get_var( ncid, varid,temp4)
366 scale_by_complement_uppersw(:) = .false.
367 where(temp4 .eq. 1) scale_by_complement_uppersw(:) = .true.
368
369 ! Close
370 status = nf90_close(ncid)
371#ifdef MPI
372 endif ! Master process
373
374 ! Other processors waiting...
375 call mpi_barrier(mpicomm, mpierr)
376
377 ! #######################################################################################
378 !
379 ! Broadcast data...
380 ! (ALL processors)
381 !
382 ! #######################################################################################
383
384 ! Real scalars
385 call mpi_bcast(press_ref_tropsw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
386 call mpi_bcast(temp_ref_psw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
387 call mpi_bcast(temp_ref_tsw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
388 call mpi_bcast(tsi_defaultsw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
389 call mpi_bcast(mg_defaultsw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
390 call mpi_bcast(sb_defaultsw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
391
392 ! Integer arrays
393 call mpi_bcast(kminor_start_lowersw, &
394 size(kminor_start_lowersw), mpi_integer, mpiroot, mpicomm, mpierr)
395 call mpi_bcast(kminor_start_uppersw, &
396 size(kminor_start_uppersw), mpi_integer, mpiroot, mpicomm, mpierr)
397 call mpi_bcast(band2gptsw, &
398 size(band2gptsw), mpi_integer, mpiroot, mpicomm, mpierr)
399 call mpi_bcast(minor_limits_gpt_lowersw, &
400 size(minor_limits_gpt_lowersw), mpi_integer, mpiroot, mpicomm, mpierr)
401 call mpi_bcast(minor_limits_gpt_uppersw, &
402 size(minor_limits_gpt_uppersw), mpi_integer, mpiroot, mpicomm, mpierr)
403 call mpi_bcast(key_speciessw, &
404 size(key_speciessw), mpi_integer, mpiroot, mpicomm, mpierr)
405
406 ! Real arrays
407 call mpi_bcast(press_refsw, &
408 size(press_refsw), mpi_double_precision, mpiroot, mpicomm, mpierr)
409 call mpi_bcast(temp_refsw, &
410 size(temp_refsw), mpi_double_precision, mpiroot, mpicomm, mpierr)
411 call mpi_bcast(solar_quietsw, &
412 size(solar_quietsw), mpi_double_precision, mpiroot, mpicomm, mpierr)
413 call mpi_bcast(solar_facularsw, &
414 size(solar_facularsw), mpi_double_precision, mpiroot, mpicomm, mpierr)
415 call mpi_bcast(solar_sunspotsw, &
416 size(solar_sunspotsw), mpi_double_precision, mpiroot, mpicomm, mpierr)
417 call mpi_bcast(band_limssw, &
418 size(band_limssw), mpi_double_precision, mpiroot, mpicomm, mpierr)
419 call mpi_bcast(vmr_refsw, &
420 size(vmr_refsw), mpi_double_precision, mpiroot, mpicomm, mpierr)
421 call mpi_bcast(kminor_lowersw, &
422 size(kminor_lowersw), mpi_double_precision, mpiroot, mpicomm, mpierr)
423 call mpi_bcast(kminor_uppersw, &
424 size(kminor_uppersw), mpi_double_precision, mpiroot, mpicomm, mpierr)
425 call mpi_bcast(rayl_lowersw, &
426 size(rayl_lowersw), mpi_double_precision, mpiroot, mpicomm, mpierr)
427 call mpi_bcast(rayl_uppersw, &
428 size(rayl_uppersw), mpi_double_precision, mpiroot, mpicomm, mpierr)
429 call mpi_bcast(kmajorsw, &
430 size(kmajorsw), mpi_double_precision, mpiroot, mpicomm, mpierr)
431
432 ! Characters
433 do ichar=1,nabsorberssw
434 call mpi_bcast(gas_namessw(ichar), &
435 len(gas_namessw(ichar)), mpi_character, mpiroot, mpicomm, mpierr)
436 enddo
437 do ichar=1,nminorabsorberssw
438 call mpi_bcast(gas_minorsw(ichar), &
439 len(gas_minorsw(ichar)), mpi_character, mpiroot, mpicomm, mpierr)
440 call mpi_bcast(identifier_minorsw(ichar), &
441 len(identifier_minorsw(ichar)), mpi_character, mpiroot, mpicomm, mpierr)
442 enddo
443 do ichar=1,nminor_absorber_intervals_lowersw
444 call mpi_bcast(minor_gases_lowersw(ichar), &
445 len(minor_gases_lowersw(ichar)), mpi_character, mpiroot, mpicomm, mpierr)
446 call mpi_bcast(scaling_gas_lowersw(ichar), &
447 len(scaling_gas_lowersw(ichar)), mpi_character, mpiroot, mpicomm, mpierr)
448 enddo
449
450 do ichar=1,nminor_absorber_intervals_uppersw
451 call mpi_bcast(minor_gases_uppersw(ichar), &
452 len(minor_gases_uppersw(ichar)), mpi_character, mpiroot, mpicomm, mpierr)
453 call mpi_bcast(scaling_gas_uppersw(ichar), &
454 len(scaling_gas_uppersw(ichar)), mpi_character, mpiroot, mpicomm, mpierr)
455 enddo
456
457 ! Logicals
458 call mpi_bcast(minor_scales_with_density_lowersw, &
459 size(minor_scales_with_density_lowersw), mpi_logical, mpiroot, mpicomm, mpierr)
460 call mpi_bcast(minor_scales_with_density_uppersw, &
461 size(minor_scales_with_density_uppersw), mpi_logical, mpiroot, mpicomm, mpierr)
462 call mpi_bcast(scale_by_complement_lowersw, &
463 size(scale_by_complement_lowersw), mpi_logical, mpiroot, mpicomm, mpierr)
464 call mpi_bcast(scale_by_complement_uppersw, &
465 size(scale_by_complement_uppersw), mpi_logical, mpiroot, mpicomm, mpierr)
466
467 call mpi_barrier(mpicomm, mpierr)
468#endif
469
470 ! #######################################################################################
471 !
472 ! Initialize RRTMGP DDT's...
473 !
474 ! #######################################################################################
475 call check_error_msg('rrtmgp_sw_gas_optics_init_gas_concs',gas_concs%init(active_gases_array))
476 call check_error_msg('rrtmgp_sw_gas_optics_init_load',sw_gas_props%load(gas_concs, &
477 gas_namessw, key_speciessw, band2gptsw, band_limssw, press_refsw, press_ref_tropsw,&
478 temp_refsw, temp_ref_psw, temp_ref_tsw, vmr_refsw, kmajorsw, kminor_lowersw, &
479 kminor_uppersw, gas_minorsw, identifier_minorsw, minor_gases_lowersw, &
480 minor_gases_uppersw, minor_limits_gpt_lowersw, minor_limits_gpt_uppersw, &
481 minor_scales_with_density_lowersw, minor_scales_with_density_uppersw, &
482 scaling_gas_lowersw, scaling_gas_uppersw, scale_by_complement_lowersw, &
483 scale_by_complement_uppersw, kminor_start_lowersw, kminor_start_uppersw, &
484 solar_quietsw, solar_facularsw, solar_sunspotsw, tsi_defaultsw, mg_defaultsw, &
485 sb_defaultsw, rayl_lowersw, rayl_uppersw))
486
487 end subroutine rrtmgp_sw_gas_optics_init
488end module rrtmgp_sw_gas_optics
489