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 (integer_varray_n_entry( &
78 aero_state%aero_sorted%size_class%inverse(b1, c1)) == 0) &
86 if (integer_varray_n_entry( &
87 aero_state%aero_sorted%size_class%inverse(b2, c2)) == 0) &
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
132 aero_state%aero_sorted%bin_grid, b1, b2, c1, c2)
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 aero_state%aero_sorted%size_class%inverse(bs, cs)), k_max, del_t, &
201 n_source_per_target, accept_factor)
203 per_particle_coag_succeeded = .false.
208 aero_state%aero_sorted%bin_grid%edges(bt))
210 aero_state%aero_sorted%bin_grid%edges(bs+1))
211 max_new_target_vol = min_target_vol + max_source_vol * n_source_per_target
213 max_target_growth_factor = max_new_target_vol / min_target_vol
215 per_particle_coag_succeeded = .false.
220 do target_unif_entry = integer_varray_n_entry( &
221 aero_state%aero_sorted%size_class%inverse(bt, ct)),1,-1
222 target_part = aero_state%aero_sorted%size_class%inverse(bt, &
223 ct)%entry(target_unif_entry)
226 target_particle = aero_state%apa%particle(target_part)
228 coag_kernel_type, bs, cs, target_particle, n_source_per_target, &
229 accept_factor, n_samp, n_coag, n_remove, source_particle)
232 target_unif_entry, source_particle, cc)
234 tot_n_samp = tot_n_samp + n_samp
235 tot_n_coag = tot_n_coag + n_coag
239 per_particle_coag_succeeded = .true.
248 b2, c1, c2, cc, bt, bs, ct, cs, correct_weight_ordering)
255 integer,
intent(in) :: b1
257 integer,
intent(in) :: b2
259 integer,
intent(in) :: c1
261 integer,
intent(in) :: c2
263 integer,
intent(in) :: cc
265 integer,
intent(out) :: bt
267 integer,
intent(out) :: bs
269 integer,
intent(out) :: ct
271 integer,
intent(out) :: cs
273 logical,
intent(out) :: correct_weight_ordering
275 real(kind=dp) :: nc1_min, nc1_max, nc2_min, nc2_max
276 logical :: monotone_increasing, monotone_decreasing
279 monotone_increasing, monotone_decreasing)
280 if (.not. monotone_decreasing)
then 281 correct_weight_ordering = .false.
286 bin_grid%edges(b1), bin_grid%edges(b1 + 1), nc1_min, nc1_max)
288 bin_grid%edges(b2), bin_grid%edges(b2 + 1), nc2_min, nc2_max)
291 correct_weight_ordering = .false.
292 if ((nc1_max < nc2_min) .and. (cc == c1))
then 297 correct_weight_ordering = .true.
299 if ((nc2_max < nc1_min) .and. (cc == c2))
then 304 correct_weight_ordering = .true.
315 integer,
intent(in) :: n_part
317 real(kind=dp),
intent(in) :: k_max
319 real(kind=dp),
intent(in) :: del_t
321 real(kind=dp),
intent(out) :: n_source_per_target
323 real(kind=dp),
intent(out) :: accept_factor
325 n_source_per_target = k_max * del_t *
real(n_part, kind=
dp)
326 accept_factor = 1d0 / k_max
334 coag_kernel_type, bs, cs, coag_particle, n_samp_mean, &
335 accept_factor, n_samp, n_coag, n_remove, source_particle)
344 integer,
intent(in) :: coag_kernel_type
346 integer,
intent(in) :: bs
348 integer,
intent(in) :: cs
352 real(kind=dp),
intent(in) :: n_samp_mean
354 real(kind=dp),
intent(in) :: accept_factor
356 integer,
intent(out) :: n_samp
358 integer,
intent(out) :: n_coag
360 integer,
intent(out) :: n_remove
364 real(kind=dp) :: prob_remove_i, prob_remove_source_max
365 real(kind=dp) :: prob_coag, prob_coag_tot, prob_coag_mean
366 real(kind=dp) :: num_conc_i, num_conc_source_min, num_conc_target, k
370 real(kind=dp) :: mean_95_conf_cv
371 integer :: n_samp_remove, n_samp_extra, n_samp_total, n_avg, i_samp
372 integer :: i_unif_entry, i_part, target_id, new_bin, ct
373 type(aero_info_t) :: aero_info
375 if (integer_varray_n_entry( &
376 aero_state%aero_sorted%size_class%inverse(bs, cs)) == 0)
then 384 coag_particle, aero_data)
385 target_id = coag_particle%id
386 ct = coag_particle%weight_class
389 aero_state%awa, cs, aero_state%aero_sorted%bin_grid%edges(bs))
390 prob_remove_source_max = num_conc_target / num_conc_source_min
391 call assert(653606684, prob_remove_source_max <= 1d0)
393 n_samp_remove =
rand_poisson(prob_remove_source_max * n_samp_mean)
394 n_samp_extra =
rand_poisson((1d0 - prob_remove_source_max) * n_samp_mean)
395 n_samp_total = n_samp_remove + n_samp_extra
407 do i_samp = 1,n_samp_total
408 if (integer_varray_n_entry( &
409 aero_state%aero_sorted%size_class%inverse(bs, cs)) == 0)
exit 410 if ((n_samp > n_samp_remove) .and. (n_avg >= 2))
then 411 vol_mean = source_particle%vol /
real(n_avg, kind=
dp)
412 where(vol_mean > 0d0) &
413 vol_cv = sqrt(vol_sq /
real(n_avg, kind=dp) - (vol_mean)**2) &
415 vol_cv_max = maxval(vol_cv)
416 mean_95_conf_cv = vol_cv_max / sqrt(
real(n_avg, kind=
dp)) &
425 aero_state%aero_sorted%size_class%inverse(bs, cs)))
426 i_part = aero_state%aero_sorted%size_class%inverse(bs, &
427 cs)%entry(i_unif_entry)
430 aero_state%apa%particle(i_part), coag_particle, cs, ct, ct, &
431 aero_data, aero_state%awa, env_state, k)
432 prob_coag = k * accept_factor
433 prob_coag_tot = prob_coag_tot + prob_coag
437 aero_state%apa%particle(i_part), source_particle)
438 vol_sq = vol_sq + aero_state%apa%particle(i_part)%vol**2
439 if (i_samp <= n_samp_remove)
then 441 aero_state%apa%particle(i_part), aero_data)
442 prob_remove_i = num_conc_target / num_conc_i
443 if (
pmc_random() < prob_remove_i / prob_remove_source_max)
then 444 n_remove = n_remove + 1
445 aero_info%id = aero_state%apa%particle(i_part)%id
446 aero_info%action = aero_info_coag
447 aero_info%other_id = target_id
460 prob_coag_mean = prob_coag_tot /
real(n_samp, kind=dp) 462 "kernel upper bound estimation is too tight")
464 source_particle%vol = source_particle%vol &
465 * (
real(n_coag, kind=dp) /
real(n_avg, kind=
dp))
473 target_unif_entry, source_particle, cc)
480 integer,
intent(in) :: bt
482 integer,
intent(in) :: ct
484 integer,
intent(in) :: target_unif_entry
488 integer,
intent(in) :: cc
490 integer :: target_part, target_id, new_bin, new_group
491 real(kind=dp) :: old_num_conc_target, new_num_conc_target
493 target_part = aero_state%aero_sorted%size_class%inverse(bt, &
494 ct)%entry(target_unif_entry)
495 target_id = aero_state%apa%particle(target_part)%id
497 aero_state%apa%particle(target_part), aero_data)
499 source_particle, aero_state%apa%particle(target_part))
500 aero_state%apa%particle(target_part)%id = target_id
509 aero_state%apa%particle(target_part), aero_data)
511 .or. (new_bin >
bin_grid_size(aero_state%aero_sorted%bin_grid)))
then 512 call die_msg(765620746,
"particle outside of bin_grid: " &
513 //
"try reducing the timestep del_t")
516 new_bin, new_group, cc)
524 aero_state%apa%particle(target_part), aero_data)
526 old_num_conc_target / new_num_conc_target, random_weight_group=.true.)
528 call assert(654300924, aero_state%apa%particle(target_part)%id &
536 subroutine per_set_coag(coag_kernel_type, k_max, env_state, aero_data, &
537 aero_state, del_t, tot_n_samp, tot_n_coag, b1, b2, c1, c2, cc)
540 integer,
intent(in) :: coag_kernel_type
542 real(kind=dp),
intent(in) :: k_max
550 real(kind=dp) :: del_t
552 integer,
intent(inout) :: tot_n_samp
554 integer,
intent(inout) :: tot_n_coag
556 integer,
intent(in) :: b1
558 integer,
intent(in) :: b2
560 integer,
intent(in) :: c1
562 integer,
intent(in) :: c2
564 integer,
intent(in) :: cc
566 real(kind=dp) :: n_samp_mean, accept_factor
567 integer :: i_samp, n_samp, n1, n2
570 n1 = integer_varray_n_entry( &
571 aero_state%aero_sorted%size_class%inverse(b1, c1))
572 n2 = integer_varray_n_entry( &
573 aero_state%aero_sorted%size_class%inverse(b2, c2))
574 call compute_n_samp(n1, n2, ((b1 == b2) .and. (c1 == c2)), k_max, del_t, &
575 n_samp_mean, n_samp, accept_factor)
576 tot_n_samp = tot_n_samp + n_samp
580 n1 = integer_varray_n_entry( &
581 aero_state%aero_sorted%size_class%inverse(b1, c1))
582 n2 = integer_varray_n_entry( &
583 aero_state%aero_sorted%size_class%inverse(b2, c2))
584 if (((n1 < 2) .and. (b1 == b2) .and. (c1 == c2)) &
585 .or. (n1 < 1) .or. (n2 < 1)) &
587 call maybe_coag_pair(env_state, aero_data, aero_state, b1, b2, c1, c2, &
588 cc, coag_kernel_type, accept_factor, did_coag)
589 if (did_coag) tot_n_coag = tot_n_coag + 1
597 subroutine compute_n_samp(ni, nj, same_bin, k_max, del_t, n_samp_mean, &
598 n_samp, accept_factor)
601 integer,
intent(in) :: ni
603 integer,
intent(in) :: nj
605 logical,
intent(in) :: same_bin
607 real(kind=dp),
intent(in) :: k_max
609 real(kind=dp),
intent(in) :: del_t
611 real(kind=dp),
intent(out) :: n_samp_mean
613 integer,
intent(out) :: n_samp
615 real(kind=dp),
intent(out) :: accept_factor
617 real(kind=dp) :: r_samp
618 real(kind=dp) :: n_possible
626 n_possible =
real(ni, kind=dp) * (
real(nj, kind=dp) - 1d0) / 2d0
628 n_possible =
real(ni, kind=dp) *
real(nj, kind=
dp)
631 r_samp = k_max * del_t
632 n_samp_mean = r_samp * n_possible
634 accept_factor = 1d0 / k_max
662 c1, c2, cc, coag_kernel_type, accept_factor, did_coag)
671 integer,
intent(in) :: b1
673 integer,
intent(in) :: b2
675 integer,
intent(in) :: c1
677 integer,
intent(in) :: c2
679 integer,
intent(in) :: cc
681 integer,
intent(in) :: coag_kernel_type
683 real(kind=dp),
intent(in) :: accept_factor
685 logical,
intent(out) :: did_coag
689 real(kind=dp) :: p, k
691 call find_rand_pair(aero_state%aero_sorted, b1, b2, c1, c2, i1, i2)
692 p1 = aero_state%aero_sorted%size_class%inverse(b1, c1)%entry(i1)
693 p2 = aero_state%aero_sorted%size_class%inverse(b2, c2)%entry(i2)
695 aero_state%apa%particle(p1), aero_state%apa%particle(p2), &
696 c1, c2, cc, aero_data, aero_state%awa, env_state, k)
697 p = k * accept_factor
699 "kernel upper bound estimation is too tight, " &
700 //
"could be caused by changing env_state")
704 call coagulate(aero_data, aero_state, p1, p2, c1, c2, cc)
720 integer,
intent(in) :: b1
722 integer,
intent(in) :: b2
724 integer,
intent(in) :: c1
726 integer,
intent(in) :: c2
728 integer,
intent(out) :: i1
730 integer,
intent(out) :: i2
733 integer_varray_n_entry(aero_sorted%size_class%inverse(b1, c1)) >= 1)
735 integer_varray_n_entry(aero_sorted%size_class%inverse(b1, c1)))
737 if ((b1 == b2) .and. (c1 == c2))
then 738 call assert(956184336, integer_varray_n_entry( &
739 aero_sorted%size_class%inverse(b2, c2)) >= 2)
741 integer_varray_n_entry(aero_sorted%size_class%inverse(b2, c2)) - 1)
743 i2 = integer_varray_n_entry(aero_sorted%size_class%inverse(b2, c2))
746 call assert(271635751, integer_varray_n_entry( &
747 aero_sorted%size_class%inverse(b2, c2)) >= 1)
749 integer_varray_n_entry(aero_sorted%size_class%inverse(b2, c2)))
760 aero_weight_array, remove_1, remove_2, create_new, id_1_lost, &
761 id_2_lost, aero_info_1, aero_info_2)
771 integer,
intent(in) :: c1
773 integer,
intent(in) :: c2
775 integer,
intent(in) :: cc
781 logical,
intent(out) :: remove_1
783 logical,
intent(out) :: remove_2
785 logical,
intent(out) :: create_new
787 logical,
intent(out) :: id_1_lost
789 logical,
intent(out) :: id_2_lost
791 type(aero_info_t),
intent(inout) :: aero_info_1
793 type(aero_info_t),
intent(inout) :: aero_info_2
795 real(kind=dp) :: r1, r2, rc, nc_min, nc1, nc2, ncc
796 real(kind=dp) :: prob_remove_1, prob_remove_2, prob_create_new
797 integer :: info_other_id, new_group
799 call assert(371947172, pt1%id /= pt2%id)
811 nc_min = min(nc1, nc2, ncc)
812 prob_remove_1 = nc_min / nc1
813 prob_remove_2 = nc_min / nc2
814 prob_create_new = nc_min / ncc
830 if (remove_1 .and. remove_2)
then 846 if (remove_1 .and. (.not. id_1_lost))
then 848 call assert(975059559, id_2_lost .eqv. remove_2)
849 elseif (remove_2 .and. (.not. id_2_lost))
then 851 call assert(246529753, id_1_lost .eqv. remove_1)
854 call assert(852038606, id_1_lost .eqv. remove_1)
855 call assert(254018921, id_2_lost .eqv. remove_2)
857 info_other_id = ptc%id
862 aero_info_1%id = pt1%id
863 aero_info_1%action = aero_info_coag
864 aero_info_1%other_id = info_other_id
866 aero_info_2%id = pt2%id
867 aero_info_2%action = aero_info_coag
868 aero_info_2%other_id = info_other_id
876 subroutine coagulate(aero_data, aero_state, p1, p2, c1, c2, cc)
883 integer,
intent(in) :: p1
885 integer,
intent(in) :: p2
887 integer,
intent(in) :: c1
889 integer,
intent(in) :: c2
891 integer,
intent(in) :: cc
895 type(aero_info_t) :: aero_info_1, aero_info_2
896 logical :: remove_1, remove_2, create_new, id_1_lost, id_2_lost
899 aero_state%apa%particle(p2), ptc, c1, c2, cc, aero_data, &
900 aero_state%awa, remove_1, remove_2, create_new, id_1_lost, &
901 id_2_lost, aero_info_1, aero_info_2)
910 id_2_lost, aero_info_2)
914 id_1_lost, aero_info_1)
919 id_1_lost, aero_info_1)
923 id_2_lost, aero_info_2)
930 allow_resort=.false.)
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...
elemental real(kind=dp) function aero_particle_radius(aero_particle, aero_data)
Total radius of the particle (m).
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...
integer function rand_poisson(mean)
Generate a Poisson-distributed random number with the given mean.
subroutine aero_state_sort(aero_state, aero_data, bin_grid, all_procs_same)
Sorts the particles if necessary.
subroutine aero_state_dup_particle(aero_state, aero_data, i_part, n_part_mean, random_weight_group)
Add copies or remove a particle, with a given mean number of resulting particles. ...
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 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 aero_weight_array_t structure and associated subroutines.
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...
subroutine compute_n_source(n_part, k_max, del_t, n_source_per_target, accept_factor)
real(kind=dp), parameter coag_accel_max_cv
Maximum allowed coefficient-of-variation due to undersampling in accelerated coagulation.
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 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).
real(kind=dp), parameter coag_accel_n_event
Minimum number of coagulation events per large particle for which accelerated coagulation is used...
The env_state_t structure and associated subroutines.
subroutine coag_target_with_source(aero_state, aero_data, bt, ct, target_unif_entry, source_particle, cc)
Coagulate a sampled source particle with a target particle.
integer function aero_sorted_n_class(aero_sorted)
Returns the number of weight classes.
integer function aero_sorted_n_bin(aero_sorted)
Returns the number of size bins.
Generic coagulation kernel.
subroutine aero_particle_set_weight(aero_particle, i_group, i_class)
Sets the aerosol particle weight group.
subroutine aero_particle_new_id(aero_particle)
Assigns a globally-unique new ID number to the particle.
The aero_sorted_t structure and assocated subroutines.
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.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
integer function coag_dest_class(aero_weight_array, aero_data, bin_grid, i_bin, j_bin, i_class, j_class)
Determine the weight class in which coagulated particles will be placed.
integer function aero_sorted_particle_in_bin(aero_sorted, aero_particle, aero_data)
Find the bin number that contains a given particle.
real(kind=dp) function pmc_random()
Returns a random number between 0 and 1.
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).
subroutine max_coag_num_conc_factor(aero_weight_array, aero_data, bin_grid, i_bin, j_bin, i_class, j_class, ij_class, f_max)
Determine the minimum and maximum number concentration factors for coagulation.
The aero_state_t structure and assocated subroutines.
Aerosol particle coagulation.
Current environment state.
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.
subroutine aero_state_remove_particle_with_info(aero_state, i_part, aero_info)
Remove the given particle and record the removal.
An array of aerosol size distribution weighting functions.
real(kind=dp) function aero_weight_array_num_conc(aero_weight_array, aero_particle, aero_data)
Compute the number concentration for a particle (m^{-3}).
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
subroutine die_msg(code, error_msg)
Error immediately.
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.
The current collection of aerosol particles.
integer function pmc_rand_int(n)
Returns a random integer between 1 and n.
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.
The aero_data_t structure and associated subroutines.
Single aerosol particle data structure.
1D grid, either logarithmic or linear.
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 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.
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.
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...
elemental real(kind=dp) function aero_particle_volume(aero_particle)
Total volume of the particle (m^3).
Sorting of particles into bins.
subroutine aero_particle_zero(aero_particle, aero_data)
Resets an aero_particle to be zero.
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...
integer, parameter dp
Kind of a double precision real number.
real(kind=dp) elemental function aero_data_vol2rad(aero_data, v)
Convert mass-equivalent volume (m^3) to geometric radius (m).
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 stats_t type and associated subroutines.
elemental integer function aero_data_n_spec(aero_data)
Return the number of aerosol species, or -1 if uninitialized.
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...
Aerosol material properties and associated data.
subroutine aero_state_add_particle(aero_state, aero_particle, aero_data, allow_resort)
Add the given particle.
Common utility subroutines.
Wrapper functions for MPI.
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.
integer function rand_binomial(n, p)
Generate a Binomial-distributed random number with the given parameters.
subroutine aero_state_remove_particle(aero_state, i_part, record_removal, aero_info)
Remove the given particle and possibly record the removal.
The aero_weight_t structure and associated subroutines.
subroutine aero_sorted_move_particle(aero_sorted, i_part, new_bin, new_group, new_class)
Move a particle to a different bin and group.
real(kind=dp), parameter max_allowable_growth_factor
Maximum allowable growth factor of a target particle volume within one timestep when using accelerate...
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...