PartMC  2.2.0
Public Member Functions | Public Attributes
pmc_coagulation Module Reference

Aerosol particle coagulation. More...

List of all members.

Public Member Functions

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 mc_coag_for_bin (coag_kernel_type, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag, i_bin, j_bin)
 Do coagulation for time del_t for the given bins.
subroutine try_per_particle_coag (coag_kernel_type, k_max, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag, i_bin, j_bin, per_particle_coag_succeeded)
 Attempt per-particle coagulation.
subroutine determine_target_and_source (aero_weight_array, bin_grid, i_bin, j_bin, target_bin, source_bin, correct_weight_ordering)
 Determine the source and target particle bin/group for per-particle coagulation, if possible.
subroutine compute_n_source (n_part, k_max, del_t, n_source_per_target, accept_factor)
subroutine sample_source_particle (aero_state, aero_data, env_state, coag_kernel_type, source_bin, coag_particle, n_samp_mean, accept_factor, n_samp, n_coag, n_remove, source_particle)
 Sample coagulation partners for a single coagulation event.
subroutine coag_target_with_source (aero_state, target_bin, target_unif_entry, source_particle)
 Coagulate a sampled source particle with a target particle.
subroutine per_set_coag (coag_kernel_type, k_max, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag, i_bin, j_bin)
 Do set-wise coagulation.
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.
subroutine maybe_coag_pair (env_state, aero_data, aero_state, b1, b2, 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.
subroutine find_rand_pair (aero_sorted, b1, b2, s1, s2)
 Given bins b1 and b2, find a random pair of particles (b1, s1) and (b2, s2) that are not the same particle particle as each other.
subroutine coagulate_weighting (particle_1, particle_2, particle_new, 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 particle_1 and particle_2 to form particle_new and compute weighting effects, including which particles should be lost and which gained.
subroutine coagulate (aero_data, aero_state, p1, p2)
 Join together particles p1 and p2, updating all particle and bin structures to reflect the change.

Public Attributes

real(kind=dp), parameter COAG_ACCEL_N_EVENT = 1d0
 Minimum number of coagulation events per large particle for which accelerated coagulation is used.
real(kind=dp), parameter COAG_ACCEL_MAX_CV = 0.1d0
 Maximum allowed coefficient-of-variation due to undersampling in accelerated coagulation.

Detailed Description

Aerosol particle coagulation.

Definition at line 9 of file coagulation.F90.


Member Function/Subroutine Documentation

subroutine pmc_coagulation::coag_target_with_source ( type(aero_state_t), intent(inout)  aero_state,
integer, intent(in)  target_bin,
integer, intent(in)  target_unif_entry,
type(aero_particle_t), intent(in)  source_particle 
)

Coagulate a sampled source particle with a target particle.

Parameters:
[in,out]aero_stateAerosol state.
[in]target_binBin of coagulating particle.
[in]target_unif_entryEntry-in-bin of coagulating particle.
[in]source_particleSampled particle to coagulate with.

Definition at line 439 of file coagulation.F90.

subroutine pmc_coagulation::coagulate ( type(aero_data_t), intent(in)  aero_data,
type(aero_state_t), intent(inout)  aero_state,
integer, intent(in)  p1,
integer, intent(in)  p2 
)

Join together particles p1 and p2, updating all particle and bin structures to reflect the change.

Parameters:
[in]aero_dataAerosol data.
[in,out]aero_stateAerosol state.
[in]p1First particle index.
[in]p2Second particle index.

Definition at line 813 of file coagulation.F90.

subroutine pmc_coagulation::coagulate_weighting ( type(aero_particle_t), intent(in)  particle_1,
type(aero_particle_t), intent(in)  particle_2,
type(aero_particle_t), intent(inout)  particle_new,
type(aero_data_t), intent(in)  aero_data,
type(aero_weight_t), dimension(:), intent(in)  aero_weight_array,
logical, intent(out)  remove_1,
logical, intent(out)  remove_2,
logical, intent(out)  create_new,
logical, intent(out)  id_1_lost,
logical, intent(out)  id_2_lost,
type(aero_info_t), intent(inout)  aero_info_1,
type(aero_info_t), intent(inout)  aero_info_2 
)

Actually coagulate particle_1 and particle_2 to form particle_new and compute weighting effects, including which particles should be lost and which gained.

Parameters:
[in]particle_1First coagulating aerosol particle.
[in]particle_2Second coagulating aerosol particle.
[in,out]particle_newCombined aerosol particle resulting from coagulation of particle_1 and particle_2.
[in]aero_dataAerosol data.
[in]aero_weight_arrayAerosol weight array.
[out]remove_1Whether to remove particle_1.
[out]remove_2Whether to remove particle_2.
[out]create_newWhether to create particle_new.
[out]id_1_lostWhether the ID of particle_1 will be lost due to coagulation.
[out]id_2_lostWhether the ID of particle_2 will be lost due to coagulation.
[in,out]aero_info_1The removal information associated with particle_1.
[in,out]aero_info_2The removal information associated with particle_2.

Definition at line 694 of file coagulation.F90.

subroutine pmc_coagulation::compute_n_samp ( integer, intent(in)  ni,
integer, intent(in)  nj,
logical, intent(in)  same_bin,
real(kind=dp), intent(in)  k_max,
real(kind=dp), intent(in)  del_t,
real(kind=dp), intent(out)  n_samp_mean,
integer, intent(out)  n_samp,
real(kind=dp), intent(out)  accept_factor 
)

Compute the number of samples required for the pair of bins.

Parameters:
[in]niNumber particles in first bin.
[in]njNumber particles in second bin.
[in]same_binWhether first bin is second bin.
[in]k_maxMaximum kernel value (s^{-1} m^3).
[in]del_tTimestep (s).
[out]n_samp_meanMean number of samples per timestep.
[out]n_sampNumber of samples per timestep.
[out]accept_factorScale factor for accept probability (1).

Definition at line 547 of file coagulation.F90.

subroutine pmc_coagulation::compute_n_source ( integer, intent(in)  n_part,
real(kind=dp), intent(in)  k_max,
real(kind=dp), intent(in)  del_t,
real(kind=dp), intent(out)  n_source_per_target,
real(kind=dp), intent(out)  accept_factor 
)
Parameters:
[in]n_partNumber of particles available as partners.
[in]k_maxMaximum coagulation kernel (s^{-1} m^3).
[in]del_tTimestep (s).
[out]n_source_per_targetMean number of source particles per target particle.
[out]accept_factorAccept factor for samples.

Definition at line 280 of file coagulation.F90.

subroutine pmc_coagulation::determine_target_and_source ( type(aero_weight_t), dimension(:), intent(in)  aero_weight_array,
type(bin_grid_t), intent(in)  bin_grid,
integer, intent(in)  i_bin,
integer, intent(in)  j_bin,
integer, intent(out)  target_bin,
integer, intent(out)  source_bin,
logical, intent(out)  correct_weight_ordering 
)

Determine the source and target particle bin/group for per-particle coagulation, if possible.

Parameters:
[in]aero_weight_arrayAero weight array.
[in]bin_gridBin grid.
[in]i_binFirst bin number.
[in]j_binSecond bin number.
[out]target_binTarget bin number.
[out]source_binSource bin number.
[out]correct_weight_orderingWhether the weight ordering is correct for per-particle coagulation.

Definition at line 228 of file coagulation.F90.

subroutine pmc_coagulation::find_rand_pair ( type(aero_sorted_t), intent(in)  aero_sorted,
integer, intent(in)  b1,
integer, intent(in)  b2,
integer, intent(out)  s1,
integer, intent(out)  s2 
)

Given bins b1 and b2, find a random pair of particles (b1, s1) and (b2, s2) that are not the same particle particle as each other.

Parameters:
[in]aero_sortedAerosol sorted data.
[in]b1Bin number of first particle.
[in]b2Bin number of second particle.
[out]s1First rand particle.
[out]s2Second rand particle.

Definition at line 660 of file coagulation.F90.

subroutine pmc_coagulation::maybe_coag_pair ( type(env_state_t), intent(in)  env_state,
type(aero_data_t), intent(in)  aero_data,
type(aero_state_t), intent(inout)  aero_state,
integer, intent(in)  b1,
integer, intent(in)  b2,
integer, intent(in)  coag_kernel_type,
real(kind=dp), intent(in)  accept_factor,
logical, intent(out)  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 probability of a coagulation will be taken as (kernel / k_max).

Parameters:
[in]env_stateEnvironment state.
[in]aero_dataAerosol data.
[in,out]aero_stateAerosol state.
[in]b1Bin of first particle.
[in]b2Bin of second particle.
[in]coag_kernel_typeCoagulation kernel type.
[in]accept_factorScale factor for accept probability (1).
[out]did_coagWhether a coagulation occured.

Definition at line 612 of file coagulation.F90.

subroutine pmc_coagulation::mc_coag ( integer, intent(in)  coag_kernel_type,
type(env_state_t), intent(in)  env_state,
type(aero_data_t), intent(in)  aero_data,
type(aero_state_t), intent(inout)  aero_state,
real(kind=dp)  del_t,
integer, intent(out)  tot_n_samp,
integer, intent(out)  tot_n_coag 
)

Do coagulation for time del_t.

Parameters:
[in]coag_kernel_typeCoagulation kernel type.
[in]env_stateEnvironment state.
[in]aero_dataAerosol data.
[in,out]aero_stateAerosol state.
del_tTimestep for coagulation.
[out]tot_n_sampTotal number of samples tested.
[out]tot_n_coagNumber of coagulation events.

Definition at line 36 of file coagulation.F90.

subroutine pmc_coagulation::mc_coag_for_bin ( integer, intent(in)  coag_kernel_type,
type(env_state_t), intent(in)  env_state,
type(aero_data_t), intent(in)  aero_data,
type(aero_state_t), intent(inout)  aero_state,
real(kind=dp)  del_t,
integer, intent(inout)  tot_n_samp,
integer, intent(inout)  tot_n_coag,
integer, intent(in)  i_bin,
integer, intent(in)  j_bin 
)

Do coagulation for time del_t for the given bins.

Parameters:
[in]coag_kernel_typeCoagulation kernel type.
[in]env_stateEnvironment state.
[in]aero_dataAerosol data.
[in,out]aero_stateAerosol state.
del_tTimestep for coagulation.
[in,out]tot_n_sampTotal number of samples tested.
[in,out]tot_n_coagNumber of coagulation events.
[in]i_binFirst bin number.
[in]j_binSecond bin number.

Definition at line 93 of file coagulation.F90.

subroutine pmc_coagulation::per_set_coag ( integer, intent(in)  coag_kernel_type,
real(kind=dp), intent(in)  k_max,
type(env_state_t), intent(in)  env_state,
type(aero_data_t), intent(in)  aero_data,
type(aero_state_t), intent(inout)  aero_state,
real(kind=dp)  del_t,
integer, intent(inout)  tot_n_samp,
integer, intent(inout)  tot_n_coag,
integer, intent(in)  i_bin,
integer, intent(in)  j_bin 
)

Do set-wise coagulation.

Parameters:
[in]coag_kernel_typeCoagulation kernel type.
[in]k_maxMaximum coagulation kernel (s^{-1} m^3).
[in]env_stateEnvironment state.
[in]aero_dataAerosol data.
[in,out]aero_stateAerosol state.
del_tTimestep for coagulation.
[in,out]tot_n_sampTotal number of samples tested.
[in,out]tot_n_coagNumber of coagulation events.
[in]i_binFirst bin number.
[in]j_binSecond bin number.

Definition at line 496 of file coagulation.F90.

subroutine pmc_coagulation::sample_source_particle ( type(aero_state_t), intent(inout)  aero_state,
type(aero_data_t), intent(in)  aero_data,
type(env_state_t), intent(in)  env_state,
integer, intent(in)  coag_kernel_type,
integer, intent(in)  source_bin,
type(aero_particle_t), intent(in)  coag_particle,
real(kind=dp), intent(in)  n_samp_mean,
real(kind=dp), intent(in)  accept_factor,
integer, intent(out)  n_samp,
integer, intent(out)  n_coag,
integer, intent(out)  n_remove,
type(aero_particle_t), intent(inout)  source_particle 
)

Sample coagulation partners for a single coagulation event.

Parameters:
[in,out]aero_stateAerosol state.
[in]aero_dataAerosol data.
[in]env_stateEnvironment state.
[in]coag_kernel_typeCoagulation kernel type.
[in]source_binBin to sample particles from.
[in]coag_particleAerosol particle that coagulation will be with.
[in]n_samp_meanMean number of samples to use.
[in]accept_factorProbability factor of accepting samples.
[out]n_sampNumber of samples used.
[out]n_coagNumber of coagulations performed.
[out]n_removeNumber of particles removed.
[in,out]source_particleSampled average coagulation partner particle.

Definition at line 302 of file coagulation.F90.

subroutine pmc_coagulation::try_per_particle_coag ( integer, intent(in)  coag_kernel_type,
real(kind=dp), intent(in)  k_max,
type(env_state_t), intent(in)  env_state,
type(aero_data_t), intent(in)  aero_data,
type(aero_state_t), intent(inout)  aero_state,
real(kind=dp)  del_t,
integer, intent(inout)  tot_n_samp,
integer, intent(inout)  tot_n_coag,
integer, intent(in)  i_bin,
integer, intent(in)  j_bin,
logical, intent(inout)  per_particle_coag_succeeded 
)

Attempt per-particle coagulation.

Parameters:
[in]coag_kernel_typeCoagulation kernel type.
[in]k_maxMaximum coagulation kernel (s^{-1} m^3).
[in]env_stateEnvironment state.
[in]aero_dataAerosol data.
[in,out]aero_stateAerosol state.
del_tTimestep for coagulation.
[in,out]tot_n_sampTotal number of samples tested.
[in,out]tot_n_coagNumber of coagulation events.
[in]i_binFirst bin number.
[in]j_binSecond bin number.
[in,out]per_particle_coag_succeededWhether we succeeded in doing per-particle coag.

Definition at line 139 of file coagulation.F90.


Member Data Documentation

real(kind=dp), parameter pmc_coagulation::COAG_ACCEL_MAX_CV = 0.1d0

Maximum allowed coefficient-of-variation due to undersampling in accelerated coagulation.

Definition at line 29 of file coagulation.F90.

real(kind=dp), parameter pmc_coagulation::COAG_ACCEL_N_EVENT = 1d0

Minimum number of coagulation events per large particle for which accelerated coagulation is used.

Definition at line 26 of file coagulation.F90.


The documentation for this module was generated from the following file: