CAMP 1.0.0
Chemistry Across Multiple Phases
rxn_wennberg_tunneling.F90
Go to the documentation of this file.
1! Copyright (C) 2021 Barcelona Supercomputing Center,
2! University of Illinois at Urbana-Champaign, and
3! National Center for Atmospheric Research
4! SPDX-License-Identifier: MIT
5
6!> \file
7!> The camp_rxn_wennberg_tunneling module.
8
9!> \page camp_rxn_wennberg_tunneling CAMP: Wennberg Tunneling Reaction
10!!
11!! Wennberg tunneling reaction rate constant equations are calculated as follows:
12!!
13!! \f[
14!! Ae^{(\frac{-B}{T})}e^{(\frac{C}{T^3})}
15!! \f]
16!!
17!! where \f$A\f$ is the pre-exponential factor
18!! (\f$(\mbox{#}\,\mbox{cm}^{-3})^{-(n-1)}\mbox{s}^{-1}\f$),
19!! and \f$B\f$ and \f$C\f$ are parameters that capture the temperature
20!! dependence as described in Wennberg et al. (2018) \cite Wennberg2018 .
21!!
22!! Input data for Wennberg tunneling equations has the following format:
23!! \code{.json}
24!! {
25!! "type" : "WENNBERG_TUNNELING",
26!! "A" : 123.45,
27!! "B" : 1200.0,
28!! "C" : 1.0e8,
29!! "time unit" : "MIN",
30!! "reactants" : {
31!! "spec1" : {},
32!! "spec2" : { "qty" : 2 },
33!! ...
34!! },
35!! "products" : {
36!! "spec3" : {},
37!! "spec4" : { "yield" : 0.65 },
38!! ...
39!! }
40!! }
41!! \endcode
42!! The key-value pairs \b reactants, and \b products are required. Reactants
43!! without a \b qty value are assumed to appear once in the reaction equation.
44!! Products without a specified \b yield are assumed to have a \b yield of
45!! 1.0.
46!!
47!! When \b A is not included, it is assumed to be 1.0, when \b B is not
48!! included, it is assumed to be 0.0 K, and when \b C is not included, it is
49!! assumed to be 0.0 \f$\mbox{K}^3\f$.
50!! The unit for time is assumed to be s, but inclusion of the optional
51!! key-value pair \b time \b unit = \b MIN can be used to indicate a rate
52!! with min as the time unit.
53
54!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
55
56!> The rxn_wennberg_tunneling_t type and associated functions.
58
62 use camp_constants, only: const
66 use camp_util, only: i_kind, dp, to_string, &
68
69 implicit none
70 private
71
72#define NUM_REACT_ this%condensed_data_int(1)
73#define NUM_PROD_ this%condensed_data_int(2)
74#define A_ this%condensed_data_real(1)
75#define B_ this%condensed_data_real(2)
76#define C_ this%condensed_data_real(3)
77#define CONV_ this%condensed_data_real(4)
78#define NUM_INT_PROP_ 2
79#define NUM_REAL_PROP_ 4
80#define NUM_ENV_PARAM_ 1
81#define REACT_(x) this%condensed_data_int(NUM_INT_PROP_ + x)
82#define PROD_(x) this%condensed_data_int(NUM_INT_PROP_ + NUM_REACT_ + x)
83#define DERIV_ID_(x) this%condensed_data_int(NUM_INT_PROP_ + NUM_REACT_ + NUM_PROD_ + x)
84#define JAC_ID_(x) this%condensed_data_int(NUM_INT_PROP_ + 2*(NUM_REACT_+NUM_PROD_) + x)
85#define YIELD_(x) this%condensed_data_real(NUM_REAL_PROP_ + x)
86
88
89 !> Generic test reaction data type
90 type, extends(rxn_data_t) :: rxn_wennberg_tunneling_t
91 contains
92 !> Reaction initialization
93 procedure :: initialize
94 !> Finalize the reaction
97
98 !> Constructor for rxn_wennberg_tunneling_t
100 procedure :: constructor
101 end interface rxn_wennberg_tunneling_t
102
103contains
104
105!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
106
107 !> Constructor for Wennberg tunneling reaction
108 function constructor() result(new_obj)
109
110 !> A new reaction instance
111 type(rxn_wennberg_tunneling_t), pointer :: new_obj
112
113 allocate(new_obj)
114 new_obj%rxn_phase = gas_rxn
115
116 end function constructor
117
118!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
119
120 !> Initialize the reaction data, validating component data and loading
121 !! any required information into the condensed data arrays for use during
122 !! solving
123 subroutine initialize(this, chem_spec_data, aero_phase, aero_rep, n_cells)
124
125 !> Reaction data
126 class(rxn_wennberg_tunneling_t), intent(inout) :: this
127 !> Chemical species data
128 type(chem_spec_data_t), intent(in) :: chem_spec_data
129 !> Aerosol phase data
130 type(aero_phase_data_ptr), intent(in) :: aero_phase(:)
131 !> Aerosol representations
132 type(aero_rep_data_ptr), pointer, intent(in) :: aero_rep(:)
133 !> Number of grid cells being solved simultaneously
134 integer(kind=i_kind), intent(in) :: n_cells
135
136 type(property_t), pointer :: spec_props, reactants, products
137 character(len=:), allocatable :: key_name, spec_name, string_val
138 integer(kind=i_kind) :: i_spec, i_qty
139
140 integer(kind=i_kind) :: temp_int
141 real(kind=dp) :: temp_real
142
143 ! Get the species involved
144 if (.not. associated(this%property_set)) call die_msg(255324828, &
145 "Missing property set needed to initialize reaction")
146 key_name = "reactants"
147 call assert_msg(362119416, &
148 this%property_set%get_property_t(key_name, reactants), &
149 "Wennberg tunneling reaction is missing reactants")
150 key_name = "products"
151 call assert_msg(409429361, &
152 this%property_set%get_property_t(key_name, products), &
153 "Wennberg tunneling reaction is missing products")
154
155 ! Count the number of reactants (including those with a qty specified)
156 call reactants%iter_reset()
157 i_spec = 0
158 do while (reactants%get_key(spec_name))
159 ! Get properties included with this reactant in the reaction data
160 call assert(463909147, reactants%get_property_t(val=spec_props))
161 key_name = "qty"
162 if (spec_props%get_int(key_name, temp_int)) i_spec = i_spec+temp_int-1
163 call reactants%iter_next()
164 i_spec = i_spec + 1
165 end do
166
167 ! Allocate space in the condensed data arrays
168 allocate(this%condensed_data_int(num_int_prop_ + &
169 (i_spec + 2) * (i_spec + products%size())))
170 allocate(this%condensed_data_real(num_real_prop_ + products%size()))
171 this%condensed_data_int(:) = int(0, kind=i_kind)
172 this%condensed_data_real(:) = real(0.0, kind=dp)
173
174 ! Save space for the environment dependent parameters
175 this%num_env_params = num_env_param_
176
177 ! Save the size of the reactant and product arrays (for reactions where
178 ! these can vary)
179 num_react_ = i_spec
180 num_prod_ = products%size()
181
182 ! Set the #/cc -> ppm conversion prefactor
183 conv_ = const%avagadro / const%univ_gas_const * 10.0d0**(-12.0d0)
184
185 ! Get reaction parameters (it might be easiest to keep these at the
186 ! beginning of the condensed data array, so they can be accessed using
187 ! compliler flags)
188 key_name = "A"
189 if (.not. this%property_set%get_real(key_name, a_)) then
190 a_ = 1.0
191 end if
192 key_name = "time unit"
193 if (this%property_set%get_string(key_name, string_val)) then
194 if (trim(string_val).eq."MIN") then
195 a_ = a_ / 60.0
196 end if
197 endif
198 key_name = "B"
199 if (.not. this%property_set%get_real(key_name, b_)) then
200 b_ = 0.0
201 end if
202 key_name = "C"
203 if (.not. this%property_set%get_real(key_name, c_)) then
204 c_ = 0.0
205 end if
206
207 ! Get the indices and chemical properties for the reactants
208 call reactants%iter_reset()
209 i_spec = 1
210 do while (reactants%get_key(spec_name))
211
212 ! Save the index of this species in the state variable array
213 react_(i_spec) = chem_spec_data%gas_state_id(spec_name)
214
215 ! Make sure the species exists
216 call assert_msg(292299004, react_(i_spec).gt.0, &
217 "Missing Wennberg tunneling reactant: "//spec_name)
218
219 ! Get properties included with this reactant in the reaction data
220 call assert(687092598, reactants%get_property_t(val=spec_props))
221 key_name = "qty"
222 if (spec_props%get_int(key_name, temp_int)) then
223 do i_qty = 1, temp_int - 1
224 react_(i_spec + i_qty) = react_(i_spec)
225 end do
226 i_spec = i_spec + temp_int - 1
227 end if
228
229 call reactants%iter_next()
230 i_spec = i_spec + 1
231 end do
232
233 ! Get the indices and chemical properties for the products
234 call products%iter_reset()
235 i_spec = 1
236 do while (products%get_key(spec_name))
237
238 ! Save the index of this species in the state variable array
239 prod_(i_spec) = chem_spec_data%gas_state_id(spec_name)
240
241 ! Make sure the species exists
242 call assert_msg(681828291, prod_(i_spec).gt.0, &
243 "Missing Wennberg tunneling product: "//spec_name)
244
245 ! Get properties included with this product in the reaction data
246 call assert(794146636, products%get_property_t(val=spec_props))
247 key_name = "yield"
248 if (spec_props%get_real(key_name, temp_real)) then
249 yield_(i_spec) = temp_real
250 else
251 yield_(i_spec) = 1.0
252 end if
253
254 call products%iter_next()
255 i_spec = i_spec + 1
256 end do
257
258 end subroutine initialize
259
260!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
261
262 !> Finalize the reaction
263 subroutine finalize(this)
264
265 !> Reaction data
266 type(rxn_wennberg_tunneling_t), intent(inout) :: this
267
268 if (associated(this%property_set)) &
269 deallocate(this%property_set)
270 if (allocated(this%condensed_data_real)) &
271 deallocate(this%condensed_data_real)
272 if (allocated(this%condensed_data_int)) &
273 deallocate(this%condensed_data_int)
274
275 end subroutine finalize
276
277!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
278
279 !> Finalize an array of reactions
280 subroutine finalize_array(this)
281
282 !> Array of reaction data
283 type(rxn_wennberg_tunneling_t), intent(inout) :: this(:)
284
285 integer(kind=i_kind) :: i_rxn
286
287 do i_rxn = 1, size(this)
288 call finalize(this(i_rxn))
289 end do
290
291 end subroutine finalize_array
292
293!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
294
Initialize the aerosol representation data, validating component data and loading any required inform...
Get the non-unique name of a chemical species by its unique name.
Interface for to_string functions.
Definition util.F90:32
The abstract aero_phase_data_t structure and associated subroutines.
subroutine finalize_array(this)
Finalize the aerosol phase data.
type(aero_phase_data_t) function, pointer constructor(phase_name, init_size)
Constructor for aero_phase_data_t.
subroutine finalize(this)
Finalize the aerosol phase data.
The abstract aero_rep_data_t structure and associated subroutines.
The camp_state_t structure and associated subroutines.
Definition camp_state.F90:9
The chem_spec_data_t structure and associated subroutines.
Physical constants.
Definition constants.F90:9
integer, parameter dp
Kind of a double precision real number.
Definition constants.F90:16
type(const_t), save const
Fixed variable for accessing the constant's values.
Definition constants.F90:77
integer, parameter i_kind
Kind of an integer.
Definition constants.F90:21
The property_t structure and associated subroutines.
Definition property.F90:9
The rxn_data_t structure and associated subroutines.
Definition rxn_data.F90:60
integer(kind=i_kind), parameter, public gas_rxn
Gas-phase reaction.
Definition rxn_data.F90:85
The rxn_wennberg_tunneling_t type and associated functions.
Common utility subroutines.
Definition util.F90:9
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition util.F90:165
subroutine die_msg(code, error_msg)
Error immediately.
Definition util.F90:196
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
Definition util.F90:130
Pointer type for building arrays.
Pointer to aero_rep_data_t extending types.
Abstract reaction data type.
Definition rxn_data.F90:99