28 real(kind=dp),
allocatable :: vol(:)
31 integer,
allocatable :: n_orig_part(:)
33 integer :: weight_group
35 integer :: weight_class
37 real(kind=dp) :: absorb_cross_sect
39 real(kind=dp) :: scatter_cross_sect
41 real(kind=dp) :: asymmetry
43 complex(kind=dc) :: refract_shell
45 complex(kind=dc) :: refract_core
47 real(kind=dp) :: core_vol
49 integer :: water_hyst_leg
53 real(kind=dp) :: least_create_time
55 real(kind=dp) :: greatest_create_time
74 call move_alloc(aero_particle_from%vol, aero_particle_to%vol)
75 call move_alloc(aero_particle_from%n_orig_part, &
76 aero_particle_to%n_orig_part)
77 aero_particle_to%weight_group = aero_particle_from%weight_group
78 aero_particle_to%weight_class = aero_particle_from%weight_class
79 aero_particle_to%absorb_cross_sect = aero_particle_from%absorb_cross_sect
80 aero_particle_to%scatter_cross_sect = &
81 aero_particle_from%scatter_cross_sect
82 aero_particle_to%asymmetry = aero_particle_from%asymmetry
83 aero_particle_to%refract_shell = aero_particle_from%refract_shell
84 aero_particle_to%refract_core = aero_particle_from%refract_core
85 aero_particle_to%core_vol = aero_particle_from%core_vol
86 aero_particle_to%water_hyst_leg = aero_particle_from%water_hyst_leg
87 aero_particle_to%id = aero_particle_from%id
88 aero_particle_to%least_create_time = aero_particle_from%least_create_time
89 aero_particle_to%greatest_create_time = &
90 aero_particle_from%greatest_create_time
105 aero_particle%vol = 0d0
108 aero_particle%n_orig_part = 0
109 aero_particle%weight_group = 0
110 aero_particle%weight_class = 0
111 aero_particle%absorb_cross_sect = 0d0
112 aero_particle%scatter_cross_sect = 0d0
113 aero_particle%asymmetry = 0d0
114 aero_particle%refract_shell = (0d0, 0d0)
115 aero_particle%refract_core = (0d0, 0d0)
116 aero_particle%core_vol = 0d0
117 aero_particle%water_hyst_leg = 0
119 aero_particle%least_create_time = 0d0
120 aero_particle%greatest_create_time = 0d0
145 real(kind=dp),
intent(in) :: create_time
147 aero_particle%least_create_time = create_time
148 aero_particle%greatest_create_time = create_time
160 real(kind=dp),
intent(in) :: vols(:)
162 aero_particle%vol = vols
174 integer,
intent(in) :: i_source
176 aero_particle%n_orig_part = 0
177 aero_particle%n_orig_part(i_source) = 1
189 integer,
intent(in),
optional :: i_group
191 integer,
intent(in),
optional :: i_class
193 if (
present(i_group)) aero_particle%weight_group = i_group
194 if (
present(i_class)) aero_particle%weight_class = i_class
221 integer,
intent(in) :: i_spec
226 * aero_data%density(i_spec)
242 aero_particle_species_masses = aero_particle%vol * aero_data%density
258 / aero_data%molec_weight)
278 aero_particle, i_spec)
283 integer,
intent(in) :: i_spec
304 if (i_spec /= aero_data%i_water)
then 306 + aero_particle%vol(i_spec)
316 aero_particle, aero_data, dry_volume)
323 logical,
intent(in) :: dry_volume
402 aero_data, env_state)
411 real(kind=dp) :: volume, mobility_radius
414 mobility_radius = fractal_vol_to_mobility_rad(aero_data%fractal, &
415 volume, env_state%temp, env_state%pressure)
439 aero_particle, aero_data, quantity)
446 real(kind=dp),
intent(in) :: quantity(:)
453 aero_data, quantity) &
469 real(kind=dp),
intent(in) :: quantity(:)
476 if (i /= aero_data%i_water)
then 477 total = total + aero_particle%vol(i) * quantity(i)
495 real(kind=dp),
intent(in) :: quantity(:)
497 call assert(420016623, aero_data%i_water > 0)
513 real(kind=dp),
intent(in) :: quantity(:)
515 call assert(223343210, aero_data%i_water > 0)
517 = aero_particle%vol(aero_data%i_water) &
518 * quantity(aero_data%i_water)
531 call assert(772012490, aero_data%i_water > 0)
533 = aero_data%molec_weight(aero_data%i_water)
550 aero_data, aero_data%molec_weight)
567 aero_data,
real(aero_data%num_ions, kind=
dp))
579 call assert(235482108, aero_data%i_water > 0)
597 aero_data, aero_data%density)
611 call assert(888636139, aero_data%i_water > 0)
613 * aero_data%density(aero_data%i_water)
629 aero_data, aero_data%density)
686 if (i_spec == aero_data%i_water)
then 688 elseif (aero_data%num_ions(i_spec) > 0)
then 689 call assert_msg(123681459, aero_data%kappa(i_spec) == 0d0, &
690 "species has nonzero num_ions and kappa: " &
691 // trim(aero_data%name(i_spec)))
692 m_a = aero_data%molec_weight(i_spec)
693 rho_a = aero_data%density(i_spec)
694 kappa(i_spec) = m_w * rho_a / (m_a * rho_w) &
695 *
real(aero_data%num_ions(i_spec), kind=
dp)
697 kappa(i_spec) = aero_data%kappa(i_spec)
710 aero_data, env_state)
719 real(kind=dp) :: kappa, diam, C, A
732 aero_data, env_state)
741 real(kind=dp) :: kappa, crit_diam, dry_diam, A
748 if (kappa < 1d-30)
then 752 / (crit_diam**3 - dry_diam**3 * (1 - kappa)) * exp(a / crit_diam)
799 aero_data, env_state)
808 integer,
parameter :: CRIT_DIAM_MAX_ITER = 100
810 real(kind=dp) :: kappa, dry_diam, A, c4, c3, c0, d, f, df, dd
816 if (kappa < 1d-30)
then 822 c4 = - 3d0 * dry_diam**3 * kappa / a
823 c3 = - dry_diam**3 * (2d0 - kappa)
824 c0 = dry_diam**6 * (1d0 - kappa)
827 d = max(sqrt(-4d0 / 3d0 * c4), (-c3)**(1d0/3d0))
828 do i_newton = 1,crit_diam_max_iter
829 f = d**6 + c4 * d**4 + c3 * d**3 + c0
830 df = 6 * d**5 + 4 * c4 * d**3 + 3 * c3 * d**2
833 if (abs(dd / d) < 1d-14)
then 838 "critical diameter Newton loop failed to converge")
840 "critical diameter Newton loop converged to invalid solution")
850 aero_particle_2, aero_particle_new)
859 call assert(203741686,
size(aero_particle_1%vol) &
860 ==
size(aero_particle_2%vol))
861 call assert(483452167,
size(aero_particle_1%n_orig_part) &
862 ==
size(aero_particle_2%n_orig_part))
863 aero_particle_new%vol = aero_particle_1%vol + aero_particle_2%vol
864 aero_particle_new%n_orig_part = aero_particle_1%n_orig_part &
865 + aero_particle_2%n_orig_part
866 aero_particle_new%weight_group = 0
867 aero_particle_new%weight_class = 0
868 aero_particle_new%absorb_cross_sect = 0d0
869 aero_particle_new%scatter_cross_sect = 0d0
870 aero_particle_new%asymmetry = 0d0
871 aero_particle_new%refract_shell = (0d0, 0d0)
872 aero_particle_new%refract_core = (0d0, 0d0)
873 aero_particle_new%core_vol = 0d0
874 if ((aero_particle_1%water_hyst_leg == 1) &
875 .and. (aero_particle_2%water_hyst_leg == 1))
then 876 aero_particle_new%water_hyst_leg = 1
878 aero_particle_new%water_hyst_leg = 0
880 aero_particle_new%id = 0
881 aero_particle_new%least_create_time = &
882 min(aero_particle_1%least_create_time, &
883 aero_particle_2%least_create_time)
884 aero_particle_new%greatest_create_time = &
885 max(aero_particle_1%greatest_create_time, &
886 aero_particle_2%greatest_create_time)
922 character,
intent(inout) :: buffer(:)
924 integer,
intent(inout) :: position
929 integer :: prev_position
931 prev_position = position
946 call assert(810223998, position - prev_position &
958 character,
intent(inout) :: buffer(:)
960 integer,
intent(inout) :: position
965 integer :: prev_position
967 prev_position = position
982 call assert(287447241, position - prev_position &
999 logical,
intent(in) :: continue_on_error
1001 if (
allocated(aero_particle%vol))
then 1003 write(0, *)
'ERROR aero_particle A:' 1004 write(0, *)
'size(aero_particle%vol)',
size(aero_particle%vol)
1005 write(0, *)
'aero_data_n_spec(aero_data)', &
1007 call assert(185878626, continue_on_error)
1010 if (
allocated(aero_particle%n_orig_part))
then 1011 if (
size(aero_particle%n_orig_part) &
1013 write(0, *)
'ERROR aero_particle A:' 1014 write(0, *)
'size(aero_particle%n_orig_part)', &
1015 size(aero_particle%n_orig_part)
1016 write(0, *)
'aero_data_n_source(aero_data)', &
1018 call assert(625490639, continue_on_error)
subroutine pmc_mpi_pack_integer_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
elemental real(kind=dp) function aero_particle_radius(aero_particle, aero_data)
Total radius of the particle (m).
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
elemental real(kind=dp) function aero_particle_mass(aero_particle, aero_data)
Total mass of the particle (kg).
integer function pmc_mpi_pack_size_integer_array(val)
Determines the number of bytes required to pack the given value.
real(kind=dp) function aero_particle_water_molec_weight(aero_data)
Returns the water molecular weight. (kg/mole)
real(kind=dp) function aero_particle_solute_radius(aero_particle, aero_data)
Returns the total solute radius (m).
real(kind=dp) function env_state_a(env_state)
Condensation parameter.
real(kind=dp) function aero_particle_mobility_diameter(aero_particle, aero_data, env_state)
Mobility diameter of the particle (m).
real(kind=dp) function aero_particle_solute_mass(aero_particle, aero_data)
Returns the total solute mass (kg).
subroutine aero_particle_set_source(aero_particle, i_source)
Sets the aerosol particle source.
real(kind=dp) function aero_particle_solute_volume(aero_particle, aero_data)
Returns the total solute volume (m^3).
integer function pmc_mpi_rank()
Returns the rank of the current process.
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 ensure_real_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine pmc_mpi_unpack_integer_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
real(kind=dp) function aero_particle_solute_molec_weight(aero_particle, aero_data)
Returns the average of the solute molecular weight (kg/mole).
subroutine ensure_integer_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
real(kind=dp) function aero_particle_crit_rel_humid(aero_particle, aero_data, env_state)
Returns the critical relative humidity (1).
real(kind=dp) function, dimension(aero_data_n_spec(aero_data)) aero_particle_species_masses(aero_particle, aero_data)
Mass of all species in the particle (kg).
The env_state_t structure and associated subroutines.
subroutine aero_particle_shift(aero_particle_from, aero_particle_to)
Shift data from one aero_particle_t to another and free the first one.
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.
subroutine pmc_mpi_unpack_aero_particle(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
real(kind=dp) function aero_particle_approx_crit_rel_humid(aero_particle, aero_data, env_state)
Returns the approximate critical relative humidity (1).
subroutine aero_particle_set_weight(aero_particle, i_group, i_class)
Sets the aerosol particle weight group.
The aero_particle_t structure and associated subroutines.
elemental real(kind=dp) function aero_particle_species_volume(aero_particle, i_spec)
Volume of a single species in the particle (m^3).
subroutine aero_particle_new_id(aero_particle)
Assigns a globally-unique new ID number to the particle.
subroutine pmc_mpi_pack_complex(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
integer, save next_id
Next unique ID number to use for a particle.
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_set_create_time(aero_particle, create_time)
Sets the creation times for the particle.
real(kind=dp) function aero_particle_water_density(aero_data)
Returns the water density (kg/m^3).
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
Current environment state.
real(kind=dp) elemental function rad2diam(r)
Convert radius (m) to diameter (m).
real(kind=dp) function aero_particle_density(aero_particle, aero_data)
Average density of the particle (kg/m^3).
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
elemental real(kind=dp) function aero_particle_diameter(aero_particle, aero_data)
Total diameter of the particle (m).
subroutine pmc_mpi_unpack_complex(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine pmc_mpi_pack_real(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.
elemental real(kind=dp) function aero_particle_dry_volume(aero_particle, aero_data)
Total dry volume of the particle (m^3).
The aero_data_t structure and associated subroutines.
Single aerosol particle data structure.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine pmc_mpi_pack_aero_particle(buffer, position, val)
Packs the given value into the buffer, advancing position.
Reading formatted text input.
real(kind=dp) function aero_particle_water_mass(aero_particle, aero_data)
Returns the water mass (kg).
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
elemental real(kind=dp) function aero_particle_volume(aero_particle)
Total volume of the particle (m^3).
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
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 aero_particle_zero(aero_particle, aero_data)
Resets an aero_particle to be zero.
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
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).
integer function pmc_mpi_pack_size_complex(val)
Determines the number of bytes required to pack the given value.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
integer function pmc_mpi_pack_size_aero_particle(val)
Determines the number of bytes required to pack the given value.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
real(kind=dp) function aero_particle_average_water_quantity(aero_particle, aero_data, quantity)
Returns the water element of quantity.
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).
elemental integer function aero_data_n_source(aero_data)
Return the number of aerosol sources, or -1 if uninitialized.
subroutine aero_particle_check(aero_particle, aero_data, continue_on_error)
Check that the particle data is consistent.
elemental integer function aero_data_n_spec(aero_data)
Return the number of aerosol species, or -1 if uninitialized.
elemental real(kind=dp) function aero_particle_moles(aero_particle, aero_data)
Total moles in the particle (1).
integer function pmc_mpi_size()
Returns the total number of processes.
Aerosol material properties and associated data.
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.
Common utility subroutines.
Wrapper functions for MPI.
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_solute_density(aero_particle, aero_data)
Returns the average of the solute densities (kg/m^3).
elemental real(kind=dp) function aero_particle_dry_diameter(aero_particle, aero_data)
Total dry diameter of the particle (m).
elemental real(kind=dp) function aero_particle_species_mass(aero_particle, i_spec, aero_data)
Mass of a single species in the particle (kg).
elemental real(kind=dp) function aero_particle_dry_radius(aero_particle, aero_data)
Total dry radius of the particle (m).