PartMC  2.3.0
coag_kernel_zero.F90
Go to the documentation of this file.
1 ! Copyright (C) 2007-2015 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_zero module.
7 
8 !> Constant kernel equal to zero.
9 !!
10 !! This is only of interest for the exact solution to the
11 !! no-coagulation, no-condensation case that can be used to test
12 !! emissions and background dilution.
14 
15  use pmc_bin_grid
16  use pmc_scenario
17  use pmc_env_state
18  use pmc_util
19  use pmc_aero_binned
20  use pmc_aero_dist
21  use pmc_aero_data
23 
24 contains
25 
26 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
27 
28  !> Zero coagulation kernel.
29  subroutine kernel_zero(aero_particle_1, aero_particle_2, &
30  aero_data, env_state, k)
31 
32  !> First particle.
33  type(aero_particle_t), intent(in) :: aero_particle_1
34  !> Second particle.
35  type(aero_particle_t), intent(in) :: aero_particle_2
36  !> Aerosol data.
37  type(aero_data_t), intent(in) :: aero_data
38  !> Environment state.
39  type(env_state_t), intent(in) :: env_state
40  !> Coagulation kernel.
41  real(kind=dp), intent(out) :: k
42 
43  k = 0d0
44 
45  end subroutine kernel_zero
46 
47 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
48 
49  !> Minimum and maximum of the zero coagulation kernel.
50  subroutine kernel_zero_minmax(v1, v2, aero_data, env_state, k_min, k_max)
51 
52  !> Volume of first particle.
53  real(kind=dp), intent(in) :: v1
54  !> Volume of second particle.
55  real(kind=dp), intent(in) :: v2
56  !> Aerosol data.
57  type(aero_data_t), intent(in) :: aero_data
58  !> Environment state.
59  type(env_state_t), intent(in) :: env_state
60  !> Coagulation kernel minimum value.
61  real(kind=dp), intent(out) :: k_min
62  !> Coagulation kernel maximum value.
63  real(kind=dp), intent(out) :: k_max
64 
65  k_min = 0d0
66  k_max = 0d0
67 
68  end subroutine kernel_zero_minmax
69 
70 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
71 
72  !> Exact solution with the zero coagulation kernel. Only useful for
73  !> testing emissions and background dilution.
74  !!
75  !! With only constant-rate emissions and dilution the number
76  !! distribution \f$n(D,t)\f$ at diameter \f$D\f$ and time \f$t\f$
77  !! satisfies:
78  !! \f[
79  !! \frac{d n(D,t)}{dt} = k_{\rm emit} n_{\rm emit}(D)
80  !! + (k_{\rm dilute} + k_{\rm loss}(D)) (n_{\rm back}(D) - n(D,t))
81  !! \f]
82  !! together with the initial condition \f$ n(D,0) = n_0(D) \f$. Here
83  !! \f$n_{\rm emit}(D)\f$ and \f$n_{\rm back}(D)\f$ are emission and
84  !! background size distributions, with corresponding rates \f$k_{\rm
85  !! emit}\f$ and \f$k_{\rm dilute}\f$. An optional loss function
86  !! \f$k_{\rm loss}(D)\f$ can be used to specify a volume-dependent
87  !! rate at which particles are lost. All values are taken at time
88  !! \f$t = 0\f$ and held constant, so there is no support for
89  !! time-varying emissions or background dilution.
90  !!
91  !! This is a family of ODEs parameterized by \f$D\f$ with
92  !! solution:
93  !! \f[
94  !! n(D,t) = n_{\infty}(D)
95  !! + (n_0(D) - n_{\infty}(D)) \exp(
96  !! -(k_{\rm dilute} + k_{\rm loss}(D)) t)
97  !! \f]
98  !! where the steady state limit is:
99  !! \f[
100  !! n_{\infty}(D) = n(D,\infty)
101  !! = n_{\rm back}(D)
102  !! + \frac{k_{\rm emit}}{k_{\rm dilute}
103  !! + k_{\rm loss}(D)} n_{\rm emit}(D)
104  !! \f]
105  subroutine soln_zero(bin_grid, aero_data, time, aero_dist_init, &
106  scenario, env_state, loss_function_type, aero_binned)
107 
108  !> Bin grid.
109  type(bin_grid_t), intent(in) :: bin_grid
110  !> Aerosol data.
111  type(aero_data_t), intent(in) :: aero_data
112  !> Current time (s).
113  real(kind=dp), intent(in) :: time
114  !> Initial distribution.
115  type(aero_dist_t), intent(in) :: aero_dist_init
116  !> Scenario.
117  type(scenario_t), intent(in) :: scenario
118  !> Environment state.
119  type(env_state_t), intent(in) :: env_state
120  !> Particle loss function type.
121  integer, intent(in) :: loss_function_type
122  !> Output state.
123  type(aero_binned_t), intent(inout) :: aero_binned
124 
125  integer :: i
126  real(kind=dp) :: emission_rate_scale, dilution_rate, p
127  real(kind=dp), allocatable :: loss_array(:)
128  type(aero_dist_t) :: emissions, background
129  type(aero_binned_t) :: background_binned, aero_binned_limit
130 
131  logical, save :: already_warned_water = .false.
132 
133  call aero_dist_allocate(emissions)
134  call aero_dist_allocate(background)
135  call aero_binned_allocate_size(background_binned, bin_grid%n_bin, &
136  aero_data%n_spec)
137 
138  call aero_dist_interp_1d(scenario%aero_emission, &
139  scenario%aero_emission_time, scenario%aero_emission_rate_scale, &
140  env_state%elapsed_time, emissions, emission_rate_scale)
141  call aero_dist_interp_1d(scenario%aero_background, &
142  scenario%aero_dilution_time, scenario%aero_dilution_rate, 0d0, &
143  background, dilution_rate)
144  call aero_binned_add_aero_dist(background_binned, bin_grid, aero_data, &
145  background)
146 
147  if (dilution_rate == 0d0 .and. &
148  loss_function_type == scenario_loss_function_invalid) then
149  call aero_binned_zero(aero_binned)
150  call aero_binned_add_aero_dist(aero_binned, bin_grid, aero_data, &
151  emissions)
152  call aero_binned_scale(aero_binned, &
153  emission_rate_scale * time / env_state%height)
154  else
155  allocate(loss_array(bin_grid%n_bin))
156 
157  if ((loss_function_type /= scenario_loss_function_zero) .or. &
158  (loss_function_type /= scenario_loss_function_invalid)) then
159  if (aero_data%name(1) /= "H2O") then
160  call warn_msg(176257943, &
161  "exact solution assumes composition is water", &
162  already_warned_water)
163  end if
164  end if
165 
166  do i = 1, bin_grid%n_bin
167  loss_array(i) = dilution_rate + scenario_loss_rate( &
168  loss_function_type, rad2vol(bin_grid%centers(i)), &
169  const%water_density, aero_data, env_state)
170  call assert_msg(181676342, loss_array(i) > 0, &
171  "non-positive loss rate")
172  loss_array(i) = 1d0 / loss_array(i)
173  end do
174 
175  ! calculate the limit steady state distribution
176  call aero_binned_allocate_size(aero_binned_limit, bin_grid%n_bin, &
177  aero_data%n_spec)
178  call aero_binned_add_aero_dist(aero_binned_limit, bin_grid, &
179  aero_data, emissions)
180  call aero_binned_scale(aero_binned_limit, &
181  emission_rate_scale / env_state%height)
182  call aero_binned_scale(background_binned, dilution_rate)
183  call aero_binned_add(aero_binned_limit, background_binned)
184  call aero_binned_scale_by_array(aero_binned_limit, loss_array)
185 
186  do i = 1, bin_grid%n_bin
187  loss_array(i) = exp(-time / loss_array(i))
188  end do
189 
190  ! calculate the current state
191  call aero_binned_zero(aero_binned)
192  call aero_binned_add_aero_dist(aero_binned, bin_grid, aero_data, &
193  aero_dist_init)
194  call aero_binned_sub(aero_binned, aero_binned_limit)
195  call aero_binned_scale_by_array(aero_binned, loss_array)
196  call aero_binned_add(aero_binned, aero_binned_limit)
197 
198  call aero_binned_deallocate(aero_binned_limit)
199  end if
200 
201  call aero_dist_deallocate(emissions)
202  call aero_dist_deallocate(background)
203  call aero_binned_deallocate(background_binned)
204 
205  end subroutine soln_zero
206 
207 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
208 
209 end module pmc_coag_kernel_zero
subroutine aero_binned_allocate_size(aero_binned, n_bin, n_spec)
Allocate an aero_binned_t of the given size.
Definition: aero_binned.F90:60
The scenario_t structure and associated subroutines.
Definition: scenario.F90:9
The aero_data_t structure and associated subroutines.
Definition: aero_data.F90:9
subroutine aero_dist_interp_1d(aero_dist_list, time_list, rate_list, time, aero_dist, rate)
Determine the current aero_dist and rate by interpolating at the current time with the lists of aero_...
Definition: aero_dist.F90:234
real(kind=dp) elemental function rad2vol(r)
Convert radius (m) to volume (m^3).
Definition: util.F90:274
The aero_dist_t structure and associated subroutines.
Definition: aero_dist.F90:18
subroutine aero_binned_sub(aero_binned, aero_binned_delta)
Subtract one aero_binned_t structure from another.
The aero_particle_t structure and associated subroutines.
The env_state_t structure and associated subroutines.
Definition: env_state.F90:9
real(kind=dp) function scenario_loss_rate(loss_function_type, vol, density, aero_data, env_state)
Evaluate a loss rate function.
Definition: scenario.F90:563
subroutine aero_binned_deallocate(aero_binned)
Free internal memory in an aero_binned_t structure.
Definition: aero_binned.F90:79
subroutine soln_zero(bin_grid, aero_data, time, aero_dist_init, scenario, env_state, loss_function_type, aero_binned)
Exact solution with the zero coagulation kernel. Only useful for testing emissions and background dil...
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
Definition: util.F90:76
subroutine aero_dist_deallocate(aero_dist)
Free all storage.
Definition: aero_dist.F90:82
subroutine kernel_zero_minmax(v1, v2, aero_data, env_state, k_min, k_max)
Minimum and maximum of the zero coagulation kernel.
A complete aerosol distribution, consisting of several modes.
Definition: aero_dist.F90:33
Common utility subroutines.
Definition: util.F90:9
subroutine aero_binned_add(aero_binned, aero_binned_delta)
Add two aero_binned_t structures together.
Current environment state.
Definition: env_state.F90:26
subroutine aero_binned_scale(aero_binned, alpha)
Scale an aero_binned_t by a real number.
Scenario data.
Definition: scenario.F90:51
Single aerosol particle data structure.
1D grid, either logarithmic or linear.
Definition: bin_grid.F90:33
The bin_grid_t structure and associated subroutines.
Definition: bin_grid.F90:9
The aero_binned_t structure and associated subroutines.
Definition: aero_binned.F90:9
subroutine warn_msg(code, warning_msg, already_warned)
Prints a warning message.
Definition: util.F90:36
subroutine kernel_zero(aero_particle_1, aero_particle_2, aero_data, env_state, k)
Zero coagulation kernel.
subroutine aero_binned_scale_by_array(aero_binned, alpha_array)
Scales an aero_binned_t element-wise by an array of reals.
Constant kernel equal to zero.
subroutine aero_binned_zero(aero_binned)
Set all internal data in an aero_binned_t structure to zero.
Definition: aero_binned.F90:92
subroutine aero_binned_add_aero_dist(aero_binned, bin_grid, aero_data, aero_dist)
Add an aero_dist_t to an aero_binned_t.
Aerosol material properties and associated data.
Definition: aero_data.F90:40
subroutine aero_dist_allocate(aero_dist)
Allocates an aero_dist.
Definition: aero_dist.F90:45
Aerosol number and volume distributions stored per bin.
Definition: aero_binned.F90:33