30 real(kind=dp),
parameter :: COAG_ACCEL_N_EVENT = 1d0
33 real(kind=dp),
parameter :: COAG_ACCEL_MAX_CV = 0.1d0
36 real(kind=dp),
parameter :: MAX_ALLOWABLE_GROWTH_FACTOR = 1.5d0
43 subroutine mc_coag(coag_kernel_type, env_state, aero_data, aero_state, &
44 del_t, tot_n_samp, tot_n_coag)
47 integer,
intent(in) :: coag_kernel_type
55 real(kind=dp) :: del_t
57 integer,
intent(out) :: tot_n_samp
59 integer,
intent(out) :: tot_n_coag
61 integer :: b1, b2, c1, c2, b2_start
64 if (.not. aero_state%aero_sorted%coag_kernel_bounds_valid)
then
66 coag_kernel_type, aero_data, env_state, &
67 aero_state%aero_sorted%coag_kernel_min, &
68 aero_state%aero_sorted%coag_kernel_max)
69 aero_state%aero_sorted%coag_kernel_bounds_valid = .true.
77 if (aero_state%aero_sorted%size_class%inverse(b1, &
86 if (aero_state%aero_sorted%size_class%inverse(b2, &
90 aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2)
102 aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2)
105 integer,
intent(in) :: coag_kernel_type
113 real(kind=dp) :: del_t
115 integer,
intent(inout) :: tot_n_samp
117 integer,
intent(inout) :: tot_n_coag
119 integer,
intent(in) :: b1
121 integer,
intent(in) :: b2
123 integer,
intent(in) :: c1
125 integer,
intent(in) :: c2
127 logical :: per_particle_coag_succeeded
129 real(kind=dp) :: f_max, k_max
131 cc =
coag_dest_class(aero_state%awa, aero_state%aero_sorted%bin_grid, b1, &
134 aero_state%aero_sorted%bin_grid, b1, b2, c1, c2, cc, f_max)
135 k_max = aero_state%aero_sorted%coag_kernel_max(b1, b2) * f_max
138 aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2, cc, &
139 per_particle_coag_succeeded)
140 if (per_particle_coag_succeeded)
return
142 call
per_set_coag(coag_kernel_type, k_max, env_state, aero_data, &
143 aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2, cc)
151 aero_data, aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2, &
152 cc, per_particle_coag_succeeded)
155 integer,
intent(in) :: coag_kernel_type
157 real(kind=dp),
intent(in) :: k_max
165 real(kind=dp) :: del_t
167 integer,
intent(inout) :: tot_n_samp
169 integer,
intent(inout) :: tot_n_coag
171 integer,
intent(in) :: b1
173 integer,
intent(in) :: b2
175 integer,
intent(in) :: c1
177 integer,
intent(in) :: c2
179 integer,
intent(in) :: cc
181 logical,
intent(inout) :: per_particle_coag_succeeded
183 logical :: correct_weight_ordering
184 integer :: target_unif_entry, target_part, n_samp, n_coag, n_remove, bt, bs
186 real(kind=dp) :: n_source_per_target, accept_factor
187 real(kind=dp) :: min_target_vol, max_source_vol, max_new_target_vol
188 real(kind=dp) :: max_target_growth_factor
192 aero_state%aero_sorted%bin_grid, b1, b2, c1, c2, cc, bt, bs, ct, cs, &
193 correct_weight_ordering)
194 if (.not. correct_weight_ordering)
then
195 per_particle_coag_succeeded = .false.
200 cs)%n_entry, k_max, del_t, n_source_per_target, accept_factor)
201 if (n_source_per_target < coag_accel_n_event)
then
202 per_particle_coag_succeeded = .false.
206 min_target_vol =
rad2vol(aero_state%aero_sorted%bin_grid%edges(bt))
207 max_source_vol =
rad2vol(aero_state%aero_sorted%bin_grid%edges(bs+1))
208 max_new_target_vol = min_target_vol + max_source_vol * n_source_per_target
210 max_target_growth_factor = max_new_target_vol / min_target_vol
211 if (max_target_growth_factor > max_allowable_growth_factor)
then
212 per_particle_coag_succeeded = .false.
220 do target_unif_entry &
221 = aero_state%aero_sorted%size_class%inverse(bt, ct)%n_entry,1,-1
222 target_part = aero_state%aero_sorted%size_class%inverse(bt, &
223 ct)%entry(target_unif_entry)
229 coag_kernel_type, bs, cs, target_particle, n_source_per_target, &
230 accept_factor, n_samp, n_coag, n_remove, source_particle)
235 tot_n_samp = tot_n_samp + n_samp
236 tot_n_coag = tot_n_coag + n_coag
243 per_particle_coag_succeeded = .true.
252 b2, c1, c2, cc, bt, bs, ct, cs, correct_weight_ordering)
259 integer,
intent(in) :: b1
261 integer,
intent(in) :: b2
263 integer,
intent(in) :: c1
265 integer,
intent(in) :: c2
267 integer,
intent(in) :: cc
269 integer,
intent(out) :: bt
271 integer,
intent(out) :: bs
273 integer,
intent(out) :: ct
275 integer,
intent(out) :: cs
277 logical,
intent(out) :: correct_weight_ordering
279 real(kind=dp) :: nc1_min, nc1_max, nc2_min, nc2_max
280 logical :: monotone_increasing, monotone_decreasing
283 monotone_increasing, monotone_decreasing)
284 if (.not. monotone_decreasing)
then
285 correct_weight_ordering = .false.
290 bin_grid%edges(b1), bin_grid%edges(b1 + 1), nc1_min, nc1_max)
292 bin_grid%edges(b2), bin_grid%edges(b2 + 1), nc2_min, nc2_max)
295 correct_weight_ordering = .false.
296 if ((nc1_max < nc2_min) .and. (cc == c1))
then
301 correct_weight_ordering = .true.
303 if ((nc2_max < nc1_min) .and. (cc == c2))
then
308 correct_weight_ordering = .true.
319 integer,
intent(in) :: n_part
321 real(kind=dp),
intent(in) :: k_max
323 real(kind=dp),
intent(in) :: del_t
325 real(kind=dp),
intent(out) :: n_source_per_target
327 real(kind=dp),
intent(out) :: accept_factor
329 n_source_per_target = k_max * del_t *
real(n_part, kind=dp)
330 accept_factor = 1d0 / k_max
338 coag_kernel_type, bs, cs, coag_particle, n_samp_mean, &
339 accept_factor, n_samp, n_coag, n_remove, source_particle)
348 integer,
intent(in) :: coag_kernel_type
350 integer,
intent(in) :: bs
352 integer,
intent(in) :: cs
356 real(kind=dp),
intent(in) :: n_samp_mean
358 real(kind=dp),
intent(in) :: accept_factor
360 integer,
intent(out) :: n_samp
362 integer,
intent(out) :: n_coag
364 integer,
intent(out) :: n_remove
368 real(kind=dp) :: prob_remove_i, prob_remove_source_max
369 real(kind=dp) :: prob_coag, prob_coag_tot, prob_coag_mean
370 real(kind=dp) :: num_conc_i, num_conc_source_min, num_conc_target, k
371 real(kind=dp) :: vol_sq(aero_data%n_spec), vol_mean(aero_data%n_spec)
372 real(kind=dp) :: vol_cv(aero_data%n_spec), vol_cv_max, mean_95_conf_cv
373 integer :: n_samp_remove, n_samp_extra, n_samp_total, n_avg, i_samp
374 integer :: i_unif_entry, i_part, target_id, new_bin, ct
376 type(aero_info_t
) :: aero_info
378 if (aero_state%aero_sorted%size_class%inverse(bs, cs)%n_entry == 0)
then
386 target_id = coag_particle%id
387 ct = coag_particle%weight_class
390 aero_state%awa, cs, aero_state%aero_sorted%bin_grid%edges(bs))
391 prob_remove_source_max = num_conc_target / num_conc_source_min
392 call
assert(653606684, prob_remove_source_max <= 1d0)
394 n_samp_remove =
rand_poisson(prob_remove_source_max * n_samp_mean)
395 n_samp_extra =
rand_poisson((1d0 - prob_remove_source_max) * n_samp_mean)
396 n_samp_total = n_samp_remove + n_samp_extra
404 aero_data%n_spec, aero_data%n_source)
409 do i_samp = 1,n_samp_total
410 if (aero_state%aero_sorted%size_class%inverse(bs, cs)%n_entry == 0)
exit
411 if ((n_samp > n_samp_remove) .and. (n_avg >= 2))
then
412 vol_mean = source_particle%vol /
real(n_avg, kind=dp)
413 where(vol_mean > 0d0) &
414 vol_cv = sqrt(vol_sq /
real(n_avg, kind=dp) - (vol_mean)**2) &
416 vol_cv_max = maxval(vol_cv)
417 mean_95_conf_cv = vol_cv_max / sqrt(
real(n_avg, kind=dp)) &
421 if (mean_95_conf_cv < coag_accel_max_cv)
exit
426 =
pmc_rand_int(aero_state%aero_sorted%size_class%inverse(bs, &
428 i_part = aero_state%aero_sorted%size_class%inverse(bs, &
429 cs)%entry(i_unif_entry)
430 i_particle => aero_state%apa%particle(i_part)
433 coag_particle, cs, ct, ct, aero_data, aero_state%awa, env_state, k)
434 prob_coag = k * accept_factor
435 prob_coag_tot = prob_coag_tot + prob_coag
440 vol_sq = vol_sq + i_particle%vol**2
441 if (i_samp <= n_samp_remove)
then
444 prob_remove_i = num_conc_target / num_conc_i
445 if (
pmc_random() < prob_remove_i / prob_remove_source_max)
then
446 n_remove = n_remove + 1
448 aero_info%id = i_particle%id
449 aero_info%action = aero_info_coag
450 aero_info%other_id = target_id
464 prob_coag_mean = prob_coag_tot /
real(n_samp, kind=dp)
466 "kernel upper bound estimation is too tight")
468 source_particle%vol = source_particle%vol &
469 * (
real(n_coag, kind=dp) /
real(n_avg, kind=dp))
482 integer,
intent(in) :: bt
484 integer,
intent(in) :: ct
486 integer,
intent(in) :: target_unif_entry
490 integer,
intent(in) :: cc
492 integer :: target_part, target_id, new_bin, new_group
493 real(kind=dp) :: old_num_conc_target, new_num_conc_target
495 target_part = aero_state%aero_sorted%size_class%inverse(bt, &
496 ct)%entry(target_unif_entry)
497 target_id = aero_state%apa%particle(target_part)%id
499 aero_state%apa%particle(target_part))
501 source_particle, aero_state%apa%particle(target_part))
502 aero_state%apa%particle(target_part)%id = target_id
510 aero_state%apa%particle(target_part))
512 .or. (new_bin > aero_state%aero_sorted%bin_grid%n_bin))
then
513 call
die_msg(765620746,
"particle outside of bin_grid: " &
514 //
"try reducing the timestep del_t")
517 new_bin, new_group, cc)
525 aero_state%apa%particle(target_part))
527 old_num_conc_target / new_num_conc_target, random_weight_group=.true.)
529 call
assert(654300924, aero_state%apa%particle(target_part)%id &
537 subroutine per_set_coag(coag_kernel_type, k_max, env_state, aero_data, &
538 aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2, cc)
541 integer,
intent(in) :: coag_kernel_type
543 real(kind=dp),
intent(in) :: k_max
551 real(kind=dp) :: del_t
553 integer,
intent(inout) :: tot_n_samp
555 integer,
intent(inout) :: tot_n_coag
557 integer,
intent(in) :: b1
559 integer,
intent(in) :: b2
561 integer,
intent(in) :: c1
563 integer,
intent(in) :: c2
565 integer,
intent(in) :: cc
567 real(kind=dp) :: n_samp_mean, accept_factor
568 integer :: i_samp, n_samp, n1, n2
571 n1 = aero_state%aero_sorted%size_class%inverse(b1, c1)%n_entry
572 n2 = aero_state%aero_sorted%size_class%inverse(b2, c2)%n_entry
573 call
compute_n_samp(n1, n2, ((b1 == b2) .and. (c1 == c2)), k_max, del_t, &
574 n_samp_mean, n_samp, accept_factor)
575 tot_n_samp = tot_n_samp + n_samp
579 n1 = aero_state%aero_sorted%size_class%inverse(b1, c1)%n_entry
580 n2 = aero_state%aero_sorted%size_class%inverse(b2, c2)%n_entry
581 if (((n1 < 2) .and. (b1 == b2) .and. (c1 == c2)) &
582 .or. (n1 < 1) .or. (n2 < 1)) &
584 call
maybe_coag_pair(env_state, aero_data, aero_state, b1, b2, c1, c2, &
585 cc, coag_kernel_type, accept_factor, did_coag)
586 if (did_coag) tot_n_coag = tot_n_coag + 1
595 n_samp, accept_factor)
598 integer,
intent(in) :: ni
600 integer,
intent(in) :: nj
602 logical,
intent(in) :: same_bin
604 real(kind=dp),
intent(in) :: k_max
606 real(kind=dp),
intent(in) :: del_t
608 real(kind=dp),
intent(out) :: n_samp_mean
610 integer,
intent(out) :: n_samp
612 real(kind=dp),
intent(out) :: accept_factor
614 real(kind=dp) :: r_samp
615 real(kind=dp) :: n_possible
623 n_possible =
real(ni, kind=dp) * (
real(nj, kind=dp) - 1d0) / 2d0
625 n_possible =
real(ni, kind=dp) *
real(nj, kind=dp)
628 r_samp = k_max * del_t
629 n_samp_mean = r_samp * n_possible
631 accept_factor = 1d0 / k_max
659 c1, c2, cc, coag_kernel_type, accept_factor, did_coag)
668 integer,
intent(in) :: b1
670 integer,
intent(in) :: b2
672 integer,
intent(in) :: c1
674 integer,
intent(in) :: c2
676 integer,
intent(in) :: cc
678 integer,
intent(in) :: coag_kernel_type
680 real(kind=dp),
intent(in) :: accept_factor
682 logical,
intent(out) :: did_coag
686 real(kind=dp) :: p, k
688 call
find_rand_pair(aero_state%aero_sorted, b1, b2, c1, c2, i1, i2)
689 p1 = aero_state%aero_sorted%size_class%inverse(b1, c1)%entry(i1)
690 p2 = aero_state%aero_sorted%size_class%inverse(b2, c2)%entry(i2)
692 aero_state%apa%particle(p1), aero_state%apa%particle(p2), &
693 c1, c2, cc, aero_data, aero_state%awa, env_state, k)
694 p = k * accept_factor
696 "kernel upper bound estimation is too tight, " &
697 //
"could be caused by changing env_state")
701 call
coagulate(aero_data, aero_state, p1, p2, c1, c2, cc)
717 integer,
intent(in) :: b1
719 integer,
intent(in) :: b2
721 integer,
intent(in) :: c1
723 integer,
intent(in) :: c2
725 integer,
intent(out) :: i1
727 integer,
intent(out) :: i2
729 call
assert(619608562, aero_sorted%size_class%inverse(b1, c1)%n_entry >= 1)
730 i1 =
pmc_rand_int(aero_sorted%size_class%inverse(b1, c1)%n_entry)
732 if ((b1 == b2) .and. (c1 == c2))
then
734 aero_sorted%size_class%inverse(b2, c2)%n_entry >= 2)
735 i2 =
pmc_rand_int(aero_sorted%size_class%inverse(b2, c2)%n_entry - 1)
737 i2 = aero_sorted%size_class%inverse(b2, c2)%n_entry
741 aero_sorted%size_class%inverse(b2, c2)%n_entry >= 1)
742 i2 =
pmc_rand_int(aero_sorted%size_class%inverse(b2, c2)%n_entry)
753 aero_weight_array, remove_1, remove_2, create_new, id_1_lost, &
754 id_2_lost, aero_info_1, aero_info_2)
764 integer,
intent(in) :: c1
766 integer,
intent(in) :: c2
768 integer,
intent(in) :: cc
774 logical,
intent(out) :: remove_1
776 logical,
intent(out) :: remove_2
778 logical,
intent(out) :: create_new
780 logical,
intent(out) :: id_1_lost
782 logical,
intent(out) :: id_2_lost
784 type(aero_info_t
),
intent(inout) :: aero_info_1
786 type(aero_info_t
),
intent(inout) :: aero_info_2
788 real(kind=dp) :: r1, r2, rc, nc_min, nc1, nc2, ncc
789 real(kind=dp) :: prob_remove_1, prob_remove_2, prob_create_new
790 integer :: info_other_id, new_group
792 call
assert(371947172, pt1%id /= pt2%id)
803 nc_min = min(nc1, nc2, ncc)
804 prob_remove_1 = nc_min / nc1
805 prob_remove_2 = nc_min / nc2
806 prob_create_new = nc_min / ncc
822 if (remove_1 .and. remove_2)
then
841 if (remove_1 .and. (.not. id_1_lost))
then
843 call
assert(975059559, id_2_lost .eqv. remove_2)
844 elseif (remove_2 .and. (.not. id_2_lost))
then
846 call
assert(246529753, id_1_lost .eqv. remove_1)
849 call
assert(852038606, id_1_lost .eqv. remove_1)
850 call
assert(254018921, id_2_lost .eqv. remove_2)
852 info_other_id = ptc%id
857 aero_info_1%id = pt1%id
858 aero_info_1%action = aero_info_coag
859 aero_info_1%other_id = info_other_id
861 aero_info_2%id = pt2%id
862 aero_info_2%action = aero_info_coag
863 aero_info_2%other_id = info_other_id
871 subroutine coagulate(aero_data, aero_state, p1, p2, c1, c2, cc)
878 integer,
intent(in) :: p1
880 integer,
intent(in) :: p2
882 integer,
intent(in) :: c1
884 integer,
intent(in) :: c2
886 integer,
intent(in) :: cc
891 type(aero_info_t
) :: aero_info_1, aero_info_2
892 logical :: remove_1, remove_2, create_new, id_1_lost, id_2_lost
898 pt1 => aero_state%apa%particle(p1)
899 pt2 => aero_state%apa%particle(p2)
902 aero_state%awa, remove_1, remove_2, create_new, id_1_lost, &
903 id_2_lost, aero_info_1, aero_info_2)
912 id_2_lost, aero_info_2)
916 id_1_lost, aero_info_1)
921 id_1_lost, aero_info_1)
925 id_2_lost, aero_info_2)
subroutine aero_particle_set_weight(aero_particle, i_group, i_class)
Sets the aerosol particle weight group.
The aero_sorted_t structure and assocated subroutines.
subroutine aero_weight_array_minmax_num_conc(aero_weight_array, i_class, radius_1, radius_2, num_conc_min, num_conc_max)
Compute the maximum and minimum number concentrations between the given radii.
subroutine sample_source_particle(aero_state, aero_data, env_state, coag_kernel_type, bs, cs, coag_particle, n_samp_mean, accept_factor, n_samp, n_coag, n_remove, source_particle)
Sample coagulation partners for a single coagulation event.
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_weight_array_check_monotonicity(aero_weight_array, monotone_increasing, monotone_decreasing)
Determine whether all weight functions in an array are monotone increasing, monotone decreasing...
subroutine coagulate(aero_data, aero_state, p1, p2, c1, c2, cc)
Join together particles p1 and p2, updating all particle and bin structures to reflect the change...
subroutine aero_info_allocate(aero_info)
Allocates and initializes.
real(kind=dp) elemental function rad2vol(r)
Convert radius (m) to volume (m^3).
integer function rand_binomial(n, p)
Generate a Binomial-distributed random number with the given parameters.
subroutine num_conc_weighted_kernel(coag_kernel_type, aero_particle_1, aero_particle_2, i_class, j_class, ij_class, aero_data, aero_weight_array, env_state, k)
Compute the kernel value with the given number concentration weighting.
The aero_weight_t structure and associated subroutines.
subroutine max_coag_num_conc_factor(aero_weight_array, bin_grid, i_bin, j_bin, i_class, j_class, ij_class, f_max)
Determine the minimum and maximum number concentration factors for coagulation.
subroutine aero_particle_deallocate(aero_particle)
Deallocates memory associated with an aero_particle_t.
subroutine mc_coag_for_bin(coag_kernel_type, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2)
Do coagulation for time del_t for the given bins.
The env_state_t structure and associated subroutines.
subroutine aero_particle_allocate_size(aero_particle, n_spec, n_source)
Allocates an aero_particle_t of the given size.
Aerosol particle coagulation.
subroutine aero_state_add_particle(aero_state, aero_particle, allow_resort)
Add the given particle.
subroutine coag_target_with_source(aero_state, bt, ct, target_unif_entry, source_particle, cc)
Coagulate a sampled source particle with a target particle.
subroutine coagulate_weighting(pt1, pt2, ptc, c1, c2, cc, aero_data, aero_weight_array, remove_1, remove_2, create_new, id_1_lost, id_2_lost, aero_info_1, aero_info_2)
Actually coagulate pt1 and pt2 to form ptc and compute weighting effects, including which particles s...
elemental real(kind=dp) function aero_particle_radius(aero_particle)
Total radius of the particle (m).
subroutine est_k_minmax_binned_unweighted(bin_grid, coag_kernel_type, aero_data, env_state, k_min, k_max)
Estimate an array of minimum and maximum kernel values. Given particles v1 in bin b1 and v2 in bin b2...
integer function aero_weight_array_rand_group(aero_weight_array, i_class, radius)
Choose a random group at the given radius, with probability inversely proportional to group weight at...
real(kind=dp) function aero_weight_array_num_conc(aero_weight_array, aero_particle)
Compute the number concentration for a particle (m^{-3}).
subroutine compute_n_source(n_part, k_max, del_t, n_source_per_target, accept_factor)
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
subroutine compute_n_samp(ni, nj, same_bin, k_max, del_t, n_samp_mean, n_samp, accept_factor)
Compute the number of samples required for the pair of bins.
Common utility subroutines.
integer function aero_sorted_n_bin(aero_sorted)
Returns the number of size bins.
subroutine find_rand_pair(aero_sorted, b1, b2, c1, c2, i1, i2)
Given (b1, c1) and (b2, c2), find a random pair of particles (b1, c1, i1) and (b2, c2, i2) that are not the same particle particle as each other.
Current environment state.
subroutine aero_state_sort(aero_state, bin_grid, all_procs_same)
Sorts the particles if necessary.
subroutine aero_state_dup_particle(aero_state, i_part, n_part_mean, random_weight_group)
Add copies or remove a particle, with a given mean number of resulting particles. ...
An array of aerosol size distribution weighting functions.
integer function pmc_rand_int(n)
Returns a random integer between 1 and n.
subroutine aero_state_remove_particle(aero_state, i_part, record_removal, aero_info)
Remove the given particle and possibly record the removal.
real(kind=dp) function student_t_95_coeff(n_sample)
Return a fairly tight upper-bound on the Student's t coefficient for the 95% confidence interval...
The aero_state_t structure and assocated subroutines.
Wrapper functions for MPI.
integer function aero_sorted_n_class(aero_sorted)
Returns the number of weight classes.
subroutine aero_sorted_move_particle(aero_sorted, i_part, new_bin, new_group, new_class)
Move a particle to a different bin and group.
subroutine try_per_particle_coag(coag_kernel_type, k_max, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2, cc, per_particle_coag_succeeded)
Attempt per-particle coagulation.
real(kind=dp) function aero_weight_array_num_conc_at_radius(aero_weight_array, i_class, radius)
Compute the total number concentration at a given radius (m^3).
The aero_weight_array_t structure and associated subroutines.
subroutine mc_coag(coag_kernel_type, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag)
Do coagulation for time del_t.
subroutine aero_particle_allocate(aero_particle)
Allocates memory in an aero_particle_t.
integer function rand_poisson(mean)
Generate a Poisson-distributed random number with the given mean.
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.
1D grid, either logarithmic or linear.
The bin_grid_t structure and associated subroutines.
subroutine aero_particle_new_id(aero_particle)
Assigns a globally-unique new ID number to the particle.
integer function coag_dest_class(aero_weight_array, bin_grid, i_bin, j_bin, i_class, j_class)
Determine the weight class in which coagulated particles will be placed.
Sorting of particles into bins.
subroutine per_set_coag(coag_kernel_type, k_max, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2, cc)
Do set-wise coagulation.
subroutine aero_particle_coagulate(aero_particle_1, aero_particle_2, aero_particle_new)
Coagulate two particles together to make a new one. The new particle will not have its ID set...
real(kind=dp) function pmc_random()
Returns a random number between 0 and 1.
real(kind=dp) elemental function vol2rad(v)
Convert volume (m^3) to radius (m).
subroutine maybe_coag_pair(env_state, aero_data, aero_state, b1, b2, c1, c2, cc, coag_kernel_type, accept_factor, did_coag)
Choose a random pair for potential coagulation and test its probability of coagulation. If it happens, do the coagulation and update all structures.
The stats_t type and associated subroutines.
elemental real(kind=dp) function aero_particle_volume(aero_particle)
Total volume of the particle (m^3).
Aerosol material properties and associated data.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine determine_target_and_source(aero_weight_array, bin_grid, b1, b2, c1, c2, cc, bt, bs, ct, cs, correct_weight_ordering)
Determine the source and target particle bin/group for per-particle coagulation, if possible...
integer function aero_sorted_particle_in_bin(aero_sorted, aero_particle)
Find the bin number that contains a given particle.
Generic coagulation kernel.
subroutine aero_particle_copy(aero_particle_from, aero_particle_to)
Copies a particle.