MFC
Exascale flow solver
Loading...
Searching...
No Matches
m_reactive_burn.fpp.f90
Go to the documentation of this file.
1# 1 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
2!>
3!! @file
4!! @brief Contains module m_reactive_burn
5
6# 1 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 1
7# 1 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 1
8# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
9# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
10# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
11# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
12# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
13# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
14
15# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
16# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
17# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
18
19# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
20# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
21# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
22
23# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
24
25# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
26
27# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
28
29# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
30
31# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
32
33# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
34
35# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
36
37# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
38! New line at end of file is required for FYPP
39# 2 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
40# 1 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 1
41# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
42# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
43# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
44# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
45# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
46# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
47
48# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
49# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
50# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
51
52# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
53# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
54# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
55
56# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
57
58# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
59
60# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
61
62# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
63
64# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
65
66# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
67
68# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
69
70# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
71! New line at end of file is required for FYPP
72# 2 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 2
73
74# 4 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
75# 5 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
76# 6 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
77# 7 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
78# 8 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
79
80# 20 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
81
82# 43 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
83
84# 48 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
85
86# 53 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
87
88# 58 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
89
90# 63 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
91
92# 68 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
93
94# 81 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
95
96# 86 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
97
98# 91 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
99
100# 96 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
101
102# 101 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
103
104# 106 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
105
106# 111 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
107
108# 116 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
109
110# 121 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
111
112# 126 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
113
114# 156 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
115
116# 197 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
117
118# 211 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
119
120# 236 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
121
122# 247 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
123
124# 249 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
125# 260 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
126
127# 310 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
128
129# 320 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
130
131# 330 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
132
133# 339 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
134
135# 356 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
136
137# 366 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
138
139# 373 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
140
141# 379 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
142
143# 385 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
144
145# 391 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
146
147# 397 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
148
149# 403 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
150! New line at end of file is required for FYPP
151# 3 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
152# 1 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 1
153# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
154# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
155# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
156# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
157# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
158# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
159
160# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
161# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
162# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
163
164# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
165# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
166# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
167
168# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
169
170# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
171
172# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
173
174# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
175
176# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
177
178# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
179
180# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
181
182# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
183! New line at end of file is required for FYPP
184# 2 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 2
185
186# 7 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
187
188# 17 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
189
190# 22 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
191
192# 27 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
193
194# 32 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
195
196# 37 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
197
198# 42 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
199
200# 47 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
201
202# 52 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
203
204# 57 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
205
206# 62 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
207
208# 73 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
209
210# 78 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
211
212# 83 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
213
214# 88 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
215
216# 103 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
217
218# 131 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
219
220# 160 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
221
222# 175 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
223
224# 193 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
225
226# 215 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
227
228# 244 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
229
230# 259 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
231
232# 269 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
233
234# 278 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
235
236# 294 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
237
238# 304 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
239
240# 311 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
241! New line at end of file is required for FYPP
242# 4 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
243
244! GPU parallel region (scalar reductions, maxval/minval)
245# 23 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
246
247! GPU parallel loop over threads (most common GPU macro)
248# 43 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
249
250! Required closing for GPU_PARALLEL_LOOP
251# 55 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
252
253! Mark routine for device compilation
254# 112 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
255
256! Declare device-resident data
257# 130 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
258
259! Inner loop within a GPU parallel region
260# 145 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
261
262! Scoped GPU data region
263# 164 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
264
265! Host code with device pointers (for MPI with GPU buffers)
266# 193 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
267
268! Allocate device memory (unscoped)
269# 207 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
270
271! Free device memory
272# 219 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
273
274! Atomic operation on device
275# 231 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
276
277! End atomic capture block
278# 242 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
279
280! Copy data between host and device
281# 254 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
282
283! Synchronization barrier
284# 266 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
285
286! Import GPU library module (openacc or omp_lib)
287# 275 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
288
289! Emit code only for AMD compiler
290# 282 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
291
292! Emit code for non-Cray compilers
293# 289 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
294
295! Emit code only for Cray compiler
296# 296 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
297
298! Emit code for non-NVIDIA compilers
299# 303 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
300
301# 305 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
302# 306 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
303! New line at end of file is required for FYPP
304# 2 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 2
305
306# 14 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
307
308! Caution: This macro requires the use of a binding script to set CUDA_VISIBLE_DEVICES, such that we have one GPU device per MPI
309! rank. That's because for both cudaMemAdvise (preferred location) and cudaMemPrefetchAsync we use location = device_id = 0. For an
310! example see misc/nvidia_uvm/bind.sh.
311# 52 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
312
313! Allocate and create GPU device memory
314# 72 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
315
316! Free GPU device memory and deallocate
317# 80 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
318
319! Cray-specific GPU pointer setup for vector fields
320# 104 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
321
322! Cray-specific GPU pointer setup for scalar fields
323# 120 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
324
325! Cray-specific GPU pointer setup for acoustic source spatials
326# 145 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
327
328# 151 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
329
330# 158 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
331! New line at end of file is required for FYPP
332# 6 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp" 2
333# 1 "/home/runner/work/MFC/MFC/src/common/include/case.fpp" 1
334! This file exists so that Fypp can be run without generating case.fpp files for
335! each target. This is useful when generating documentation, for example. This
336! should also let MFC be built with CMake directly, without invoking mfc.sh.
337
338! For pre-process.
339# 8 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
340
341! For moving immersed boundaries in simulation
342# 12 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
343# 7 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp" 2
344
345!> @brief Condensed-phase reactive burn: a pressure-driven programmed-burn source that converts a "reactant" fluid into a "product"
346!! fluid on the multi-fluid model (num_fluids=2, chemistry='F'). The two fluids share the same stiffened-gas EOS (gamma, pi_inf) and
347!! differ only in their reference energy qv, so the reactant->product conversion releases (qv_reactant - qv_product) per unit mass
348!! through the mixture EOS with no explicit energy source. Because the two fluids are mechanically identical the volume-fraction
349!! swap is exact (the product volume fraction is the reaction progress), making this a reactive-Euler/ZND detonation model expressed
350!! through the diffuse-interface framework. A shock raises the pressure above rburn%pign, the reactant burns, and the energy release
351!! sustains the shock -- a self-propagating condensed-phase detonation.
352module m_reactive_burn
353
354 use m_global_parameters
355 use m_variables_conversion, only: s_phase_temperature, s_compute_mixture_coefficients, f_pressure
356
357 implicit none
358
359 private; public :: s_compute_reactive_burn, s_reactive_burn_substep
360
361contains
362
363 !> Programmed-burn rate dlambda/dt for one cell state. Both the RHS source and the operator-split integrator call this, so the
364 !! rate law is stated once.
365 !! @param pres Mixture pressure
366 !! @param lambda Reaction progress, i.e. the product volume fraction
367 !! @param alpha_rho_react Reactant partial density, for the optional Arrhenius factor
368 !! @param alpha_react Reactant volume fraction, for the optional Arrhenius factor
369 !! @param rate dlambda/dt; zero below the ignition pressure and once the reactant is spent
370 subroutine s_burn_rate(pres, lambda, alpha_rho_react, alpha_react, rate)
371
372
373# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
374#ifdef _CRAYFTN
375# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
376#if MFC_OpenACC
377# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
378!$acc routine seq
379# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
380#elif MFC_OpenMP
381# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
382
383# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
384
385# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
386!$omp declare target device_type(any)
387# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
388#else
389# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
390!DIR$ INLINEALWAYS s_burn_rate
391# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
392#endif
393# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
394#elif MFC_OpenACC
395# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
396!$acc routine seq
397# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
398#elif MFC_OpenMP
399# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
400
401# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
402
403# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
404!$omp declare target device_type(any)
405# 35 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
406#endif
407
408 real(wp), intent(in) :: pres, lambda, alpha_rho_react, alpha_react
409 real(wp), intent(out) :: rate
410 real(wp) :: drive, T_r
411
412 ! Pressure-driven programmed burn: fires only behind the shock (p > rburn%pign).
413 drive = (pres - rburn%pign)/rburn%pref
414 if (drive > 0._wp .and. lambda < 1._wp) then
415 rate = rburn%k*(1._wp - lambda)*drive**rburn%n
416 ! Optional Arrhenius dependence on the reactant's phasic temperature, from the reactant's own
417 ! EOS. rburn%ta = 0, the default, leaves the pure pressure-driven rate unchanged.
418 if (rburn%ta > 0._wp) then
419 call s_phase_temperature(alpha_rho_react/max(alpha_react, sgm_eps), pres, 1, t_r)
420 rate = rate*exp(-rburn%ta/t_r)
421 end if
422 else
423 rate = 0._wp
424 end if
425
426 end subroutine s_burn_rate
427
428 !> Add the programmed-burn reaction source to the continuity and volume-fraction RHS.
429 !! @param rhs_vf Right-hand-side accumulator (inout)
430 !! @param q_cons_vf Conserved variables (partial densities live here)
431 !! @param q_prim_vf Primitive variables (pressure and volume fractions live here)
432 !! @param bounds Interior cell bounds
433 subroutine s_compute_reactive_burn(rhs_vf, q_cons_vf, q_prim_vf, bounds)
434
435 type(scalar_field), dimension(sys_size), intent(inout) :: rhs_vf
436 type(scalar_field), dimension(sys_size), intent(in) :: q_cons_vf, q_prim_vf
437 type(int_bounds_info), dimension(1:3), intent(in) :: bounds
438 integer :: x, y, z
439 real(wp) :: rho, pres, lambda, rate, mdot
440 real(wp) :: alpha_rho_react, alpha_react
441
442
443# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
444
445# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
446#if defined(MFC_OpenACC)
447# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
448!$acc parallel loop collapse(3) gang vector default(present) private(rho, pres, lambda, rate, mdot, alpha_rho_react, alpha_react) copyin(bounds)
449# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
450#elif defined(MFC_OpenMP)
451# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
452
453# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
454
455# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
456
457# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
458!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) &
459# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
460!$omp& private(rho, pres, lambda, rate, mdot, alpha_rho_react, alpha_react) map(to:bounds)
461# 71 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
462#endif
463 do z = bounds(3)%beg, bounds(3)%end
464 do y = bounds(2)%beg, bounds(2)%end
465 do x = bounds(1)%beg, bounds(1)%end
466 ! reactant is fluid 1, product is fluid 2
467 rho = q_cons_vf(eqn_idx%cont%beg)%sf(x, y, z) + q_cons_vf(eqn_idx%cont%beg + 1)%sf(x, y, z)
468 pres = q_prim_vf(eqn_idx%E)%sf(x, y, z)
469 lambda = q_prim_vf(eqn_idx%adv%beg + 1)%sf(x, y, z) ! reaction progress = product volume fraction
470
471 alpha_rho_react = q_cons_vf(eqn_idx%cont%beg)%sf(x, y, z)
472 alpha_react = q_prim_vf(eqn_idx%adv%beg)%sf(x, y, z)
473 call s_burn_rate(pres, lambda, alpha_rho_react, alpha_react, rate)
474 if (rate > 0._wp) then
475 mdot = rho*rate ! mass reactant -> product
476
477 ! continuity: reactant loses mass, product gains it
478 rhs_vf(eqn_idx%cont%beg)%sf(x, y, z) = rhs_vf(eqn_idx%cont%beg)%sf(x, y, z) - mdot
479 rhs_vf(eqn_idx%cont%beg + 1)%sf(x, y, z) = rhs_vf(eqn_idx%cont%beg + 1)%sf(x, y, z) + mdot
480
481 ! volume fraction: exact swap (fluids share the EOS), so d(alpha)/dt = +/- rate
482 rhs_vf(eqn_idx%adv%beg)%sf(x, y, z) = rhs_vf(eqn_idx%adv%beg)%sf(x, y, z) - rate
483 rhs_vf(eqn_idx%adv%beg + 1)%sf(x, y, z) = rhs_vf(eqn_idx%adv%beg + 1)%sf(x, y, z) + rate
484 end if
485 end do
486 end do
487 end do
488
489# 97 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
490#if defined(MFC_OpenACC)
491# 97 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
492!$acc end parallel loop
493# 97 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
494#elif defined(MFC_OpenMP)
495# 97 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
496
497# 97 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
498!$omp end target teams loop
499# 97 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
500#endif
501
502 end subroutine s_compute_reactive_burn
503
504 !> Operator-split alternative to s_compute_reactive_burn, used when rburn%substeps > 0. The flow is frozen and the burn ODE is
505 !! integrated over one time step in equal sub-steps, so the reaction time scale is decoupled from the acoustic CFL. The mixture
506 !! pressure is re-evaluated from the frozen internal energy each sub-step, which is what carries the rate's own feedback: the
507 !! coefficients move as the reactant becomes product.
508 !! @param q_cons_vf Conserved variables, updated in place
509 !! @param dtime Time step to integrate across
510 !! @param bounds Interior cell bounds
511 subroutine s_reactive_burn_substep(q_cons_vf, dtime, bounds)
512
513 type(scalar_field), dimension(sys_size), intent(inout) :: q_cons_vf
514 real(wp), intent(in) :: dtime
515 type(int_bounds_info), dimension(1:3), intent(in) :: bounds
516 integer :: x, y, z, i, sub
517 real(wp) :: rho, pres, lambda, rate
518 real(wp) :: dt_sub, e_int, gamma_mix, pi_inf_mix, qv_mix
519 real(wp) :: rho_mix, dlambda, dmass
520 real(wp) :: alpha_rho_react, alpha_react
521
522 ! Bounded by num_fluids_max, not num_fluids: under case optimization num_fluids is a compile-time
523 ! constant that can be 1, and the reactant/product indices below are literal. The validator holds
524 ! reactive_burn to two fluids, so only the first two entries are ever used.
525 real(wp), dimension(num_fluids_max) :: alpha_rho, alpha
526
527 dt_sub = dtime/real(rburn%substeps, wp)
528
529
530# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
531
532# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
533#if defined(MFC_OpenACC)
534# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
535!$acc parallel loop collapse(3) gang vector default(present) &
536# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
537!$acc& private(alpha_rho, alpha, rho, pres, lambda, rate, e_int, gamma_mix, pi_inf_mix, qv_mix, rho_mix, dlambda, dmass, alpha_rho_react, alpha_react, i, sub) copyin(bounds, dt_sub)
538# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
539#elif defined(MFC_OpenMP)
540# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
541
542# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
543
544# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
545
546# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
547!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) &
548# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
549!$omp& private(alpha_rho, alpha, rho, pres, lambda, rate, e_int, gamma_mix, pi_inf_mix, qv_mix, rho_mix, dlambda, dmass, alpha_rho_react, alpha_react, i, sub) map(to:bounds, dt_sub)
550# 126 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
551#endif
552# 128 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
553 do z = bounds(3)%beg, bounds(3)%end
554 do y = bounds(2)%beg, bounds(2)%end
555 do x = bounds(1)%beg, bounds(1)%end
556
557# 131 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
558#if defined(MFC_OpenACC)
559# 131 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
560!$acc loop seq
561# 131 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
562#elif defined(MFC_OpenMP)
563# 131 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
564
565# 131 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
566#endif
567 do i = 1, num_fluids
568 alpha_rho(i) = q_cons_vf(i + eqn_idx%cont%beg - 1)%sf(x, y, z)
569 alpha(i) = q_cons_vf(i + eqn_idx%adv%beg - 1)%sf(x, y, z)
570 end do
571 rho = alpha_rho(1) + alpha_rho(2)
572
573 ! Internal energy per unit volume is what the burn conserves: it moves mass and volume
574 ! between the phases, and the qv difference surfaces as pressure through the mixture EOS.
575 e_int = q_cons_vf(eqn_idx%E)%sf(x, y, z)
576
577# 141 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
578#if defined(MFC_OpenACC)
579# 141 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
580!$acc loop seq
581# 141 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
582#elif defined(MFC_OpenMP)
583# 141 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
584
585# 141 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
586#endif
587 do i = eqn_idx%mom%beg, eqn_idx%mom%end
588 e_int = e_int - 0.5_wp*q_cons_vf(i)%sf(x, y, z)**2/rho
589 end do
590
591
592# 146 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
593#if defined(MFC_OpenACC)
594# 146 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
595!$acc loop seq
596# 146 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
597#elif defined(MFC_OpenMP)
598# 146 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
599
600# 146 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
601#endif
602 do sub = 1, rburn%substeps
603 lambda = alpha(2)
604 call s_compute_mixture_coefficients(alpha_rho, alpha, rho_mix, gamma_mix, pi_inf_mix, qv_mix)
605 pres = f_pressure(e_int, gamma_mix, pi_inf_mix, qv_mix)
606 alpha_rho_react = alpha_rho(1)
607 alpha_react = alpha(1)
608 call s_burn_rate(pres, lambda, alpha_rho_react, alpha_react, rate)
609 if (rate <= 0._wp) exit
610 ! A sub-step longer than the reaction time would carry the progress variable past one.
611 ! Stop it there and hand over the reactant's remaining mass in the same sub-step: capping
612 ! the two independently strands mass at zero volume, and the EOS divides one by the other.
613 if (rate*dt_sub >= 1._wp - lambda) then
614 dlambda = 1._wp - lambda
615 dmass = alpha_rho(1)
616 else
617 dlambda = rate*dt_sub
618 dmass = min(rho*dlambda, alpha_rho(1))
619 end if
620 alpha(1) = alpha(1) - dlambda
621 alpha(2) = alpha(2) + dlambda
622 alpha_rho(1) = alpha_rho(1) - dmass
623 alpha_rho(2) = alpha_rho(2) + dmass
624 end do
625
626
627# 171 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
628#if defined(MFC_OpenACC)
629# 171 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
630!$acc loop seq
631# 171 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
632#elif defined(MFC_OpenMP)
633# 171 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
634
635# 171 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
636#endif
637 do i = 1, num_fluids
638 q_cons_vf(i + eqn_idx%cont%beg - 1)%sf(x, y, z) = alpha_rho(i)
639 q_cons_vf(i + eqn_idx%adv%beg - 1)%sf(x, y, z) = alpha(i)
640 end do
641 end do
642 end do
643 end do
644
645# 179 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
646#if defined(MFC_OpenACC)
647# 179 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
648!$acc end parallel loop
649# 179 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
650#elif defined(MFC_OpenMP)
651# 179 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
652
653# 179 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
654!$omp end target teams loop
655# 179 "/home/runner/work/MFC/MFC/src/simulation/m_reactive_burn.fpp"
656#endif
657
658 end subroutine s_reactive_burn_substep
659
660end module m_reactive_burn
q_cons_vf
type(scalar_field), dimension(sys_size), intent(inout) q_cons_vf
Definition m_phase_change.fpp.f90:1005
m_global_parameters
Global parameters for the computational domain, fluid properties, and simulation algorithm configurat...
Definition m_global_parameters.fpp.f90:346
m_reactive_burn
Condensed-phase reactive burn: a pressure-driven programmed-burn source that converts a "reactant" fl...
Definition m_reactive_burn.fpp.f90:352
m_reactive_burn::s_compute_reactive_burn
subroutine, public s_compute_reactive_burn(rhs_vf, q_cons_vf, q_prim_vf, bounds)
Add the programmed-burn reaction source to the continuity and volume-fraction RHS.
Definition m_reactive_burn.fpp.f90:434
m_reactive_burn::s_burn_rate
subroutine s_burn_rate(pres, lambda, alpha_rho_react, alpha_react, rate)
Programmed-burn rate dlambda/dt for one cell state. Both the RHS source and the operator-split integr...
Definition m_reactive_burn.fpp.f90:371
m_reactive_burn::s_reactive_burn_substep
subroutine, public s_reactive_burn_substep(q_cons_vf, dtime, bounds)
Operator-split alternative to s_compute_reactive_burn, used when rburnsubsteps > 0....
Definition m_reactive_burn.fpp.f90:512
m_variables_conversion
Conservative-to-primitive variable conversion, mixture property evaluation, and pressure computation.
Definition m_variables_conversion.fpp.f90:346
m_variables_conversion::s_phase_temperature
subroutine, public s_phase_temperature(rho, pres, i, t)
Temperature of phase i at (rho, p): the stiffened-gas relation, or T_ref(rho) + (e - e_ref)/c_v.
Definition m_variables_conversion.fpp.f90:3631
m_variables_conversion::s_compute_mixture_coefficients
subroutine, public s_compute_mixture_coefficients(alpha_rho_k, alpha_k, rho_k, gamma_k, pi_inf_k, qv_k)
Mixture coefficients of one state. Under bubbles_euler with num_fluids == 1 the sole advection slot a...
Definition m_variables_conversion.fpp.f90:2683
m_variables_conversion::f_pressure
real(wp) function, public f_pressure(e_int, gamma, pi_inf, qv)
Pressure of a stiffened gas from its internal energy density - the inverse of s_compute_energy....
Definition m_variables_conversion.fpp.f90:3882