56 character(len=AERO_MODE_NAME_LEN) :: name
60 real(kind=dp) :: char_radius
62 real(kind=dp) :: log10_std_dev_radius
64 real(kind=dp),
allocatable :: sample_radius(:)
66 real(kind=dp),
allocatable :: sample_num_conc(:)
68 real(kind=dp) :: num_conc
70 real(kind=dp),
allocatable :: vol_frac(:)
73 real(kind=dp),
allocatable :: vol_frac_std(:)
86 integer,
intent(in) :: type
120 call die_msg(719625922,
"unknown aero_mode type: " &
130 log10_sigma_g, bin_grid, num_conc)
133 real(kind=dp),
intent(in) :: total_num_conc
135 real(kind=dp),
intent(in) :: geom_mean_radius
137 real(kind=dp),
intent(in) :: log10_sigma_g
141 real(kind=dp),
intent(out) :: num_conc(
bin_grid_size(bin_grid))
146 num_conc(k) = total_num_conc / (sqrt(2d0 *
const%pi) &
147 * log10_sigma_g) * dexp(-(dlog10(bin_grid%centers(k)) &
148 - dlog10(geom_mean_radius))**2d0 &
149 / (2d0 * log10_sigma_g**2d0)) / dlog(10d0)
164 geom_mean_radius, log10_sigma_g, bin_grid, aero_data, vol_conc)
167 real(kind=dp),
intent(in) :: total_num_conc
169 real(kind=dp),
intent(in) :: geom_mean_radius
171 real(kind=dp),
intent(in) :: log10_sigma_g
177 real(kind=dp),
intent(out) :: vol_conc(
bin_grid_size(bin_grid))
182 log10_sigma_g, bin_grid, num_conc)
194 subroutine num_conc_exp(total_num_conc, radius_at_mean_vol, &
195 bin_grid, aero_data, num_conc)
198 real(kind=dp),
intent(in) :: total_num_conc
200 real(kind=dp),
intent(in) :: radius_at_mean_vol
206 real(kind=dp),
intent(out) :: num_conc(
bin_grid_size(bin_grid))
209 real(kind=dp) :: mean_vol, num_conc_vol
213 num_conc_vol = total_num_conc / mean_vol &
216 call vol_to_lnr(bin_grid%centers(k), num_conc_vol, num_conc(k))
224 subroutine vol_conc_exp(total_num_conc, radius_at_mean_vol, &
225 bin_grid, aero_data, vol_conc)
228 real(kind=dp),
intent(in) :: total_num_conc
230 real(kind=dp),
intent(in) :: radius_at_mean_vol
236 real(kind=dp),
intent(out) :: vol_conc(
bin_grid_size(bin_grid))
241 bin_grid, aero_data, num_conc)
250 subroutine num_conc_mono(total_num_conc, radius, bin_grid, num_conc)
253 real(kind=dp),
intent(in) :: total_num_conc
255 real(kind=dp),
intent(in) :: radius
259 real(kind=dp),
intent(out) :: num_conc(
bin_grid_size(bin_grid))
266 call warn_msg(825666877,
"monodisperse radius outside of bin_grid")
268 num_conc(k) = total_num_conc / bin_grid%widths(k)
277 bin_grid, aero_data, vol_conc)
280 real(kind=dp),
intent(in) :: total_num_conc
282 real(kind=dp),
intent(in) :: radius
288 real(kind=dp),
intent(out) :: vol_conc(
bin_grid_size(bin_grid))
295 call warn_msg(420930707,
"monodisperse radius outside of bin_grid")
297 vol_conc(k) = total_num_conc / bin_grid%widths(k) &
310 real(kind=dp),
intent(in) :: sample_radius(:)
312 real(kind=dp),
intent(in) :: sample_num_conc(:)
316 real(kind=dp),
intent(out) :: num_conc(
bin_grid_size(bin_grid))
318 integer :: i_sample, n_sample, i_lower, i_upper, i_bin
319 real(kind=dp) :: r_lower, r_upper
320 real(kind=dp) :: r_bin_lower, r_bin_upper, r1, r2, ratio
322 n_sample =
size(sample_num_conc)
323 call assert(188766208,
size(sample_radius) == n_sample + 1)
324 call assert(295384037, n_sample >= 1)
327 do i_sample = 1,n_sample
328 r_lower = sample_radius(i_sample)
329 r_upper = sample_radius(i_sample + 1)
332 if (i_upper < 1) cycle
334 i_lower = max(1, i_lower)
336 do i_bin = i_lower,i_upper
337 r_bin_lower = bin_grid%edges(i_bin)
338 r_bin_upper = bin_grid%edges(i_bin + 1)
339 r1 = max(r_lower, r_bin_lower)
340 r2 = min(r_upper, r_bin_upper)
341 ratio = (log(r2) - log(r1)) / (log(r_upper) - log(r_lower))
342 num_conc(i_bin) = num_conc(i_bin) + ratio &
343 * sample_num_conc(i_sample) / bin_grid%widths(i_bin)
353 bin_grid, aero_data, vol_conc)
356 real(kind=dp),
intent(in) :: sample_radius(:)
358 real(kind=dp),
intent(in) :: sample_num_conc(:)
364 real(kind=dp),
intent(out) :: vol_conc(
bin_grid_size(bin_grid))
386 real(kind=dp),
intent(out) :: num_conc(
bin_grid_size(bin_grid))
390 aero_mode%log10_std_dev_radius, bin_grid, num_conc)
393 aero_mode%char_radius, bin_grid, aero_data, num_conc)
395 call num_conc_mono(aero_mode%num_conc, aero_mode%char_radius, &
399 aero_mode%sample_num_conc, bin_grid, num_conc)
401 call die_msg(223903246,
"unknown aero_mode type: " &
421 real(kind=dp),
intent(out) :: vol_conc(
bin_grid_size(bin_grid), &
422 aero_data_n_spec(aero_data))
429 aero_mode%char_radius, aero_mode%log10_std_dev_radius, &
430 bin_grid, aero_data, vol_conc_total)
433 aero_mode%char_radius, bin_grid, aero_data, vol_conc_total)
436 aero_mode%char_radius, bin_grid, aero_data, vol_conc_total)
439 aero_mode%sample_num_conc, bin_grid, aero_data, vol_conc_total)
441 call die_msg(314169653,
"Unknown aero_mode type: " &
444 call assert_msg(756593082, sum(aero_mode%vol_frac_std) == 0d0, &
445 "cannot convert species fractions with non-zero standard deviation " &
446 //
"to binned distributions")
447 do i_spec = 1,aero_data_n_spec(aero_data)
448 vol_conc(:,i_spec) = vol_conc_total * aero_mode%vol_frac(i_spec)
464 real(kind=dp),
intent(out) :: weighted_num_conc(:)
467 real(kind=dp) :: x0, x1
471 size(weighted_num_conc) ==
size(aero_mode%sample_num_conc))
474 weighted_num_conc = aero_mode%sample_num_conc
477 do i_sample = 1,
size(aero_mode%sample_num_conc)
478 x0 = log(aero_mode%sample_radius(i_sample))
479 x1 = log(aero_mode%sample_radius(i_sample + 1))
480 weighted_num_conc(i_sample) = aero_mode%sample_num_conc(i_sample) &
481 / aero_weight%exponent * (exp(- aero_weight%exponent * x0) &
482 - exp(- aero_weight%exponent * x1)) / (x1 - x0)
485 call die_msg(576124393,
"unknown aero_weight type: " &
501 real(kind=dp) :: x_mean_prime
502 real(kind=dp),
allocatable :: weighted_num_conc(:)
510 x_mean_prime = log10(aero_mode%char_radius) &
511 - aero_weight%exponent * aero_mode%log10_std_dev_radius**2 &
514 * exp((x_mean_prime**2 - log10(aero_mode%char_radius)**2) &
515 / (2d0 * aero_mode%log10_std_dev_radius**2))
517 call die_msg(466668240,
"unknown aero_weight type: " &
525 "cannot use exponential modes with weighting")
530 aero_mode%char_radius)
532 allocate(weighted_num_conc(
size(aero_mode%sample_num_conc)))
536 deallocate(weighted_num_conc)
538 call die_msg(901140225,
"unknown aero_mode type: " &
557 real(kind=dp),
intent(out) :: radius
559 real(kind=dp) :: x_mean_prime, x0, x1, x, r, inv_nc0, inv_nc1, inv_nc
561 real(kind=dp),
allocatable :: weighted_num_conc(:)
565 x_mean_prime = log10(aero_mode%char_radius)
568 x_mean_prime = log10(aero_mode%char_radius) &
569 - aero_weight%exponent * aero_mode%log10_std_dev_radius**2 &
572 call die_msg(517376844,
"unknown aero_weight type: " &
576 aero_mode%log10_std_dev_radius)
578 allocate(weighted_num_conc(
size(aero_mode%sample_num_conc)))
582 deallocate(weighted_num_conc)
586 .and. (aero_weight%exponent == 0d0)))
then 587 x0 = log(aero_mode%sample_radius(i_sample))
588 x1 = log(aero_mode%sample_radius(i_sample + 1))
590 x = (1d0 - r) * x0 + r * x1
595 aero_mode%sample_radius(i_sample))
597 aero_mode%sample_radius(i_sample + 1))
599 inv_nc = (1d0 - r) * inv_nc0 + r * inv_nc1
602 call die_msg(769131141,
"unknown aero_weight type: " &
612 "cannot use exponential modes with weighting")
614 call die_msg(301787712,
"unknown aero_weight type: " &
618 radius = aero_mode%char_radius
620 call die_msg(749122931,
"Unknown aero_mode type: " &
634 real(kind=dp),
intent(in) :: total_vol
636 real(kind=dp),
intent(out) :: vols(size(aero_mode%vol_frac))
638 integer,
parameter :: AERO_MODE_MAX_SAMPLE_LOOPS = 1000000
641 real(kind=dp) :: offset
646 do i_sample = 1,aero_mode_max_sample_loops
651 offset = 1d0 - sum(vols)
652 vols = vols + offset * aero_mode%vol_frac
653 if (minval(vols) >= 0d0)
exit 655 if (i_sample == aero_mode_max_sample_loops)
then 656 call die_msg(549015143,
"Unable to sample non-negative volumes for " &
657 //
"mode: " // trim(aero_mode%name))
659 vols = vols / sum(vols) * total_vol
666 subroutine spec_file_read_vol_frac(file, aero_data, vol_frac, vol_frac_std)
673 real(kind=dp),
allocatable,
intent(inout) :: vol_frac(:)
675 real(kind=dp),
allocatable,
intent(inout) :: vol_frac_std(:)
677 integer :: n_species, species, i
678 character(len=SPEC_LINE_MAX_VAR_LEN),
allocatable :: species_name(:)
679 real(kind=dp),
allocatable :: species_data(:,:)
680 real(kind=dp) :: tot_vol_frac
728 n_species =
size(species_data, 1)
729 if (n_species < 1)
then 730 call die_msg(628123166,
'file ' // trim(file%name) &
731 //
' must contain at least one line of data')
733 if ((
size(species_data, 2) /= 1) .and. (
size(species_data, 2) /= 2))
then 734 call die_msg(427666881,
'each line in file ' // trim(file%name) &
735 //
' must contain exactly one or two data values')
739 if (
allocated(vol_frac))
deallocate(vol_frac)
740 if (
allocated(vol_frac_std))
deallocate(vol_frac_std)
747 if (species == 0)
then 748 call die_msg(775942501,
'unknown species ' // trim(species_name(i)) &
749 //
' in file ' // trim(file%name))
751 vol_frac(species) = species_data(i, 1)
752 if (
size(species_data, 2) == 2)
then 753 vol_frac_std(species) = species_data(i, 2)
758 vol_frac = vol_frac / aero_data%density
759 vol_frac_std = vol_frac_std / aero_data%density
762 tot_vol_frac = sum(vol_frac)
763 if ((minval(vol_frac) < 0d0) .or. (tot_vol_frac <= 0d0))
then 764 call die_msg(356648030,
'fractions in ' // trim(file%name) &
765 //
' are not positive')
767 if (minval(vol_frac_std) < 0d0)
then 768 call die_msg(676576501,
'standard deviations in ' // trim(file%name) &
769 //
' are not positive')
771 vol_frac = vol_frac / tot_vol_frac
772 vol_frac_std = vol_frac_std / tot_vol_frac
774 end subroutine spec_file_read_vol_frac
779 subroutine spec_file_read_size_dist(file, sample_radius, sample_num_conc)
784 real(kind=dp),
allocatable,
intent(inout) :: sample_radius(:)
786 real(kind=dp),
allocatable,
intent(inout) :: sample_num_conc(:)
788 character(len=SPEC_LINE_MAX_VAR_LEN),
allocatable :: names(:)
789 real(kind=dp),
allocatable :: data(:,:)
790 integer :: n_sample, i_sample
827 'must contain a line starting with "diam"')
829 n_sample =
size(
data,2) - 1
831 'must have at least two diam values')
833 if (
allocated(sample_radius))
deallocate(sample_radius)
834 allocate(sample_radius(n_sample + 1))
836 do i_sample = 1,n_sample
838 sample_radius(i_sample) < sample_radius(i_sample + 1), &
839 'diam values must be strictly increasing')
844 'must contain a line starting with "num_conc"')
848 'must have one fewer num_conc than diam values')
850 if (
allocated(sample_num_conc))
deallocate(sample_num_conc)
851 allocate(sample_num_conc(n_sample))
852 sample_num_conc =
data(1,:)
853 do i_sample = 1,n_sample
855 sample_num_conc(i_sample) >= 0d0, &
856 'num_conc values must be non-negative')
859 end subroutine spec_file_read_size_dist
865 subroutine spec_file_read_aero_mode(file, aero_data, aero_mode, eof)
876 character(len=SPEC_LINE_MAX_VAR_LEN) :: tmp_str, mode_type, diam_type_str
877 character(len=SPEC_LINE_MAX_VAR_LEN) :: mass_frac_filename
878 character(len=SPEC_LINE_MAX_VAR_LEN) :: size_dist_filename
879 type(spec_line_t) :: line
880 type(
spec_file_t) :: mass_frac_file, size_dist_file
881 real(kind=dp) :: diam, temp, pressure
882 integer :: diam_type, i_radius
883 real(kind=dp) :: log10_r0_mob, log10_r1_mob, log10_r2_mob, &
884 r0_mob, r2_mob, r0_geom, r2_geom, log10_r0_geom, log10_r2_geom
983 tmp_str = line%data(1)
989 call spec_file_read_vol_frac(mass_frac_file, aero_data, &
990 aero_mode%vol_frac, aero_mode%vol_frac_std)
994 if (trim(diam_type_str) ==
'geometric')
then 996 elseif (trim(diam_type_str) ==
'mobility')
then 1002 "Unknown diam_type: " // trim(diam_type_str))
1006 aero_mode%sample_radius = [
real(kind=dp) :: ]
1007 aero_mode%sample_num_conc = [
real(kind=dp) :: ]
1008 if (trim(mode_type) ==
'log_normal')
then 1013 aero_mode%log10_std_dev_radius)
1015 aero_mode%char_radius =
diam2rad(diam)
1017 aero_mode%char_radius &
1026 log10_r1_mob = log10(
diam2rad(diam))
1027 log10_r0_mob = log10_r1_mob - aero_mode%log10_std_dev_radius
1028 log10_r2_mob = log10_r1_mob + aero_mode%log10_std_dev_radius
1029 r0_mob = 10**log10_r0_mob
1030 r2_mob = 10**log10_r2_mob
1032 r0_mob, temp, pressure)
1034 r2_mob, temp, pressure)
1035 log10_r0_geom = log10(r0_geom)
1036 log10_r2_geom = log10(r2_geom)
1037 aero_mode%log10_std_dev_radius &
1038 = (log10_r2_geom - log10_r0_geom) / 2d0
1040 call die_msg(532966100,
"Diameter type not handled: " &
1043 elseif (trim(mode_type) ==
'exp')
then 1048 aero_mode%char_radius =
diam2rad(diam)
1050 aero_mode%char_radius &
1054 call die_msg(585104460,
"Diameter type not handled: " &
1057 elseif (trim(mode_type) ==
'mono')
then 1062 aero_mode%char_radius =
diam2rad(diam)
1064 aero_mode%char_radius &
1068 call die_msg(902864269,
"Diameter type not handled: " &
1071 elseif (trim(mode_type) ==
'sampled')
then 1075 call spec_file_read_size_dist(size_dist_file, &
1076 aero_mode%sample_radius, aero_mode%sample_num_conc)
1081 do i_radius = 1,
size(aero_mode%sample_radius)
1082 aero_mode%sample_radius(i_radius) &
1084 aero_mode%sample_radius(i_radius), temp, pressure)
1087 call die_msg(239088838,
"Diameter type not handled: " &
1092 "Unknown distribution mode type: " // trim(mode_type))
1096 end subroutine spec_file_read_aero_mode
1126 character,
intent(inout) :: buffer(:)
1128 integer,
intent(inout) :: position
1133 integer :: prev_position
1135 prev_position = position
1158 character,
intent(inout) :: buffer(:)
1160 integer,
intent(inout) :: position
1165 integer :: prev_position
1167 prev_position = position
real(kind=dp) function aero_weight_radius_at_num_conc(aero_weight, num_conc)
Compute the radius at a given number concentration (m). Inverse of aero_weight_num_conc_at_radius().
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
integer, parameter aero_mode_type_len
Maximum length of a mode type.
Random number generators.
An input file with extra data for printing messages.
real(kind=dp) function aero_mode_number(aero_mode, aero_weight)
Return the total number of computational particles for an aero_mode.
subroutine num_conc_log_normal(total_num_conc, geom_mean_radius, log10_sigma_g, bin_grid, num_conc)
Compute a log-normal distribution.
integer function pmc_mpi_pack_size_string(val)
Determines the number of bytes required to pack the given value.
subroutine aero_mode_vol_conc(aero_mode, bin_grid, aero_data, vol_conc)
Return the binned per-species volume concentration for an aero_mode.
subroutine pmc_mpi_unpack_string(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine spec_file_close(file)
Close a spec file.
The bin_grid_t structure and associated subroutines.
real(kind=dp) elemental function aero_data_rad2vol(aero_data, r)
Convert geometric radius (m) to mass-equivalent volume (m^3).
subroutine num_conc_exp(total_num_conc, radius_at_mean_vol, bin_grid, aero_data, num_conc)
Exponential distribution in volume.
integer, parameter aero_weight_type_power
Type code for power function weighting.
subroutine warn_msg(code, warning_msg, already_warned)
Prints a warning message.
subroutine num_conc_sampled(sample_radius, sample_num_conc, bin_grid, num_conc)
Sampled distribution, not normalized.
integer, parameter aero_mode_diam_type_geometric
Type code for geometric diameter.
subroutine pmc_mpi_pack_aero_mode(buffer, position, val)
Packs the given value into the buffer, advancing position.
real(kind=dp) function rand_normal(mean, stddev)
Generates a normally distributed random number with the given mean and standard deviation.
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
subroutine spec_file_check_line_length(file, line, length)
Check that the length of the line data is as given.
subroutine vol_conc_log_normal(total_num_conc, geom_mean_radius, log10_sigma_g, bin_grid, aero_data, vol_conc)
Compute a log-normal distribution in volume.
integer, parameter aero_mode_diam_type_mobility
Type code for mobility equivalent diameter.
real(kind=dp) function aero_mode_total_num_conc(aero_mode)
Returns the total number concentration of a mode. (#/m^3)
An aerosol size distribution weighting function.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
integer function aero_data_source_by_name(aero_data, name)
Returns the number of the source in aero_data with the given name, or adds the source if it doesn't e...
real(kind=dp) function pmc_random()
Returns a random number between 0 and 1.
subroutine spec_file_read_real(file, name, var)
Read a real number from a spec file that must have the given name.
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
An aerosol size distribution mode.
character(len=aero_mode_type_len) function aero_mode_type_to_string(type)
Return a string representation of a kernel type.
subroutine die_msg(code, error_msg)
Error immediately.
integer function sample_cts_pdf(pdf)
Sample the given continuous probability density function.
integer, parameter aero_mode_type_invalid
Type code for an undefined or invalid mode.
elemental integer function bin_grid_size(bin_grid)
Return the number of bins in the grid, or -1 if the bin grid is not allocated.
subroutine vol_conc_exp(total_num_conc, radius_at_mean_vol, bin_grid, aero_data, vol_conc)
Exponential distribution in volume.
subroutine pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine aero_mode_sample_vols(aero_mode, total_vol, vols)
Return an array of volumes randomly sampled from the volume fractions.
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
integer, parameter aero_mode_type_sampled
Type code for a sampled mode.
The aero_data_t structure and associated subroutines.
real(kind=dp) elemental function diam2rad(d)
Convert diameter (m) to radius (m).
1D grid, either logarithmic or linear.
subroutine aero_mode_weighted_sampled_num_conc(aero_mode, aero_weight, weighted_num_conc)
Compute weighted sampled number concentrations.
subroutine spec_file_read_line(file, line, eof)
Read a spec_line from the spec_file.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine rand_normal_array_1d(mean, stddev, val)
Generates a vector of normally distributed random numbers with the given means and standard deviation...
Reading formatted text input.
integer function bin_grid_find(bin_grid, val)
Find the bin number that contains a given value.
subroutine spec_file_check_name(file, name, read_name)
Checks that the read_name is the same as name.
subroutine pmc_mpi_pack_string(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine spec_file_open(filename, file)
Open a spec file for reading.
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
integer, parameter aero_mode_type_log_normal
Type code for a log-normal mode.
subroutine pmc_mpi_unpack_aero_mode(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
integer function pmc_mpi_pack_size_aero_mode(val)
Determines the number of bytes required to pack the given value.
integer, parameter aero_mode_type_mono
Type code for a mono-disperse mode.
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine spec_file_read_real_named_array(file, max_lines, names, vals)
Read an array of named lines with real data. All lines must have the same number of data elements...
subroutine num_conc_mono(total_num_conc, radius, bin_grid, num_conc)
Mono-disperse distribution. Normalized so that sum(num_conc(k) * log_width) = 1.
integer, parameter aero_mode_name_len
Maximum length of a mode name.
subroutine vol_to_lnr(r, f_vol, f_lnr)
Convert a concentration f(vol)d(vol) to f(ln(r))d(ln(r)) where vol = 4/3 pi r^3.
real(kind=dp) elemental function aero_data_vol2rad(aero_data, v)
Convert mass-equivalent volume (m^3) to geometric radius (m).
subroutine spec_file_check_line_name(file, line, name)
Check that the name of the line data is as given.
subroutine aero_mode_sample_radius(aero_mode, aero_data, aero_weight, radius)
Return a radius randomly sampled from the mode distribution.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
subroutine spec_file_assert_msg(code, file, condition_ok, msg)
Exit with an error message containing filename and line number if condition_ok is ...
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine spec_file_read_string(file, name, var)
Read a string from a spec file that must have a given name.
integer, parameter aero_mode_type_exp
Type code for an exponential mode.
type(const_t), save const
Fixed variable for accessing the constant's values.
integer, parameter aero_weight_type_none
Type code for no (or flat) weighting.
elemental integer function aero_data_n_spec(aero_data)
Return the number of aerosol species, or -1 if uninitialized.
integer, parameter aero_weight_type_mfa
Type code for MFA weighting.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
subroutine vol_conc_sampled(sample_radius, sample_num_conc, bin_grid, aero_data, vol_conc)
Sampled distribution in volume.
subroutine vol_conc_mono(total_num_conc, radius, bin_grid, aero_data, vol_conc)
Mono-disperse distribution in volume.
real(kind=dp) function aero_weight_num_conc_at_radius(aero_weight, radius)
Compute the number concentration at a given radius (m^{-3}).
Aerosol material properties and associated data.
subroutine aero_mode_num_conc(aero_mode, bin_grid, aero_data, num_conc)
Return the binned number concentration for an aero_mode.
Common utility subroutines.
Wrapper functions for MPI.
real(kind=dp) function aero_data_mobility_rad_to_geometric_rad(aero_data, mobility_rad, temp, pressure)
Convert mobility equivalent radius (m) to geometric radius (m^3).
The aero_weight_t structure and associated subroutines.
integer function aero_data_spec_by_name(aero_data, name)
Returns the number of the species in aero_data with the given name, or returns 0 if there is no such ...
The aero_mode_t structure and associated subroutines.
integer, parameter aero_mode_diam_type_invalid
Type code for an undefined for invalid diameter type.