CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
rrtmgp_lw_cloud_optics.F90
1
9 use machine, only: kind_phys
10 use mo_rte_kind, only: wl
11 use mo_cloud_optics, only: ty_cloud_optics
12 use rrtmgp_lw_gas_optics, only: lw_gas_props
13 use radiation_tools, only: check_error_msg
14 use netcdf
15#ifdef MPI
16 use mpi_f08
17#endif
18
19 implicit none
20
21 type(ty_cloud_optics) :: lw_cloud_props
22 integer :: &
23 nrghice_fromfilelw, nbandlw, nsize_liqlw, nsize_icelw, nsizereglw, &
24 ncoeff_extlw, ncoeff_ssa_glw, nboundlw, npairslw
25 real(kind_phys) :: &
26 radliq_faclw, & !< Factor for calculating LUT interpolation indices for liquid
27 radice_faclw
28 real(kind_phys), dimension(:,:), allocatable :: &
29 lut_extliqlw, & !< LUT shortwave liquid extinction coefficient
30 lut_ssaliqlw, & !< LUT shortwave liquid single scattering albedo
31 lut_asyliqlw, & !< LUT shortwave liquid asymmetry parameter
32 band_limscldlw
33 real(kind_phys), dimension(:,:,:), allocatable :: &
34 lut_exticelw, & !< LUT shortwave ice extinction coefficient
35 lut_ssaicelw, & !< LUT shortwave ice single scattering albedo
36 lut_asyicelw
37 real(kind_phys), dimension(:), allocatable :: &
38 pade_sizereg_extliqlw, & !< Particle size regime boundaries for shortwave liquid extinction
39
40 pade_sizereg_ssaliqlw, &
42 pade_sizereg_asyliqlw, &
44 pade_sizereg_exticelw, &
46 pade_sizereg_ssaicelw, &
48 pade_sizereg_asyicelw
50 real(kind_phys), dimension(:,:,:), allocatable :: &
51 pade_extliqlw, & !< PADE coefficients for shortwave liquid extinction
52 pade_ssaliqlw, & !< PADE coefficients for shortwave liquid single scattering albedo
53 pade_asyliqlw
54 real(kind_phys), dimension(:,:,:,:), allocatable :: &
55 pade_exticelw, & !< PADE coefficients for shortwave ice extinction
56 pade_ssaicelw, & !< PADE coefficients for shortwave ice single scattering albedo
57 pade_asyicelw
58
59 ! Parameters used for rain and snow(+groupel) RRTMGP cloud-optics
60 real(kind_phys), parameter :: &
61 absrain = 0.33e-3, &
62 abssnow0 = 1.5, &
63 abssnow1 = 2.34e-3
64 real(kind_phys) :: &
65 radliq_lwrlw, & !< Liquid particle size lower bound for LUT interpolation
66 radliq_uprlw, & !< Liquid particle size upper bound for LUT interpolation
67 radice_lwrlw, & !< Ice particle size upper bound for LUT interpolation
68 radice_uprlw
69
70contains
71
72 ! ######################################################################################
73 ! SUBROUTINE rrtmgp_lw_cloud_optics_init()
74 ! ######################################################################################
76 subroutine rrtmgp_lw_cloud_optics_init(rrtmgp_root_dir, rrtmgp_lw_file_clouds, &
77 doGP_cldoptics_PADE, doGP_cldoptics_LUT, nrghice, mpicomm, mpirank, mpiroot, &
78 errmsg, errflg)
79
80 ! Inputs
81 character(len=128),intent(in) :: &
82 rrtmgp_root_dir, & !< RTE-RRTMGP root directory
83 rrtmgp_lw_file_clouds
84
85 logical, intent(in) :: &
86 doGP_cldoptics_PADE,& !< Use RRTMGP cloud-optics: PADE approximation?
87 doGP_cldoptics_LUT
88 integer, intent(inout) :: &
89 nrghice
90 type(mpi_comm), intent(in) :: &
91 mpicomm
92 integer, intent(in) :: &
93 mpirank, & !< Current MPI rank
94 mpiroot
95
96 ! Outputs
97 character(len=*), intent(out) :: &
98 errmsg
99 integer, intent(out) :: &
100 errflg
101
102 ! Local variables
103 integer :: dimID,varID,status,ncid,mpierr
104 character(len=264) :: lw_cloud_props_file
105
106 ! Initialize
107 errmsg = ''
108 errflg = 0
109
110 ! Filenames are set in the physics_nml
111 lw_cloud_props_file = trim(rrtmgp_root_dir)//trim(rrtmgp_lw_file_clouds)
112
113 ! #######################################################################################
114 !
115 ! Read dimensions for longwave cloud-optics fields...
116 ! (ONLY master processor(0), if MPI enabled)
117 !
118 ! #######################################################################################
119#ifdef MPI
120 if (mpirank .eq. mpiroot) then
121#endif
122 write (*,*) 'Reading RRTMGP longwave cloud-optics metadata ... '
123
124 ! Open file
125 status = nf90_open(trim(lw_cloud_props_file), nf90_nowrite, ncid)
126
127 ! Read dimensions
128 status = nf90_inq_dimid(ncid, 'nband', dimid)
129 status = nf90_inquire_dimension(ncid, dimid, len=nbandlw)
130 status = nf90_inq_dimid(ncid, 'nrghice', dimid)
131 status = nf90_inquire_dimension(ncid, dimid, len=nrghice_fromfilelw)
132 status = nf90_inq_dimid(ncid, 'nsize_liq', dimid)
133 status = nf90_inquire_dimension(ncid, dimid, len=nsize_liqlw)
134 status = nf90_inq_dimid(ncid, 'nsize_ice', dimid)
135 status = nf90_inquire_dimension(ncid, dimid, len=nsize_icelw)
136 status = nf90_inq_dimid(ncid, 'nsizereg', dimid)
137 status = nf90_inquire_dimension(ncid, dimid, len=nsizereglw)
138 status = nf90_inq_dimid(ncid, 'ncoeff_ext', dimid)
139 status = nf90_inquire_dimension(ncid, dimid, len=ncoeff_extlw)
140 status = nf90_inq_dimid(ncid, 'ncoeff_ssa_g', dimid)
141 status = nf90_inquire_dimension(ncid, dimid, len=ncoeff_ssa_glw)
142 status = nf90_inq_dimid(ncid, 'nbound', dimid)
143 status = nf90_inquire_dimension(ncid, dimid, len=nboundlw)
144 status = nf90_inq_dimid(ncid, 'pair', dimid)
145 status = nf90_inquire_dimension(ncid, dimid, len=npairslw)
146
147#ifdef MPI
148 endif ! On master processor
149
150 ! Other processors waiting...
151 call mpi_barrier(mpicomm, mpierr)
152
153 ! #######################################################################################
154 !
155 ! Broadcast dimensions...
156 ! (ALL processors)
157 !
158 ! #######################################################################################
159 call mpi_bcast(nbandlw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
160 call mpi_bcast(nsize_liqlw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
161 call mpi_bcast(nsize_icelw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
162 call mpi_bcast(nsizereglw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
163 call mpi_bcast(ncoeff_extlw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
164 call mpi_bcast(ncoeff_ssa_glw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
165 call mpi_bcast(nboundlw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
166 call mpi_bcast(npairslw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
167#endif
168
169 ! Has the number of ice-roughnesses to use been provided from the namelist?
170 ! If so, override nrghice from cloud-optics file
171 if (nrghice .ne. 0) nrghice_fromfilelw = nrghice
172#ifdef MPI
173 call mpi_bcast(nrghice_fromfilelw, 1, mpi_integer, mpiroot, mpicomm, mpierr)
174#endif
175
176 ! #######################################################################################
177 !
178 ! Allocate space for arrays...
179 ! (ALL processors)
180 !
181 ! #######################################################################################
182 if (dogp_cldoptics_lut) then
183 allocate(lut_extliqlw(nsize_liqlw, nbandlw))
184 allocate(lut_ssaliqlw(nsize_liqlw, nbandlw))
185 allocate(lut_asyliqlw(nsize_liqlw, nbandlw))
186 allocate(lut_exticelw(nsize_icelw, nbandlw, nrghice_fromfilelw))
187 allocate(lut_ssaicelw(nsize_icelw, nbandlw, nrghice_fromfilelw))
188 allocate(lut_asyicelw(nsize_icelw, nbandlw, nrghice_fromfilelw))
189 endif
190 if (dogp_cldoptics_pade) then
191 allocate(pade_extliqlw(nbandlw, nsizereglw, ncoeff_extlw ))
192 allocate(pade_ssaliqlw(nbandlw, nsizereglw, ncoeff_ssa_glw))
193 allocate(pade_asyliqlw(nbandlw, nsizereglw, ncoeff_ssa_glw))
194 allocate(pade_exticelw(nbandlw, nsizereglw, ncoeff_extlw, nrghice_fromfilelw))
195 allocate(pade_ssaicelw(nbandlw, nsizereglw, ncoeff_ssa_glw, nrghice_fromfilelw))
196 allocate(pade_asyicelw(nbandlw, nsizereglw, ncoeff_ssa_glw, nrghice_fromfilelw))
197 allocate(pade_sizereg_extliqlw(nboundlw))
198 allocate(pade_sizereg_ssaliqlw(nboundlw))
199 allocate(pade_sizereg_asyliqlw(nboundlw))
200 allocate(pade_sizereg_exticelw(nboundlw))
201 allocate(pade_sizereg_ssaicelw(nboundlw))
202 allocate(pade_sizereg_asyicelw(nboundlw))
203 endif
204 allocate(band_limscldlw(2,nbandlw))
205
206 ! #######################################################################################
207 !
208 ! Read in data ...
209 ! (ONLY master processor(0), if MPI enabled)
210 !
211 ! #######################################################################################
212#ifdef MPI
213 if (mpirank .eq. mpiroot) then
214#endif
215 ! Read in fields from file
216 if (dogp_cldoptics_lut) then
217 write (*,*) 'Reading RRTMGP longwave cloud data (LUT) ... '
218 status = nf90_inq_varid(ncid,'radliq_lwr',varid)
219 status = nf90_get_var(ncid,varid,radliq_lwrlw)
220 status = nf90_inq_varid(ncid,'radliq_upr',varid)
221 status = nf90_get_var(ncid,varid,radliq_uprlw)
222 status = nf90_inq_varid(ncid,'radliq_fac',varid)
223 status = nf90_get_var(ncid,varid,radliq_faclw)
224 status = nf90_inq_varid(ncid,'radice_lwr',varid)
225 status = nf90_get_var(ncid,varid,radice_lwrlw)
226 status = nf90_inq_varid(ncid,'radice_upr',varid)
227 status = nf90_get_var(ncid,varid,radice_uprlw)
228 status = nf90_inq_varid(ncid,'radice_fac',varid)
229 status = nf90_get_var(ncid,varid,radice_faclw)
230 status = nf90_inq_varid(ncid,'lut_extliq',varid)
231 status = nf90_get_var(ncid,varid,lut_extliqlw)
232 status = nf90_inq_varid(ncid,'lut_ssaliq',varid)
233 status = nf90_get_var(ncid,varid,lut_ssaliqlw)
234 status = nf90_inq_varid(ncid,'lut_asyliq',varid)
235 status = nf90_get_var(ncid,varid,lut_asyliqlw)
236 status = nf90_inq_varid(ncid,'lut_extice',varid)
237 status = nf90_get_var(ncid,varid,lut_exticelw)
238 status = nf90_inq_varid(ncid,'lut_ssaice',varid)
239 status = nf90_get_var(ncid,varid,lut_ssaicelw)
240 status = nf90_inq_varid(ncid,'lut_asyice',varid)
241 status = nf90_get_var(ncid,varid,lut_asyicelw)
242 status = nf90_inq_varid(ncid,'bnd_limits_wavenumber',varid)
243 status = nf90_get_var(ncid,varid,band_limscldlw)
244 endif
245 if (dogp_cldoptics_pade) then
246 write (*,*) 'Reading RRTMGP longwave cloud data (PADE) ... '
247 status = nf90_inq_varid(ncid,'radliq_lwr',varid)
248 status = nf90_get_var(ncid,varid,radliq_lwrlw)
249 status = nf90_inq_varid(ncid,'radliq_upr',varid)
250 status = nf90_get_var(ncid,varid,radliq_uprlw)
251 status = nf90_inq_varid(ncid,'radliq_fac',varid)
252 status = nf90_get_var(ncid,varid,radliq_faclw)
253 status = nf90_inq_varid(ncid,'radice_lwr',varid)
254 status = nf90_get_var(ncid,varid,radice_lwrlw)
255 status = nf90_inq_varid(ncid,'radice_upr',varid)
256 status = nf90_get_var(ncid,varid,radice_uprlw)
257 status = nf90_inq_varid(ncid,'radice_fac',varid)
258 status = nf90_get_var(ncid,varid,radice_faclw)
259 status = nf90_inq_varid(ncid,'pade_extliq',varid)
260 status = nf90_get_var(ncid,varid,pade_extliqlw)
261 status = nf90_inq_varid(ncid,'pade_ssaliq',varid)
262 status = nf90_get_var(ncid,varid,pade_ssaliqlw)
263 status = nf90_inq_varid(ncid,'pade_asyliq',varid)
264 status = nf90_get_var(ncid,varid,pade_asyliqlw)
265 status = nf90_inq_varid(ncid,'pade_extice',varid)
266 status = nf90_get_var(ncid,varid,pade_exticelw)
267 status = nf90_inq_varid(ncid,'pade_ssaice',varid)
268 status = nf90_get_var(ncid,varid,pade_ssaicelw)
269 status = nf90_inq_varid(ncid,'pade_asyice',varid)
270 status = nf90_get_var(ncid,varid,pade_asyicelw)
271 status = nf90_inq_varid(ncid,'pade_sizreg_extliq',varid)
272 status = nf90_get_var(ncid,varid,pade_sizereg_extliqlw)
273 status = nf90_inq_varid(ncid,'pade_sizreg_ssaliq',varid)
274 status = nf90_get_var(ncid,varid,pade_sizereg_ssaliqlw)
275 status = nf90_inq_varid(ncid,'pade_sizreg_asyliq',varid)
276 status = nf90_get_var(ncid,varid,pade_sizereg_asyliqlw)
277 status = nf90_inq_varid(ncid,'pade_sizreg_extice',varid)
278 status = nf90_get_var(ncid,varid,pade_sizereg_exticelw)
279 status = nf90_inq_varid(ncid,'pade_sizreg_ssaice',varid)
280 status = nf90_get_var(ncid,varid,pade_sizereg_ssaicelw)
281 status = nf90_inq_varid(ncid,'pade_sizreg_asyice',varid)
282 status = nf90_get_var(ncid,varid,pade_sizereg_asyicelw)
283 status = nf90_inq_varid(ncid,'bnd_limits_wavenumber',varid)
284 status = nf90_get_var(ncid,varid,band_limscldlw)
285 endif
286
287 ! Close file
288 status = nf90_close(ncid)
289#ifdef MPI
290 endif ! Master process
291
292 ! Other processors waiting...
293 call mpi_barrier(mpicomm, mpierr)
294
295 ! #######################################################################################
296 !
297 ! Broadcast data...
298 ! (ALL processors)
299 !
300 ! #######################################################################################
301
302 ! Real scalars
303 call mpi_bcast(radliq_faclw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
304 call mpi_bcast(radice_faclw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
305 call mpi_bcast(radliq_lwrlw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
306 call mpi_bcast(radliq_uprlw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
307 call mpi_bcast(radice_lwrlw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
308 call mpi_bcast(radice_uprlw, 1, mpi_double_precision, mpiroot, mpicomm, mpierr)
309
310 ! Real arrays
311 call mpi_bcast(band_limscldlw, size(band_limscldlw), &
312 mpi_double_precision, mpiroot, mpicomm, mpierr)
313 if (dogp_cldoptics_lut) then
314 call mpi_bcast(lut_extliqlw, size(lut_extliqlw), &
315 mpi_double_precision, mpiroot, mpicomm, mpierr)
316 call mpi_bcast(lut_ssaliqlw, size(lut_ssaliqlw), &
317 mpi_double_precision, mpiroot, mpicomm, mpierr)
318 call mpi_bcast(lut_asyliqlw, size(lut_asyliqlw), &
319 mpi_double_precision, mpiroot, mpicomm, mpierr)
320 call mpi_bcast(lut_exticelw, size(lut_exticelw), &
321 mpi_double_precision, mpiroot, mpicomm, mpierr)
322 call mpi_bcast(lut_ssaicelw, size(lut_ssaicelw), &
323 mpi_double_precision, mpiroot, mpicomm, mpierr)
324 call mpi_bcast(lut_asyicelw, size(lut_asyicelw), &
325 mpi_double_precision, mpiroot, mpicomm, mpierr)
326 endif
327 if (dogp_cldoptics_pade) then
328 call mpi_bcast(pade_extliqlw, size(pade_extliqlw), &
329 mpi_double_precision, mpiroot, mpicomm, mpierr)
330 call mpi_bcast(pade_ssaliqlw, size(pade_ssaliqlw), &
331 mpi_double_precision, mpiroot, mpicomm, mpierr)
332 call mpi_bcast(pade_asyliqlw, size(pade_asyliqlw), &
333 mpi_double_precision, mpiroot, mpicomm, mpierr)
334 call mpi_bcast(pade_exticelw, size(pade_exticelw), &
335 mpi_double_precision, mpiroot, mpicomm, mpierr)
336 call mpi_bcast(pade_ssaicelw, size(pade_ssaicelw), &
337 mpi_double_precision, mpiroot, mpicomm, mpierr)
338 call mpi_bcast(pade_asyicelw, size(pade_asyicelw), &
339 mpi_double_precision, mpiroot, mpicomm, mpierr)
340 call mpi_bcast(pade_sizereg_extliqlw, size(pade_sizereg_extliqlw), &
341 mpi_double_precision, mpiroot, mpicomm, mpierr)
342 call mpi_bcast(pade_sizereg_ssaliqlw, size(pade_sizereg_ssaliqlw), &
343 mpi_double_precision, mpiroot, mpicomm, mpierr)
344 call mpi_bcast(pade_sizereg_asyliqlw, size(pade_sizereg_asyliqlw), &
345 mpi_double_precision, mpiroot, mpicomm, mpierr)
346 call mpi_bcast(pade_sizereg_exticelw, size(pade_sizereg_exticelw), &
347 mpi_double_precision, mpiroot, mpicomm, mpierr)
348 call mpi_bcast(pade_sizereg_ssaicelw, size(pade_sizereg_ssaicelw), &
349 mpi_double_precision, mpiroot, mpicomm, mpierr)
350 call mpi_bcast(pade_sizereg_asyicelw, size(pade_sizereg_asyicelw), &
351 mpi_double_precision, mpiroot, mpicomm, mpierr)
352 endif
353#endif
354
355 ! #######################################################################################
356 !
357 ! Initialize RRTMGP DDT's...
358 !
359 ! #######################################################################################
360 if (dogp_cldoptics_lut) then
361 call check_error_msg('lw_cloud_optics_init',lw_cloud_props%load(band_limscldlw, &
362 radliq_lwrlw, radliq_uprlw, radliq_faclw, radice_lwrlw, radice_uprlw, &
363 radice_faclw, lut_extliqlw, lut_ssaliqlw, lut_asyliqlw, lut_exticelw, &
364 lut_ssaicelw, lut_asyicelw))
365 endif
366
367 if (dogp_cldoptics_pade) then
368 call check_error_msg('lw_cloud_optics_init', lw_cloud_props%load(band_limscldlw, &
369 pade_extliqlw, pade_ssaliqlw, pade_asyliqlw, pade_exticelw, pade_ssaicelw, &
370 pade_asyicelw, pade_sizereg_extliqlw, pade_sizereg_ssaliqlw, &
371 pade_sizereg_asyliqlw, pade_sizereg_exticelw, pade_sizereg_ssaicelw, &
372 pade_sizereg_asyicelw))
373 endif
374
375 call check_error_msg('lw_cloud_optics_init',lw_cloud_props%set_ice_roughness(nrghice))
376
377 end subroutine rrtmgp_lw_cloud_optics_init
378end module rrtmgp_lw_cloud_optics