24 integer,
parameter :: AERO_MODE_NAME_LEN = 300
26 integer,
parameter :: AERO_MODE_TYPE_LEN = 20
29 integer,
parameter :: AERO_MODE_TYPE_INVALID = 0
31 integer,
parameter :: AERO_MODE_TYPE_LOG_NORMAL = 1
33 integer,
parameter :: AERO_MODE_TYPE_EXP = 2
35 integer,
parameter :: AERO_MODE_TYPE_MONO = 3
37 integer,
parameter :: AERO_MODE_TYPE_SAMPLED = 4
49 character(len=AERO_MODE_NAME_LEN) :: name
53 real(kind=dp) :: char_radius
55 real(kind=dp) :: log10_std_dev_radius
57 real(kind=dp),
pointer :: sample_radius(:)
59 real(kind=dp),
pointer :: sample_num_conc(:)
61 real(kind=dp) :: num_conc
63 real(kind=dp),
pointer :: vol_frac(:)
65 real(kind=dp),
pointer :: vol_frac_std(:)
78 integer,
intent(in) :: type
80 if (
type == aero_mode_type_invalid)
then
82 elseif (
type == aero_mode_type_log_normal)
then
84 elseif (
type == aero_mode_type_exp)
then
86 elseif (
type == aero_mode_type_mono)
then
88 elseif (
type == aero_mode_type_sampled)
then
104 allocate(aero_mode%vol_frac(0))
105 allocate(aero_mode%vol_frac_std(0))
106 allocate(aero_mode%sample_radius(0))
107 allocate(aero_mode%sample_num_conc(0))
119 integer,
intent(in) :: n_spec
121 allocate(aero_mode%vol_frac(n_spec))
122 allocate(aero_mode%vol_frac_std(n_spec))
123 allocate(aero_mode%sample_radius(0))
124 allocate(aero_mode%sample_num_conc(0))
136 deallocate(aero_mode%vol_frac)
137 deallocate(aero_mode%vol_frac_std)
138 deallocate(aero_mode%sample_radius)
139 deallocate(aero_mode%sample_num_conc)
155 aero_mode_to%name = aero_mode_from%name
156 aero_mode_to%type = aero_mode_from%type
157 aero_mode_to%char_radius = aero_mode_from%char_radius
158 aero_mode_to%log10_std_dev_radius = aero_mode_from%log10_std_dev_radius
159 aero_mode_to%num_conc = aero_mode_from%num_conc
160 aero_mode_to%vol_frac = aero_mode_from%vol_frac
161 aero_mode_to%vol_frac_std = aero_mode_from%vol_frac_std
162 aero_mode_to%source = aero_mode_from%source
163 deallocate(aero_mode_to%sample_radius)
164 deallocate(aero_mode_to%sample_num_conc)
165 allocate(aero_mode_to%sample_radius(
size(aero_mode_from%sample_radius)))
166 allocate(aero_mode_to%sample_num_conc( &
167 size(aero_mode_from%sample_num_conc)))
168 aero_mode_to%sample_radius = aero_mode_from%sample_radius
169 aero_mode_to%sample_num_conc = aero_mode_from%sample_num_conc
182 if ((aero_mode%type == aero_mode_type_log_normal) &
183 .or. (aero_mode%type == aero_mode_type_exp) &
184 .or. (aero_mode%type == aero_mode_type_mono))
then
186 elseif (aero_mode%type == aero_mode_type_sampled)
then
189 call
die_msg(719625922,
"unknown aero_mode type: " &
199 log10_sigma_g, bin_grid, num_conc)
202 real(kind=dp),
intent(in) :: total_num_conc
204 real(kind=dp),
intent(in) :: geom_mean_radius
206 real(kind=dp),
intent(in) :: log10_sigma_g
210 real(kind=dp),
intent(out) :: num_conc(bin_grid%n_bin)
214 do k = 1,bin_grid%n_bin
215 num_conc(k) = total_num_conc / (sqrt(2d0 * const%pi) &
216 * log10_sigma_g) * dexp(-(dlog10(bin_grid%centers(k)) &
217 - dlog10(geom_mean_radius))**2d0 &
218 / (2d0 * log10_sigma_g**2d0)) / dlog(10d0)
233 log10_sigma_g, bin_grid, vol_conc)
236 real(kind=dp),
intent(in) :: total_num_conc
238 real(kind=dp),
intent(in) :: geom_mean_radius
240 real(kind=dp),
intent(in) :: log10_sigma_g
244 real(kind=dp),
intent(out) :: vol_conc(bin_grid%n_bin)
246 real(kind=dp) :: num_conc(bin_grid%n_bin)
249 log10_sigma_g, bin_grid, num_conc)
250 vol_conc = num_conc *
rad2vol(bin_grid%centers)
260 subroutine num_conc_exp(total_num_conc, radius_at_mean_vol, bin_grid, &
264 real(kind=dp),
intent(in) :: total_num_conc
266 real(kind=dp),
intent(in) :: radius_at_mean_vol
270 real(kind=dp),
intent(out) :: num_conc(bin_grid%n_bin)
273 real(kind=dp) :: mean_vol, num_conc_vol
275 mean_vol =
rad2vol(radius_at_mean_vol)
276 do k = 1,bin_grid%n_bin
277 num_conc_vol = total_num_conc / mean_vol &
278 * exp(-(
rad2vol(bin_grid%centers(k)) / mean_vol))
279 call
vol_to_lnr(bin_grid%centers(k), num_conc_vol, num_conc(k))
287 subroutine vol_conc_exp(total_num_conc, radius_at_mean_vol, bin_grid, &
291 real(kind=dp),
intent(in) :: total_num_conc
293 real(kind=dp),
intent(in) :: radius_at_mean_vol
297 real(kind=dp),
intent(out) :: vol_conc(bin_grid%n_bin)
299 real(kind=dp) :: num_conc(bin_grid%n_bin)
301 call
num_conc_exp(total_num_conc, radius_at_mean_vol, bin_grid, num_conc)
302 vol_conc = num_conc *
rad2vol(bin_grid%centers)
313 real(kind=dp),
intent(in) :: total_num_conc
315 real(kind=dp),
intent(in) :: radius
319 real(kind=dp),
intent(out) :: num_conc(bin_grid%n_bin)
325 if ((k < 1) .or. (k > bin_grid%n_bin))
then
326 call
warn_msg(825666877,
"monodisperse radius outside of bin_grid")
328 num_conc(k) = total_num_conc / bin_grid%widths(k)
339 real(kind=dp),
intent(in) :: total_num_conc
341 real(kind=dp),
intent(in) :: radius
345 real(kind=dp),
intent(out) :: vol_conc(bin_grid%n_bin)
351 if ((k < 1) .or. (k > bin_grid%n_bin))
then
352 call
warn_msg(420930707,
"monodisperse radius outside of bin_grid")
354 vol_conc(k) = total_num_conc / bin_grid%widths(k) &
367 real(kind=dp),
intent(in) :: sample_radius(:)
369 real(kind=dp),
intent(in) :: sample_num_conc(:)
373 real(kind=dp),
intent(out) :: num_conc(bin_grid%n_bin)
375 integer :: i_sample, n_sample, i_lower, i_upper, i_bin
376 real(kind=dp) :: r_lower, r_upper
377 real(kind=dp) :: r_bin_lower, r_bin_upper, r1, r2, ratio
379 n_sample =
size(sample_num_conc)
380 call
assert(188766208,
size(sample_radius) == n_sample + 1)
381 call
assert(295384037, n_sample >= 1)
384 do i_sample = 1,n_sample
385 r_lower = sample_radius(i_sample)
386 r_upper = sample_radius(i_sample + 1)
389 if (i_upper < 1) cycle
390 if (i_lower > bin_grid%n_bin) cycle
391 i_lower = max(1, i_lower)
392 i_upper = min(bin_grid%n_bin, i_upper)
393 do i_bin = i_lower,i_upper
394 r_bin_lower = bin_grid%edges(i_bin)
395 r_bin_upper = bin_grid%edges(i_bin + 1)
396 r1 = max(r_lower, r_bin_lower)
397 r2 = min(r_upper, r_bin_upper)
398 ratio = (log(r2) - log(r1)) / (log(r_upper) - log(r_lower))
399 num_conc(i_bin) = num_conc(i_bin) + ratio &
400 * sample_num_conc(i_sample) / bin_grid%widths(i_bin)
413 real(kind=dp),
intent(in) :: sample_radius(:)
415 real(kind=dp),
intent(in) :: sample_num_conc(:)
419 real(kind=dp),
intent(out) :: vol_conc(bin_grid%n_bin)
421 real(kind=dp) :: num_conc(bin_grid%n_bin)
424 vol_conc = num_conc *
rad2vol(bin_grid%centers)
441 real(kind=dp),
intent(out) :: num_conc(bin_grid%n_bin)
443 if (aero_mode%type == aero_mode_type_log_normal)
then
445 aero_mode%log10_std_dev_radius, bin_grid, num_conc)
446 elseif (aero_mode%type == aero_mode_type_exp)
then
447 call
num_conc_exp(aero_mode%num_conc, aero_mode%char_radius, &
449 elseif (aero_mode%type == aero_mode_type_mono)
then
450 call
num_conc_mono(aero_mode%num_conc, aero_mode%char_radius, &
452 elseif (aero_mode%type == aero_mode_type_sampled)
then
454 aero_mode%sample_num_conc, bin_grid, num_conc)
456 call
die_msg(223903246,
"unknown aero_mode type: " &
476 real(kind=dp),
intent(out) :: vol_conc(bin_grid%n_bin, aero_data%n_spec)
479 real(kind=dp) :: vol_conc_total(bin_grid%n_bin)
481 if (aero_mode%type == aero_mode_type_log_normal)
then
483 aero_mode%log10_std_dev_radius, bin_grid, vol_conc_total)
484 elseif (aero_mode%type == aero_mode_type_exp)
then
485 call
vol_conc_exp(aero_mode%num_conc, aero_mode%char_radius, &
486 bin_grid, vol_conc_total)
487 elseif (aero_mode%type == aero_mode_type_mono)
then
488 call
vol_conc_mono(aero_mode%num_conc, aero_mode%char_radius, &
489 bin_grid, vol_conc_total)
490 elseif (aero_mode%type == aero_mode_type_sampled)
then
492 aero_mode%sample_num_conc, bin_grid, vol_conc_total)
494 call
die_msg(314169653,
"Unknown aero_mode type: " &
497 call
assert_msg(756593082, sum(aero_mode%vol_frac_std) == 0d0, &
498 "cannot convert species fractions with non-zero standard deviation " &
499 //
"to binned distributions")
500 do i_spec = 1,aero_data%n_spec
501 vol_conc(:,i_spec) = vol_conc_total * aero_mode%vol_frac(i_spec)
517 real(kind=dp),
intent(out) :: weighted_num_conc(:)
520 real(kind=dp) :: x0, x1
522 call
assert(256667423, aero_mode%type == aero_mode_type_sampled)
524 size(weighted_num_conc) ==
size(aero_mode%sample_num_conc))
526 if (aero_weight%type == aero_weight_type_none)
then
527 weighted_num_conc = aero_mode%sample_num_conc
528 elseif ((aero_weight%type == aero_weight_type_power) &
529 .or. (aero_weight%type == aero_weight_type_mfa))
then
530 do i_sample = 1,
size(aero_mode%sample_num_conc)
531 x0 = log(aero_mode%sample_radius(i_sample))
532 x1 = log(aero_mode%sample_radius(i_sample + 1))
533 weighted_num_conc(i_sample) = aero_mode%sample_num_conc(i_sample) &
534 / aero_weight%exponent * (exp(- aero_weight%exponent * x0) &
535 - exp(- aero_weight%exponent * x1)) / (x1 - x0)
538 call
die_msg(576124393,
"unknown aero_weight type: " &
554 real(kind=dp) :: x_mean_prime
555 real(kind=dp),
allocatable :: weighted_num_conc(:)
558 if (aero_mode%type == aero_mode_type_log_normal)
then
559 if (aero_weight%type == aero_weight_type_none)
then
561 elseif ((aero_weight%type == aero_weight_type_power) &
562 .or. (aero_weight%type == aero_weight_type_mfa))
then
563 x_mean_prime = log10(aero_mode%char_radius) &
564 - aero_weight%exponent * aero_mode%log10_std_dev_radius**2 &
567 * exp((x_mean_prime**2 - log10(aero_mode%char_radius)**2) &
568 / (2d0 * aero_mode%log10_std_dev_radius**2))
570 call
die_msg(466668240,
"unknown aero_weight type: " &
573 elseif (aero_mode%type == aero_mode_type_exp)
then
574 if (aero_weight%type == aero_weight_type_none)
then
578 "cannot use exponential modes with weighting")
580 elseif (aero_mode%type == aero_mode_type_mono)
then
583 aero_mode%char_radius)
584 elseif (aero_mode%type == aero_mode_type_sampled)
then
585 allocate(weighted_num_conc(
size(aero_mode%sample_num_conc)))
589 deallocate(weighted_num_conc)
591 call
die_msg(901140225,
"unknown aero_mode type: " &
607 real(kind=dp),
intent(out) :: radius
609 real(kind=dp) :: x_mean_prime, x0, x1, x, r, inv_nc0, inv_nc1, inv_nc
611 real(kind=dp),
allocatable :: weighted_num_conc(:)
613 if (aero_mode%type == aero_mode_type_log_normal)
then
614 if (aero_weight%type == aero_weight_type_none)
then
615 x_mean_prime = log10(aero_mode%char_radius)
616 elseif ((aero_weight%type == aero_weight_type_power) &
617 .or. (aero_weight%type == aero_weight_type_mfa))
then
618 x_mean_prime = log10(aero_mode%char_radius) &
619 - aero_weight%exponent * aero_mode%log10_std_dev_radius**2 &
622 call
die_msg(517376844,
"unknown aero_weight type: " &
626 aero_mode%log10_std_dev_radius)
627 elseif (aero_mode%type == aero_mode_type_sampled)
then
628 allocate(weighted_num_conc(
size(aero_mode%sample_num_conc)))
632 deallocate(weighted_num_conc)
633 if ((aero_weight%type == aero_weight_type_none) &
634 .or. (((aero_weight%type == aero_weight_type_power) &
635 .or. (aero_weight%type == aero_weight_type_mfa)) &
636 .and. (aero_weight%exponent == 0d0)))
then
637 x0 = log(aero_mode%sample_radius(i_sample))
638 x1 = log(aero_mode%sample_radius(i_sample + 1))
640 x = (1d0 - r) * x0 + r * x1
642 elseif ((aero_weight%type == aero_weight_type_power) &
643 .or. (aero_weight%type == aero_weight_type_mfa))
then
645 aero_mode%sample_radius(i_sample))
647 aero_mode%sample_radius(i_sample + 1))
649 inv_nc = (1d0 - r) * inv_nc0 + r * inv_nc1
652 call
die_msg(769131141,
"unknown aero_weight type: " &
655 elseif (aero_mode%type == aero_mode_type_exp)
then
656 if (aero_weight%type == aero_weight_type_none)
then
659 elseif ((aero_weight%type == aero_weight_type_power) &
660 .or. (aero_weight%type == aero_weight_type_mfa))
then
662 "cannot use exponential modes with weighting")
664 call
die_msg(301787712,
"unknown aero_weight type: " &
667 elseif (aero_mode%type == aero_mode_type_mono)
then
668 radius = aero_mode%char_radius
670 call
die_msg(749122931,
"Unknown aero_mode type: " &
684 real(kind=dp),
intent(in) :: total_vol
686 real(kind=dp),
intent(out) :: vols(size(aero_mode%vol_frac))
689 vols = max(vols, 0d0)
690 vols = vols / sum(vols) * total_vol
697 subroutine spec_file_read_vol_frac(file, aero_data, vol_frac, vol_frac_std)
704 real(kind=dp),
intent(inout) :: vol_frac(:)
706 real(kind=dp),
intent(inout) :: vol_frac_std(:)
708 integer :: n_species, species, i
709 character(len=SPEC_LINE_MAX_VAR_LEN),
pointer :: species_name(:)
710 real(kind=dp),
pointer :: species_data(:,:)
711 real(kind=dp) :: tot_vol_frac
755 allocate(species_name(0))
756 allocate(species_data(0,0))
761 n_species =
size(species_data, 1)
762 if (n_species < 1)
then
763 call
die_msg(628123166,
'file ' // trim(file%name) &
764 //
' must contain at least one line of data')
766 if ((
size(species_data, 2) /= 1) .and. (
size(species_data, 2) /= 2))
then
767 call
die_msg(427666881,
'each line in file ' // trim(file%name) &
768 //
' must contain exactly one or two data values')
776 if (species == 0)
then
777 call
die_msg(775942501,
'unknown species ' // trim(species_name(i)) &
778 //
' in file ' // trim(file%name))
780 vol_frac(species) = species_data(i, 1)
781 if (
size(species_data, 2) == 2)
then
782 vol_frac_std(species) = species_data(i, 2)
785 deallocate(species_name)
786 deallocate(species_data)
789 vol_frac = vol_frac / aero_data%density
790 vol_frac_std = vol_frac_std / aero_data%density
793 tot_vol_frac = sum(vol_frac)
794 if ((minval(vol_frac) < 0d0) .or. (tot_vol_frac <= 0d0))
then
795 call
die_msg(356648030,
'fractions in ' // trim(file%name) &
796 //
' are not positive')
798 if (minval(vol_frac_std) < 0d0)
then
799 call
die_msg(676576501,
'standard deviations in ' // trim(file%name) &
800 //
' are not positive')
802 vol_frac = vol_frac / tot_vol_frac
803 vol_frac_std = vol_frac_std / tot_vol_frac
805 end subroutine spec_file_read_vol_frac
810 subroutine spec_file_read_size_dist(file, sample_radius, sample_num_conc)
815 real(kind=dp),
pointer :: sample_radius(:)
817 real(kind=dp),
pointer :: sample_num_conc(:)
819 character(len=SPEC_LINE_MAX_VAR_LEN),
pointer :: names(:)
820 real(kind=dp),
pointer :: data(:,:)
821 integer :: n_sample, i_sample
860 n_sample =
size(
data,2) - 1
862 'must have at least two diam values')
864 deallocate(sample_radius)
865 allocate(sample_radius(n_sample + 1))
867 do i_sample = 1,n_sample
869 sample_radius(i_sample) < sample_radius(i_sample + 1), &
870 'diam values must be strictly increasing')
877 'must have one fewer num_conc than diam values')
879 deallocate(sample_num_conc)
880 allocate(sample_num_conc(n_sample))
881 sample_num_conc =
data(1,:)
882 do i_sample = 1,n_sample
884 sample_num_conc(i_sample) >= 0d0, &
885 'num_conc values must be non-negative')
888 end subroutine spec_file_read_size_dist
894 subroutine spec_file_read_aero_mode(file, aero_data, aero_mode, eof)
905 character(len=SPEC_LINE_MAX_VAR_LEN) :: tmp_str, mode_type
906 character(len=SPEC_LINE_MAX_VAR_LEN) :: mass_frac_filename
907 character(len=SPEC_LINE_MAX_VAR_LEN) :: size_dist_filename
908 type(spec_line_t
) :: line
909 type(spec_file_t) :: mass_frac_file, size_dist_file
910 real(kind=dp) :: diam
1001 tmp_str = line%data(1)
1002 aero_mode%name = tmp_str(1:aero_mode_name_len)
1006 call spec_file_read_vol_frac(mass_frac_file, aero_data, &
1007 aero_mode%vol_frac, aero_mode%vol_frac_std)
1010 if (trim(mode_type) ==
'log_normal')
then
1011 aero_mode%type = aero_mode_type_log_normal
1014 aero_mode%char_radius =
diam2rad(diam)
1016 aero_mode%log10_std_dev_radius)
1017 elseif (trim(mode_type) ==
'exp')
then
1018 aero_mode%type = aero_mode_type_exp
1021 aero_mode%char_radius =
diam2rad(diam)
1022 elseif (trim(mode_type) ==
'mono')
then
1023 aero_mode%type = aero_mode_type_mono
1026 aero_mode%char_radius =
diam2rad(diam)
1027 elseif (trim(mode_type) ==
'sampled')
then
1028 aero_mode%type = aero_mode_type_sampled
1031 call spec_file_read_size_dist(size_dist_file, &
1032 aero_mode%sample_radius, aero_mode%sample_num_conc)
1036 "Unknown distribution mode type: " // trim(mode_type))
1041 end subroutine spec_file_read_aero_mode
1071 character,
intent(inout) :: buffer(:)
1073 integer,
intent(inout) :: position
1078 integer :: prev_position
1080 prev_position = position
1103 character,
intent(inout) :: buffer(:)
1105 integer,
intent(inout) :: position
1110 integer :: prev_position
1112 prev_position = position
real(kind=dp) function rand_normal(mean, stddev)
Generates a normally distributed random number with the given mean and standard deviation.
real(kind=dp) function aero_weight_num_conc_at_radius(aero_weight, radius)
Compute the number concentration at a given radius (m^{-3}).
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
subroutine num_conc_sampled(sample_radius, sample_num_conc, bin_grid, num_conc)
Sampled distribution, not normalized.
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
subroutine aero_mode_weighted_sampled_num_conc(aero_mode, aero_weight, weighted_num_conc)
Compute weighted sampled number concentrations.
An input file with extra data for printing messages.
subroutine spec_file_check_line_name(file, line, name)
Check that the name of the line data is as given.
The aero_data_t structure and associated subroutines.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
subroutine die_msg(code, error_msg)
Error immediately.
subroutine aero_mode_sample_vols(aero_mode, total_vol, vols)
Return an array of volumes randomly sampled from the volume fractions.
subroutine pmc_mpi_unpack_aero_mode(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
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) elemental function rad2vol(r)
Convert radius (m) to volume (m^3).
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine num_conc_log_normal(total_num_conc, geom_mean_radius, log10_sigma_g, bin_grid, num_conc)
Compute a log-normal distribution.
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.
The aero_weight_t structure and associated subroutines.
subroutine spec_file_open(filename, file)
Open a spec file for reading.
subroutine num_conc_exp(total_num_conc, radius_at_mean_vol, bin_grid, num_conc)
Exponential distribution in volume.
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 ...
subroutine spec_file_assert_msg(code, file, condition_ok, msg)
Exit with an error message containing filename and line number if condition_ok is ...
integer function sample_cts_pdf(pdf)
Sample the given continuous probability density function.
integer function bin_grid_find(bin_grid, val)
Find the bin number that contains a given value.
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
subroutine vol_conc_log_normal(total_num_conc, geom_mean_radius, log10_sigma_g, bin_grid, vol_conc)
Compute a log-normal distribution in volume.
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_pack_string(buffer, position, val)
Packs the given value into the buffer, advancing position.
integer function pmc_mpi_pack_size_string(val)
Determines the number of bytes required to pack the given value.
An aerosol size distribution weighting function.
Random number generators.
integer function pmc_mpi_pack_size_aero_mode(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.
real(kind=dp) function aero_mode_number(aero_mode, aero_weight)
Return the total number of computational particles for an aero_mode.
subroutine aero_mode_num_conc(aero_mode, bin_grid, aero_data, num_conc)
Return the binned number concentration for an aero_mode.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine aero_mode_allocate(aero_mode)
Allocates an aero_mode.
real(kind=dp) function aero_mode_total_num_conc(aero_mode)
Returns the total number concentration of a mode. (#/m^3)
real(kind=dp) elemental function diam2rad(d)
Convert diameter (m) to radius (m).
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 spec_file_check_name(file, name, read_name)
Checks that the read_name is the same as name.
Common utility subroutines.
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
character(len=aero_mode_type_len) function aero_mode_type_to_string(type)
Return a string representation of a kernel type.
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.
An aerosol size distribution mode.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine aero_mode_copy(aero_mode_from, aero_mode_to)
Copy an aero_mode.
Wrapper functions for MPI.
subroutine spec_file_read_string(file, name, var)
Read a string from a spec file that must have a given name.
subroutine rand_normal_array_1d(mean, stddev, val)
Generates a vector of normally distributed random numbers with the given means and standard deviation...
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
subroutine aero_mode_sample_radius(aero_mode, aero_weight, radius)
Return a radius randomly sampled from the mode distribution.
The aero_mode_t structure and associated subroutines.
1D grid, either logarithmic or linear.
The bin_grid_t structure and associated subroutines.
subroutine pmc_mpi_unpack_string(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Reading formatted text input.
subroutine spec_file_close(file)
Close a spec file.
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 spec_line_allocate(spec_line)
Allocates memory for a spec_line.
subroutine warn_msg(code, warning_msg, already_warned)
Prints a warning message.
subroutine spec_file_read_line(file, line, eof)
Read a spec_line from the spec_file.
subroutine pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
real(kind=dp) function pmc_random()
Returns a random number between 0 and 1.
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine spec_line_deallocate(spec_line)
Frees all storage.
real(kind=dp) elemental function vol2rad(v)
Convert volume (m^3) to radius (m).
subroutine pmc_mpi_pack_aero_mode(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine aero_mode_allocate_size(aero_mode, n_spec)
Allocates an aero_mode of the given size.
subroutine spec_file_read_real(file, name, var)
Read a real number from a spec file that must have the given name.
Aerosol material properties and associated data.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine vol_conc_mono(total_num_conc, radius, bin_grid, vol_conc)
Mono-disperse distribution in volume.
subroutine vol_conc_exp(total_num_conc, radius_at_mean_vol, bin_grid, vol_conc)
Exponential distribution in volume.
subroutine aero_mode_deallocate(aero_mode)
Free all storage.
subroutine vol_conc_sampled(sample_radius, sample_num_conc, bin_grid, vol_conc)
Sampled distribution in volume.