PartMC  2.3.0
coag_kernel.F90
Go to the documentation of this file.
1 ! Copyright (C) 2005-2012 Nicole Riemer and Matthew West
2 ! Licensed under the GNU General Public License version 2 or (at your
3 ! option) any later version. See the file COPYING for details.
4 
5 !> \file
6 !> The pmc_coag_kernel module.
7 
8 !> Generic coagulation kernel.
10 
11  use pmc_env_state
12  use pmc_bin_grid
14  use pmc_aero_data
15  use pmc_aero_weight
22 
23  !> Maximum length of a mode type.
24  integer, parameter :: COAG_KERNEL_TYPE_LEN = 20
25 
26  !> Type code for an undefined or invalid kernel.
27  integer, parameter :: COAG_KERNEL_TYPE_INVALID = 0
28  !> Type code for a sedimentation kernel.
29  integer, parameter :: COAG_KERNEL_TYPE_SEDI = 1
30  !> Type code for an additive kernel.
31  integer, parameter :: COAG_KERNEL_TYPE_ADDITIVE = 2
32  !> Type code for a constant kernel.
33  integer, parameter :: COAG_KERNEL_TYPE_CONSTANT = 3
34  !> Type code for a Brownian kernel.
35  integer, parameter :: COAG_KERNEL_TYPE_BROWN = 4
36  !> Type code for a zero kernel.
37  integer, parameter :: COAG_KERNEL_TYPE_ZERO = 5
38 
39 contains
40 
41 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
42 
43  !> Return a string representation of a kernel type.
44  character(len=COAG_KERNEL_TYPE_LEN) function coag_kernel_type_to_string( &
45  coag_kernel_type)
46 
47  !> Coagulation kernel type.
48  integer, intent(in) :: coag_kernel_type
49 
50  if (coag_kernel_type == coag_kernel_type_invalid) then
51  coag_kernel_type_to_string = "invalid"
52  elseif (coag_kernel_type == coag_kernel_type_sedi) then
54  elseif (coag_kernel_type == coag_kernel_type_additive) then
55  coag_kernel_type_to_string = "additive"
56  elseif (coag_kernel_type == coag_kernel_type_constant) then
57  coag_kernel_type_to_string = "constant"
58  elseif (coag_kernel_type == coag_kernel_type_brown) then
60  elseif (coag_kernel_type == coag_kernel_type_zero) then
62  else
63  coag_kernel_type_to_string = "unknown"
64  end if
65 
66  end function coag_kernel_type_to_string
67 
68 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
69 
70  !> Evalulate a coagulation kernel function.
71  subroutine kernel(coag_kernel_type, aero_particle_1, aero_particle_2, &
72  aero_data, env_state, k)
73 
74  !> Coagulation kernel type.
75  integer, intent(in) :: coag_kernel_type
76  !> First particle.
77  type(aero_particle_t), intent(in) :: aero_particle_1
78  !> Second particle.
79  type(aero_particle_t), intent(in) :: aero_particle_2
80  !> Aerosol data.
81  type(aero_data_t), intent(in) :: aero_data
82  !> Environment state.
83  type(env_state_t), intent(in) :: env_state
84  !> Kernel k(a,b) (m^3/s).
85  real(kind=dp), intent(out) :: k
86 
87  if (coag_kernel_type == coag_kernel_type_sedi) then
88  call kernel_sedi(aero_particle_1, aero_particle_2, &
89  aero_data, env_state, k)
90  elseif (coag_kernel_type == coag_kernel_type_additive) then
91  call kernel_additive(aero_particle_1, aero_particle_2, &
92  aero_data, env_state, k)
93  elseif (coag_kernel_type == coag_kernel_type_constant) then
94  call kernel_constant(aero_particle_1, aero_particle_2, &
95  aero_data, env_state, k)
96  elseif (coag_kernel_type == coag_kernel_type_brown) then
97  call kernel_brown(aero_particle_1, aero_particle_2, &
98  aero_data, env_state, k)
99  elseif (coag_kernel_type == coag_kernel_type_zero) then
100  call kernel_zero(aero_particle_1, aero_particle_2, &
101  aero_data, env_state, k)
102  else
103  call die_msg(200724934, "Unknown kernel type: " &
104  // trim(integer_to_string(coag_kernel_type)))
105  end if
106 
107  end subroutine kernel
108 
109 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
110 
111  !> Compute the minimum and maximum coagulation kernel.
112  subroutine kernel_minmax(coag_kernel_type, v1, v2, aero_data, env_state, &
113  k_min, k_max)
114 
115  !> Coagulation kernel type.
116  integer, intent(in) :: coag_kernel_type
117  !> Volume of first particle (m^3).
118  real(kind=dp), intent(in) :: v1
119  !> Volume of second particle (m^3).
120  real(kind=dp), intent(in) :: v2
121  !> Aerosol data.
122  type(aero_data_t), intent(in) :: aero_data
123  !> Environment state.
124  type(env_state_t), intent(in) :: env_state
125  !> Minimum kernel value (m^3/s).
126  real(kind=dp), intent(out) :: k_min
127  !> Maximum kernel value (m^3/s).
128  real(kind=dp), intent(out) :: k_max
129 
130  if (coag_kernel_type == coag_kernel_type_sedi) then
131  call kernel_sedi_minmax(v1, v2, aero_data, env_state, k_min, k_max)
132  elseif (coag_kernel_type == coag_kernel_type_additive) then
133  call kernel_additive_minmax(v1, v2, aero_data, env_state, k_min, k_max)
134  elseif (coag_kernel_type == coag_kernel_type_constant) then
135  call kernel_constant_minmax(v1, v2, aero_data, env_state, k_min, k_max)
136  elseif (coag_kernel_type == coag_kernel_type_brown) then
137  call kernel_brown_minmax(v1, v2, aero_data, env_state, k_min, k_max)
138  elseif (coag_kernel_type == coag_kernel_type_zero) then
139  call kernel_zero_minmax(v1, v2, aero_data, env_state, k_min, k_max)
140  else
141  call die_msg(330498208, "Unknown kernel type: " &
142  // trim(integer_to_string(coag_kernel_type)))
143  end if
144 
145  end subroutine kernel_minmax
146 
147 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
148 
149  !> Compute the kernel value with the given number concentration
150  !> weighting.
151  subroutine num_conc_weighted_kernel(coag_kernel_type, aero_particle_1, &
152  aero_particle_2, i_class, j_class, ij_class, aero_data, &
153  aero_weight_array, env_state, k)
154 
155  !> Coagulation kernel type.
156  integer, intent(in) :: coag_kernel_type
157  !> First particle.
158  type(aero_particle_t), intent(in) :: aero_particle_1
159  !> Second particle.
160  type(aero_particle_t), intent(in) :: aero_particle_2
161  !> Weight class of first particle.
162  integer, intent(in) :: i_class
163  !> Weight class of second particle.
164  integer, intent(in) :: j_class
165  !> Weight class of combined particle.
166  integer, intent(in) :: ij_class
167  !> Aerosol data.
168  type(aero_data_t), intent(in) :: aero_data
169  !> Aerosol weight array.
170  type(aero_weight_array_t), intent(in) :: aero_weight_array
171  !> Environment state.
172  type(env_state_t), intent(in) :: env_state
173  !> Coagulation kernel.
174  real(kind=dp), intent(out) :: k
175 
176  real(kind=dp) :: unweighted_k, i_r, j_r
177 
178  call kernel(coag_kernel_type, aero_particle_1, aero_particle_2, &
179  aero_data, env_state, unweighted_k)
180  i_r = aero_particle_radius(aero_particle_1)
181  j_r = aero_particle_radius(aero_particle_2)
182  k = unweighted_k * coag_num_conc_factor(aero_weight_array, i_r, j_r, &
183  i_class, j_class, ij_class)
184 
185  end subroutine num_conc_weighted_kernel
186 
187 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
188 
189  !> Computes an array of kernel values for each bin pair. k(i,j) is
190  !> the kernel value at the centers of bins i and j. This assumes the
191  !> kernel is only a function of the particle volumes.
192  subroutine bin_kernel(n_bin, bin_r, aero_data, coag_kernel_type, &
193  env_state, k)
194 
195  !> Number of bins.
196  integer, intent(in) :: n_bin
197  !> Radii of particles in bins (m).
198  real(kind=dp), intent(in) :: bin_r(n_bin)
199  !> Aerosol data.
200  type(aero_data_t), intent(in) :: aero_data
201  !> Coagulation kernel type.
202  integer, intent(in) :: coag_kernel_type
203  !> Environment state.
204  type(env_state_t), intent(in) :: env_state
205  !> Kernel values.
206  real(kind=dp), intent(out) :: k(n_bin,n_bin)
207 
208  integer :: i, j
209  type(aero_particle_t) :: aero_particle_1, aero_particle_2
210 
211  call aero_particle_allocate_size(aero_particle_1, aero_data%n_spec, &
212  aero_data%n_source)
213  call aero_particle_allocate_size(aero_particle_2, aero_data%n_spec, &
214  aero_data%n_source)
215  do i = 1,n_bin
216  do j = 1,n_bin
217  aero_particle_1%vol(1) = rad2vol(bin_r(i))
218  aero_particle_2%vol(1) = rad2vol(bin_r(j))
219  call kernel(coag_kernel_type, aero_particle_1, aero_particle_2, &
220  aero_data, env_state, k(i,j))
221  end do
222  end do
223  call aero_particle_deallocate(aero_particle_1)
224  call aero_particle_deallocate(aero_particle_2)
225 
226  end subroutine bin_kernel
227 
228 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
229 
230  !> Estimate an array of minimum and maximum kernel values. Given
231  !> particles v1 in bin b1 and v2 in bin b2, it is probably true that
232  !> <tt>k_min(b1,b2) <= kernel(v1,v2) <= k_max(b1,b2)</tt>.
233  subroutine est_k_minmax_binned_unweighted(bin_grid, coag_kernel_type, &
234  aero_data, env_state, k_min, k_max)
235 
236  !> Bin_grid.
237  type(bin_grid_t), intent(in) :: bin_grid
238  !> Coagulation kernel type.
239  integer, intent(in) :: coag_kernel_type
240  !> Aerosol data.
241  type(aero_data_t), intent(in) :: aero_data
242  !> Environment state.
243  type(env_state_t), intent(in) :: env_state
244  !> Minimum kernel vals.
245  real(kind=dp), intent(out) :: k_min(bin_grid%n_bin,bin_grid%n_bin)
246  !> Maximum kernel vals.
247  real(kind=dp), intent(out) :: k_max(bin_grid%n_bin,bin_grid%n_bin)
248 
249  integer i, j
250 
251  do i = 1,bin_grid%n_bin
252  do j = 1,bin_grid%n_bin
253  call est_k_minmax_for_bin_unweighted(bin_grid, coag_kernel_type, &
254  i, j, aero_data, env_state, k_min(i,j), k_max(i,j))
255  end do
256  end do
257 
258  end subroutine est_k_minmax_binned_unweighted
259 
260 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
261 
262  !> Samples within bins b1 and b2 to find the minimum and maximum
263  !> value of the kernel between particles from the two bins.
264  subroutine est_k_minmax_for_bin_unweighted(bin_grid, coag_kernel_type, &
265  b1, b2, aero_data, env_state, k_min, k_max)
266 
267  !> Bin_grid.
268  type(bin_grid_t), intent(in) :: bin_grid
269  !> Coagulation kernel type.
270  integer, intent(in) :: coag_kernel_type
271  !> First bin.
272  integer, intent(in) :: b1
273  !> Second bin.
274  integer, intent(in) :: b2
275  !> Aerosol data.
276  type(aero_data_t), intent(in) :: aero_data
277  !> Environment state.
278  type(env_state_t), intent(in) :: env_state
279  !> Minimum kernel value.
280  real(kind=dp), intent(out) :: k_min
281  !> Maximum kernel value.
282  real(kind=dp), intent(out) :: k_max
283 
284  !> Number of sample points per bin.
285  integer, parameter :: n_sample = 3
286  !> Over-estimation scale factor parameter.
287  real(kind=dp), parameter :: over_scale = 2d0
288 
289  real(kind=dp) :: v1, v2, v1_high, v1_low, v2_high, v2_low
290  real(kind=dp) :: new_k_min, new_k_max
291  integer :: i, j
292 
293  ! v1_low < bin_v(b1) < v1_high
294  v1_low = rad2vol(bin_grid%edges(b1))
295  v1_high = rad2vol(bin_grid%edges(b1 + 1))
296 
297  ! v2_low < bin_v(b2) < v2_high
298  v2_low = rad2vol(bin_grid%edges(b2))
299  v2_high = rad2vol(bin_grid%edges(b2 + 1))
300 
301  do i = 1,n_sample
302  do j = 1,n_sample
303  v1 = interp_linear_disc(v1_low, v1_high, n_sample, i)
304  v2 = interp_linear_disc(v2_low, v2_high, n_sample, j)
305  call kernel_minmax(coag_kernel_type, v1, v2, aero_data, &
306  env_state, new_k_min, new_k_max)
307  if ((i == 1) .and. (j == 1)) then
308  k_min = new_k_min
309  k_max = new_k_max
310  else
311  k_min = min(k_min, new_k_min)
312  k_max = max(k_max, new_k_max)
313  end if
314  end do
315  end do
316 
317  k_max = k_max * over_scale
318 
319  end subroutine est_k_minmax_for_bin_unweighted
320 
321 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
322 
323  !> Coagulation scale factor due to number concentrations.
324  real(kind=dp) function coag_num_conc_factor(aero_weight_array, &
325  i_r, j_r, i_class, j_class, ij_class)
326 
327  !> Aerosol weight array.
328  type(aero_weight_array_t), intent(in) :: aero_weight_array
329  !> Radius of first particle.
330  real(kind=dp), intent(in) :: i_r
331  !> Radius of second particle.
332  real(kind=dp), intent(in) :: j_r
333  !> Weight class of first particle.
334  integer, intent(in) :: i_class
335  !> Weight class of second particle.
336  integer, intent(in) :: j_class
337  !> Weight class of combined particle.
338  integer, intent(in) :: ij_class
339 
340  real(kind=dp) :: ij_r, i_nc, j_nc, ij_nc, nc_min
341 
342  ij_r = vol2rad(rad2vol(i_r) + rad2vol(j_r))
343  i_nc = aero_weight_array_num_conc_at_radius(aero_weight_array, i_class, &
344  i_r)
345  j_nc = aero_weight_array_num_conc_at_radius(aero_weight_array, j_class, &
346  j_r)
347  ij_nc = aero_weight_array_num_conc_at_radius(aero_weight_array, ij_class, &
348  ij_r)
349  nc_min = min(i_nc, j_nc, ij_nc)
350  coag_num_conc_factor = i_nc * j_nc / nc_min
351 
352  end function coag_num_conc_factor
353 
354 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
355 
356  !> Determine the weight class in which coagulated particles will be placed.
357  integer function coag_dest_class(aero_weight_array, bin_grid, i_bin, j_bin, &
358  i_class, j_class)
359 
360  !> Aerosol weight array.
361  type(aero_weight_array_t), intent(in) :: aero_weight_array
362  !> Bin grid.
363  type(bin_grid_t), intent(in) :: bin_grid
364  !> First bin number.
365  integer, intent(in) :: i_bin
366  !> Second bin number.
367  integer, intent(in) :: j_bin
368  !> Weight class of first particle.
369  integer, intent(in) :: i_class
370  !> Weight class of second particle.
371  integer, intent(in) :: j_class
372 
373  real(kind=dp) :: i_r, j_r, ij_r, ij_nc_i, ij_nc_j
374 
375  i_r = bin_grid%centers(i_bin)
376  j_r = bin_grid%centers(i_bin)
377  ij_r = vol2rad(rad2vol(i_r) + rad2vol(j_r))
378  ij_nc_i = aero_weight_array_num_conc_at_radius(aero_weight_array, &
379  i_class, ij_r)
380  ij_nc_j = aero_weight_array_num_conc_at_radius(aero_weight_array, &
381  j_class, ij_r)
382  if (ij_nc_i < ij_nc_j) then
383  coag_dest_class = i_class
384  else
385  coag_dest_class = j_class
386  end if
387 
388  end function coag_dest_class
389 
390 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
391 
392  !> Determine the minimum and maximum number concentration factors
393  !> for coagulation.
394  subroutine max_coag_num_conc_factor(aero_weight_array, bin_grid, &
395  i_bin, j_bin, i_class, j_class, ij_class, f_max)
396 
397  !> Aerosol weight array.
398  type(aero_weight_array_t), intent(in) :: aero_weight_array
399  !> Bin grid.
400  type(bin_grid_t), intent(in) :: bin_grid
401  !> First bin number.
402  integer, intent(in) :: i_bin
403  !> Second bin number.
404  integer, intent(in) :: j_bin
405  !> Weight class of first particle.
406  integer, intent(in) :: i_class
407  !> Weight class of second particle.
408  integer, intent(in) :: j_class
409  !> Weight class of coagulated particle.
410  integer, intent(in) :: ij_class
411  !> Maximum coagulation factor.
412  real(kind=dp), intent(out) :: f_max
413 
414  integer, parameter :: n_sample = 5
415 
416  real(kind=dp) :: i_r_min, i_r_max, j_r_min, j_r_max, i_r, j_r, f
417  integer :: i_sample, j_sample
418 
419  i_r_min = bin_grid%edges(i_bin)
420  i_r_max = bin_grid%edges(i_bin + 1)
421  j_r_min = bin_grid%edges(j_bin)
422  j_r_max = bin_grid%edges(j_bin + 1)
423 
424  f_max = 0d0
425  do i_sample = 1,n_sample
426  do j_sample = 1,n_sample
427  i_r = interp_linear_disc(i_r_min, i_r_max, n_sample, i_sample)
428  j_r = interp_linear_disc(j_r_min, j_r_max, n_sample, j_sample)
429  f = coag_num_conc_factor(aero_weight_array, i_r, j_r, i_class, &
430  j_class, ij_class)
431  f_max = max(f_max, f)
432  end do
433  end do
434 
435  end subroutine max_coag_num_conc_factor
436 
437 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
438 
439  !> Read the specification for a kernel type from a spec file and
440  !> generate it.
441  subroutine spec_file_read_coag_kernel_type(file, coag_kernel_type)
442 
443  !> Spec file.
444  type(spec_file_t), intent(inout) :: file
445  !> Kernel type.
446  integer, intent(out) :: coag_kernel_type
447 
448  character(len=SPEC_LINE_MAX_VAR_LEN) :: kernel_name
449 
450  !> \page input_format_coag_kernel Input File Format: Coagulation Kernel
451  !!
452  !! The coagulation kernel is specified by the parameter:
453  !! - \b coag_kernel (string): the type of coagulation kernel ---
454  !! must be one of: \c sedi for the gravitational sedimentation
455  !! kernel; \c additive for the additive kernel; \c constant
456  !! for the constant kernel; \c brown for the Brownian kernel,
457  !! or \c zero for no coagulation
458  !!
459  !! See also:
460  !! - \ref spec_file_format --- the input file text format
461 
462  call spec_file_read_string(file, 'coag_kernel', kernel_name)
463  if (trim(kernel_name) == 'sedi') then
464  coag_kernel_type = coag_kernel_type_sedi
465  elseif (trim(kernel_name) == 'additive') then
466  coag_kernel_type = coag_kernel_type_additive
467  elseif (trim(kernel_name) == 'constant') then
468  coag_kernel_type = coag_kernel_type_constant
469  elseif (trim(kernel_name) == 'brown') then
470  coag_kernel_type = coag_kernel_type_brown
471  elseif (trim(kernel_name) == 'zero') then
472  coag_kernel_type = coag_kernel_type_zero
473  else
474  call spec_file_die_msg(920761229, file, &
475  "Unknown coagulation kernel type: " // trim(kernel_name))
476  end if
477 
478  end subroutine spec_file_read_coag_kernel_type
479 
480 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
481 
482 end module pmc_coag_kernel
An input file with extra data for printing messages.
Definition: spec_file.F90:59
The aero_data_t structure and associated subroutines.
Definition: aero_data.F90:9
subroutine die_msg(code, error_msg)
Error immediately.
Definition: util.F90:133
Additive coagulation kernel.
subroutine kernel_additive_minmax(v1, v2, aero_data, env_state, k_min, k_max)
Minimum and maximum values of the additive kernel.
real(kind=dp) elemental function rad2vol(r)
Convert radius (m) to volume (m^3).
Definition: util.F90:274
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.
Definition: aero_weight.F90:9
subroutine kernel_additive(aero_particle_1, aero_particle_2, aero_data, env_state, k)
Additive coagulation kernel.
The aero_particle_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.
The env_state_t structure and associated subroutines.
Definition: env_state.F90:9
subroutine aero_particle_allocate_size(aero_particle, n_spec, n_source)
Allocates an aero_particle_t of the given size.
subroutine kernel_constant(aero_particle_1, aero_particle_2, aero_data, env_state, k)
Constant coagulation kernel.
subroutine kernel_minmax(coag_kernel_type, v1, v2, aero_data, env_state, k_min, k_max)
Compute the minimum and maximum coagulation kernel.
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...
subroutine kernel_zero_minmax(v1, v2, aero_data, env_state, k_min, k_max)
Minimum and maximum of the zero coagulation kernel.
real(kind=dp) function interp_linear_disc(x_1, x_n, n, i)
Linear interpolation over discrete indices.
Definition: util.F90:644
character(len=coag_kernel_type_len) function coag_kernel_type_to_string(coag_kernel_type)
Return a string representation of a kernel type.
Definition: coag_kernel.F90:44
Current environment state.
Definition: env_state.F90:26
subroutine est_k_minmax_for_bin_unweighted(bin_grid, coag_kernel_type, b1, b2, aero_data, env_state, k_min, k_max)
Samples within bins b1 and b2 to find the minimum and maximum value of the kernel between particles f...
An array of aerosol size distribution weighting functions.
Gravitational sedimentation coagulation kernel.
subroutine kernel_brown_minmax(v1, v2, aero_data, env_state, k_min, k_max)
Compute the minimum and maximum Brownian coagulation kernel.
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).
Constant coagulation kernel.
subroutine spec_file_read_string(file, name, var)
Read a string from a spec file that must have a given name.
Definition: spec_file.F90:625
The aero_weight_array_t structure and associated subroutines.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
Definition: util.F90:743
Single aerosol particle data structure.
1D grid, either logarithmic or linear.
Definition: bin_grid.F90:33
Brownian coagulation kernel.
The bin_grid_t structure and associated subroutines.
Definition: bin_grid.F90:9
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.
subroutine kernel(coag_kernel_type, aero_particle_1, aero_particle_2, aero_data, env_state, k)
Evalulate a coagulation kernel function.
Definition: coag_kernel.F90:71
subroutine kernel_zero(aero_particle_1, aero_particle_2, aero_data, env_state, k)
Zero coagulation kernel.
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
Definition: spec_file.F90:73
Constant kernel equal to zero.
subroutine kernel_sedi(aero_particle_1, aero_particle_2, aero_data, env_state, k)
Sedimentation coagulation kernel.
subroutine bin_kernel(n_bin, bin_r, aero_data, coag_kernel_type, env_state, k)
Computes an array of kernel values for each bin pair. k(i,j) is the kernel value at the centers of bi...
real(kind=dp) elemental function vol2rad(v)
Convert volume (m^3) to radius (m).
Definition: util.F90:238
subroutine kernel_sedi_minmax(v1, v2, aero_data, env_state, k_min, k_max)
Minimum and maximum values of the sedimentation coagulation.
subroutine kernel_constant_minmax(v1, v2, aero_data, env_state, k_min, k_max)
Minimum and maximum values of the constant coagulation kernel.
subroutine kernel_brown(aero_particle_1, aero_particle_2, aero_data, env_state, k)
Compute the Brownian coagulation kernel.
Aerosol material properties and associated data.
Definition: aero_data.F90:40
real(kind=dp) function coag_num_conc_factor(aero_weight_array, i_r, j_r, i_class, j_class, ij_class)
Coagulation scale factor due to number concentrations.
Generic coagulation kernel.
Definition: coag_kernel.F90:9