25 integer,
parameter :: SCENARIO_LOSS_FUNCTION_INVALID = 0
27 integer,
parameter :: SCENARIO_LOSS_FUNCTION_ZERO = 1
29 integer,
parameter :: SCENARIO_LOSS_FUNCTION_CONSTANT = 2
31 integer,
parameter :: SCENARIO_LOSS_FUNCTION_VOLUME = 3
33 integer,
parameter :: SCENARIO_LOSS_FUNCTION_DRY_DEP = 4
38 real(kind=dp),
parameter :: alg_threshold = 1.0d0
53 real(kind=dp),
pointer :: temp_time(:)
55 real(kind=dp),
pointer :: temp(:)
58 real(kind=dp),
pointer :: pressure_time(:)
60 real(kind=dp),
pointer :: pressure(:)
63 real(kind=dp),
pointer :: height_time(:)
65 real(kind=dp),
pointer :: height(:)
68 real(kind=dp),
pointer :: gas_emission_time(:)
70 real(kind=dp),
pointer :: gas_emission_rate_scale(:)
75 real(kind=dp),
pointer :: gas_dilution_time(:)
77 real(kind=dp),
pointer :: gas_dilution_rate(:)
82 real(kind=dp),
pointer :: aero_emission_time(:)
84 real(kind=dp),
pointer :: aero_emission_rate_scale(:)
89 real(kind=dp),
pointer :: aero_dilution_time(:)
91 real(kind=dp),
pointer :: aero_dilution_rate(:)
106 allocate(scenario%temp_time(0))
107 allocate(scenario%temp(0))
109 allocate(scenario%pressure_time(0))
110 allocate(scenario%pressure(0))
112 allocate(scenario%height_time(0))
113 allocate(scenario%height(0))
115 allocate(scenario%gas_emission_time(0))
116 allocate(scenario%gas_emission_rate_scale(0))
117 allocate(scenario%gas_emission(0))
119 allocate(scenario%gas_dilution_time(0))
120 allocate(scenario%gas_dilution_rate(0))
121 allocate(scenario%gas_background(0))
123 allocate(scenario%aero_emission_time(0))
124 allocate(scenario%aero_emission_rate_scale(0))
125 allocate(scenario%aero_emission(0))
127 allocate(scenario%aero_dilution_time(0))
128 allocate(scenario%aero_dilution_rate(0))
129 allocate(scenario%aero_background(0))
143 deallocate(scenario%temp_time)
144 deallocate(scenario%temp)
146 deallocate(scenario%pressure_time)
147 deallocate(scenario%pressure)
149 deallocate(scenario%height_time)
150 deallocate(scenario%height)
152 do i = 1,
size(scenario%gas_emission)
155 deallocate(scenario%gas_emission_time)
156 deallocate(scenario%gas_emission_rate_scale)
157 deallocate(scenario%gas_emission)
159 do i = 1,
size(scenario%gas_background)
162 deallocate(scenario%gas_dilution_time)
163 deallocate(scenario%gas_dilution_rate)
164 deallocate(scenario%gas_background)
166 do i = 1,
size(scenario%aero_emission)
169 deallocate(scenario%aero_emission_time)
170 deallocate(scenario%aero_emission_rate_scale)
171 deallocate(scenario%aero_emission)
173 do i = 1,
size(scenario%aero_background)
176 deallocate(scenario%aero_dilution_time)
177 deallocate(scenario%aero_dilution_rate)
178 deallocate(scenario%aero_background)
190 type(scenario_t),
intent(inout) :: scenario_to
196 allocate(scenario_to%temp_time( &
197 size(scenario_from%temp_time)))
198 scenario_to%temp_time = scenario_from%temp_time
199 allocate(scenario_to%temp( &
200 size(scenario_from%temp)))
201 scenario_to%temp = scenario_from%temp
203 allocate(scenario_to%pressure_time( &
204 size(scenario_from%pressure_time)))
205 scenario_to%pressure_time = scenario_from%pressure_time
206 allocate(scenario_to%pressure( &
207 size(scenario_from%pressure)))
208 scenario_to%pressure = scenario_from%pressure
210 allocate(scenario_to%height_time( &
211 size(scenario_from%height_time)))
212 scenario_to%height_time = scenario_from%height_time
213 allocate(scenario_to%height( &
214 size(scenario_from%height)))
215 scenario_to%height = scenario_from%height
217 allocate(scenario_to%gas_emission_time( &
218 size(scenario_from%gas_emission_time)))
219 scenario_to%gas_emission_time = scenario_from%gas_emission_time
220 allocate(scenario_to%gas_emission_rate_scale( &
221 size(scenario_from%gas_emission_rate_scale)))
222 scenario_to%gas_emission_rate_scale = scenario_from%gas_emission_rate_scale
223 allocate(scenario_to%gas_emission( &
224 size(scenario_from%gas_emission)))
225 do i = 1,
size(scenario_from%gas_emission)
228 scenario_to%gas_emission(i))
231 allocate(scenario_to%gas_dilution_time( &
232 size(scenario_from%gas_dilution_time)))
233 scenario_to%gas_dilution_time = scenario_from%gas_dilution_time
234 allocate(scenario_to%gas_dilution_rate( &
235 size(scenario_from%gas_dilution_rate)))
236 scenario_to%gas_dilution_rate = scenario_from%gas_dilution_rate
237 allocate(scenario_to%gas_background( &
238 size(scenario_from%gas_background)))
239 do i = 1,
size(scenario_from%gas_background)
242 scenario_to%gas_background(i))
245 allocate(scenario_to%aero_emission_time( &
246 size(scenario_from%aero_emission_time)))
247 scenario_to%aero_emission_time = scenario_from%aero_emission_time
248 allocate(scenario_to%aero_emission_rate_scale( &
249 size(scenario_from%aero_emission_rate_scale)))
250 scenario_to%aero_emission_rate_scale &
251 = scenario_from%aero_emission_rate_scale
252 allocate(scenario_to%aero_emission( &
253 size(scenario_from%aero_emission)))
254 do i = 1,
size(scenario_from%aero_emission)
257 scenario_to%aero_emission(i))
260 allocate(scenario_to%aero_dilution_time( &
261 size(scenario_from%aero_dilution_time)))
262 scenario_to%aero_dilution_time = scenario_from%aero_dilution_time
263 allocate(scenario_to%aero_dilution_rate( &
264 size(scenario_from%aero_dilution_rate)))
265 scenario_to%aero_dilution_rate = scenario_from%aero_dilution_rate
266 allocate(scenario_to%aero_background( &
267 size(scenario_from%aero_background)))
268 do i = 1,
size(scenario_from%aero_background)
271 scenario_to%aero_background(i))
287 real(kind=dp),
intent(in) :: time
289 env_state%temp =
interp_1d(scenario%temp_time, scenario%temp, time)
290 env_state%pressure =
interp_1d(scenario%pressure_time, &
291 scenario%pressure, time)
292 env_state%height =
interp_1d(scenario%height_time, scenario%height, time)
293 env_state%elapsed_time = time
308 real(kind=dp),
intent(in) :: time
311 real(kind=dp) :: pmv_old, pmv_new
313 real(kind=dp) :: pressure_old
315 real(kind=dp) :: temp_old
321 pressure_old = env_state%pressure
322 temp_old = env_state%temp
324 env_state%temp =
interp_1d(scenario%temp_time, scenario%temp, time)
325 env_state%pressure =
interp_1d(scenario%pressure_time, &
326 scenario%pressure, time)
328 pmv_new = pmv_old * env_state%pressure / pressure_old
331 env_state%height =
interp_1d(scenario%height_time, scenario%height, time)
332 env_state%elapsed_time = time
376 old_env_state, gas_data, gas_state)
381 real(kind=dp),
intent(in) :: delta_t
391 real(kind=dp) :: emission_rate_scale, dilution_rate, p
397 scenario%gas_emission_time, scenario%gas_emission_rate_scale, &
398 env_state%elapsed_time, emissions, emission_rate_scale)
400 p = emission_rate_scale * delta_t / env_state%height
407 scenario%gas_dilution_time, scenario%gas_dilution_rate, &
408 env_state%elapsed_time, background, dilution_rate)
409 p = exp(- dilution_rate * delta_t)
410 if (env_state%height > old_env_state%height)
then
411 p = p * old_env_state%height / env_state%height
429 old_env_state, aero_data, aero_state, n_emit, n_dil_in, n_dil_out, &
430 allow_doubling, allow_halving)
435 real(kind=dp),
intent(in) :: delta_t
445 integer,
intent(out) :: n_emit
447 integer,
intent(out) :: n_dil_in
449 integer,
intent(out) :: n_dil_out
451 logical,
intent(in) :: allow_doubling
453 logical,
intent(in) :: allow_halving
455 real(kind=dp) :: emission_rate_scale, dilution_rate, p
462 scenario%aero_emission_time, scenario%aero_emission_rate_scale, &
463 env_state%elapsed_time, emissions, emission_rate_scale)
464 p = emission_rate_scale * delta_t / env_state%height
466 emissions, p, env_state%elapsed_time, allow_doubling, allow_halving, &
473 scenario%aero_dilution_time, scenario%aero_dilution_rate, &
474 env_state%elapsed_time, background, dilution_rate)
475 p = exp(- dilution_rate * delta_t)
476 if (env_state%height > old_env_state%height)
then
477 p = p * old_env_state%height / env_state%height
482 1d0 - p, aero_info_dilution)
487 background, 1d0 - p, env_state%elapsed_time, allow_doubling, &
488 allow_halving, n_dil_in)
493 old_env_state%temp * env_state%pressure &
494 / (env_state%temp * old_env_state%pressure))
505 old_env_state, bin_grid, aero_data, aero_binned)
510 real(kind=dp),
intent(in) :: delta_t
520 type(aero_binned_t
),
intent(inout) :: aero_binned
522 real(kind=dp) :: emission_rate_scale, dilution_rate, p
524 type(aero_binned_t
) :: emissions_binned, background_binned
531 scenario%aero_emission_time, scenario%aero_emission_rate_scale, &
532 env_state%elapsed_time, emissions, emission_rate_scale)
535 p = emission_rate_scale * delta_t / env_state%height
545 scenario%aero_dilution_time, scenario%aero_dilution_rate, &
546 env_state%elapsed_time, background, dilution_rate)
549 p = exp(- dilution_rate * delta_t)
550 if (env_state%height > old_env_state%height)
then
551 p = p * old_env_state%height / env_state%height
564 aero_data, env_state)
567 integer,
intent(in) :: loss_function_type
569 real(kind=dp),
intent(in) :: vol
571 real(kind=dp),
intent(in) :: density
578 if (loss_function_type == scenario_loss_function_invalid)
then
580 else if (loss_function_type == scenario_loss_function_zero)
then
582 else if (loss_function_type == scenario_loss_function_constant)
then
584 else if (loss_function_type == scenario_loss_function_volume)
then
586 else if (loss_function_type == scenario_loss_function_dry_dep)
then
589 call
die_msg(201594391,
"Unknown loss function id: " &
603 real(kind=dp),
intent(in) :: vol
605 real(kind=dp),
intent(in) :: density
612 real(kind=dp) :: density_air
613 real(kind=dp) :: visc_d, visc_k
614 real(kind=dp) :: gas_speed, gas_mean_free_path
615 real(kind=dp) :: knud, cunning
616 real(kind=dp) :: grav
617 real(kind=dp) :: r_s, r_a
618 real(kind=dp) :: alpha, beta, gamma, a, eps_0
619 real(kind=dp) :: diff_p
620 real(kind=dp) :: von_karman
621 real(kind=dp) :: st, sc, u_star
622 real(kind=dp) :: e_b, e_im, e_in, r1
623 real(kind=dp) :: u_mean, z_ref, z_rough
638 density_air = (const%air_molec_weight * env_state%pressure) &
639 / (const%univ_gas_const * env_state%temp)
641 visc_d = 1.8325d-5 * (416.16 / (env_state%temp + 120.0d0)) &
642 * (env_state%temp / 296.16)**1.5d0
644 visc_k = visc_d / density_air
647 sqrt((8.0d0 * const%boltzmann * env_state%temp * const%avagadro) / &
648 (const%pi * const%air_molec_weight))
650 gas_mean_free_path = (2.0d0 * visc_d) / (density_air * gas_speed)
652 knud = (2.0d0 * gas_mean_free_path) / d_p
654 cunning = 1.0d0 + knud * (1.257d0 + 0.4d0 * exp(-1.1d0 / knud))
658 v_s = (density * d_p**2.0d0 * grav * cunning) / (18.0d0 * visc_d)
662 u_star = .4d0 * u_mean / log(z_ref / z_rough)
663 r_a = (1.0d0 / (.4d0 * u_star)) * log(z_ref / z_rough)
665 diff_p = (const%boltzmann * env_state%temp * cunning) / &
666 (3.d0 * const%pi * visc_d * d_p)
672 e_in = .5d0 * (d_p / a)**2.0d0
675 st = (v_s * u_star) / (grav * a)
676 e_im = (st / (alpha + st))**beta
683 r_s = 1.0d0 / (eps_0 * u_star * (e_b + e_in + e_im) * r1)
686 v_d = v_s + (1.0d0 / (r_a + r_s + r_a * r_s * v_s))
697 aero_data, env_state)
700 integer,
intent(in) :: loss_function_type
702 real(kind=dp),
intent(in) :: vol
709 integer,
parameter :: n_sample = 3
711 real(kind=dp) :: d, d_min, d_max, loss
714 d_min = minval(aero_data%density)
715 d_max = maxval(aero_data%density)
735 aero_data, env_state, loss_max)
738 integer,
intent(in) :: loss_function_type
746 real(kind=dp),
intent(out) :: loss_max(bin_grid%n_bin)
749 integer,
parameter :: n_sample = 3
751 real(kind=dp),
parameter :: over_scale = 2d0
753 real(kind=dp) :: v_low, v_high, vol, r, r_max
756 do b = 1,bin_grid%n_bin
757 v_low =
rad2vol(bin_grid%edges(b))
758 v_high =
rad2vol(bin_grid%edges(b + 1))
764 r_max = max(r_max, r)
766 loss_max(b) = r_max * over_scale
782 aero_state, env_state)
785 integer,
intent(in) :: loss_function_type
787 real(kind=dp),
intent(in) :: delta_t
795 integer :: c, b, s, i_part
796 real(kind=dp) :: over_rate, over_prob, rand_real, rand_geom
798 if (loss_function_type == scenario_loss_function_zero .or. &
799 loss_function_type == scenario_loss_function_invalid)
return
801 if (alg_threshold <= 0d0)
then
803 do i_part = aero_state%apa%n_part, 1, -1
805 aero_data, aero_state, env_state, i_part, 1d0)
812 if (.not. aero_state%aero_sorted%removal_rate_bounds_valid)
then
814 aero_state%aero_sorted%bin_grid, aero_data, env_state, &
815 aero_state%aero_sorted%removal_rate_max)
816 aero_state%aero_sorted%removal_rate_bounds_valid = .true.
821 over_rate = aero_state%aero_sorted%removal_rate_max(b)
822 if (delta_t * over_rate <= 0d0) cycle
823 over_prob = 1d0 - exp(-delta_t * over_rate)
824 if (over_prob >= alg_threshold)
then
826 do s = aero_state%aero_sorted%size_class%inverse(b, c)%n_entry, &
829 aero_state%aero_sorted%size_class%inverse(b, c)%entry(s)
831 delta_t, aero_data, aero_state, env_state, i_part, 1d0)
835 s = aero_state%aero_sorted%size_class%inverse(b, c)%n_entry + 1
838 if (rand_real <= 0d0)
exit
839 rand_geom = -log(rand_real) / (delta_t * over_rate) + 1d0
840 if (rand_geom >=
real(s, kind=dp)) exit
841 s = s - floor(rand_geom)
847 aero_state%aero_sorted%size_class%inverse(b, c)%entry(s)
849 delta_t, aero_data, aero_state, env_state, i_part, &
868 aero_data, aero_state, env_state, i_part, over_prob)
870 integer,
intent(in) :: loss_function_type
872 real(kind=dp),
intent(in) :: delta_t
880 integer,
intent(in) :: i_part
882 real(kind=dp),
intent(in) :: over_prob
884 real(kind=dp) :: prob, rate, vol, density
886 type(aero_info_t
) :: aero_info
888 aero_particle => aero_state%apa%particle(i_part)
893 prob = 1d0 - exp(-delta_t * rate)
895 "particle loss upper bound estimation is too tight: " &
901 aero_info%id = aero_particle%id
902 aero_info%action = aero_info_dilution
903 aero_info%other_id = 0
918 integer,
intent(in) :: aero_mode_type
931 subroutine spec_file_read_scenario(file, gas_data, aero_data, scenario)
942 character(len=PMC_MAX_FILENAME_LEN) :: sub_filename
981 scenario%temp_time, scenario%temp)
988 scenario%pressure_time, scenario%pressure)
995 scenario%height_time, scenario%height)
1001 call spec_file_read_gas_states_times_rates(sub_file, gas_data, &
1002 scenario%gas_emission_time, scenario%gas_emission_rate_scale, &
1003 scenario%gas_emission)
1009 call spec_file_read_gas_states_times_rates(sub_file, gas_data, &
1010 scenario%gas_dilution_time, scenario%gas_dilution_rate, &
1011 scenario%gas_background)
1017 call spec_file_read_aero_dists_times_rates(sub_file, aero_data, &
1018 scenario%aero_emission_time, scenario%aero_emission_rate_scale, &
1019 scenario%aero_emission)
1025 call spec_file_read_aero_dists_times_rates(sub_file, aero_data, &
1026 scenario%aero_dilution_time, scenario%aero_dilution_rate, &
1027 scenario%aero_background)
1030 end subroutine spec_file_read_scenario
1141 integer :: total_size, i, n
1158 do i = 1,
size(val%gas_emission)
1159 total_size = total_size &
1162 do i = 1,
size(val%gas_background)
1163 total_size = total_size &
1166 do i = 1,
size(val%aero_emission)
1167 total_size = total_size &
1170 do i = 1,
size(val%aero_background)
1171 total_size = total_size &
1185 character,
intent(inout) :: buffer(:)
1187 integer,
intent(inout) :: position
1192 integer :: prev_position, i
1194 prev_position = position
1207 val%aero_emission_rate_scale)
1210 do i = 1,
size(val%gas_emission)
1213 do i = 1,
size(val%gas_background)
1216 do i = 1,
size(val%aero_emission)
1219 do i = 1,
size(val%aero_background)
1234 character,
intent(inout) :: buffer(:)
1236 integer,
intent(inout) :: position
1241 integer :: prev_position, i
1245 prev_position = position
1254 val%gas_emission_rate_scale)
1259 val%aero_emission_rate_scale)
1262 deallocate(val%gas_emission)
1263 deallocate(val%gas_background)
1264 deallocate(val%aero_emission)
1265 deallocate(val%aero_background)
1266 allocate(val%gas_emission(
size(val%gas_emission_time)))
1267 allocate(val%gas_background(
size(val%gas_dilution_time)))
1268 allocate(val%aero_emission(
size(val%aero_emission_time)))
1269 allocate(val%aero_background(
size(val%aero_dilution_time)))
1270 do i = 1,
size(val%gas_emission)
1274 do i = 1,
size(val%gas_background)
1278 do i = 1,
size(val%aero_emission)
1282 do i = 1,
size(val%aero_background)
1296 subroutine spec_file_read_loss_function_type(file, loss_function_type)
1301 integer,
intent(out) :: loss_function_type
1303 character(len=SPEC_LINE_MAX_VAR_LEN) :: function_name
1319 if (trim(function_name) ==
'zero')
then
1320 loss_function_type = scenario_loss_function_zero
1321 else if (trim(function_name) ==
'constant')
then
1322 loss_function_type = scenario_loss_function_constant
1323 else if (trim(function_name) ==
'volume')
then
1324 loss_function_type = scenario_loss_function_volume
1325 else if (trim(function_name) ==
'drydep')
then
1326 loss_function_type = scenario_loss_function_dry_dep
1329 "Unknown loss function type: " // trim(function_name))
1332 end subroutine spec_file_read_loss_function_type
real(kind=dp) elemental function vol2diam(v)
Convert volume (m^3) to diameter (m).
subroutine aero_binned_allocate_size(aero_binned, n_bin, n_spec)
Allocate an aero_binned_t of the given size.
The scenario_t structure and associated subroutines.
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
An input file with extra data for printing messages.
subroutine scenario_update_aero_state(scenario, delta_t, env_state, old_env_state, aero_data, aero_state, n_emit, n_dil_in, n_dil_out, allow_doubling, allow_halving)
Do emissions and background dilution for a particle aerosol distribution.
The aero_data_t structure and associated subroutines.
subroutine die_msg(code, error_msg)
Error immediately.
subroutine aero_info_deallocate(aero_info)
Deallocates.
subroutine aero_dist_interp_1d(aero_dist_list, time_list, rate_list, time, aero_dist, rate)
Determine the current aero_dist and rate by interpolating at the current time with the lists of aero_...
integer function pmc_mpi_pack_size_gas_state(val)
Determines the number of bytes required to pack the given value.
subroutine aero_info_allocate(aero_info)
Allocates and initializes.
real(kind=dp) elemental function rad2vol(r)
Convert radius (m) to volume (m^3).
subroutine scenario_particle_loss(loss_function_type, delta_t, aero_data, aero_state, env_state)
Performs stochastic particle loss for one time-step. If a particle i_part has a scenario_loss_rate() ...
integer function pmc_mpi_pack_size_aero_dist(val)
Determines the number of bytes required to pack the given value.
The aero_dist_t structure and associated subroutines.
subroutine spec_file_open(filename, file)
Open a spec file for reading.
The env_state_t structure and associated subroutines.
subroutine aero_state_allocate(aero_state)
Allocates aerosol arrays.
subroutine scenario_init_env_state(scenario, env_state, time)
Initialize the time-dependent contents of the environment. Thereafter scenario_update_env_state() sho...
subroutine aero_weight_array_scale(aero_weight_array, factor)
Scale the weights by the given factor, so new_weight = old_weight * factor.
real(kind=dp) function scenario_loss_rate(loss_function_type, vol, density, aero_data, env_state)
Evaluate a loss rate function.
subroutine spec_file_read_timed_real_array(file, name, times, vals)
Read an a time-indexed array of real data.
subroutine aero_binned_deallocate(aero_binned)
Free internal memory in an aero_binned_t structure.
subroutine pmc_mpi_pack_gas_state(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine scenario_loss_rate_bin_max(loss_function_type, bin_grid, aero_data, env_state, loss_max)
Compute an upper bound on the maximum kernel value for each bin. Value over_scale is multiplied to th...
subroutine aero_dist_deallocate(aero_dist)
Free all storage.
subroutine scenario_update_env_state(scenario, env_state, time)
Update time-dependent contents of the environment. scenario_init_env_state() should have been called ...
integer function pmc_mpi_pack_size_scenario(val)
Determines the number of bytes required to pack the given value.
subroutine gas_state_add_scaled(gas_state, gas_state_delta, alpha)
Adds the given gas_state_delta scaled by alpha.
real(kind=dp) function interp_linear_disc(x_1, x_n, n, i)
Linear interpolation over discrete indices.
elemental logical function scenario_contains_aero_mode_type(scenario, aero_mode_type)
Whether any of the contained aerosol modes are of the given type.
subroutine aero_state_add_aero_dist_sample(aero_state, aero_data, aero_dist, sample_prop, create_time, allow_doubling, allow_halving, n_part_add)
Generates a Poisson sample of an aero_dist, adding to aero_state, with the given sample proportion...
The gas_data_t structure and associated subroutines.
subroutine scenario_copy(scenario_from, scenario_to)
Copy structure.
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
real(kind=dp) function aero_particle_density(aero_particle, aero_data)
Average density of the particle (kg/m^3).
subroutine scenario_allocate(scenario)
Allocate an scenario.
A complete aerosol distribution, consisting of several modes.
Common utility subroutines.
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
integer function aero_sorted_n_bin(aero_sorted)
Returns the number of size bins.
Current environment state.
subroutine aero_state_sort(aero_state, bin_grid, all_procs_same)
Sorts the particles if necessary.
subroutine aero_binned_scale(aero_binned, alpha)
Scale an aero_binned_t by a real number.
subroutine scenario_try_single_particle_loss(loss_function_type, delta_t, aero_data, aero_state, env_state, i_part, over_prob)
Test a candidate particle to see if it should be removed, and remove if necessary. Particle is removed with probability (1d0 - exp(-delta_t*rate))/over_prob, where rate is the loss function evaluated for the given particle.
The aero_state_t structure and assocated subroutines.
Wrapper functions for MPI.
real(kind=dp) function scenario_loss_rate_max(loss_function_type, vol, aero_data, env_state)
Compute and return the max loss rate function for a given volume.
integer function aero_sorted_n_class(aero_sorted)
Returns the number of weight classes.
subroutine gas_state_deallocate(gas_state)
Free all storage.
subroutine spec_file_read_string(file, name, var)
Read a string from a spec file that must have a given name.
The gas_state_t structure and associated subroutines.
subroutine gas_state_allocate(gas_state)
Allocate storage for gas species.
integer function aero_state_total_particles(aero_state, i_group, i_class)
Returns the total number of particles in an aerosol distribution.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
The current collection of aerosol particles.
subroutine aero_state_remove_particle_with_info(aero_state, i_part, aero_info)
Remove the given particle and record the removal.
Single aerosol particle data structure.
subroutine scenario_update_gas_state(scenario, delta_t, env_state, old_env_state, gas_data, gas_state)
Do gas emissions and background dilution.
1D grid, either logarithmic or linear.
subroutine gas_state_scale(gas_state, alpha)
Scale a gas state.
subroutine gas_state_ensure_nonnegative(gas_state)
Set any negative values to zero.
subroutine aero_state_deallocate(aero_state)
Deallocates a previously allocated aerosol.
Reading formatted text input.
subroutine spec_file_close(file)
Close a spec file.
subroutine aero_dist_copy(aero_dist_from, aero_dist_to)
Copy an aero_dist.
elemental logical function aero_dist_contains_aero_mode_type(aero_dist, aero_mode_type)
Whether any of the modes are of the given type.
real(kind=dp) function scenario_loss_rate_dry_dep(vol, density, env_state)
Compute and return the dry deposition rate for a given particle. All equations used here are written ...
subroutine aero_state_sample_particles(aero_state_from, aero_state_to, sample_prob, removal_action)
Generates a random sample by removing particles from aero_state_from and adding them to aero_state_to...
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
subroutine pmc_mpi_unpack_aero_dist(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
real(kind=dp) function pmc_random()
Returns a random number between 0 and 1.
Current state of the gas mixing ratios in the system.
subroutine scenario_update_aero_binned(scenario, delta_t, env_state, old_env_state, bin_grid, aero_data, aero_binned)
Do emissions and background dilution from the environment for a binned aerosol distribution.
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
real(kind=dp) function interp_1d(x_vals, y_vals, x)
1D linear interpolation.
subroutine pmc_mpi_pack_scenario(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine pmc_mpi_pack_aero_dist(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine gas_state_allocate_size(gas_state, n_spec)
Allocate storage for gas species of the given size.
subroutine gas_state_interp_1d(gas_state_list, time_list, rate_list, time, gas_state, rate)
Determine the current gas_state and rate by interpolating at the current time with the lists of gas_s...
subroutine aero_binned_add_aero_dist(aero_binned, bin_grid, aero_data, aero_dist)
Add an aero_dist_t to an aero_binned_t.
elemental real(kind=dp) function aero_particle_volume(aero_particle)
Total volume of the particle (m^3).
real(kind=dp) function env_state_sat_vapor_pressure(env_state)
Computes the current saturation vapor pressure (Pa).
character(len=pmc_util_convert_string_len) function real_to_string(val)
Convert a real to a string format.
subroutine gas_state_mole_dens_to_ppb(gas_state, env_state)
Convert (mol m^{-3}) to (ppb).
Aerosol material properties and associated data.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine pmc_mpi_unpack_scenario(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine aero_binned_add_scaled(aero_binned, aero_binned_delta, alpha)
Add a scaled aero_binned_t structure to an existing one.
subroutine aero_dist_allocate(aero_dist)
Allocates an aero_dist.
subroutine scenario_deallocate(scenario)
Free all storage.
subroutine gas_state_copy(from_state, to_state)
Copy to an already allocated to_state.
subroutine pmc_mpi_unpack_gas_state(buffer, position, val)
Unpacks the given value from the buffer, advancing position.