Line data Source code
1 : module modal_aerosol_properties_mod
2 : use shr_kind_mod, only: r8 => shr_kind_r8
3 : use physconst, only: pi
4 : use aerosol_properties_mod, only: aerosol_properties, aero_name_len
5 : use rad_constituents, only: rad_cnst_get_info, rad_cnst_get_mode_props, rad_cnst_get_aer_props
6 :
7 : implicit none
8 :
9 : private
10 :
11 : public :: modal_aerosol_properties
12 :
13 : type, extends(aerosol_properties) :: modal_aerosol_properties
14 : private
15 : real(r8), allocatable :: exp45logsig_(:)
16 : real(r8), allocatable :: voltonumblo_(:)
17 : real(r8), allocatable :: voltonumbhi_(:)
18 : integer, allocatable :: sulfate_mode_ndxs_(:)
19 : integer, allocatable :: dust_mode_ndxs_(:)
20 : integer, allocatable :: ssalt_mode_ndxs_(:)
21 : integer, allocatable :: ammon_mode_ndxs_(:)
22 : integer, allocatable :: nitrate_mode_ndxs_(:)
23 : integer, allocatable :: msa_mode_ndxs_(:)
24 : integer, allocatable :: bcarbon_mode_ndxs_(:,:)
25 : integer, allocatable :: porganic_mode_ndxs_(:,:)
26 : integer, allocatable :: sorganic_mode_ndxs_(:,:)
27 : integer :: num_soa_ = 0
28 : integer :: num_poa_ = 0
29 : integer :: num_bc_ = 0
30 : contains
31 : procedure :: number_transported
32 : procedure :: get
33 : procedure :: amcube
34 : procedure :: actfracs
35 : procedure :: num_names
36 : procedure :: mmr_names
37 : procedure :: amb_num_name
38 : procedure :: amb_mmr_name
39 : procedure :: species_type
40 : procedure :: icenuc_updates_num
41 : procedure :: icenuc_updates_mmr
42 : procedure :: apply_number_limits
43 : procedure :: hetfrz_species
44 : procedure :: optics_params
45 : procedure :: nbins_rlist
46 : procedure :: nspecies_per_bin_rlist
47 : procedure :: alogsig_rlist
48 : procedure :: soluble
49 : procedure :: min_mass_mean_rad
50 : procedure :: bin_name
51 : procedure :: scav_diam
52 : procedure :: resuspension_resize
53 : procedure :: rebin_bulk_fluxes
54 : procedure :: hydrophilic
55 :
56 : final :: destructor
57 : end type modal_aerosol_properties
58 :
59 : interface modal_aerosol_properties
60 : procedure :: constructor
61 : end interface modal_aerosol_properties
62 :
63 : logical, parameter :: debug = .false.
64 :
65 : contains
66 :
67 : !------------------------------------------------------------------------------
68 : !------------------------------------------------------------------------------
69 3072 : function constructor() result(newobj)
70 :
71 : type(modal_aerosol_properties), pointer :: newobj
72 :
73 : integer :: l, m, nmodes, ncnst_tot, mm
74 : real(r8) :: dgnumlo
75 : real(r8) :: dgnumhi
76 3072 : integer,allocatable :: nspecies(:)
77 3072 : real(r8),allocatable :: sigmag(:)
78 3072 : real(r8),allocatable :: alogsig(:)
79 3072 : real(r8),allocatable :: f1(:)
80 3072 : real(r8),allocatable :: f2(:)
81 : integer :: ierr
82 :
83 : character(len=aero_name_len) :: spectype
84 :
85 : integer :: npoa, nsoa, nbc
86 :
87 3072 : allocate(newobj,stat=ierr)
88 3072 : if( ierr /= 0 ) then
89 0 : nullify(newobj)
90 0 : return
91 : end if
92 :
93 3072 : call rad_cnst_get_info(0, nmodes=nmodes)
94 :
95 9216 : allocate(nspecies(nmodes),stat=ierr)
96 3072 : if( ierr /= 0 ) then
97 0 : nullify(newobj)
98 : return
99 : end if
100 9216 : allocate(alogsig(nmodes),stat=ierr)
101 3072 : if( ierr /= 0 ) then
102 0 : nullify(newobj)
103 : return
104 : end if
105 6144 : allocate( f1(nmodes),stat=ierr )
106 3072 : if( ierr /= 0 ) then
107 0 : nullify(newobj)
108 : return
109 : end if
110 6144 : allocate( f2(nmodes),stat=ierr )
111 3072 : if( ierr /= 0 ) then
112 0 : nullify(newobj)
113 : return
114 : end if
115 :
116 6144 : allocate(sigmag(nmodes),stat=ierr)
117 3072 : if( ierr /= 0 ) then
118 0 : nullify(newobj)
119 : return
120 : end if
121 6144 : allocate(newobj%exp45logsig_(nmodes),stat=ierr)
122 3072 : if( ierr /= 0 ) then
123 0 : nullify(newobj)
124 : return
125 : end if
126 6144 : allocate(newobj%voltonumblo_(nmodes),stat=ierr)
127 3072 : if( ierr /= 0 ) then
128 0 : nullify(newobj)
129 : return
130 : end if
131 6144 : allocate(newobj%voltonumbhi_(nmodes),stat=ierr)
132 3072 : if( ierr /= 0 ) then
133 0 : nullify(newobj)
134 : return
135 : end if
136 :
137 3072 : ncnst_tot = 0
138 :
139 18432 : do m = 1, nmodes
140 15360 : call rad_cnst_get_info(0, m, nspec=nspecies(m))
141 :
142 15360 : ncnst_tot = ncnst_tot + nspecies(m) + 1
143 :
144 : call rad_cnst_get_mode_props(0, m, sigmag=sigmag(m), &
145 15360 : dgnumhi=dgnumhi, dgnumlo=dgnumlo )
146 :
147 15360 : alogsig(m) = log(sigmag(m))
148 :
149 15360 : newobj%exp45logsig_(m) = exp(4.5_r8*alogsig(m)*alogsig(m))
150 :
151 15360 : f1(m) = 0.5_r8*exp(2.5_r8*alogsig(m)*alogsig(m))
152 15360 : f2(m) = 1._r8 + 0.25_r8*alogsig(m)
153 :
154 0 : newobj%voltonumblo_(m) = 1._r8 / ( (pi/6._r8)* &
155 15360 : (dgnumlo**3._r8)*exp(4.5_r8*alogsig(m)**2._r8) )
156 0 : newobj%voltonumbhi_(m) = 1._r8 / ( (pi/6._r8)* &
157 18432 : (dgnumhi**3._r8)*exp(4.5_r8*alogsig(m)**2._r8) )
158 :
159 : end do
160 :
161 3072 : call newobj%initialize(nmodes,ncnst_tot,nspecies,nspecies,alogsig,f1,f2,ierr)
162 :
163 3072 : npoa = 0
164 3072 : nsoa = 0
165 3072 : nbc = 0
166 :
167 3072 : m = 1
168 21504 : do l = 1,newobj%nspecies(m)
169 18432 : mm = newobj%indexer(m,l)
170 18432 : call newobj%species_type(m, l, spectype)
171 39936 : select case ( trim(spectype) )
172 : case('p-organic')
173 3072 : npoa = npoa + 1
174 : case('s-organic')
175 3072 : nsoa = nsoa + 1
176 : case('black-c')
177 36864 : nbc = nbc + 1
178 : end select
179 : end do
180 :
181 3072 : newobj%num_soa_ = nsoa
182 3072 : newobj%num_poa_ = npoa
183 3072 : newobj%num_bc_ = nbc
184 :
185 9216 : allocate(newobj%sulfate_mode_ndxs_(newobj%nbins()),stat=ierr)
186 3072 : if( ierr /= 0 ) then
187 0 : nullify(newobj)
188 : return
189 : end if
190 9216 : allocate(newobj%dust_mode_ndxs_(newobj%nbins()),stat=ierr)
191 3072 : if( ierr /= 0 ) then
192 0 : nullify(newobj)
193 : return
194 : end if
195 9216 : allocate(newobj%ssalt_mode_ndxs_(newobj%nbins()),stat=ierr)
196 3072 : if( ierr /= 0 ) then
197 0 : nullify(newobj)
198 : return
199 : end if
200 9216 : allocate(newobj%ammon_mode_ndxs_(newobj%nbins()),stat=ierr)
201 3072 : if( ierr /= 0 ) then
202 0 : nullify(newobj)
203 : return
204 : end if
205 9216 : allocate(newobj%nitrate_mode_ndxs_(newobj%nbins()),stat=ierr)
206 3072 : if( ierr /= 0 ) then
207 0 : nullify(newobj)
208 : return
209 : end if
210 9216 : allocate(newobj%msa_mode_ndxs_(newobj%nbins()),stat=ierr)
211 3072 : if( ierr /= 0 ) then
212 0 : nullify(newobj)
213 : return
214 : end if
215 :
216 18432 : newobj%sulfate_mode_ndxs_ = 0
217 18432 : newobj%dust_mode_ndxs_ = 0
218 18432 : newobj%ssalt_mode_ndxs_ = 0
219 18432 : newobj%ammon_mode_ndxs_ = 0
220 18432 : newobj%nitrate_mode_ndxs_ = 0
221 18432 : newobj%msa_mode_ndxs_ = 0
222 :
223 12288 : allocate(newobj%porganic_mode_ndxs_(newobj%nbins(),npoa),stat=ierr)
224 3072 : if( ierr /= 0 ) then
225 0 : nullify(newobj)
226 : return
227 : end if
228 12288 : allocate(newobj%sorganic_mode_ndxs_(newobj%nbins(),nsoa),stat=ierr)
229 3072 : if( ierr /= 0 ) then
230 0 : nullify(newobj)
231 : return
232 : end if
233 12288 : allocate(newobj%bcarbon_mode_ndxs_(newobj%nbins(),nbc),stat=ierr)
234 3072 : if( ierr /= 0 ) then
235 0 : nullify(newobj)
236 : return
237 : end if
238 :
239 21504 : newobj%porganic_mode_ndxs_ = 0._r8
240 21504 : newobj%sorganic_mode_ndxs_ = 0._r8
241 21504 : newobj%bcarbon_mode_ndxs_ = 0._r8
242 :
243 18432 : do m = 1,newobj%nbins()
244 15360 : npoa = 0
245 15360 : nsoa = 0
246 15360 : nbc = 0
247 :
248 67584 : do l = 1,newobj%nspecies(m)
249 49152 : mm = newobj%indexer(m,l)
250 49152 : call newobj%species_type(m, l, spectype)
251 :
252 113664 : select case ( trim(spectype) )
253 : case('sulfate')
254 12288 : newobj%sulfate_mode_ndxs_(m) = mm
255 : case('dust')
256 9216 : newobj%dust_mode_ndxs_(m) = mm
257 : case('nitrate')
258 0 : newobj%nitrate_mode_ndxs_(m) = mm
259 : case('ammonium')
260 0 : newobj%ammon_mode_ndxs_(m) = mm
261 : case('seasalt')
262 9216 : newobj%ssalt_mode_ndxs_(m) = mm
263 : case('msa')
264 0 : newobj%msa_mode_ndxs_(m) = mm
265 : case('p-organic')
266 6144 : npoa = npoa + 1
267 6144 : newobj%porganic_mode_ndxs_(m,npoa) = mm
268 : case('s-organic')
269 6144 : nsoa = nsoa + 1
270 6144 : newobj%sorganic_mode_ndxs_(m,nsoa) = mm
271 : case('black-c')
272 6144 : nbc = nbc + 1
273 98304 : newobj%bcarbon_mode_ndxs_(m,nbc) = mm
274 : end select
275 :
276 : end do
277 : end do
278 :
279 3072 : if( ierr /= 0 ) then
280 0 : nullify(newobj)
281 : return
282 : end if
283 3072 : deallocate(nspecies)
284 3072 : deallocate(alogsig)
285 3072 : deallocate(sigmag)
286 3072 : deallocate(f1)
287 3072 : deallocate(f2)
288 :
289 3072 : end function constructor
290 :
291 : !------------------------------------------------------------------------------
292 : !------------------------------------------------------------------------------
293 768 : subroutine destructor(self)
294 : type(modal_aerosol_properties), intent(inout) :: self
295 :
296 768 : if (allocated(self%exp45logsig_)) then
297 768 : deallocate(self%exp45logsig_)
298 : end if
299 768 : if (allocated(self%voltonumblo_)) then
300 768 : deallocate(self%voltonumblo_)
301 : end if
302 768 : if (allocated(self%voltonumbhi_)) then
303 768 : deallocate(self%voltonumbhi_)
304 : end if
305 :
306 768 : if (allocated(self%sulfate_mode_ndxs_)) then
307 768 : deallocate(self%sulfate_mode_ndxs_)
308 : end if
309 768 : if (allocated(self%dust_mode_ndxs_)) then
310 768 : deallocate(self%dust_mode_ndxs_)
311 : end if
312 768 : if (allocated(self%ssalt_mode_ndxs_)) then
313 768 : deallocate(self%ssalt_mode_ndxs_)
314 : end if
315 768 : if (allocated(self%ammon_mode_ndxs_)) then
316 768 : deallocate(self%ammon_mode_ndxs_)
317 : end if
318 768 : if (allocated(self%nitrate_mode_ndxs_)) then
319 768 : deallocate(self%nitrate_mode_ndxs_)
320 : end if
321 768 : if (allocated(self%msa_mode_ndxs_)) then
322 768 : deallocate(self%msa_mode_ndxs_)
323 : end if
324 768 : if (allocated(self%porganic_mode_ndxs_)) then
325 768 : deallocate(self%porganic_mode_ndxs_)
326 : end if
327 768 : if (allocated(self%sorganic_mode_ndxs_)) then
328 768 : deallocate(self%sorganic_mode_ndxs_)
329 : end if
330 768 : if (allocated(self%bcarbon_mode_ndxs_)) then
331 768 : deallocate(self%bcarbon_mode_ndxs_)
332 : end if
333 :
334 768 : call self%final()
335 :
336 768 : end subroutine destructor
337 :
338 : !------------------------------------------------------------------------------
339 : ! returns number of transported aerosol constituents
340 : !------------------------------------------------------------------------------
341 768 : integer function number_transported(self)
342 : class(modal_aerosol_properties), intent(in) :: self
343 : ! to be implemented later
344 768 : number_transported = -1
345 768 : end function number_transported
346 :
347 : !------------------------------------------------------------------------
348 : ! returns aerosol properties:
349 : ! density
350 : ! hygroscopicity
351 : ! species type
352 : ! species name
353 : ! short wave species refractive indices
354 : ! long wave species refractive indices
355 : ! species morphology
356 : !------------------------------------------------------------------------
357 272910432 : subroutine get(self, bin_ndx, species_ndx, list_ndx, density, hygro, &
358 : spectype, specname, specmorph, refindex_sw, refindex_lw)
359 :
360 : class(modal_aerosol_properties), intent(in) :: self
361 : integer, intent(in) :: bin_ndx ! bin index
362 : integer, intent(in) :: species_ndx ! species index
363 : integer, optional, intent(in) :: list_ndx ! climate or a diagnostic list number
364 : real(r8), optional, intent(out) :: density ! density (kg/m3)
365 : real(r8), optional, intent(out) :: hygro ! hygroscopicity
366 : character(len=*), optional, intent(out) :: spectype ! species type
367 : character(len=*), optional, intent(out) :: specname ! species name
368 : character(len=*), optional, intent(out) :: specmorph ! species morphology
369 : complex(r8), pointer, optional, intent(out) :: refindex_sw(:) ! short wave species refractive indices
370 : complex(r8), pointer, optional, intent(out) :: refindex_lw(:) ! long wave species refractive indices
371 :
372 : integer :: ilist
373 :
374 272910432 : if (present(list_ndx)) then
375 210898944 : ilist = list_ndx
376 : else
377 62011488 : ilist = 0
378 : end if
379 :
380 : call rad_cnst_get_aer_props(ilist, bin_ndx, species_ndx, &
381 : density_aer=density, hygro_aer=hygro, spectype=spectype, &
382 780083232 : refindex_aer_sw=refindex_sw, refindex_aer_lw=refindex_lw)
383 :
384 272910432 : if (present(specname)) then
385 8884224 : call rad_cnst_get_info(ilist, bin_ndx, species_ndx, spec_name=specname)
386 : end if
387 :
388 272910432 : if (present(specmorph)) then
389 0 : specmorph = 'UNKNOWN'
390 : end if
391 :
392 272910432 : end subroutine get
393 :
394 : !------------------------------------------------------------------------
395 : ! returns optics type and table parameters
396 : !------------------------------------------------------------------------
397 46080 : subroutine optics_params(self, list_ndx, bin_ndx, opticstype, extpsw, abspsw, asmpsw, absplw, &
398 : refrtabsw, refitabsw, refrtablw, refitablw, ncoef, prefr, prefi, sw_hygro_ext_wtp, &
399 : sw_hygro_ssa_wtp, sw_hygro_asm_wtp, lw_hygro_ext_wtp, wgtpct, nwtp, &
400 : sw_hygro_coreshell_ext, sw_hygro_coreshell_ssa, sw_hygro_coreshell_asm, lw_hygro_coreshell_ext, &
401 : corefrac, bcdust, kap, relh, nfrac, nbcdust, nkap, nrelh )
402 :
403 : class(modal_aerosol_properties), intent(in) :: self
404 : integer, intent(in) :: bin_ndx ! bin index
405 : integer, intent(in) :: list_ndx ! rad climate/diags list
406 :
407 : character(len=*), optional, intent(out) :: opticstype
408 :
409 : ! refactive index table parameters
410 : real(r8), optional, pointer :: extpsw(:,:,:,:) ! short wave specific extinction
411 : real(r8), optional, pointer :: abspsw(:,:,:,:) ! short wave specific absorption
412 : real(r8), optional, pointer :: asmpsw(:,:,:,:) ! short wave asymmetry factor
413 : real(r8), optional, pointer :: absplw(:,:,:,:) ! long wave specific absorption
414 : real(r8), optional, pointer :: refrtabsw(:,:) ! table of short wave real refractive indices for aerosols
415 : real(r8), optional, pointer :: refitabsw(:,:) ! table of short wave imaginary refractive indices for aerosols
416 : real(r8), optional, pointer :: refrtablw(:,:) ! table of long wave real refractive indices for aerosols
417 : real(r8), optional, pointer :: refitablw(:,:) ! table of long wave imaginary refractive indices for aerosols
418 : integer, optional, intent(out) :: ncoef ! number of chebychev polynomials
419 : integer, optional, intent(out) :: prefr ! number of real refractive indices in table
420 : integer, optional, intent(out) :: prefi ! number of imaginary refractive indices in table
421 :
422 : ! hygrowghtpct table parameters
423 : real(r8), optional, pointer :: sw_hygro_ext_wtp(:,:) ! short wave extinction table
424 : real(r8), optional, pointer :: sw_hygro_ssa_wtp(:,:) ! short wave single-scatter albedo table
425 : real(r8), optional, pointer :: sw_hygro_asm_wtp(:,:) ! short wave asymmetry table
426 : real(r8), optional, pointer :: lw_hygro_ext_wtp(:,:) ! long wave absorption table
427 : real(r8), optional, pointer :: wgtpct(:) ! weight precent of H2SO4/H2O solution
428 : integer, optional, intent(out) :: nwtp ! number of weight precent values
429 :
430 : ! hygrocoreshell table parameters
431 : real(r8), optional, pointer :: sw_hygro_coreshell_ext(:,:,:,:,:) ! short wave extinction table
432 : real(r8), optional, pointer :: sw_hygro_coreshell_ssa(:,:,:,:,:) ! short wave single-scatter albedo table
433 : real(r8), optional, pointer :: sw_hygro_coreshell_asm(:,:,:,:,:) ! short wave asymmetry table
434 : real(r8), optional, pointer :: lw_hygro_coreshell_ext(:,:,:,:,:) ! long wave absorption table
435 : real(r8), optional, pointer :: corefrac(:) ! core fraction dimension values
436 : real(r8), optional, pointer :: bcdust(:) ! bc/(bc + dust) fraction dimension values
437 : real(r8), optional, pointer :: kap(:) ! hygroscopicity dimension values
438 : real(r8), optional, pointer :: relh(:) ! relative humidity dimension values
439 : integer, optional, intent(out) :: nfrac ! core fraction dimension size
440 : integer, optional, intent(out) :: nbcdust ! bc/(bc + dust) fraction dimension size
441 : integer, optional, intent(out) :: nkap ! hygroscopicity dimension size
442 : integer, optional, intent(out) :: nrelh ! relative humidity dimension size
443 :
444 : ! refactive index table parameters
445 : call rad_cnst_get_mode_props(list_ndx, bin_ndx, &
446 : opticstype=opticstype, &
447 : extpsw=extpsw, &
448 : abspsw=abspsw, &
449 : asmpsw=asmpsw, &
450 : absplw=absplw, &
451 : refrtabsw=refrtabsw, &
452 : refitabsw=refitabsw, &
453 : refrtablw=refrtablw, &
454 : refitablw=refitablw, &
455 : ncoef=ncoef, &
456 : prefr=prefr, &
457 253440 : prefi=prefi)
458 :
459 : ! hygrowghtpct table parameters
460 46080 : if (present(sw_hygro_ext_wtp)) then
461 0 : nullify(sw_hygro_ext_wtp)
462 : end if
463 46080 : if (present(sw_hygro_ssa_wtp)) then
464 0 : nullify(sw_hygro_ssa_wtp)
465 : end if
466 46080 : if (present(sw_hygro_asm_wtp)) then
467 0 : nullify(sw_hygro_asm_wtp)
468 : end if
469 46080 : if (present(lw_hygro_ext_wtp)) then
470 0 : nullify(lw_hygro_ext_wtp)
471 : end if
472 46080 : if (present(wgtpct)) then
473 0 : nullify(wgtpct)
474 : end if
475 46080 : if (present(nwtp)) then
476 0 : nwtp = -1
477 : end if
478 :
479 : ! hygrocoreshell table parameters
480 46080 : if (present(sw_hygro_coreshell_ext)) then
481 0 : nullify(sw_hygro_coreshell_ext)
482 : end if
483 46080 : if (present(sw_hygro_coreshell_ssa)) then
484 0 : nullify(sw_hygro_coreshell_ssa)
485 : end if
486 46080 : if (present(sw_hygro_coreshell_asm)) then
487 0 : nullify(sw_hygro_coreshell_asm)
488 : end if
489 46080 : if (present(lw_hygro_coreshell_ext)) then
490 0 : nullify(lw_hygro_coreshell_ext)
491 : end if
492 46080 : if (present(corefrac)) then
493 0 : nullify(corefrac)
494 : end if
495 46080 : if (present(bcdust)) then
496 0 : nullify(bcdust)
497 : end if
498 46080 : if (present(kap)) then
499 0 : nullify(kap)
500 : end if
501 46080 : if (present(relh)) then
502 0 : nullify(relh)
503 : end if
504 46080 : if (present(nfrac)) then
505 0 : nfrac = -1
506 : end if
507 46080 : if (present(nbcdust)) then
508 0 : nbcdust = -1
509 : end if
510 46080 : if (present(nkap)) then
511 0 : nkap = -1
512 : end if
513 46080 : if (present(nrelh)) then
514 0 : nrelh = -1
515 : end if
516 :
517 46080 : end subroutine optics_params
518 :
519 : !------------------------------------------------------------------------------
520 : ! returns radius^3 (m3) of a given bin number
521 : !------------------------------------------------------------------------------
522 57836378 : pure elemental real(r8) function amcube(self, bin_ndx, volconc, numconc)
523 :
524 : class(modal_aerosol_properties), intent(in) :: self
525 : integer, intent(in) :: bin_ndx ! bin number
526 : real(r8), intent(in) :: volconc ! volume conc (m3/m3)
527 : real(r8), intent(in) :: numconc ! number conc (1/m3)
528 :
529 57836378 : amcube = (3._r8*volconc/(4._r8*pi*self%exp45logsig_(bin_ndx)*numconc))
530 :
531 57836378 : end function amcube
532 :
533 : !------------------------------------------------------------------------------
534 : ! returns mass and number activation fractions
535 : !------------------------------------------------------------------------------
536 8514030 : subroutine actfracs(self, bin_ndx, smc, smax, fn, fm )
537 : use shr_spfn_mod, only: erf => shr_spfn_erf
538 : class(modal_aerosol_properties), intent(in) :: self
539 : integer, intent(in) :: bin_ndx ! bin index
540 : real(r8),intent(in) :: smc ! critical supersaturation for particles of bin radius
541 : real(r8),intent(in) :: smax ! maximum supersaturation for multiple competing aerosols
542 : real(r8),intent(out) :: fn ! activation fraction for aerosol number
543 : real(r8),intent(out) :: fm ! activation fraction for aerosol mass
544 :
545 : real(r8) :: x,y
546 : real(r8), parameter :: twothird = 2._r8/3._r8
547 : real(r8), parameter :: sq2 = sqrt(2._r8)
548 :
549 8514030 : x=twothird*(log(smc)-log(smax))/(sq2*self%alogsig(bin_ndx))
550 8514030 : y=x-1.5_r8*sq2*self%alogsig(bin_ndx)
551 :
552 8514030 : fn = 0.5_r8*(1._r8-erf(x))
553 8514030 : fm = 0.5_r8*(1._r8-erf(y))
554 :
555 8514030 : end subroutine actfracs
556 :
557 : !------------------------------------------------------------------------
558 : ! returns constituents names of aerosol number mixing ratios
559 : !------------------------------------------------------------------------
560 3139200 : subroutine num_names(self, bin_ndx, name_a, name_c)
561 : class(modal_aerosol_properties), intent(in) :: self
562 : integer, intent(in) :: bin_ndx ! bin number
563 : character(len=*), intent(out) :: name_a ! constituent name of ambient aerosol number dens
564 : character(len=*), intent(out) :: name_c ! constituent name of cloud-borne aerosol number dens
565 :
566 3139200 : call rad_cnst_get_info(0,bin_ndx, num_name=name_a, num_name_cw=name_c)
567 3139200 : end subroutine num_names
568 :
569 : !------------------------------------------------------------------------
570 : ! returns constituents names of aerosol mass mixing ratios
571 : !------------------------------------------------------------------------
572 1198080 : subroutine mmr_names(self, bin_ndx, species_ndx, name_a, name_c)
573 : class(modal_aerosol_properties), intent(in) :: self
574 : integer, intent(in) :: bin_ndx ! bin number
575 : integer, intent(in) :: species_ndx ! species number
576 : character(len=*), intent(out) :: name_a ! constituent name of ambient aerosol MMR
577 : character(len=*), intent(out) :: name_c ! constituent name of cloud-borne aerosol MMR
578 :
579 1198080 : call rad_cnst_get_info(0, bin_ndx, species_ndx, spec_name=name_a, spec_name_cw=name_c)
580 1198080 : end subroutine mmr_names
581 :
582 : !------------------------------------------------------------------------
583 : ! returns constituent name of ambient aerosol number mixing ratios
584 : !------------------------------------------------------------------------
585 768 : subroutine amb_num_name(self, bin_ndx, name)
586 : class(modal_aerosol_properties), intent(in) :: self
587 : integer, intent(in) :: bin_ndx ! bin number
588 : character(len=*), intent(out) :: name ! constituent name of ambient aerosol number dens
589 :
590 768 : call rad_cnst_get_info(0,bin_ndx, num_name=name)
591 :
592 768 : end subroutine amb_num_name
593 :
594 : !------------------------------------------------------------------------
595 : ! returns constituent name of ambient aerosol mass mixing ratios
596 : !------------------------------------------------------------------------
597 768 : subroutine amb_mmr_name(self, bin_ndx, species_ndx, name)
598 : class(modal_aerosol_properties), intent(in) :: self
599 : integer, intent(in) :: bin_ndx ! bin number
600 : integer, intent(in) :: species_ndx ! species number
601 : character(len=*), intent(out) :: name ! constituent name of ambient aerosol MMR
602 :
603 768 : call rad_cnst_get_info(0, bin_ndx, species_ndx, spec_name=name)
604 :
605 768 : end subroutine amb_mmr_name
606 :
607 : !------------------------------------------------------------------------
608 : ! returns species type
609 : !------------------------------------------------------------------------
610 75654342 : subroutine species_type(self, bin_ndx, species_ndx, spectype)
611 : class(modal_aerosol_properties), intent(in) :: self
612 : integer, intent(in) :: bin_ndx ! bin number
613 : integer, intent(in) :: species_ndx ! species number
614 : character(len=*), intent(out) :: spectype ! species type
615 :
616 75654342 : call rad_cnst_get_info(0, bin_ndx, species_ndx, spec_type=spectype)
617 :
618 75654342 : end subroutine species_type
619 :
620 : !------------------------------------------------------------------------------
621 : ! returns TRUE if Ice Nucleation tendencies are applied to given aerosol bin number
622 : !------------------------------------------------------------------------------
623 51262770 : function icenuc_updates_num(self, bin_ndx) result(res)
624 : class(modal_aerosol_properties), intent(in) :: self
625 : integer, intent(in) :: bin_ndx ! bin number
626 :
627 : logical :: res
628 :
629 : character(len=aero_name_len) :: spectype
630 : character(len=aero_name_len) :: modetype
631 : integer :: spc_ndx
632 :
633 51262770 : res = .false.
634 :
635 51262770 : call rad_cnst_get_info(0, bin_ndx, mode_type=modetype)
636 51262770 : if (.not.(modetype=='coarse' .or. modetype=='coarse_dust')) then
637 : return
638 : end if
639 :
640 41010216 : do spc_ndx = 1, self%nspecies(bin_ndx)
641 30757662 : call self%species_type( bin_ndx, spc_ndx, spectype)
642 41010216 : if (spectype=='dust') res = .true.
643 : end do
644 :
645 10252554 : end function icenuc_updates_num
646 :
647 : !------------------------------------------------------------------------------
648 : ! returns TRUE if Ice Nucleation tendencies are applied to a given species within a bin
649 : !------------------------------------------------------------------------------
650 30758430 : function icenuc_updates_mmr(self, bin_ndx, species_ndx) result(res)
651 : class(modal_aerosol_properties), intent(in) :: self
652 : integer, intent(in) :: bin_ndx ! bin number
653 : integer, intent(in) :: species_ndx ! species number
654 :
655 : logical :: res
656 :
657 : character(len=32) :: spectype
658 : character(len=32) :: modetype
659 :
660 30758430 : res = .false.
661 :
662 30758430 : if (species_ndx>0) then
663 :
664 30757662 : call rad_cnst_get_info(0, bin_ndx, mode_type=modetype)
665 30757662 : if (.not.(modetype=='coarse' .or. modetype=='coarse_dust')) then
666 : return
667 : end if
668 :
669 30757662 : call self%species_type( bin_ndx, species_ndx, spectype)
670 30757662 : if (spectype=='dust') res = .true.
671 : end if
672 :
673 : end function icenuc_updates_mmr
674 :
675 : !------------------------------------------------------------------------------
676 : ! apply max / min to number concentration
677 : !------------------------------------------------------------------------------
678 13715310 : subroutine apply_number_limits( self, naerosol, vaerosol, istart, istop, m )
679 : class(modal_aerosol_properties), intent(in) :: self
680 : real(r8), intent(inout) :: naerosol(:) ! number conc (1/m3)
681 : real(r8), intent(in) :: vaerosol(:) ! volume conc (m3/m3)
682 : integer, intent(in) :: istart ! start column index (1 <= istart <= istop <= pcols)
683 : integer, intent(in) :: istop ! stop column index
684 : integer, intent(in) :: m ! mode or bin index
685 :
686 : integer :: i
687 :
688 : ! adjust number so that dgnumlo < dgnum < dgnumhi
689 : ! -- the diameter falls within the lower and upper limits which are
690 : ! represented by voltonumhi and voltonumblo values, respectively
691 84644700 : do i = istart, istop
692 70929390 : naerosol(i) = max(naerosol(i), vaerosol(i)*self%voltonumbhi_(m))
693 84644700 : naerosol(i) = min(naerosol(i), vaerosol(i)*self%voltonumblo_(m))
694 : end do
695 :
696 13715310 : end subroutine apply_number_limits
697 :
698 : !------------------------------------------------------------------------------
699 : ! returns TRUE if species `spc_ndx` in aerosol subset `bin_ndx` contributes to
700 : ! the particles' ability to act as heterogeneous freezing nuclei
701 : !------------------------------------------------------------------------------
702 24576 : function hetfrz_species(self, bin_ndx, spc_ndx) result(res)
703 : class(modal_aerosol_properties), intent(in) :: self
704 : integer, intent(in) :: bin_ndx ! bin number
705 : integer, intent(in) :: spc_ndx ! species number
706 :
707 : logical :: res
708 :
709 : character(len=aero_name_len) :: mode_name, species_type
710 :
711 24576 : res = .false.
712 :
713 24576 : call rad_cnst_get_info(0, bin_ndx, mode_type=mode_name)
714 :
715 24576 : if ((trim(mode_name)/='aitken')) then
716 :
717 18432 : call self%species_type(bin_ndx, spc_ndx, species_type)
718 :
719 18432 : if ((trim(species_type)=='black-c').or.(trim(species_type)=='dust')) then
720 :
721 6144 : res = .true.
722 :
723 : end if
724 :
725 : end if
726 :
727 24576 : end function hetfrz_species
728 :
729 : !------------------------------------------------------------------------------
730 : ! returns TRUE if soluble
731 : !------------------------------------------------------------------------------
732 96768 : logical function soluble(self,bin_ndx)
733 : class(modal_aerosol_properties), intent(in) :: self
734 : integer, intent(in) :: bin_ndx ! bin number
735 :
736 : character(len=aero_name_len) :: mode_name
737 :
738 96768 : call rad_cnst_get_info(0, bin_ndx, mode_type=mode_name)
739 :
740 96768 : soluble = trim(mode_name)/='primary_carbon'
741 :
742 96768 : end function soluble
743 :
744 : !------------------------------------------------------------------------------
745 : ! returns minimum mass mean radius (meters)
746 : !------------------------------------------------------------------------------
747 96768 : function min_mass_mean_rad(self,bin_ndx,species_ndx) result(minrad)
748 : class(modal_aerosol_properties), intent(in) :: self
749 : integer, intent(in) :: bin_ndx ! bin number
750 : integer, intent(in) :: species_ndx ! species number
751 :
752 : real(r8) :: minrad ! meters
753 :
754 : integer :: nmodes
755 : character(len=aero_name_len) :: species_type, mode_type
756 :
757 96768 : call self%species_type(bin_ndx, species_ndx, spectype=species_type)
758 193536 : select case ( trim(species_type) )
759 : case('dust')
760 48384 : call rad_cnst_get_info(0, bin_ndx, mode_type=mode_type)
761 145152 : select case ( trim(mode_type) )
762 : case ('accum','fine_dust')
763 24192 : minrad = 0.258e-6_r8
764 : case ('coarse','coarse_dust')
765 24192 : minrad = 1.576e-6_r8
766 : case default
767 96768 : minrad = -huge(1._r8)
768 : end select
769 : case('black-c')
770 48384 : call rad_cnst_get_info(0, nmodes=nmodes)
771 48384 : if (nmodes==3) then
772 : minrad = 0.04e-6_r8
773 : else
774 48384 : minrad = 0.067e-6_r8 ! from emission size
775 : endif
776 : case default
777 193536 : minrad = -huge(1._r8)
778 : end select
779 :
780 96768 : end function min_mass_mean_rad
781 :
782 : !------------------------------------------------------------------------------
783 : ! returns the total number of bins for a given radiation list index
784 : !------------------------------------------------------------------------------
785 4608 : function nbins_rlist(self, list_ndx) result(res)
786 : class(modal_aerosol_properties), intent(in) :: self
787 : integer, intent(in) :: list_ndx ! radiation list number
788 :
789 : integer :: res
790 :
791 4608 : call rad_cnst_get_info(list_ndx, nmodes=res)
792 :
793 4608 : end function nbins_rlist
794 :
795 : !------------------------------------------------------------------------------
796 : ! returns number of species in a bin for a given radiation list index
797 : !------------------------------------------------------------------------------
798 62933760 : function nspecies_per_bin_rlist(self, list_ndx, bin_ndx) result(res)
799 : class(modal_aerosol_properties), intent(in) :: self
800 : integer, intent(in) :: list_ndx ! radiation list number
801 : integer, intent(in) :: bin_ndx ! bin number
802 :
803 : integer :: res
804 :
805 62933760 : call rad_cnst_get_info(list_ndx, bin_ndx, nspec=res)
806 :
807 62933760 : end function nspecies_per_bin_rlist
808 :
809 : !------------------------------------------------------------------------------
810 : ! returns the natural log of geometric standard deviation of the number
811 : ! distribution for radiation list number and aerosol bin
812 : !------------------------------------------------------------------------------
813 23040 : function alogsig_rlist(self, list_ndx, bin_ndx) result(res)
814 : class(modal_aerosol_properties), intent(in) :: self
815 : integer, intent(in) :: list_ndx ! radiation list number
816 : integer, intent(in) :: bin_ndx ! bin number
817 :
818 : real(r8) :: res
819 :
820 : real(r8) :: sig
821 :
822 23040 : call rad_cnst_get_mode_props(list_ndx, bin_ndx, sigmag=sig)
823 23040 : res = log(sig)
824 :
825 23040 : end function alogsig_rlist
826 :
827 : !------------------------------------------------------------------------------
828 : ! returns name for a given radiation list number and aerosol bin
829 : !------------------------------------------------------------------------------
830 15360 : function bin_name(self, list_ndx, bin_ndx) result(name)
831 : class(modal_aerosol_properties), intent(in) :: self
832 : integer, intent(in) :: list_ndx ! radiation list number
833 : integer, intent(in) :: bin_ndx ! bin number
834 :
835 : character(len=32) name
836 :
837 15360 : call rad_cnst_get_info(list_ndx, bin_ndx, mode_type=name)
838 :
839 15360 : end function bin_name
840 :
841 : !------------------------------------------------------------------------------
842 : ! returns scavenging diameter (cm) for a given aerosol bin number
843 : !------------------------------------------------------------------------------
844 11256265 : function scav_diam(self, bin_ndx) result(diam)
845 : use modal_aero_data, only: dgnum_amode
846 :
847 : class(modal_aerosol_properties), intent(in) :: self
848 : integer, intent(in) :: bin_ndx ! bin number
849 :
850 : real(r8) :: diam
851 :
852 11256265 : diam = dgnum_amode(bin_ndx)
853 :
854 11256265 : end function scav_diam
855 :
856 : !------------------------------------------------------------------------------
857 : ! adjust aerosol concentration tendencies to create larger sizes of aerosols
858 : ! during resuspension
859 : !------------------------------------------------------------------------------
860 851719 : subroutine resuspension_resize(self, dcondt)
861 :
862 11256265 : use modal_aero_data, only: mode_size_order
863 :
864 : class(modal_aerosol_properties), intent(in) :: self
865 : real(r8), intent(inout) :: dcondt(:)
866 :
867 : integer :: i
868 : character(len=4) :: spcstr
869 :
870 851719 : call accumulate_to_larger_mode( 'SO4', self%sulfate_mode_ndxs_, dcondt )
871 851719 : call accumulate_to_larger_mode( 'DUST',self%dust_mode_ndxs_,dcondt )
872 851719 : call accumulate_to_larger_mode( 'NACL',self%ssalt_mode_ndxs_,dcondt )
873 851719 : call accumulate_to_larger_mode( 'MSA', self%msa_mode_ndxs_, dcondt )
874 851719 : call accumulate_to_larger_mode( 'NH4', self%ammon_mode_ndxs_, dcondt )
875 851719 : call accumulate_to_larger_mode( 'NO3', self%nitrate_mode_ndxs_, dcondt )
876 :
877 851719 : spcstr = ' '
878 1703438 : do i = 1,self%num_soa_
879 851719 : write(spcstr,'(i4)') i
880 1703438 : call accumulate_to_larger_mode( 'SOA'//adjustl(spcstr), self%sorganic_mode_ndxs_(:,i), dcondt )
881 : enddo
882 851719 : spcstr = ' '
883 1703438 : do i = 1,self%num_poa_
884 851719 : write(spcstr,'(i4)') i
885 1703438 : call accumulate_to_larger_mode( 'POM'//adjustl(spcstr), self%porganic_mode_ndxs_(:,i), dcondt )
886 : enddo
887 851719 : spcstr = ' '
888 2555157 : do i = 1,self%num_bc_
889 851719 : write(spcstr,'(i4)') i
890 1703438 : call accumulate_to_larger_mode( 'BC'//adjustl(spcstr), self%bcarbon_mode_ndxs_(:,i), dcondt )
891 : enddo
892 :
893 : contains
894 :
895 : !------------------------------------------------------------------------------
896 7665471 : subroutine accumulate_to_larger_mode( spc_name, lptr, prevap )
897 :
898 851719 : use cam_logfile, only: iulog
899 : use spmd_utils, only: masterproc
900 :
901 : character(len=*), intent(in) :: spc_name
902 : integer, intent(in) :: lptr(:)
903 : real(r8), intent(inout) :: prevap(:)
904 :
905 : integer :: m,n, nl,ns
906 :
907 : logical, parameter :: debug = .false.
908 :
909 : ! find constituent index of the largest mode for the species
910 22996413 : loop1: do m = 1,self%nbins()-1
911 20441256 : nl = lptr(mode_size_order(m))
912 22996413 : if (nl>0) exit loop1
913 : end do loop1
914 :
915 7665471 : if (.not. nl>0) return
916 :
917 : ! accumulate the smaller modes into the largest mode
918 20441256 : do n = m+1,self%nbins()
919 15330942 : ns = lptr(mode_size_order(n))
920 20441256 : if (ns>0) then
921 8517190 : prevap(nl) = prevap(nl) + prevap(ns)
922 8517190 : prevap(ns) = 0._r8
923 : if (masterproc .and. debug) then
924 : write(iulog,'(a,i3,a,i3)') trim(spc_name)//' mode number accumulate ',ns,'->',nl
925 : endif
926 : endif
927 : end do
928 :
929 : end subroutine accumulate_to_larger_mode
930 : !------------------------------------------------------------------------------
931 :
932 : end subroutine resuspension_resize
933 :
934 : !------------------------------------------------------------------------------
935 : ! returns bulk deposition fluxes of the specified species type
936 : ! rebinned to specified diameter limits
937 : !------------------------------------------------------------------------------
938 552960 : subroutine rebin_bulk_fluxes(self, bulk_type, dep_fluxes, diam_edges, bulk_fluxes, &
939 : error_code, error_string)
940 : use infnan, only: nan, assignment(=)
941 :
942 : class(modal_aerosol_properties), intent(in) :: self
943 : character(len=*),intent(in) :: bulk_type ! aerosol type to rebin
944 : real(r8), intent(in) :: dep_fluxes(:) ! kg/m2
945 : real(r8), intent(in) :: diam_edges(:) ! meters
946 : real(r8), intent(out) :: bulk_fluxes(:) ! kg/m2
947 : integer, intent(out) :: error_code ! error code (0 if no error)
948 : character(len=*), intent(out) :: error_string ! error string
949 :
950 : real(r8) :: dns_dst ! kg/m3
951 1105920 : real(r8) :: sigma_g, vmd, tmp, massfrac_bin(size(bulk_fluxes))
952 : real(r8) :: Ntype, Mtype, Mtotal, Ntot
953 : integer :: k,l,m,mm, nbulk
954 : logical :: has_type, type_not_found
955 :
956 : character(len=aero_name_len) :: spectype
957 : character(len=aero_name_len) :: modetype
958 :
959 : real(r8), parameter :: sqrtwo = sqrt(2._r8)
960 : real(r8), parameter :: onethrd = 1._r8/3._r8
961 :
962 552960 : error_code = 0
963 552960 : error_string = ' '
964 :
965 552960 : type_not_found = .true.
966 :
967 552960 : nbulk = size(bulk_fluxes)
968 :
969 2764800 : bulk_fluxes(:) = 0.0_r8
970 :
971 3317760 : do m = 1,self%nbins()
972 2764800 : Mtype = 0._r8
973 2764800 : Mtotal = 0._r8
974 2764800 : mm = self%indexer(m,0)
975 2764800 : Ntot = dep_fluxes(mm) ! #/m2
976 :
977 2764800 : has_type = .false.
978 :
979 11612160 : do l = 1,self%nspecies(m)
980 8847360 : mm = self%indexer(m,l)
981 8847360 : call self%get(m,l, spectype=spectype, density=dns_dst) ! kg/m3
982 8847360 : if (spectype==bulk_type) then
983 1658880 : Mtype = dep_fluxes(mm) ! kg/m2
984 1658880 : has_type = .true.
985 1658880 : type_not_found = .false.
986 : end if
987 20459520 : Mtotal = Mtotal + dep_fluxes(mm) ! kg/m2
988 : end do
989 3317760 : mode_has_type: if (has_type) then
990 1658880 : call rad_cnst_get_info(0, m, mode_type=modetype)
991 1658880 : if (Ntot>1.e-40_r8 .and. Mtype>1.e-40_r8 .and. Mtotal>1.e-40_r8) then
992 :
993 1551183 : call rad_cnst_get_mode_props(0, m, sigmag=sigma_g)
994 1551183 : tmp = sqrtwo*log(sigma_g)
995 :
996 : ! type number concentration
997 1551183 : Ntype = Ntot * Mtype/Mtotal ! #/m2
998 :
999 : ! volume median diameter (meters)
1000 1551183 : vmd = (6._r8*Mtype/(pi*Ntype*dns_dst))**onethrd * exp(1.5_r8*(log(sigma_g))**2)
1001 :
1002 7755915 : massfrac_bin = 0._r8
1003 :
1004 7755915 : do k = 1,nbulk
1005 12409464 : massfrac_bin(k) = 0.5_r8*( erf((log(diam_edges(k+1)/vmd))/tmp) &
1006 18614196 : - erf((log(diam_edges(k )/vmd))/tmp) )
1007 7755915 : bulk_fluxes(k) = bulk_fluxes(k) + massfrac_bin(k) * Mtype
1008 : end do
1009 :
1010 : if (debug) then
1011 : if (abs(1._r8-sum(massfrac_bin)) > 1.e-6_r8) then
1012 : write(*,*) 'rebin_bulk_fluxes WARNING mode-num, massfrac_bin, sum(massfrac_bin) = ', &
1013 : m, massfrac_bin, sum(massfrac_bin)
1014 : end if
1015 : end if
1016 :
1017 : end if
1018 : end if mode_has_type
1019 : end do
1020 :
1021 552960 : if (type_not_found) then
1022 0 : bulk_fluxes(:) = nan
1023 0 : error_code = 1
1024 0 : write(error_string,*) 'aerosol_properties::rebin_bulk_fluxes ERROR : ',trim(bulk_type),' not found'
1025 : end if
1026 :
1027 552960 : end subroutine rebin_bulk_fluxes
1028 :
1029 : !------------------------------------------------------------------------------
1030 : ! Returns TRUE if bin is hydrophilic, otherwise FALSE
1031 : !------------------------------------------------------------------------------
1032 23040 : logical function hydrophilic(self, bin_ndx)
1033 : class(modal_aerosol_properties), intent(in) :: self
1034 : integer, intent(in) :: bin_ndx ! bin number
1035 :
1036 : character(len=aero_name_len) :: modetype
1037 :
1038 23040 : call rad_cnst_get_info(0, bin_ndx, mode_type=modetype)
1039 :
1040 23040 : hydrophilic = (trim(modetype) == 'accum')
1041 :
1042 23040 : end function hydrophilic
1043 :
1044 2304 : end module modal_aerosol_properties_mod
|
