23 integer,
parameter :: AERO_WEIGHT_TYPE_INVALID = 0
25 integer,
parameter :: AERO_WEIGHT_TYPE_NONE = 1
27 integer,
parameter :: AERO_WEIGHT_TYPE_POWER = 2
29 integer,
parameter :: AERO_WEIGHT_TYPE_MFA = 3
36 real(kind=dp) :: magnitude
38 real(kind=dp) :: exponent
85 aero_weight%type = aero_weight_type_invalid
86 aero_weight%magnitude = 1d0
87 aero_weight%exponent = 0d0
99 aero_weight%magnitude = 1d0
113 aero_weight_to%type = aero_weight_from%type
114 aero_weight_to%magnitude = aero_weight_from%magnitude
115 aero_weight_to%exponent = aero_weight_from%exponent
128 real(kind=dp),
intent(in) :: factor
130 aero_weight%magnitude = aero_weight%magnitude * factor
144 aero_weight%magnitude =
harmonic_mean(aero_weight%magnitude, &
145 aero_weight_delta%magnitude)
155 aero_weight_to, sample_prop, overwrite_to)
162 real(kind=dp),
intent(in) :: sample_prop
164 logical,
intent(in),
optional :: overwrite_to
166 real(kind=dp) :: magnitude_transfer
167 logical :: use_overwrite_to
169 magnitude_transfer = aero_weight_from%magnitude / sample_prop
170 aero_weight_from%magnitude =
harmonic_mean(aero_weight_from%magnitude, &
171 - magnitude_transfer)
172 use_overwrite_to = .false.
173 if (present(overwrite_to))
then
174 if (overwrite_to)
then
175 use_overwrite_to = .true.
178 if (use_overwrite_to)
then
179 aero_weight_to%magnitude = magnitude_transfer
181 aero_weight_to%magnitude =
harmonic_mean(aero_weight_to%magnitude, &
195 real(kind=dp),
intent(in) :: radius
198 if (aero_weight%type == aero_weight_type_none)
then
200 elseif ((aero_weight%type == aero_weight_type_power) &
201 .or. (aero_weight%type == aero_weight_type_mfa))
then
203 * radius**aero_weight%exponent
205 call
die_msg(700421478,
"unknown aero_weight type: " &
235 real(kind=dp),
intent(in) :: num_conc
238 if (aero_weight%type == aero_weight_type_none)
then
239 call
die_msg(545688537,
"cannot invert weight type 'none'")
240 elseif ((aero_weight%type == aero_weight_type_power) &
241 .or. (aero_weight%type == aero_weight_type_mfa))
then
242 call
assert_msg(902242996, aero_weight%exponent /= 0d0, &
243 "cannot invert weight with zero exponent")
245 = exp(log(num_conc / aero_weight%magnitude) / aero_weight%exponent)
247 call
die_msg(307766845,
"unknown aero_weight type: " &
262 call
assert(269952052, (aero_weight%type == aero_weight_type_none) &
263 .or. (aero_weight%type == aero_weight_type_power) &
264 .or. (aero_weight%type == aero_weight_type_mfa))
265 if (aero_weight%type == aero_weight_type_none)
then
266 call
assert(853998284, aero_weight%exponent == 0d0)
268 if (aero_weight%type == aero_weight_type_mfa)
then
269 call
assert(829651126, aero_weight%exponent == -3d0)
279 monotone_increasing, monotone_decreasing)
284 logical,
intent(out) :: monotone_increasing
286 logical,
intent(out) :: monotone_decreasing
289 call
assert(610698264, (aero_weight%type == aero_weight_type_none) &
290 .or. (aero_weight%type == aero_weight_type_power) &
291 .or. (aero_weight%type == aero_weight_type_mfa))
292 monotone_increasing = .true.
293 monotone_decreasing = .true.
294 if (aero_weight%type == aero_weight_type_power)
then
295 if (aero_weight%exponent < 0d0)
then
296 monotone_increasing = .false.
298 if (aero_weight%exponent > 0d0)
then
299 monotone_decreasing = .false.
302 if (aero_weight%type == aero_weight_type_mfa)
then
303 monotone_increasing = .false.
311 subroutine spec_file_read_aero_weight(file, aero_weight)
318 character(len=SPEC_LINE_MAX_VAR_LEN) :: weight_type
353 aero_weight%magnitude = 0d0
354 if (trim(weight_type) ==
'none')
then
355 aero_weight%type = aero_weight_type_none
356 elseif (trim(weight_type) ==
'power')
then
357 aero_weight%type = aero_weight_type_power
359 elseif (trim(weight_type) ==
'mfa')
then
360 aero_weight%type = aero_weight_type_mfa
361 aero_weight%exponent = -3d0
364 // trim(weight_type))
367 end subroutine spec_file_read_aero_weight
390 character,
intent(inout) :: buffer(:)
392 integer,
intent(inout) :: position
397 integer :: prev_position
399 prev_position = position
415 character,
intent(inout) :: buffer(:)
417 integer,
intent(inout) :: position
422 integer :: prev_position
424 prev_position = position
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 pmc_mpi_unpack_aero_weight(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
An input file with extra data for printing messages.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
elemental subroutine aero_weight_normalize(aero_weight)
Sets the aero_weight to a non-zero normalized value.
subroutine die_msg(code, error_msg)
Error immediately.
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
elemental subroutine aero_weight_copy(aero_weight_from, aero_weight_to)
Copy an aero_weight.
The aero_weight_t structure and associated subroutines.
The aero_particle_t structure and associated subroutines.
elemental subroutine aero_weight_allocate(aero_weight)
Allocates an aero_weight.
elemental subroutine aero_weight_shift(aero_weight_from, aero_weight_to, sample_prop, overwrite_to)
Adjust source and destination weights to reflect moving sample_prop proportion of particles from aero...
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
subroutine aero_weight_check_valid_exponent(aero_weight)
Ensures that a weight function exponent is valid.
elemental subroutine aero_weight_scale(aero_weight, factor)
Scale the weight by the given fraction, so new_weight = old_weight * factor.
elemental real(kind=dp) function aero_particle_radius(aero_particle)
Total radius of the particle (m).
An aerosol size distribution weighting function.
Random number generators.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
integer function pmc_mpi_pack_size_aero_weight(val)
Determines the number of bytes required to pack the given value.
elemental subroutine aero_weight_combine(aero_weight, aero_weight_delta)
Combine aero_weight_delta into aero_weight with a harmonic mean.
subroutine pmc_mpi_pack_aero_weight(buffer, position, val)
Packs the given value into the buffer, advancing position.
Common utility subroutines.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
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.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
Single aerosol particle data structure.
Reading formatted text input.
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().
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
subroutine pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine aero_weight_check_monotonicity(aero_weight, monotone_increasing, monotone_decreasing)
Determine whether a weight function is monotone increasing, monotone decreasing, or neither...
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
elemental real(kind=dp) function harmonic_mean(x1, x2)
Compute the harmonic mean of two numbers.
elemental subroutine aero_weight_allocate_size(aero_weight)
Allocates an aero_weight of the given size.
real(kind=dp) function aero_weight_num_conc(aero_weight, aero_particle)
Compute the number concentration for a particle (m^{-3}).
elemental subroutine aero_weight_zero(aero_weight)
Zeros the contents of the aero_weight.
subroutine spec_file_read_real(file, name, var)
Read a real number from a spec file that must have the given name.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
elemental subroutine aero_weight_deallocate(aero_weight)
Free all storage.