PartMC  2.3.0
aero_particle.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_aero_particle module.
7 
8 !> The aero_particle_t structure and associated subroutines.
10 
11  use pmc_util
12  use pmc_aero_data
13  use pmc_spec_file
14  use pmc_env_state
15  use pmc_mpi
16 #ifdef PMC_USE_MPI
17  use mpi
18 #endif
19 
20  !> Single aerosol particle data structure.
21  !!
22  !! The \c vol array stores the total volumes of the different
23  !! species that make up the particle. This array must have length
24  !! equal to aero_data%%n_spec, so that \c vol(i) is the volume (in
25  !! m^3) of the i'th aerosol species.
27  !> Constituent species volumes [length aero_data%%n_spec] (m^3).
28  real(kind=dp), pointer :: vol(:)
29  !> Number of original particles from each source that coagulated
30  !> to form this one [length aero_data%%n_source].
31  integer, pointer :: n_orig_part(:)
32  !> Weighting function group number (see \c aero_weight_array_t).
33  integer :: weight_group
34  !> Weighting function class number (see \c aero_weight_array_t).
35  integer :: weight_class
36  !> Absorption cross-section (m^2).
37  real(kind=dp) :: absorb_cross_sect
38  !> Scattering cross-section (m^2).
39  real(kind=dp) :: scatter_cross_sect
40  !> Asymmetry parameter (1).
41  real(kind=dp) :: asymmetry
42  !> Refractive index of the shell (1).
43  complex(kind=dc) :: refract_shell
44  !> Refractive index of the core (1).
45  complex(kind=dc) :: refract_core
46  !> Volume of the core (m^3).
47  real(kind=dp) :: core_vol
48  !> Water hysteresis curve section (0 = lower, 1 = upper)
49  integer :: water_hyst_leg
50  !> Unique ID number.
51  integer :: id
52  !> First time a constituent was created (s).
53  real(kind=dp) :: least_create_time
54  !> Last time a constituent was created (s).
55  real(kind=dp) :: greatest_create_time
56  end type aero_particle_t
57 
58  !> Next unique ID number to use for a particle.
59  integer, save :: next_id = 1
60 
61 contains
62 
63 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
64 
65  !> Allocates memory in an aero_particle_t.
66  subroutine aero_particle_allocate(aero_particle)
67 
68  !> Particle to init.
69  type(aero_particle_t), intent(out) :: aero_particle
70 
71  allocate(aero_particle%vol(0))
72  allocate(aero_particle%n_orig_part(0))
73  call aero_particle_zero(aero_particle)
74 
75  end subroutine aero_particle_allocate
76 
77 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
78 
79  !> Allocates an aero_particle_t of the given size.
80  subroutine aero_particle_allocate_size(aero_particle, n_spec, n_source)
81 
82  !> Particle to init.
83  type(aero_particle_t), intent(out) :: aero_particle
84  !> Number of species.
85  integer, intent(in) :: n_spec
86  !> Number of sources.
87  integer, intent(in) :: n_source
88 
89  allocate(aero_particle%vol(n_spec))
90  allocate(aero_particle%n_orig_part(n_source))
91  call aero_particle_zero(aero_particle)
92 
93  end subroutine aero_particle_allocate_size
94 
95 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
96 
97  !> Deallocates memory associated with an aero_particle_t.
98  subroutine aero_particle_deallocate(aero_particle)
99 
100  !> Particle to free.
101  type(aero_particle_t), intent(inout) :: aero_particle
102 
103  deallocate(aero_particle%vol)
104  deallocate(aero_particle%n_orig_part)
105 
106  end subroutine aero_particle_deallocate
107 
108 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
109 
110  !> Copies a particle.
111  subroutine aero_particle_copy(aero_particle_from, aero_particle_to)
112 
113  !> Reference particle.
114  type(aero_particle_t), intent(in) :: aero_particle_from
115  !> Destination particle (already alloced on entry).
116  type(aero_particle_t), intent(inout) :: aero_particle_to
117 
118  integer :: n_spec, n_source
119 
120  n_spec = size(aero_particle_from%vol)
121  n_source = size(aero_particle_from%n_orig_part)
122  if ((n_spec /= size(aero_particle_to%vol)) &
123  .or. (n_source /= size(aero_particle_to%n_orig_part))) then
124  call aero_particle_deallocate(aero_particle_to)
125  call aero_particle_allocate_size(aero_particle_to, n_spec, n_source)
126  end if
127  call assert(651178226, size(aero_particle_from%vol) &
128  == size(aero_particle_to%vol))
129  call assert(105980551, size(aero_particle_from%n_orig_part) &
130  == size(aero_particle_to%n_orig_part))
131  aero_particle_to%vol = aero_particle_from%vol
132  aero_particle_to%n_orig_part = aero_particle_from%n_orig_part
133  aero_particle_to%weight_group = aero_particle_from%weight_group
134  aero_particle_to%weight_class = aero_particle_from%weight_class
135  aero_particle_to%absorb_cross_sect = aero_particle_from%absorb_cross_sect
136  aero_particle_to%scatter_cross_sect = &
137  aero_particle_from%scatter_cross_sect
138  aero_particle_to%asymmetry = aero_particle_from%asymmetry
139  aero_particle_to%refract_shell = aero_particle_from%refract_shell
140  aero_particle_to%refract_core = aero_particle_from%refract_core
141  aero_particle_to%core_vol = aero_particle_from%core_vol
142  aero_particle_to%water_hyst_leg = aero_particle_from%water_hyst_leg
143  aero_particle_to%id = aero_particle_from%id
144  aero_particle_to%least_create_time = aero_particle_from%least_create_time
145  aero_particle_to%greatest_create_time = &
146  aero_particle_from%greatest_create_time
147 
148  end subroutine aero_particle_copy
149 
150 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
151 
152  !> Shift data from one aero_particle_t to another and free the first
153  !> one.
154  !!
155  !! This is roughly equivalent to aero_particle_copy(from, to)
156  !! followed by aero_particle_deallocate(from), but faster and with
157  !! different memory allocation requirements.
158  subroutine aero_particle_shift(aero_particle_from, aero_particle_to)
159 
160  !> Reference particle (will be deallocated on return).
161  type(aero_particle_t), intent(inout) :: aero_particle_from
162  !> Destination particle (not allocated on entry).
163  type(aero_particle_t), intent(inout) :: aero_particle_to
164 
165  aero_particle_to%vol => aero_particle_from%vol
166  nullify(aero_particle_from%vol)
167  aero_particle_to%n_orig_part => aero_particle_from%n_orig_part
168  nullify(aero_particle_from%n_orig_part)
169  aero_particle_to%weight_group = aero_particle_from%weight_group
170  aero_particle_to%weight_class = aero_particle_from%weight_class
171  aero_particle_to%absorb_cross_sect = aero_particle_from%absorb_cross_sect
172  aero_particle_to%scatter_cross_sect = &
173  aero_particle_from%scatter_cross_sect
174  aero_particle_to%asymmetry = aero_particle_from%asymmetry
175  aero_particle_to%refract_shell = aero_particle_from%refract_shell
176  aero_particle_to%refract_core = aero_particle_from%refract_core
177  aero_particle_to%core_vol = aero_particle_from%core_vol
178  aero_particle_to%water_hyst_leg = aero_particle_from%water_hyst_leg
179  aero_particle_to%id = aero_particle_from%id
180  aero_particle_to%least_create_time = aero_particle_from%least_create_time
181  aero_particle_to%greatest_create_time = &
182  aero_particle_from%greatest_create_time
183 
184  end subroutine aero_particle_shift
185 
186 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
187 
188  !> Resets an aero_particle to be zero.
189  subroutine aero_particle_zero(aero_particle)
190 
191  !> Particle to zero.
192  type(aero_particle_t), intent(inout) :: aero_particle
193 
194  aero_particle%vol = 0d0
195  aero_particle%n_orig_part = 0
196  aero_particle%weight_group = 0
197  aero_particle%weight_class = 0
198  aero_particle%absorb_cross_sect = 0d0
199  aero_particle%scatter_cross_sect = 0d0
200  aero_particle%asymmetry = 0d0
201  aero_particle%refract_shell = (0d0, 0d0)
202  aero_particle%refract_core = (0d0, 0d0)
203  aero_particle%core_vol = 0d0
204  aero_particle%water_hyst_leg = 0
205  aero_particle%id = 0
206  aero_particle%least_create_time = 0d0
207  aero_particle%greatest_create_time = 0d0
208 
209  end subroutine aero_particle_zero
210 
211 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
212 
213  !> Assigns a globally-unique new ID number to the particle.
214  subroutine aero_particle_new_id(aero_particle)
215 
216  !> Particle to set ID for.
217  type(aero_particle_t), intent(inout) :: aero_particle
218 
219  aero_particle%id = (next_id - 1) * pmc_mpi_size() + pmc_mpi_rank() + 1
220  next_id = next_id + 1
221 
222  end subroutine aero_particle_new_id
223 
224 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
225 
226  !> Sets the creation times for the particle.
227  subroutine aero_particle_set_create_time(aero_particle, create_time)
228 
229  !> Particle to set time for.
230  type(aero_particle_t), intent(inout) :: aero_particle
231  !> Creation time.
232  real(kind=dp), intent(in) :: create_time
233 
234  aero_particle%least_create_time = create_time
235  aero_particle%greatest_create_time = create_time
236 
237  end subroutine aero_particle_set_create_time
238 
239 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
240 
241  !> Sets the aerosol particle volumes.
242  subroutine aero_particle_set_vols(aero_particle, vols)
243 
244  !> Particle.
245  type(aero_particle_t), intent(inout) :: aero_particle
246  !> New volumes.
247  real(kind=dp), intent(in) :: vols(size(aero_particle%vol))
248 
249  aero_particle%vol = vols
250 
251  end subroutine aero_particle_set_vols
252 
253 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
254 
255  !> Sets the aerosol particle source.
256  subroutine aero_particle_set_source(aero_particle, i_source)
257 
258  !> Particle.
259  type(aero_particle_t), intent(inout) :: aero_particle
260  !> Source number for the particle.
261  integer, intent(in) :: i_source
262 
263  aero_particle%n_orig_part(i_source) = 1
264 
265  end subroutine aero_particle_set_source
266 
267 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
268 
269  !> Sets the aerosol particle weight group.
270  subroutine aero_particle_set_weight(aero_particle, i_group, i_class)
271 
272  !> Particle.
273  type(aero_particle_t), intent(inout) :: aero_particle
274  !> Weight group number for the particle.
275  integer, intent(in), optional :: i_group
276  !> Weight class number for the particle.
277  integer, intent(in), optional :: i_class
278 
279  if (present(i_group)) aero_particle%weight_group = i_group
280  if (present(i_class)) aero_particle%weight_class = i_class
281 
282  end subroutine aero_particle_set_weight
283 
284 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
285 
286  !> Total mass of the particle (kg).
287  elemental real(kind=dp) function aero_particle_mass(aero_particle, aero_data)
288 
289  !> Particle.
290  type(aero_particle_t), intent(in) :: aero_particle
291  !> Aerosol data.
292  type(aero_data_t), intent(in) :: aero_data
293 
294  aero_particle_mass = sum(aero_particle%vol * aero_data%density)
295 
296  end function aero_particle_mass
297 
298 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
299 
300  !> Mass of a single species in the particle (kg).
301  elemental real(kind=dp) function aero_particle_species_mass(aero_particle, &
302  i_spec, aero_data)
303 
304  !> Particle.
305  type(aero_particle_t), intent(in) :: aero_particle
306  !> Species number to find mass of.
307  integer, intent(in) :: i_spec
308  !> Aerosol data.
309  type(aero_data_t), intent(in) :: aero_data
310 
311  aero_particle_species_mass = aero_particle%vol(i_spec) &
312  * aero_data%density(i_spec)
313 
314  end function aero_particle_species_mass
315 
316 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
317 
318  !> Mass of all species in the particle (kg).
319  function aero_particle_species_masses(aero_particle, aero_data)
320 
321  !> Particle.
322  type(aero_particle_t), intent(in) :: aero_particle
323  !> Aerosol data.
324  type(aero_data_t), intent(in) :: aero_data
325  !> Return value.
326  real(kind=dp) :: aero_particle_species_masses(aero_data%n_spec)
327 
328  aero_particle_species_masses = aero_particle%vol * aero_data%density
329 
330  end function aero_particle_species_masses
331 
332 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
333 
334  !> Total moles in the particle (1).
335  elemental real(kind=dp) function aero_particle_moles(aero_particle, &
336  aero_data)
337 
338  !> Particle.
339  type(aero_particle_t), intent(in) :: aero_particle
340  !> Aerosol data.
341  type(aero_data_t), intent(in) :: aero_data
342 
343  aero_particle_moles = sum(aero_particle%vol * aero_data%density &
344  / aero_data%molec_weight)
345 
346  end function aero_particle_moles
347 
348 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
349 
350  !> Total volume of the particle (m^3).
351  elemental real(kind=dp) function aero_particle_volume(aero_particle)
352 
353  !> Particle.
354  type(aero_particle_t), intent(in) :: aero_particle
355 
356  aero_particle_volume = sum(aero_particle%vol)
357 
358  end function aero_particle_volume
359 
360 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
361 
362  !> Total dry volume of the particle (m^3).
363  elemental real(kind=dp) function aero_particle_dry_volume(aero_particle, &
364  aero_data)
365 
366  !> Particle.
367  type(aero_particle_t), intent(in) :: aero_particle
368  !> Aerosol data.
369  type(aero_data_t), intent(in) :: aero_data
370 
371  integer :: i_spec
372 
374  do i_spec = 1,aero_data%n_spec
375  if (i_spec /= aero_data%i_water) then
377  + aero_particle%vol(i_spec)
378  end if
379  end do
380 
381  end function aero_particle_dry_volume
382 
383 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
384 
385  !> Total volume (dry or wet) of the particle (m^3).
386  elemental real(kind=dp) function aero_particle_volume_maybe_dry( &
387  aero_particle, aero_data, dry_volume)
388 
389  !> Particle.
390  type(aero_particle_t), intent(in) :: aero_particle
391  !> Aerosol data.
392  type(aero_data_t), intent(in) :: aero_data
393  !> Whether the desired volume is dry (otherwise wet).
394  logical, intent(in) :: dry_volume
395 
396  if (dry_volume) then
398  = aero_particle_dry_volume(aero_particle, aero_data)
399  else
401  end if
402 
403  end function aero_particle_volume_maybe_dry
404 
405 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
406 
407  !> Total radius of the particle (m).
408  elemental real(kind=dp) function aero_particle_radius(aero_particle)
409 
410  !> Particle.
411  type(aero_particle_t), intent(in) :: aero_particle
412 
414 
415  end function aero_particle_radius
416 
417 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
418 
419  !> Total dry radius of the particle (m).
420  elemental real(kind=dp) function aero_particle_dry_radius(aero_particle, &
421  aero_data)
422 
423  !> Particle.
424  type(aero_particle_t), intent(in) :: aero_particle
425  !> Aerosol data.
426  type(aero_data_t), intent(in) :: aero_data
427 
429  aero_particle_dry_volume(aero_particle, aero_data))
430 
431  end function aero_particle_dry_radius
432 
433 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
434 
435  !> Total diameter of the particle (m).
436  elemental real(kind=dp) function aero_particle_diameter(aero_particle)
437 
438  !> Particle.
439  type(aero_particle_t), intent(in) :: aero_particle
440 
441  aero_particle_diameter = 2d0 * aero_particle_radius(aero_particle)
442 
443  end function aero_particle_diameter
444 
445 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
446 
447  !> Total dry diameter of the particle (m).
448  elemental real(kind=dp) function aero_particle_dry_diameter(aero_particle, &
449  aero_data)
450 
451  !> Particle.
452  type(aero_particle_t), intent(in) :: aero_particle
453  !> Aerosol data.
454  type(aero_data_t), intent(in) :: aero_data
455 
457  = 2d0 * aero_particle_dry_radius(aero_particle, aero_data)
458 
459  end function aero_particle_dry_diameter
460 
461 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
462 
463  !> Average density of the particle (kg/m^3).
464  real(kind=dp) function aero_particle_density(aero_particle, aero_data)
465 
466  !> Particle.
467  type(aero_particle_t), intent(in) :: aero_particle
468  !> Aerosol data.
469  type(aero_data_t), intent(in) :: aero_data
470 
471  aero_particle_density = aero_particle_mass(aero_particle, aero_data) &
472  / aero_particle_volume(aero_particle)
473 
474  end function aero_particle_density
475 
476 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
477 
478  !> Returns the volume-average of the non-water elements of quantity.
479  real(kind=dp) function aero_particle_average_solute_quantity( &
480  aero_particle, aero_data, quantity)
481 
482  !> Aerosol particle.
483  type(aero_particle_t), intent(in) :: aero_particle
484  !> Aerosol data.
485  type(aero_data_t), intent(in) :: aero_data
486  !> Quantity to average.
487  real(kind=dp), intent(in) :: quantity(:)
488 
489  real(kind=dp) :: ones(aero_data%n_spec)
490 
491  ones = 1d0
493  aero_particle_total_solute_quantity(aero_particle, &
494  aero_data, quantity) &
495  / aero_particle_total_solute_quantity(aero_particle, aero_data, ones)
496 
498 
499 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
500 
501  !> Returns the volume-total of the non-water elements of quantity.
502  real(kind=dp) function aero_particle_total_solute_quantity(aero_particle, &
503  aero_data, quantity)
504 
505  !> Aerosol particle.
506  type(aero_particle_t), intent(in) :: aero_particle
507  !> Aerosol data.
508  type(aero_data_t), intent(in) :: aero_data
509  !> Quantity to total.
510  real(kind=dp), intent(in) :: quantity(:)
511 
512  real(kind=dp) total
513  integer i
514 
515  total = 0d0
516  do i = 1,aero_data%n_spec
517  if (i /= aero_data%i_water) then
518  total = total + aero_particle%vol(i) * quantity(i)
519  end if
520  end do
522 
524 
525 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
526 
527  !> Returns the water element of quantity.
528  real(kind=dp) function aero_particle_average_water_quantity(aero_particle, &
529  aero_data, quantity)
530 
531  !> Aerosol particle.
532  type(aero_particle_t), intent(in) :: aero_particle
533  !> Aerosol data.
534  type(aero_data_t), intent(in) :: aero_data
535  !> Quantity to average.
536  real(kind=dp), intent(in) :: quantity(:)
537 
538  call assert(420016623, aero_data%i_water > 0)
539  aero_particle_average_water_quantity = quantity(aero_data%i_water)
540 
542 
543 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
544 
545  !> Returns the volume-total of the water element of quantity.
546  real(kind=dp) function aero_particle_total_water_quantity(aero_particle, &
547  aero_data, quantity)
548 
549  !> Aerosol particle.
550  type(aero_particle_t), intent(in) :: aero_particle
551  !> Aerosol data.
552  type(aero_data_t), intent(in) :: aero_data
553  !> Quantity to total.
554  real(kind=dp), intent(in) :: quantity(:)
555 
556  call assert(223343210, aero_data%i_water > 0)
558  = aero_particle%vol(aero_data%i_water) &
559  * quantity(aero_data%i_water)
560 
562 
563 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
564 
565  !> Returns the water molecular weight.
566  !> (kg/mole)
567  real(kind=dp) function aero_particle_water_molec_weight(aero_data)
568 
569  !> Aerosol data.
570  type(aero_data_t), intent(in) :: aero_data
571 
572  call assert(772012490, aero_data%i_water > 0)
574  = aero_data%molec_weight(aero_data%i_water)
575 
577 
578 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
579 
580  !> Returns the average of the solute molecular weight (kg/mole).
581  real(kind=dp) function aero_particle_solute_molec_weight(aero_particle, &
582  aero_data)
583 
584  !> Aerosol particle.
585  type(aero_particle_t), intent(in) :: aero_particle
586  !> Aerosol data.
587  type(aero_data_t), intent(in) :: aero_data
588 
590  = aero_particle_average_solute_quantity(aero_particle, &
591  aero_data, aero_data%molec_weight)
592 
594 
595 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
596 
597  !> Returns the average of the solute ion number (1).
598  real(kind=dp) function aero_particle_solute_num_ions(aero_particle, &
599  aero_data)
600 
601  !> Aerosol data.
602  type(aero_data_t), intent(in) :: aero_data
603  !> Aerosol particle.
604  type(aero_particle_t), intent(in) :: aero_particle
605 
607  = aero_particle_average_solute_quantity(aero_particle, &
608  aero_data, real(aero_data%num_ions, kind=dp))
609 
610  end function aero_particle_solute_num_ions
611 
612 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
613 
614  !> Returns the water density (kg/m^3).
615  real(kind=dp) function aero_particle_water_density(aero_data)
616 
617  !> Aerosol data.
618  type(aero_data_t), intent(in) :: aero_data
619 
620  call assert(235482108, aero_data%i_water > 0)
621  aero_particle_water_density = aero_data%density(aero_data%i_water)
622 
623  end function aero_particle_water_density
624 
625 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
626 
627  !> Returns the average of the solute densities (kg/m^3).
628  real(kind=dp) function aero_particle_solute_density(aero_particle, &
629  aero_data)
630 
631  !> Aerosol data.
632  type(aero_data_t), intent(in) :: aero_data
633  !> Aerosol particle.
634  type(aero_particle_t), intent(in) :: aero_particle
635 
637  = aero_particle_average_solute_quantity(aero_particle, &
638  aero_data, aero_data%density)
639 
640  end function aero_particle_solute_density
641 
642 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
643 
644  !> Returns the water mass (kg).
645  real(kind=dp) function aero_particle_water_mass(aero_particle, aero_data)
646 
647  !> Aerosol data.
648  type(aero_data_t), intent(in) :: aero_data
649  !> Aerosol particle.
650  type(aero_particle_t), intent(in) :: aero_particle
651 
652  call assert(888636139, aero_data%i_water > 0)
653  aero_particle_water_mass = aero_particle%vol(aero_data%i_water) &
654  * aero_data%density(aero_data%i_water)
655 
656  end function aero_particle_water_mass
657 
658 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
659 
660  !> Returns the total solute mass (kg).
661  real(kind=dp) function aero_particle_solute_mass(aero_particle, aero_data)
662 
663  !> Aerosol data.
664  type(aero_data_t), intent(in) :: aero_data
665  !> Aerosol particle.
666  type(aero_particle_t), intent(in) :: aero_particle
667 
669  = aero_particle_total_solute_quantity(aero_particle, &
670  aero_data, aero_data%density)
671 
672  end function aero_particle_solute_mass
673 
674 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
675 
676  !> Returns the total solute volume (m^3).
677  real(kind=dp) function aero_particle_solute_volume(aero_particle, &
678  aero_data)
679 
680  !> Aerosol data.
681  type(aero_data_t), intent(in) :: aero_data
682  !> Aerosol particle.
683  type(aero_particle_t), intent(in) :: aero_particle
684 
685  real(kind=dp) :: ones(aero_data%n_spec)
686 
687  ones = 1d0
689  = aero_particle_total_solute_quantity(aero_particle, &
690  aero_data, ones)
691 
692  end function aero_particle_solute_volume
693 
694 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
695 
696  !> Returns the total solute radius (m).
697  real(kind=dp) function aero_particle_solute_radius(aero_particle, &
698  aero_data)
699 
700  !> Aerosol data.
701  type(aero_data_t), intent(in) :: aero_data
702  !> Aerosol particle.
703  type(aero_particle_t), intent(in) :: aero_particle
704 
706  = vol2rad(aero_particle_solute_volume(aero_particle, aero_data))
707 
708  end function aero_particle_solute_radius
709 
710 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
711 
712  !> Returns the average of the solute kappas (1).
713  real(kind=dp) function aero_particle_solute_kappa(aero_particle, aero_data)
714 
715  !> Aerosol data.
716  type(aero_data_t), intent(in) :: aero_data
717  !> Aerosol particle.
718  type(aero_particle_t), intent(in) :: aero_particle
719 
720  real(kind=dp) :: kappa(aero_data%n_spec), m_w, rho_w, m_a, rho_a
721  integer :: i_spec
722 
723  m_w = aero_particle_water_molec_weight(aero_data)
724  rho_w = aero_particle_water_density(aero_data)
725  do i_spec = 1,aero_data%n_spec
726  if (i_spec == aero_data%i_water) then
727  kappa(i_spec) = 0d0
728  elseif (aero_data%num_ions(i_spec) > 0) then
729  call assert_msg(123681459, aero_data%kappa(i_spec) == 0d0, &
730  "species has nonzero num_ions and kappa: " &
731  // trim(aero_data%name(i_spec)))
732  m_a = aero_data%molec_weight(i_spec)
733  rho_a = aero_data%density(i_spec)
734  kappa(i_spec) = m_w * rho_a / (m_a * rho_w) &
735  * real(aero_data%num_ions(i_spec), kind=dp)
736  else
737  kappa(i_spec) = aero_data%kappa(i_spec)
738  end if
739  end do
741  = aero_particle_average_solute_quantity(aero_particle, &
742  aero_data, kappa)
743 
744  end function aero_particle_solute_kappa
745 
746 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
747 
748  !> Returns the approximate critical relative humidity (1).
749  real(kind=dp) function aero_particle_approx_crit_rel_humid(aero_particle, &
750  aero_data, env_state)
751 
752  !> Aerosol particle.
753  type(aero_particle_t), intent(in) :: aero_particle
754  !> Aerosol data.
755  type(aero_data_t), intent(in) :: aero_data
756  !> Environment state.
757  type(env_state_t), intent(in) :: env_state
758 
759  real(kind=dp) :: kappa, diam, c, a
760 
761  kappa = aero_particle_solute_kappa(aero_particle, aero_data)
762  c = sqrt(4d0 * env_state_a(env_state)**3 / 27d0)
763  diam = aero_particle_diameter(aero_particle)
764  aero_particle_approx_crit_rel_humid = c / sqrt(kappa * diam**3) + 1d0
765 
767 
768 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
769 
770  !> Returns the critical relative humidity (1).
771  real(kind=dp) function aero_particle_crit_rel_humid(aero_particle, &
772  aero_data, env_state)
773 
774  !> Aerosol particle.
775  type(aero_particle_t), intent(in) :: aero_particle
776  !> Aerosol data.
777  type(aero_data_t), intent(in) :: aero_data
778  !> Environment state.
779  type(env_state_t), intent(in) :: env_state
780 
781  real(kind=dp) :: kappa, crit_diam, dry_diam, a
782 
783  a = env_state_a(env_state)
784  dry_diam = aero_particle_dry_diameter(aero_particle, aero_data)
785  crit_diam = aero_particle_crit_diameter(aero_particle, aero_data, &
786  env_state)
787  kappa = aero_particle_solute_kappa(aero_particle, aero_data)
788  if (kappa < 1d-30) then
789  aero_particle_crit_rel_humid = exp(a / crit_diam)
790  else
791  aero_particle_crit_rel_humid = (crit_diam**3 - dry_diam**3) &
792  / (crit_diam**3 - dry_diam**3 * (1 - kappa)) * exp(a / crit_diam)
793  end if
794 
795  end function aero_particle_crit_rel_humid
796 
797 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
798 
799  !> Returns the critical diameter (m).
800  !!
801  !! The method is as follows. We need to solve the polynomial
802  !! equation \f$f(D)\f$ given in Eq. (7) of
803  !!
804  !! N. Riemer, M. West, R. Zaveri, D. Easter, "Estimating black
805  !! carbon aging time-scales with a particle-resolved aerosol
806  !! model", Aerosol Science 41, 143-158, 2010.
807  !!
808  !! This is the equation:
809  !! \f[
810  !! f(D) = D^6 + c_4 D^4 + c_3 D^3 + c_0,
811  !! \f]
812  !! where \f$c_4 < 0\f$, \f$c_3 < 0\f$, and \f$c_0 > 0\f$. There is
813  !! unique solution for \f$D > D_{\rm dry}\f$, as shown in the above
814  !! paper.
815  !!
816  !! The polynomial has first derivative which factors like
817  !! \f[
818  !! f'(D) = 6 D^5 + 4 c_4 D^3 + 3 c_3 D^2
819  !! = (3 D^2 + 4 c_4) D^3 + (D^3 + c_3) 3 D^2.
820  !! \f]
821  !! The first term is positive for \f$D > (-4 c_4 / 3)^{1/2}\f$ and
822  !! the second is positive for \f$D > (-c_3)^{1/3}\f$. If we take
823  !! \f[
824  !! D_0 = max((-4 c_4 / 3)^{1/2}, (-c_3)^{1/3})
825  !! \f]
826  !! then we have that \f$f'(D) > 0\f$ for \f$D > D_0\f$. Similarly,
827  !! \f[
828  !! f''(D) = 30 D^4 + 12 c_4 D^2 + 6 c_3 D
829  !! = (5 D^2 + 4 c_4) 3 D^2 + (5 D^3 + 2 c_3) 3 D,
830  !! \f]
831  !! so \f$f''(D) > 0\f$ for \f$D > D_0\f$ (as this ensures that \f$D
832  !! > (-4 c_4 / 5)^{1/2}\f$ and \f$D > (-2 c_3 / 5)^{1/3}\f$).
833  !!
834  !! Thus for $D > D_0$ we have that the first and second derivatives
835  !! of $f(D)$ are positive, so Newton's method starting from
836  !! \f$D_0\f$ will converge quadratically. This is the scheme used
837  !! here.
838  real(kind=dp) function aero_particle_crit_diameter(aero_particle, &
839  aero_data, env_state)
840 
841  !> Aerosol particle.
842  type(aero_particle_t), intent(in) :: aero_particle
843  !> Aerosol data.
844  type(aero_data_t), intent(in) :: aero_data
845  !> Environment state.
846  type(env_state_t), intent(in) :: env_state
847 
848  integer, parameter :: crit_diam_max_iter = 100
849 
850  real(kind=dp) :: kappa, dry_diam, a, c4, c3, c0, d, f, df, dd
851  integer :: i_newton
852 
853  a = env_state_a(env_state)
854  dry_diam = aero_particle_dry_diameter(aero_particle, aero_data)
855  kappa = aero_particle_solute_kappa(aero_particle, aero_data)
856  if (kappa < 1d-30) then
857  ! bail out early for hydrophobic particles
858  aero_particle_crit_diameter = dry_diam
859  return
860  end if
861 
862  c4 = - 3d0 * dry_diam**3 * kappa / a
863  c3 = - dry_diam**3 * (2d0 - kappa)
864  c0 = dry_diam**6 * (1d0 - kappa)
865 
866  ! Newton's method for f(d) = d^6 + c4 d^4 + c3 d^3 + c0
867  d = max(sqrt(-4d0 / 3d0 * c4), (-c3)**(1d0/3d0))
868  do i_newton = 1,crit_diam_max_iter
869  f = d**6 + c4 * d**4 + c3 * d**3 + c0
870  df = 6 * d**5 + 4 * c4 * d**3 + 3 * c3 * d**2
871  dd = f / df
872  d = d - dd
873  if (abs(dd / d) < 1d-14) then
874  exit
875  end if
876  end do
877  call warn_assert_msg(408545686, i_newton < crit_diam_max_iter, &
878  "critical diameter Newton loop failed to converge")
879  call warn_assert_msg(353290871, d >= dry_diam, &
880  "critical diameter Newton loop converged to invalid solution")
882 
883  end function aero_particle_crit_diameter
884 
885 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
886 
887  !> Coagulate two particles together to make a new one. The new
888  !> particle will not have its ID set.
889  subroutine aero_particle_coagulate(aero_particle_1, &
890  aero_particle_2, aero_particle_new)
891 
892  !> First particle.
893  type(aero_particle_t), intent(in) :: aero_particle_1
894  !> Second particle.
895  type(aero_particle_t), intent(in) :: aero_particle_2
896  !> Result particle.
897  type(aero_particle_t), intent(inout) :: aero_particle_new
898 
899  call assert(203741686, size(aero_particle_1%vol) &
900  == size(aero_particle_new%vol))
901  call assert(586181003, size(aero_particle_2%vol) &
902  == size(aero_particle_new%vol))
903  call assert(483452167, size(aero_particle_1%n_orig_part) &
904  == size(aero_particle_new%n_orig_part))
905  call assert(126911035, size(aero_particle_2%n_orig_part) &
906  == size(aero_particle_new%n_orig_part))
907  aero_particle_new%vol = aero_particle_1%vol + aero_particle_2%vol
908  aero_particle_new%n_orig_part = aero_particle_1%n_orig_part &
909  + aero_particle_2%n_orig_part
910  aero_particle_new%weight_group = 0
911  aero_particle_new%weight_class = 0
912  aero_particle_new%absorb_cross_sect = 0d0
913  aero_particle_new%scatter_cross_sect = 0d0
914  aero_particle_new%asymmetry = 0d0
915  aero_particle_new%refract_shell = (0d0, 0d0)
916  aero_particle_new%refract_core = (0d0, 0d0)
917  aero_particle_new%core_vol = 0d0
918  if ((aero_particle_1%water_hyst_leg == 1) &
919  .and. (aero_particle_2%water_hyst_leg == 1)) then
920  aero_particle_new%water_hyst_leg = 1
921  else
922  aero_particle_new%water_hyst_leg = 0
923  end if
924  aero_particle_new%id = 0
925  aero_particle_new%least_create_time = &
926  min(aero_particle_1%least_create_time, &
927  aero_particle_2%least_create_time)
928  aero_particle_new%greatest_create_time = &
929  max(aero_particle_1%greatest_create_time, &
930  aero_particle_2%greatest_create_time)
931 
932  end subroutine aero_particle_coagulate
933 
934 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
935 
936  !> Determines the number of bytes required to pack the given value.
937  integer function pmc_mpi_pack_size_aero_particle(val)
938 
939  !> Value to pack.
940  type(aero_particle_t), intent(in) :: val
941 
944  + pmc_mpi_pack_size_integer_array(val%n_orig_part) &
945  + pmc_mpi_pack_size_integer(val%weight_group) &
946  + pmc_mpi_pack_size_integer(val%weight_class) &
947  + pmc_mpi_pack_size_real(val%absorb_cross_sect) &
948  + pmc_mpi_pack_size_real(val%scatter_cross_sect) &
949  + pmc_mpi_pack_size_real(val%asymmetry) &
950  + pmc_mpi_pack_size_complex(val%refract_shell) &
951  + pmc_mpi_pack_size_complex(val%refract_core) &
952  + pmc_mpi_pack_size_real(val%core_vol) &
953  + pmc_mpi_pack_size_integer(val%water_hyst_leg) &
954  + pmc_mpi_pack_size_integer(val%id) &
955  + pmc_mpi_pack_size_real(val%least_create_time) &
956  + pmc_mpi_pack_size_real(val%greatest_create_time)
957 
959 
960 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
961 
962  !> Packs the given value into the buffer, advancing position.
963  subroutine pmc_mpi_pack_aero_particle(buffer, position, val)
964 
965  !> Memory buffer.
966  character, intent(inout) :: buffer(:)
967  !> Current buffer position.
968  integer, intent(inout) :: position
969  !> Value to pack.
970  type(aero_particle_t), intent(in) :: val
971 
972 #ifdef PMC_USE_MPI
973  integer :: prev_position
974 
975  prev_position = position
976  call pmc_mpi_pack_real_array(buffer, position, val%vol)
977  call pmc_mpi_pack_integer_array(buffer, position, val%n_orig_part)
978  call pmc_mpi_pack_integer(buffer, position, val%weight_group)
979  call pmc_mpi_pack_integer(buffer, position, val%weight_class)
980  call pmc_mpi_pack_real(buffer, position, val%absorb_cross_sect)
981  call pmc_mpi_pack_real(buffer, position, val%scatter_cross_sect)
982  call pmc_mpi_pack_real(buffer, position, val%asymmetry)
983  call pmc_mpi_pack_complex(buffer, position, val%refract_shell)
984  call pmc_mpi_pack_complex(buffer, position, val%refract_core)
985  call pmc_mpi_pack_real(buffer, position, val%core_vol)
986  call pmc_mpi_pack_integer(buffer, position, val%water_hyst_leg)
987  call pmc_mpi_pack_integer(buffer, position, val%id)
988  call pmc_mpi_pack_real(buffer, position, val%least_create_time)
989  call pmc_mpi_pack_real(buffer, position, val%greatest_create_time)
990  call assert(810223998, position - prev_position &
992 #endif
993 
994  end subroutine pmc_mpi_pack_aero_particle
995 
996 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
997 
998  !> Unpacks the given value from the buffer, advancing position.
999  subroutine pmc_mpi_unpack_aero_particle(buffer, position, val)
1000 
1001  !> Memory buffer.
1002  character, intent(inout) :: buffer(:)
1003  !> Current buffer position.
1004  integer, intent(inout) :: position
1005  !> Value to pack.
1006  type(aero_particle_t), intent(inout) :: val
1007 
1008 #ifdef PMC_USE_MPI
1009  integer :: prev_position
1010 
1011  prev_position = position
1012  call pmc_mpi_unpack_real_array(buffer, position, val%vol)
1013  call pmc_mpi_unpack_integer_array(buffer, position, val%n_orig_part)
1014  call pmc_mpi_unpack_integer(buffer, position, val%weight_group)
1015  call pmc_mpi_unpack_integer(buffer, position, val%weight_class)
1016  call pmc_mpi_unpack_real(buffer, position, val%absorb_cross_sect)
1017  call pmc_mpi_unpack_real(buffer, position, val%scatter_cross_sect)
1018  call pmc_mpi_unpack_real(buffer, position, val%asymmetry)
1019  call pmc_mpi_unpack_complex(buffer, position, val%refract_shell)
1020  call pmc_mpi_unpack_complex(buffer, position, val%refract_core)
1021  call pmc_mpi_unpack_real(buffer, position, val%core_vol)
1022  call pmc_mpi_unpack_integer(buffer, position, val%water_hyst_leg)
1023  call pmc_mpi_unpack_integer(buffer, position, val%id)
1024  call pmc_mpi_unpack_real(buffer, position, val%least_create_time)
1025  call pmc_mpi_unpack_real(buffer, position, val%greatest_create_time)
1026  call assert(287447241, position - prev_position &
1028 #endif
1029 
1030  end subroutine pmc_mpi_unpack_aero_particle
1031 
1032 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1033 
1034 end module pmc_aero_particle
subroutine aero_particle_set_weight(aero_particle, i_group, i_class)
Sets the aerosol particle weight group.
subroutine pmc_mpi_pack_aero_particle(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine aero_particle_set_vols(aero_particle, vols)
Sets the aerosol particle volumes.
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:367
real(kind=dp) function aero_particle_total_solute_quantity(aero_particle, aero_data, quantity)
Returns the volume-total of the non-water elements of quantity.
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:471
elemental real(kind=dp) function aero_particle_moles(aero_particle, aero_data)
Total moles in the particle (1).
The aero_data_t structure and associated subroutines.
Definition: aero_data.F90:9
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:347
subroutine pmc_mpi_unpack_aero_particle(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
real(kind=dp) function aero_particle_solute_kappa(aero_particle, aero_data)
Returns the average of the solute kappas (1).
real(kind=dp) function aero_particle_crit_diameter(aero_particle, aero_data, env_state)
Returns the critical diameter (m).
subroutine aero_particle_zero(aero_particle)
Resets an aero_particle to be zero.
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:796
The aero_particle_t structure and associated subroutines.
subroutine aero_particle_deallocate(aero_particle)
Deallocates memory associated with an aero_particle_t.
elemental real(kind=dp) function aero_particle_species_mass(aero_particle, i_spec, aero_data)
Mass of a single species in the particle (kg).
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 aero_particle_set_create_time(aero_particle, create_time)
Sets the creation times for the particle.
real(kind=dp) function aero_particle_solute_mass(aero_particle, aero_data)
Returns the total solute mass (kg).
subroutine aero_particle_shift(aero_particle_from, aero_particle_to)
Shift data from one aero_particle_t to another and free the first one.
elemental real(kind=dp) function aero_particle_dry_radius(aero_particle, aero_data)
Total dry radius of the particle (m).
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
Definition: util.F90:76
elemental real(kind=dp) function aero_particle_diameter(aero_particle)
Total diameter of the particle (m).
subroutine pmc_mpi_unpack_complex(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:874
real(kind=dp) function aero_particle_water_molec_weight(aero_data)
Returns the water molecular weight. (kg/mole)
subroutine aero_particle_set_source(aero_particle, i_source)
Sets the aerosol particle source.
elemental real(kind=dp) function aero_particle_radius(aero_particle)
Total radius of the particle (m).
real(kind=dp) function env_state_a(env_state)
Condensation parameter.
Definition: env_state.F90:232
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:534
real(kind=dp) function aero_particle_water_density(aero_data)
Returns the water density (kg/m^3).
integer function pmc_mpi_size()
Returns the total number of processes.
Definition: mpi.F90:133
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
Definition: util.F90:58
real(kind=dp) function aero_particle_density(aero_particle, aero_data)
Average density of the particle (kg/m^3).
Common utility subroutines.
Definition: util.F90:9
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:927
Current environment state.
Definition: env_state.F90:26
real(kind=dp) function aero_particle_solute_density(aero_particle, aero_data)
Returns the average of the solute densities (kg/m^3).
real(kind=dp) function aero_particle_total_water_quantity(aero_particle, aero_data, quantity)
Returns the volume-total of the water element of quantity.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:771
integer function pmc_mpi_rank()
Returns the rank of the current process.
Definition: mpi.F90:116
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine pmc_mpi_unpack_integer_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:899
subroutine pmc_mpi_pack_complex(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:636
subroutine aero_particle_allocate(aero_particle)
Allocates memory in an aero_particle_t.
integer function pmc_mpi_pack_size_complex(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:429
real(kind=dp) function aero_particle_solute_num_ions(aero_particle, aero_data)
Returns the average of the solute ion number (1).
real(kind=dp) function aero_particle_crit_rel_humid(aero_particle, aero_data, env_state)
Returns the critical relative humidity (1).
Single aerosol particle data structure.
subroutine aero_particle_new_id(aero_particle)
Assigns a globally-unique new ID number to the particle.
elemental real(kind=dp) function aero_particle_dry_volume(aero_particle, aero_data)
Total dry volume of the particle (m^3).
Reading formatted text input.
Definition: spec_file.F90:43
real(kind=dp) function aero_particle_solute_molec_weight(aero_particle, aero_data)
Returns the average of the solute molecular weight (kg/mole).
subroutine pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:559
real(kind=dp) function aero_particle_solute_volume(aero_particle, aero_data)
Returns the total solute volume (m^3).
integer function pmc_mpi_pack_size_integer_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:449
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 aero_particle_average_water_quantity(aero_particle, aero_data, quantity)
Returns the water element of quantity.
real(kind=dp) function aero_particle_average_solute_quantity(aero_particle, aero_data, quantity)
Returns the volume-average of the non-water elements of quantity.
real(kind=dp) function, dimension(aero_data%n_spec) aero_particle_species_masses(aero_particle, aero_data)
Mass of all species in the particle (kg).
elemental real(kind=dp) function aero_particle_mass(aero_particle, aero_data)
Total mass of the particle (kg).
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:688
real(kind=dp) function aero_particle_approx_crit_rel_humid(aero_particle, aero_data, env_state)
Returns the approximate critical relative humidity (1).
real(kind=dp) elemental function vol2rad(v)
Convert volume (m^3) to radius (m).
Definition: util.F90:238
elemental real(kind=dp) function aero_particle_volume(aero_particle)
Total volume of the particle (m^3).
real(kind=dp) function aero_particle_solute_radius(aero_particle, aero_data)
Returns the total solute radius (m).
integer function pmc_mpi_pack_size_aero_particle(val)
Determines the number of bytes required to pack the given value.
Aerosol material properties and associated data.
Definition: aero_data.F90:40
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:102
elemental real(kind=dp) function aero_particle_dry_diameter(aero_particle, aero_data)
Total dry diameter of the particle (m).
subroutine pmc_mpi_pack_integer_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:661
subroutine aero_particle_copy(aero_particle_from, aero_particle_to)
Copies a particle.
elemental real(kind=dp) function aero_particle_volume_maybe_dry(aero_particle, aero_data, dry_volume)
Total volume (dry or wet) of the particle (m^3).
real(kind=dp) function aero_particle_water_mass(aero_particle, aero_data)
Returns the water mass (kg).