Skip to content

Commit 2db8695

Browse files
mabruzzomabruzzo
committed
shift to using rxn_rate_bufs for tracking collisional reaction rates
1 parent ab58434 commit 2db8695

File tree

7 files changed

+734
-748
lines changed

7 files changed

+734
-748
lines changed

‎src/clib/chemistry_solver_funcs.hpp‎

Lines changed: 648 additions & 647 deletions
Large diffs are not rendered by default.

‎src/clib/lookup_cool_rates1d.hpp‎

Lines changed: 17 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -128,7 +128,7 @@ inline void interpolate_h2_heating_terms_(
128128
/// interpolate the "standard" collisional reaction rates for each each index
129129
/// in the index-range that is also selected by the given itmask
130130
///
131-
/// @param[out] kcol_buf A struct containing the buffers that are filled by
131+
/// @param[out] rxn_rate_buf A struct containing the buffers that are filled by
132132
/// this function
133133
/// @param[in] idx_range Specifies the current index-range
134134
/// @param[in] kcol_rate_tables The 1D tables of "standard" collisional
@@ -141,7 +141,7 @@ inline void interpolate_h2_heating_terms_(
141141
/// @param[in] logTlininterp_buf Specifies the information related to the
142142
/// position in the logT interpolations (for a number of chemistry zones)
143143
inline void interpolate_kcol_rate_tables_(
144-
grackle::impl::CollisionalRxnRateCollection kcol_buf, IndexRange idx_range,
144+
FullRxnRateBuf rxn_rate_buf, IndexRange idx_range,
145145
grackle::impl::CollisionalRxnRateCollection kcol_rate_tables,
146146
int* kcol_lut_indices, int n_rates, const gr_mask_type* itmask,
147147
grackle::impl::LogTLinInterpScratchBuf logTlininterp_buf) {
@@ -153,17 +153,18 @@ inline void interpolate_kcol_rate_tables_(
153153
// -> to maximize performance, we may want to change order of values in the
154154
// tables and possibly change transpose the data layout. But, we may want
155155
// to get GPU support working before we head down this path
156+
double* const* kcol_buf = FullRxnRateBuf_kcol_bufs(&rxn_rate_buf);
156157

157158
for (int i = idx_range.i_start; i < idx_range.i_stop; i++) {
158159
if (itmask[i] != MASK_FALSE) {
159160
for (int counter = 0; counter < n_rates; counter++) {
160161
int idx = kcol_lut_indices[counter];
161162
const double* table = kcol_rate_tables.data[idx];
162163

163-
kcol_buf.data[idx][i] = table[logTlininterp_buf.indixe[i] - 1] +
164-
(table[logTlininterp_buf.indixe[i]] -
165-
table[logTlininterp_buf.indixe[i] - 1]) *
166-
logTlininterp_buf.tdef[i];
164+
kcol_buf[idx][i] = table[logTlininterp_buf.indixe[i] - 1] +
165+
(table[logTlininterp_buf.indixe[i]] -
166+
table[logTlininterp_buf.indixe[i] - 1]) *
167+
logTlininterp_buf.tdef[i];
167168
}
168169
}
169170
}
@@ -172,7 +173,7 @@ inline void interpolate_kcol_rate_tables_(
172173
/// interpolate collisional reaction rates for each each index in the
173174
/// index-range that is also selected by the given itmask
174175
///
175-
/// @param[out] kcol_buf A struct containing the buffers that are filled by
176+
/// @param[out] rxn_rate_buf A struct containing the buffers that are filled by
176177
/// this function
177178
/// @param[in] idx_range Specifies the current index-range
178179
/// @param[in] tgas1d specifies the gas temperatures for the `idx_range`
@@ -184,17 +185,16 @@ inline void interpolate_kcol_rate_tables_(
184185
/// @param[in] logTlininterp_buf Specifies the information related to the
185186
/// position in the logT interpolations (for a number of chemistry zones)
186187
inline void interpolate_collisional_rxn_rates_(
187-
grackle::impl::CollisionalRxnRateCollection kcol_buf, IndexRange idx_range,
188-
const double* tgas1d, const gr_mask_type* itmask, double dom,
189-
chemistry_data* my_chemistry, grackle_field_data* my_fields,
190-
chemistry_data_storage* my_rates,
188+
FullRxnRateBuf rxn_rate_buf, IndexRange idx_range, const double* tgas1d,
189+
const gr_mask_type* itmask, double dom, chemistry_data* my_chemistry,
190+
grackle_field_data* my_fields, chemistry_data_storage* my_rates,
191191
grackle::impl::LogTLinInterpScratchBuf logTlininterp_buf) {
192192
// There are 2 parts to this function
193193
// ----------------------------------
194194

195195
// Part 1: Handle interpolations for most collisional rxn rates
196196
interpolate_kcol_rate_tables_(
197-
kcol_buf, idx_range, *(my_rates->opaque_storage->kcol_rate_tables),
197+
rxn_rate_buf, idx_range, *(my_rates->opaque_storage->kcol_rate_tables),
198198
my_rates->opaque_storage->used_kcol_rate_indices,
199199
my_rates->opaque_storage->n_kcol_rate_indices, itmask, logTlininterp_buf);
200200

@@ -229,11 +229,12 @@ inline void interpolate_collisional_rxn_rates_(
229229

230230
#define USE_DENSITY_DEPENDENT_H2_DISSOCIATION_RATE
231231
#ifdef USE_DENSITY_DEPENDENT_H2_DISSOCIATION_RATE
232+
double* const* kcol_buf = FullRxnRateBuf_kcol_bufs(&rxn_rate_buf);
232233
if (my_chemistry->three_body_rate == 0) {
233234
for (int i = idx_range.i_start; i < idx_range.i_stop; i++) {
234235
if (itmask[i] != MASK_FALSE) {
235236
double nh = std::fmin(HI(i, idx_range.j, idx_range.k) * dom, 1.0e9);
236-
kcol_buf.data[CollisionalRxnLUT::k13][i] = tiny8;
237+
kcol_buf[CollisionalRxnLUT::k13][i] = tiny8;
237238
if (tgas1d[i] >= 500. && tgas1d[i] < 1.0e6) {
238239
// define a local buffer and fill it with values
239240
double k13dd[14];
@@ -254,7 +255,7 @@ inline void interpolate_collisional_rxn_rates_(
254255
k13dd[10] / (1. + std::pow((nh / k13dd[12]), k13dd[13]));
255256
k13_DT = std::fmax(std::pow(10., k13_DT), tiny8);
256257
//
257-
kcol_buf.data[CollisionalRxnLUT::k13][i] = k13_DT + k13_CID;
258+
kcol_buf[CollisionalRxnLUT::k13][i] = k13_DT + k13_CID;
258259
}
259260
}
260261
}
@@ -783,7 +784,6 @@ inline void lookup_cool_rates1d(
783784
grackle::impl::GrainSpeciesCollection grain_growth_rates,
784785
grackle::impl::GrainSpeciesCollection grain_temperatures,
785786
grackle::impl::LogTLinInterpScratchBuf logTlininterp_buf,
786-
grackle::impl::CollisionalRxnRateCollection kcol_buf,
787787
grackle::impl::PhotoRxnRateCollection kshield_buf,
788788
FullRxnRateBuf rxn_rate_buf,
789789
grackle::impl::ChemHeatingRates chemheatrates_buf,
@@ -826,8 +826,8 @@ inline void lookup_cool_rates1d(
826826
}
827827

828828
// interpolate all collisional reaction rates
829-
interpolate_collisional_rxn_rates_(kcol_buf, idx_range, tgas1d, itmask, dom,
830-
my_chemistry, my_fields, my_rates,
829+
interpolate_collisional_rxn_rates_(rxn_rate_buf, idx_range, tgas1d, itmask,
830+
dom, my_chemistry, my_fields, my_rates,
831831
logTlininterp_buf);
832832

833833
// interpolate terms used to compute H2 formation heating terms.

‎src/clib/rate_timestep_g.cpp‎

Lines changed: 53 additions & 53 deletions
Original file line numberDiff line numberDiff line change
@@ -31,7 +31,6 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
3131
const gr_mask_type* itmask, double* edot, double chunit,
3232
double dom, chemistry_data* my_chemistry,
3333
grackle_field_data* my_fields, IndexRange idx_range,
34-
grackle::impl::CollisionalRxnRateCollection kcr_buf,
3534
grackle::impl::PhotoRxnRateCollection kshield_buf,
3635
grackle::impl::ChemHeatingRates chemheatrates_buf,
3736
FullRxnRateBuf rxn_rate_buf) {
@@ -116,6 +115,7 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
116115
my_fields->grid_dimension[1], my_fields->grid_dimension[2]);
117116

118117
// locals
118+
const double* const* kcol_buf = FullRxnRateBuf_kcol_bufs(&rxn_rate_buf);
119119

120120
int i;
121121
double atten;
@@ -128,28 +128,28 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
128128
// Compute the electron density rate-of-change
129129

130130
dedot[i] =
131-
+kcr_buf.data[CollisionalRxnLUT::k1][i] *
131+
+kcol_buf[CollisionalRxnLUT::k1][i] *
132132
HI(i, idx_range.j, idx_range.k) *
133133
de(i, idx_range.j, idx_range.k) +
134-
kcr_buf.data[CollisionalRxnLUT::k3][i] *
134+
kcol_buf[CollisionalRxnLUT::k3][i] *
135135
HeI(i, idx_range.j, idx_range.k) *
136136
de(i, idx_range.j, idx_range.k) / 4. +
137-
kcr_buf.data[CollisionalRxnLUT::k5][i] *
137+
kcol_buf[CollisionalRxnLUT::k5][i] *
138138
HeII(i, idx_range.j, idx_range.k) *
139139
de(i, idx_range.j, idx_range.k) / 4. -
140-
kcr_buf.data[CollisionalRxnLUT::k2][i] *
140+
kcol_buf[CollisionalRxnLUT::k2][i] *
141141
HII(i, idx_range.j, idx_range.k) *
142142
de(i, idx_range.j, idx_range.k) -
143-
kcr_buf.data[CollisionalRxnLUT::k4][i] *
143+
kcol_buf[CollisionalRxnLUT::k4][i] *
144144
HeII(i, idx_range.j, idx_range.k) *
145145
de(i, idx_range.j, idx_range.k) / 4. -
146-
kcr_buf.data[CollisionalRxnLUT::k6][i] *
146+
kcol_buf[CollisionalRxnLUT::k6][i] *
147147
HeIII(i, idx_range.j, idx_range.k) *
148148
de(i, idx_range.j, idx_range.k) / 4. +
149-
kcr_buf.data[CollisionalRxnLUT::k57][i] *
149+
kcol_buf[CollisionalRxnLUT::k57][i] *
150150
HI(i, idx_range.j, idx_range.k) *
151151
HI(i, idx_range.j, idx_range.k) +
152-
kcr_buf.data[CollisionalRxnLUT::k58][i] *
152+
kcol_buf[CollisionalRxnLUT::k58][i] *
153153
HI(i, idx_range.j, idx_range.k) *
154154
HeI(i, idx_range.j, idx_range.k) / 4. +
155155
(kshield_buf.k24[i] * HI(i, idx_range.j, idx_range.k) +
@@ -158,16 +158,16 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
158158

159159
// Compute the HI density rate-of-change
160160

161-
HIdot[i] = -kcr_buf.data[CollisionalRxnLUT::k1][i] *
161+
HIdot[i] = -kcol_buf[CollisionalRxnLUT::k1][i] *
162162
HI(i, idx_range.j, idx_range.k) *
163163
de(i, idx_range.j, idx_range.k) +
164-
kcr_buf.data[CollisionalRxnLUT::k2][i] *
164+
kcol_buf[CollisionalRxnLUT::k2][i] *
165165
HII(i, idx_range.j, idx_range.k) *
166166
de(i, idx_range.j, idx_range.k) -
167-
kcr_buf.data[CollisionalRxnLUT::k57][i] *
167+
kcol_buf[CollisionalRxnLUT::k57][i] *
168168
HI(i, idx_range.j, idx_range.k) *
169169
HI(i, idx_range.j, idx_range.k) -
170-
kcr_buf.data[CollisionalRxnLUT::k58][i] *
170+
kcol_buf[CollisionalRxnLUT::k58][i] *
171171
HI(i, idx_range.j, idx_range.k) *
172172
HeI(i, idx_range.j, idx_range.k) / 4. -
173173
kshield_buf.k24[i] * HI(i, idx_range.j, idx_range.k);
@@ -179,55 +179,55 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
179179
for (i = idx_range.i_start; i <= idx_range.i_end; i++) {
180180
if (itmask[i] != MASK_FALSE) {
181181
HIdot[i] =
182-
-kcr_buf.data[CollisionalRxnLUT::k1][i] *
182+
-kcol_buf[CollisionalRxnLUT::k1][i] *
183183
de(i, idx_range.j, idx_range.k) *
184184
HI(i, idx_range.j, idx_range.k) -
185-
kcr_buf.data[CollisionalRxnLUT::k7][i] *
185+
kcol_buf[CollisionalRxnLUT::k7][i] *
186186
de(i, idx_range.j, idx_range.k) *
187187
HI(i, idx_range.j, idx_range.k) -
188-
kcr_buf.data[CollisionalRxnLUT::k8][i] *
188+
kcol_buf[CollisionalRxnLUT::k8][i] *
189189
HM(i, idx_range.j, idx_range.k) *
190190
HI(i, idx_range.j, idx_range.k) -
191-
kcr_buf.data[CollisionalRxnLUT::k9][i] *
191+
kcol_buf[CollisionalRxnLUT::k9][i] *
192192
HII(i, idx_range.j, idx_range.k) *
193193
HI(i, idx_range.j, idx_range.k) -
194-
kcr_buf.data[CollisionalRxnLUT::k10][i] *
194+
kcol_buf[CollisionalRxnLUT::k10][i] *
195195
H2II(i, idx_range.j, idx_range.k) *
196196
HI(i, idx_range.j, idx_range.k) / 2. -
197-
2. * kcr_buf.data[CollisionalRxnLUT::k22][i] *
197+
2. * kcol_buf[CollisionalRxnLUT::k22][i] *
198198
std::pow(HI(i, idx_range.j, idx_range.k), 2) *
199199
HI(i, idx_range.j, idx_range.k) +
200-
kcr_buf.data[CollisionalRxnLUT::k2][i] *
200+
kcol_buf[CollisionalRxnLUT::k2][i] *
201201
HII(i, idx_range.j, idx_range.k) *
202202
de(i, idx_range.j, idx_range.k) +
203-
2. * kcr_buf.data[CollisionalRxnLUT::k13][i] *
203+
2. * kcol_buf[CollisionalRxnLUT::k13][i] *
204204
HI(i, idx_range.j, idx_range.k) *
205205
H2I(i, idx_range.j, idx_range.k) / 2. +
206-
kcr_buf.data[CollisionalRxnLUT::k11][i] *
206+
kcol_buf[CollisionalRxnLUT::k11][i] *
207207
HII(i, idx_range.j, idx_range.k) *
208208
H2I(i, idx_range.j, idx_range.k) / 2. +
209-
2. * kcr_buf.data[CollisionalRxnLUT::k12][i] *
209+
2. * kcol_buf[CollisionalRxnLUT::k12][i] *
210210
de(i, idx_range.j, idx_range.k) *
211211
H2I(i, idx_range.j, idx_range.k) / 2. +
212-
kcr_buf.data[CollisionalRxnLUT::k14][i] *
212+
kcol_buf[CollisionalRxnLUT::k14][i] *
213213
HM(i, idx_range.j, idx_range.k) *
214214
de(i, idx_range.j, idx_range.k) +
215-
kcr_buf.data[CollisionalRxnLUT::k15][i] *
215+
kcol_buf[CollisionalRxnLUT::k15][i] *
216216
HM(i, idx_range.j, idx_range.k) *
217217
HI(i, idx_range.j, idx_range.k) +
218-
2. * kcr_buf.data[CollisionalRxnLUT::k16][i] *
218+
2. * kcol_buf[CollisionalRxnLUT::k16][i] *
219219
HM(i, idx_range.j, idx_range.k) *
220220
HII(i, idx_range.j, idx_range.k) +
221-
2. * kcr_buf.data[CollisionalRxnLUT::k18][i] *
221+
2. * kcol_buf[CollisionalRxnLUT::k18][i] *
222222
H2II(i, idx_range.j, idx_range.k) *
223223
de(i, idx_range.j, idx_range.k) / 2. +
224-
kcr_buf.data[CollisionalRxnLUT::k19][i] *
224+
kcol_buf[CollisionalRxnLUT::k19][i] *
225225
H2II(i, idx_range.j, idx_range.k) *
226226
HM(i, idx_range.j, idx_range.k) / 2. -
227-
kcr_buf.data[CollisionalRxnLUT::k57][i] *
227+
kcol_buf[CollisionalRxnLUT::k57][i] *
228228
HI(i, idx_range.j, idx_range.k) *
229229
HI(i, idx_range.j, idx_range.k) -
230-
kcr_buf.data[CollisionalRxnLUT::k58][i] *
230+
kcol_buf[CollisionalRxnLUT::k58][i] *
231231
HI(i, idx_range.j, idx_range.k) *
232232
HeI(i, idx_range.j, idx_range.k) / 4. -
233233
kshield_buf.k24[i] * HI(i, idx_range.j, idx_range.k) +
@@ -246,46 +246,46 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
246246
// Compute the electron density rate-of-change
247247

248248
dedot[i] =
249-
+kcr_buf.data[CollisionalRxnLUT::k1][i] *
249+
+kcol_buf[CollisionalRxnLUT::k1][i] *
250250
HI(i, idx_range.j, idx_range.k) *
251251
de(i, idx_range.j, idx_range.k) +
252-
kcr_buf.data[CollisionalRxnLUT::k3][i] *
252+
kcol_buf[CollisionalRxnLUT::k3][i] *
253253
HeI(i, idx_range.j, idx_range.k) *
254254
de(i, idx_range.j, idx_range.k) / 4. +
255-
kcr_buf.data[CollisionalRxnLUT::k5][i] *
255+
kcol_buf[CollisionalRxnLUT::k5][i] *
256256
HeII(i, idx_range.j, idx_range.k) *
257257
de(i, idx_range.j, idx_range.k) / 4. +
258-
kcr_buf.data[CollisionalRxnLUT::k8][i] *
258+
kcol_buf[CollisionalRxnLUT::k8][i] *
259259
HM(i, idx_range.j, idx_range.k) *
260260
HI(i, idx_range.j, idx_range.k) +
261-
kcr_buf.data[CollisionalRxnLUT::k15][i] *
261+
kcol_buf[CollisionalRxnLUT::k15][i] *
262262
HM(i, idx_range.j, idx_range.k) *
263263
HI(i, idx_range.j, idx_range.k) +
264-
kcr_buf.data[CollisionalRxnLUT::k17][i] *
264+
kcol_buf[CollisionalRxnLUT::k17][i] *
265265
HM(i, idx_range.j, idx_range.k) *
266266
HII(i, idx_range.j, idx_range.k) +
267-
kcr_buf.data[CollisionalRxnLUT::k14][i] *
267+
kcol_buf[CollisionalRxnLUT::k14][i] *
268268
HM(i, idx_range.j, idx_range.k) *
269269
de(i, idx_range.j, idx_range.k) -
270-
kcr_buf.data[CollisionalRxnLUT::k2][i] *
270+
kcol_buf[CollisionalRxnLUT::k2][i] *
271271
HII(i, idx_range.j, idx_range.k) *
272272
de(i, idx_range.j, idx_range.k) -
273-
kcr_buf.data[CollisionalRxnLUT::k4][i] *
273+
kcol_buf[CollisionalRxnLUT::k4][i] *
274274
HeII(i, idx_range.j, idx_range.k) *
275275
de(i, idx_range.j, idx_range.k) / 4. -
276-
kcr_buf.data[CollisionalRxnLUT::k6][i] *
276+
kcol_buf[CollisionalRxnLUT::k6][i] *
277277
HeIII(i, idx_range.j, idx_range.k) *
278278
de(i, idx_range.j, idx_range.k) / 4. -
279-
kcr_buf.data[CollisionalRxnLUT::k7][i] *
279+
kcol_buf[CollisionalRxnLUT::k7][i] *
280280
HI(i, idx_range.j, idx_range.k) *
281281
de(i, idx_range.j, idx_range.k) -
282-
kcr_buf.data[CollisionalRxnLUT::k18][i] *
282+
kcol_buf[CollisionalRxnLUT::k18][i] *
283283
H2II(i, idx_range.j, idx_range.k) *
284284
de(i, idx_range.j, idx_range.k) / 2. +
285-
kcr_buf.data[CollisionalRxnLUT::k57][i] *
285+
kcol_buf[CollisionalRxnLUT::k57][i] *
286286
HI(i, idx_range.j, idx_range.k) *
287287
HI(i, idx_range.j, idx_range.k) +
288-
kcr_buf.data[CollisionalRxnLUT::k58][i] *
288+
kcol_buf[CollisionalRxnLUT::k58][i] *
289289
HI(i, idx_range.j, idx_range.k) *
290290
HeI(i, idx_range.j, idx_range.k) / 4. +
291291
(kshield_buf.k24[i] * HI(i, idx_range.j, idx_range.k) +
@@ -295,22 +295,22 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
295295
// HII, HeII, HeIII recombination heating
296296

297297
edot[i] = edot[i] -
298-
chunit * (13.6 * (kcr_buf.data[CollisionalRxnLUT::k1][i] *
298+
chunit * (13.6 * (kcol_buf[CollisionalRxnLUT::k1][i] *
299299
HI(i, idx_range.j, idx_range.k) *
300300
de(i, idx_range.j, idx_range.k) -
301-
kcr_buf.data[CollisionalRxnLUT::k2][i] *
301+
kcol_buf[CollisionalRxnLUT::k2][i] *
302302
HII(i, idx_range.j, idx_range.k) *
303303
de(i, idx_range.j, idx_range.k)) +
304-
24.6 * (kcr_buf.data[CollisionalRxnLUT::k3][i] *
304+
24.6 * (kcol_buf[CollisionalRxnLUT::k3][i] *
305305
HeI(i, idx_range.j, idx_range.k) *
306306
de(i, idx_range.j, idx_range.k) / 4. -
307-
kcr_buf.data[CollisionalRxnLUT::k4][i] *
307+
kcol_buf[CollisionalRxnLUT::k4][i] *
308308
HeII(i, idx_range.j, idx_range.k) *
309309
de(i, idx_range.j, idx_range.k) / 4.) +
310-
79.0 * (kcr_buf.data[CollisionalRxnLUT::k5][i] *
310+
79.0 * (kcol_buf[CollisionalRxnLUT::k5][i] *
311311
HeII(i, idx_range.j, idx_range.k) *
312312
de(i, idx_range.j, idx_range.k) / 4. -
313-
kcr_buf.data[CollisionalRxnLUT::k6][i] *
313+
kcol_buf[CollisionalRxnLUT::k6][i] *
314314
HeIII(i, idx_range.j, idx_range.k) *
315315
de(i, idx_range.j, idx_range.k) / 4.));
316316

@@ -329,12 +329,12 @@ void rate_timestep_g(double* dedot, double* HIdot, gr_mask_type anydust,
329329
// need to apply it outside the delta calculation.
330330

331331
H2delta[i] = HI(i, idx_range.j, idx_range.k) *
332-
((3.53 * kcr_buf.data[CollisionalRxnLUT::k8][i] *
332+
((3.53 * kcol_buf[CollisionalRxnLUT::k8][i] *
333333
HM(i, idx_range.j, idx_range.k) +
334-
4.48 * kcr_buf.data[CollisionalRxnLUT::k22][i] *
334+
4.48 * kcol_buf[CollisionalRxnLUT::k22][i] *
335335
std::pow(HI(i, idx_range.j, idx_range.k), 2.)) *
336336
h2heatfac[i] -
337-
4.48 * kcr_buf.data[CollisionalRxnLUT::k13][i] *
337+
4.48 * kcol_buf[CollisionalRxnLUT::k13][i] *
338338
H2I(i, idx_range.j, idx_range.k) / 2.);
339339
// ! corrected by GC 202002
340340

0 commit comments

Comments
 (0)