MFC
Exascale flow solver
Loading...
Searching...
No Matches
m_riemann_solver_lf.fpp.f90
Go to the documentation of this file.
1# 1 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2!>
3!! @file
4!! @brief Contains module m_riemann_solver_lf
5
6!> @brief Lax-Friedrichs (Rusanov) approximate Riemann solver
7# 1 "/home/runner/work/MFC/MFC/src/common/include/case.fpp" 1
8! This file exists so that Fypp can be run without generating case.fpp files for
9! each target. This is useful when generating documentation, for example. This
10! should also let MFC be built with CMake directly, without invoking mfc.sh.
11
12! For pre-process.
13# 8 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
14
15! For moving immersed boundaries in simulation
16# 12 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
17# 7 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp" 2
18# 1 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 1
19# 1 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 1
20# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
21# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
22# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
23# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
24# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
25# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
26
27# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
28# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
29# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
30
31# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
32
33# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
34
35# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
36
37# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
38
39# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
40
41# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
42
43# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
44
45# 167 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
46! New line at end of file is required for FYPP
47# 2 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
48# 1 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 1
49# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
50# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
51# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
52# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
53# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
54# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
55
56# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
57# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
58# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
59
60# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
61
62# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
63
64# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
65
66# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
67
68# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
69
70# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
71
72# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
73
74# 167 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
75! New line at end of file is required for FYPP
76# 2 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 2
77
78# 4 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
79# 5 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
80# 6 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
81# 7 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
82# 8 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
83
84# 20 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
85
86# 43 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
87
88# 48 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
89
90# 53 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
91
92# 58 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
93
94# 63 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
95
96# 68 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
97
98# 76 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
99
100# 81 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
101
102# 86 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
103
104# 91 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
105
106# 96 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
107
108# 101 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
109
110# 106 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
111
112# 111 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
113
114# 116 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
115
116# 121 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
117
118# 151 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
119
120# 192 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
121
122# 206 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
123
124# 231 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
125
126# 242 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
127
128# 244 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
129# 255 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
130
131# 284 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
132
133# 294 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
134
135# 304 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
136
137# 313 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
138
139# 330 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
140
141# 340 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
142
143# 347 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
144
145# 353 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
146
147# 359 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
148
149# 365 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
150
151# 371 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
152
153# 377 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
154! New line at end of file is required for FYPP
155# 3 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
156# 1 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 1
157# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
158# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
159# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
160# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
161# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
162# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
163
164# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
165# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
166# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
167
168# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
169
170# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
171
172# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
173
174# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
175
176# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
177
178# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
179
180# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
181
182# 167 "/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# 55 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
312
313! Allocate and create GPU device memory
314# 75 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
315
316! Free GPU device memory and deallocate
317# 83 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
318
319! Cray-specific GPU pointer setup for vector fields
320# 107 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
321
322! Cray-specific GPU pointer setup for scalar fields
323# 123 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
324
325! Cray-specific GPU pointer setup for acoustic source spatials
326# 148 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
327
328# 154 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
329
330# 161 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
331! New line at end of file is required for FYPP
332# 8 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp" 2
333# 1 "/home/runner/work/MFC/MFC/src/simulation/include/inline_riemann.fpp" 1
334# 13 "/home/runner/work/MFC/MFC/src/simulation/include/inline_riemann.fpp"
335
336# 60 "/home/runner/work/MFC/MFC/src/simulation/include/inline_riemann.fpp"
337
338# 70 "/home/runner/work/MFC/MFC/src/simulation/include/inline_riemann.fpp"
339
340# 94 "/home/runner/work/MFC/MFC/src/simulation/include/inline_riemann.fpp"
341
342# 109 "/home/runner/work/MFC/MFC/src/simulation/include/inline_riemann.fpp"
343
344# 116 "/home/runner/work/MFC/MFC/src/simulation/include/inline_riemann.fpp"
345# 9 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp" 2
346
348
353 use m_thermochem, only: gas_constant, get_mixture_molecular_weight, get_mixture_specific_heat_cv_mass, &
354 & get_mixture_energy_mass, get_species_specific_heats_r, get_mixture_specific_heat_cp_mass, molecular_weights
356
357 implicit none
358
359contains
360
361 !> Lax-Friedrichs (Rusanov) approximate Riemann solver
362 subroutine s_lf_riemann_solver(qL_prim_rsx_vf, dqL_prim_dx_vf, dqL_prim_dy_vf, dqL_prim_dz_vf, qL_prim_vf, qR_prim_rsx_vf, &
363 & dqR_prim_dx_vf, dqR_prim_dy_vf, dqR_prim_dz_vf, qR_prim_vf, q_prim_vf, flux_vf, flux_src_vf, &
364 & flux_gsrc_vf, norm_dir, ix, iy, iz)
365
366 real(wp), dimension(idwbuff(1)%beg:,idwbuff(2)%beg:,idwbuff(3)%beg:,1:), intent(inout) :: qL_prim_rsx_vf, qR_prim_rsx_vf
367 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf
368 type(scalar_field), allocatable, dimension(:), intent(inout) :: qL_prim_vf, qR_prim_vf
369 type(scalar_field), allocatable, dimension(:), intent(inout) :: dqL_prim_dx_vf, dqR_prim_dx_vf, dqL_prim_dy_vf, &
370 & dqR_prim_dy_vf, dqL_prim_dz_vf, dqR_prim_dz_vf
371
372 ! Intercell fluxes
373 type(scalar_field), dimension(sys_size), intent(inout) :: flux_vf, flux_src_vf, flux_gsrc_vf
374 integer, intent(in) :: norm_dir
375 type(int_bounds_info), intent(in) :: ix, iy, iz
376
377# 49 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
378 real(wp), dimension(num_fluids) :: alpha_rho_L, alpha_rho_R
379 real(wp), dimension(num_vels) :: vel_L, vel_R
380 real(wp), dimension(num_fluids) :: alpha_L, alpha_R
381 real(wp), dimension(num_species) :: Ys_L, Ys_R
382 real(wp), dimension(num_species) :: Cp_iL, Cp_iR, Xs_L, Xs_R, Gamma_iL, Gamma_iR
383 real(wp), dimension(num_species) :: Yi_avg, Phi_avg, h_iL, h_iR, h_avg_2
384 !> Averaged velocity gradient tensor `d(vel_i)/d(coord_j)`.
385 real(wp), dimension(num_dims, num_dims) :: vel_grad_L, vel_grad_R
386# 58 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
387 real(wp) :: rho_L, rho_R
388 real(wp) :: pres_L, pres_R
389 real(wp) :: E_L, E_R
390 real(wp) :: H_L, H_R
391 real(wp) :: Cp_avg, Cv_avg, T_avg, eps, c_sum_Yi_Phi
392 real(wp) :: T_L, T_R
393 real(wp) :: Y_L, Y_R
394 real(wp) :: MW_L, MW_R
395 real(wp) :: R_gas_L, R_gas_R
396 real(wp) :: Cp_L, Cp_R
397 real(wp) :: Cv_L, Cv_R
398 real(wp) :: Gamm_L, Gamm_R
399 real(wp) :: gamma_L, gamma_R
400 real(wp) :: pi_inf_L, pi_inf_R
401 real(wp) :: qv_L, qv_R
402 real(wp) :: c_L, c_R
403 real(wp), dimension(2) :: Re_L, Re_R
404 real(wp) :: rho_avg
405 real(wp) :: H_avg
406 real(wp) :: gamma_avg
407 real(wp) :: c_avg
408 real(wp) :: s_L, s_R, s_M, s_P, s_S
409 real(wp) :: xi_M, xi_P
410 real(wp) :: ptilde_L, ptilde_R
411 real(wp) :: vel_L_rms, vel_R_rms, vel_avg_rms
412 real(wp) :: vel_L_tmp, vel_R_tmp
413 real(wp) :: Ms_L, Ms_R, pres_SL, pres_SR
414 real(wp) :: alpha_L_sum, alpha_R_sum
415 real(wp) :: zcoef, pcorr !< low Mach number correction
416 integer :: i, j, k, l !< Generic loop iterators
417 integer :: Re_size_loc1, Re_size_loc2 !< host copies of Re_size; amdflang reads the declare-target original stale cross-TU
418 integer, dimension(3) :: idx_right_phys !< Physical (j,k,l) indices for right state.
419 ! Populating the buffers of the left and right Riemann problem states variables, based on the choice of boundary conditions
420
421 call s_populate_riemann_states_variables_buffers(ql_prim_rsx_vf, dql_prim_dx_vf, dql_prim_dy_vf, dql_prim_dz_vf, &
422 & qr_prim_rsx_vf, dqr_prim_dx_vf, dqr_prim_dy_vf, dqr_prim_dz_vf, norm_dir, ix, iy, iz)
423
424 ! Reshaping inputted data based on dimensional splitting direction
425 call s_initialize_riemann_solver(flux_src_vf, norm_dir)
426 re_size_loc1 = re_size(1); re_size_loc2 = re_size(2)
427# 102 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
428# 103 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
429# 104 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
430 if (norm_dir == 1) then
431
432# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
433
434# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
435#if defined(MFC_OpenACC)
436# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
437!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, l, alpha_rho_L, alpha_rho_R, vel_L, vel_R, alpha_L, alpha_R, Re_L, Re_R, rho_avg, h_avg, gamma_avg, s_L, s_R, s_S, Ys_L, &
438# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
439!$acc& Ys_R, Cp_iL, Cp_iR, Xs_L, Xs_R, Gamma_iL, Gamma_iR, Yi_avg, Phi_avg, h_iL, h_iR, h_avg_2, pcorr, zcoef, vel_grad_L, vel_grad_R, idx_right_phys, vel_L_rms, vel_R_rms, vel_avg_rms, vel_L_tmp, &
440# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
441!$acc& vel_R_tmp, Ms_L, Ms_R, pres_SL, pres_SR, alpha_L_sum, alpha_R_sum, c_avg, pres_L, pres_R, rho_L, rho_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, qv_L, qv_R, c_L, c_R, E_L, E_R, H_L, H_R, ptilde_L, &
442# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
443!$acc& ptilde_R, s_M, s_P, xi_M, xi_P, Cp_avg, Cv_avg, T_avg, eps, c_sum_Yi_Phi, Cp_L, Cp_R, Cv_L, Cv_R, R_gas_L, R_gas_R, MW_L, MW_R, T_L, T_R, Y_L, Y_R) firstprivate(Re_size_loc1, Re_size_loc2)
444# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
445#elif defined(MFC_OpenMP)
446# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
447
448# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
449
450# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
451
452# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
453!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(i, j, k, l, &
454# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
455!$omp& alpha_rho_L, alpha_rho_R, vel_L, vel_R, alpha_L, alpha_R, Re_L, Re_R, rho_avg, h_avg, gamma_avg, s_L, s_R, s_S, Ys_L, Ys_R, Cp_iL, Cp_iR, Xs_L, Xs_R, Gamma_iL, Gamma_iR, Yi_avg, Phi_avg, h_iL, &
456# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
457!$omp& h_iR, h_avg_2, pcorr, zcoef, vel_grad_L, vel_grad_R, idx_right_phys, vel_L_rms, vel_R_rms, vel_avg_rms, vel_L_tmp, vel_R_tmp, Ms_L, Ms_R, pres_SL, pres_SR, alpha_L_sum, alpha_R_sum, c_avg, &
458# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
459!$omp& pres_L, pres_R, rho_L, rho_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, qv_L, qv_R, c_L, c_R, E_L, E_R, H_L, H_R, ptilde_L, ptilde_R, s_M, s_P, xi_M, xi_P, Cp_avg, Cv_avg, T_avg, eps, c_sum_Yi_Phi, &
460# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
461!$omp& Cp_L, Cp_R, Cv_L, Cv_R, R_gas_L, R_gas_R, MW_L, MW_R, T_L, T_R, Y_L, Y_R) firstprivate(Re_size_loc1, Re_size_loc2)
462# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
463#endif
464# 113 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
465 do l = is3%beg, is3%end
466 do k = is2%beg, is2%end
467 do j = is1%beg, is1%end
468
469# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
470#if defined(MFC_OpenACC)
471# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
472!$acc loop seq
473# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
474#elif defined(MFC_OpenMP)
475# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
476
477# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
478#endif
479 do i = 1, eqn_idx%cont%end
480 alpha_rho_l(i) = ql_prim_rsx_vf(j, k, l, i)
481 alpha_rho_r(i) = qr_prim_rsx_vf(j + 1, k, l, i)
482 end do
483
484 vel_l_rms = 0._wp; vel_r_rms = 0._wp
485
486
487# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
488#if defined(MFC_OpenACC)
489# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
490!$acc loop seq
491# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
492#elif defined(MFC_OpenMP)
493# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
494
495# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
496#endif
497 do i = 1, num_vels
498 vel_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%cont%end + i)
499 vel_r(i) = qr_prim_rsx_vf(j + 1, k, l, eqn_idx%cont%end + i)
500 vel_l_rms = vel_l_rms + vel_l(i)**2._wp
501 vel_r_rms = vel_r_rms + vel_r(i)**2._wp
502 end do
503
504
505# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
506#if defined(MFC_OpenACC)
507# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
508!$acc loop seq
509# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
510#elif defined(MFC_OpenMP)
511# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
512
513# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
514#endif
515 do i = 1, num_fluids
516 alpha_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%E + i)
517 alpha_r(i) = qr_prim_rsx_vf(j + 1, k, l, eqn_idx%E + i)
518 end do
519
520 pres_l = ql_prim_rsx_vf(j, k, l, eqn_idx%E)
521 pres_r = qr_prim_rsx_vf(j + 1, k, l, eqn_idx%E)
522
523 rho_l = 0._wp
524 gamma_l = 0._wp
525 pi_inf_l = 0._wp
526 qv_l = 0._wp
527
528 rho_r = 0._wp
529 gamma_r = 0._wp
530 pi_inf_r = 0._wp
531 qv_r = 0._wp
532
533 alpha_l_sum = 0._wp
534 alpha_r_sum = 0._wp
535
536 if (mpp_lim) then
537
538# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
539#if defined(MFC_OpenACC)
540# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
541!$acc loop seq
542# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
543#elif defined(MFC_OpenMP)
544# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
545
546# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
547#endif
548 do i = 1, num_fluids
549 alpha_rho_l(i) = max(0._wp, alpha_rho_l(i))
550 alpha_l(i) = min(max(0._wp, alpha_l(i)), 1._wp)
551 alpha_l_sum = alpha_l_sum + alpha_l(i)
552 alpha_rho_r(i) = max(0._wp, alpha_rho_r(i))
553 alpha_r(i) = min(max(0._wp, alpha_r(i)), 1._wp)
554 alpha_r_sum = alpha_r_sum + alpha_r(i)
555 end do
556
557 alpha_l = alpha_l/max(alpha_l_sum, sgm_eps)
558 alpha_r = alpha_r/max(alpha_r_sum, sgm_eps)
559 end if
560
561 call s_accumulate_mixture_properties(num_fluids, alpha_rho_l, alpha_l, rho_l, gamma_l, pi_inf_l, qv_l)
562 call s_accumulate_mixture_properties(num_fluids, alpha_rho_r, alpha_r, rho_r, gamma_r, pi_inf_r, qv_r)
563
564 if (viscous) then
565 call s_compute_interface_reynolds(alpha_l, re_l, re_size_loc1, re_size_loc2)
566 call s_compute_interface_reynolds(alpha_r, re_r, re_size_loc1, re_size_loc2)
567 end if
568
569 if (chemistry) then
570
571# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
572#if defined(MFC_OpenACC)
573# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
574!$acc loop seq
575# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
576#elif defined(MFC_OpenMP)
577# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
578
579# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
580#endif
581 do i = eqn_idx%species%beg, eqn_idx%species%end
582 ys_l(i - eqn_idx%species%beg + 1) = ql_prim_rsx_vf(j, k, l, i)
583 ys_r(i - eqn_idx%species%beg + 1) = qr_prim_rsx_vf(j + 1, k, l, i)
584 end do
585
586 call get_mixture_molecular_weight(ys_l, mw_l)
587 call get_mixture_molecular_weight(ys_r, mw_r)
588
589 xs_l(:) = ys_l(:)*mw_l/molecular_weights(:)
590 xs_r(:) = ys_r(:)*mw_r/molecular_weights(:)
591
592 r_gas_l = gas_constant/mw_l
593 r_gas_r = gas_constant/mw_r
594 t_l = pres_l/rho_l/r_gas_l
595 t_r = pres_r/rho_r/r_gas_r
596
597 call get_species_specific_heats_r(t_l, cp_il)
598 call get_species_specific_heats_r(t_r, cp_ir)
599
600 if (chem_params%gamma_method == 1) then
601 ! gamma_method = 1: Ref. Section 2.3.1 Formulation of doi:10.7907/ZKW8-ES97.
602 gamma_il = cp_il/(cp_il - 1.0_wp)
603 gamma_ir = cp_ir/(cp_ir - 1.0_wp)
604
605 gamma_l = sum(xs_l(:)/(gamma_il(:) - 1.0_wp))
606 gamma_r = sum(xs_r(:)/(gamma_ir(:) - 1.0_wp))
607 else if (chem_params%gamma_method == 2) then
608 ! gamma_method = 2: c_p / c_v where c_p, c_v are specific heats.
609 call get_mixture_specific_heat_cp_mass(t_l, ys_l, cp_l)
610 call get_mixture_specific_heat_cp_mass(t_r, ys_r, cp_r)
611 call get_mixture_specific_heat_cv_mass(t_l, ys_l, cv_l)
612 call get_mixture_specific_heat_cv_mass(t_r, ys_r, cv_r)
613
614 gamm_l = cp_l/cv_l
615 gamma_l = 1.0_wp/(gamm_l - 1.0_wp)
616 gamm_r = cp_r/cv_r
617 gamma_r = 1.0_wp/(gamm_r - 1.0_wp)
618 end if
619
620 call get_mixture_energy_mass(t_l, ys_l, e_l)
621 call get_mixture_energy_mass(t_r, ys_r, e_r)
622
623 e_l = rho_l*e_l + 5.e-1*rho_l*vel_l_rms
624 e_r = rho_r*e_r + 5.e-1*rho_r*vel_r_rms
625 h_l = (e_l + pres_l)/rho_l
626 h_r = (e_r + pres_r)/rho_r
627 else
628 e_l = gamma_l*pres_l + pi_inf_l + 5.e-1*rho_l*vel_l_rms + qv_l
629 e_r = gamma_r*pres_r + pi_inf_r + 5.e-1*rho_r*vel_r_rms + qv_r
630 h_l = (e_l + pres_l)/rho_l
631 h_r = (e_r + pres_r)/rho_r
632 end if
633
634 call s_compute_speed_of_sound(pres_l, rho_l, gamma_l, pi_inf_l, h_l, alpha_l, vel_l_rms, 0._wp, c_l, &
635 & qv_l)
636
637 call s_compute_speed_of_sound(pres_r, rho_r, gamma_r, pi_inf_r, h_r, alpha_r, vel_r_rms, 0._wp, c_r, &
638 & qv_r)
639
640 s_l = 0._wp; s_r = 0._wp
641
642
643# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
644#if defined(MFC_OpenACC)
645# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
646!$acc loop seq
647# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
648#elif defined(MFC_OpenMP)
649# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
650
651# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
652#endif
653 do i = 1, num_dims
654 s_l = s_l + vel_l(i)**2._wp
655 s_r = s_r + vel_r(i)**2._wp
656 end do
657
658 s_l = sqrt(s_l)
659 s_r = sqrt(s_r)
660
661 s_p = max(s_l, s_r) + max(c_l, c_r)
662 s_m = -s_p
663
664 s_l = s_m
665 s_r = s_p
666
667 ! Low Mach correction
668 if (low_mach == 1) then
669 if (riemann_solver == riemann_solver_hll .or. riemann_solver == riemann_solver_lax_friedrichs) then
670# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
671 zcoef = min(1._wp, max(vel_l_rms**5.e-1_wp/c_l, vel_r_rms**5.e-1_wp/c_r))
672# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
673 pcorr = 0._wp
674# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
675
676# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
677 if (low_mach == 1) then
678# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
679 pcorr = -(s_p - s_m)*(rho_l + rho_r)/8._wp*(zcoef - 1._wp)
680# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
681 end if
682# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
683 else if (riemann_solver == riemann_solver_hllc) then
684# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
685 zcoef = min(1._wp, max(vel_l_rms**5.e-1_wp/c_l, vel_r_rms**5.e-1_wp/c_r))
686# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
687 pcorr = 0._wp
688# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
689
690# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
691 if (low_mach == 1) then
692# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
693 pcorr = rho_l*rho_r*(s_l - vel_l(dir_idx(1)))*(s_r - vel_r(dir_idx(1)))*(vel_r(dir_idx(1)) - vel_l(dir_idx(1))) &
694# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
695 & /(rho_r*(s_r - vel_r(dir_idx(1))) - rho_l*(s_l - vel_l(dir_idx(1))))*(zcoef - 1._wp)
696# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
697 else if (low_mach == 2) then
698# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
699 vel_l_tmp = 5.e-1_wp*((vel_l(dir_idx(1)) + vel_r(dir_idx(1))) + zcoef*(vel_l(dir_idx(1)) - vel_r(dir_idx(1))))
700# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
701 vel_r_tmp = 5.e-1_wp*((vel_l(dir_idx(1)) + vel_r(dir_idx(1))) + zcoef*(vel_r(dir_idx(1)) - vel_l(dir_idx(1))))
702# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
703 vel_l(dir_idx(1)) = vel_l_tmp
704# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
705 vel_r(dir_idx(1)) = vel_r_tmp
706# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
707 end if
708# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
709 end if
710 else
711 pcorr = 0._wp
712 end if
713
714 ! Mass
715
716# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
717#if defined(MFC_OpenACC)
718# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
719!$acc loop seq
720# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
721#elif defined(MFC_OpenMP)
722# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
723
724# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
725#endif
726 do i = 1, eqn_idx%cont%end
727 flux_rsx_vf(j, k, l, &
728 & i) = (s_m*alpha_rho_r(i)*vel_r(norm_dir) - s_p*alpha_rho_l(i)*vel_l(norm_dir) &
729 & + s_m*s_p*(alpha_rho_l(i) - alpha_rho_r(i)))/(s_m - s_p)
730 end do
731
732 ! Momentum
733 if (bubbles_euler) then
734
735# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
736#if defined(MFC_OpenACC)
737# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
738!$acc loop seq
739# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
740#elif defined(MFC_OpenMP)
741# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
742
743# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
744#endif
745 do i = 1, num_vels
746 flux_rsx_vf(j, k, l, &
747 & eqn_idx%cont%end + dir_idx(i)) = (s_m*(rho_r*vel_r(dir_idx(1))*vel_r(dir_idx(i)) &
748 & + dir_flg(dir_idx(i))*(pres_r - ptilde_r)) - s_p*(rho_l*vel_l(dir_idx(1)) &
749 & *vel_l(dir_idx(i)) + dir_flg(dir_idx(i))*(pres_l - ptilde_l)) &
750 & + s_m*s_p*(rho_l*vel_l(dir_idx(i)) - rho_r*vel_r(dir_idx(i))))/(s_m - s_p) &
751 & + (s_m/s_l)*(s_p/s_r)*pcorr*(vel_r(dir_idx(i)) - vel_l(dir_idx(i)))
752 end do
753 else
754
755# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
756#if defined(MFC_OpenACC)
757# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
758!$acc loop seq
759# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
760#elif defined(MFC_OpenMP)
761# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
762
763# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
764#endif
765 do i = 1, num_vels
766 flux_rsx_vf(j, k, l, &
767 & eqn_idx%cont%end + dir_idx(i)) = (s_m*(rho_r*vel_r(dir_idx(1))*vel_r(dir_idx(i)) &
768 & + dir_flg(dir_idx(i))*pres_r) - s_p*(rho_l*vel_l(dir_idx(1))*vel_l(dir_idx(i)) &
769 & + dir_flg(dir_idx(i))*pres_l) + s_m*s_p*(rho_l*vel_l(dir_idx(i)) &
770 & - rho_r*vel_r(dir_idx(i))))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r) &
771 & *pcorr*(vel_r(dir_idx(i)) - vel_l(dir_idx(i)))
772 end do
773 end if
774
775 ! Energy
776 if (bubbles_euler) then
777 flux_rsx_vf(j, k, l, &
778 & eqn_idx%E) = (s_m*vel_r(dir_idx(1))*(e_r + pres_r - ptilde_r) - s_p*vel_l(dir_idx(1) &
779 & )*(e_l + pres_l - ptilde_l) + s_m*s_p*(e_l - e_r))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r) &
780 & *pcorr*(vel_r_rms - vel_l_rms)/2._wp
781 else
782 flux_rsx_vf(j, k, l, &
783 & eqn_idx%E) = (s_m*vel_r(dir_idx(1))*(e_r + pres_r) - s_p*vel_l(dir_idx(1))*(e_l &
784 & + pres_l) + s_m*s_p*(e_l - e_r))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r)*pcorr*(vel_r_rms &
785 & - vel_l_rms)/2._wp
786 end if
787
788 ! Advection flux and source: interface velocity for volume fraction transport
789
790# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
791#if defined(MFC_OpenACC)
792# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
793!$acc loop seq
794# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
795#elif defined(MFC_OpenMP)
796# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
797
798# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
799#endif
800 do i = eqn_idx%adv%beg, eqn_idx%adv%end
801 flux_rsx_vf(j, k, l, i) = (ql_prim_rsx_vf(j, k, l, i) - qr_prim_rsx_vf(j + 1, k, l, &
802 & i))*s_m*s_p/(s_m - s_p)
803 flux_src_rsx_vf(j, k, l, i) = (s_m*qr_prim_rsx_vf(j + 1, k, l, &
804 & i) - s_p*ql_prim_rsx_vf(j, k, l, i))/(s_m - s_p)
805 end do
806
807 if (bubbles_euler) then
808 ! From HLLC: Kills mass transport @ bubble gas density
809 if (num_fluids > 1) then
810 flux_rsx_vf(j, k, l, eqn_idx%cont%end) = 0._wp
811 end if
812 end if
813
814 if (chemistry) then
815
816# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
817#if defined(MFC_OpenACC)
818# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
819!$acc loop seq
820# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
821#elif defined(MFC_OpenMP)
822# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
823
824# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
825#endif
826 do i = eqn_idx%species%beg, eqn_idx%species%end
827 y_l = ql_prim_rsx_vf(j, k, l, i)
828 y_r = qr_prim_rsx_vf(j + 1, k, l, i)
829
830 flux_rsx_vf(j, k, l, &
831 & i) = (s_m*y_r*rho_r*vel_r(dir_idx(1)) - s_p*y_l*rho_l*vel_l(dir_idx(1)) &
832 & + s_m*s_p*(y_l*rho_l - y_r*rho_r))/(s_m - s_p)
833 flux_src_rsx_vf(j, k, l, i) = 0._wp
834 end do
835 end if
836
837# 352 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
838 end do
839 end do
840 end do
841
842# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
843#if defined(MFC_OpenACC)
844# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
845!$acc end parallel loop
846# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
847#elif defined(MFC_OpenMP)
848# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
849
850# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
851!$omp end target teams loop
852# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
853#endif
854 end if
855# 102 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
856# 103 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
857# 104 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
858 if (norm_dir == 2) then
859
860# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
861
862# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
863#if defined(MFC_OpenACC)
864# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
865!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, l, alpha_rho_L, alpha_rho_R, vel_L, vel_R, alpha_L, alpha_R, Re_L, Re_R, rho_avg, h_avg, gamma_avg, s_L, s_R, s_S, Ys_L, &
866# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
867!$acc& Ys_R, Cp_iL, Cp_iR, Xs_L, Xs_R, Gamma_iL, Gamma_iR, Yi_avg, Phi_avg, h_iL, h_iR, h_avg_2, pcorr, zcoef, vel_grad_L, vel_grad_R, idx_right_phys, vel_L_rms, vel_R_rms, vel_avg_rms, vel_L_tmp, &
868# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
869!$acc& vel_R_tmp, Ms_L, Ms_R, pres_SL, pres_SR, alpha_L_sum, alpha_R_sum, c_avg, pres_L, pres_R, rho_L, rho_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, qv_L, qv_R, c_L, c_R, E_L, E_R, H_L, H_R, ptilde_L, &
870# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
871!$acc& ptilde_R, s_M, s_P, xi_M, xi_P, Cp_avg, Cv_avg, T_avg, eps, c_sum_Yi_Phi, Cp_L, Cp_R, Cv_L, Cv_R, R_gas_L, R_gas_R, MW_L, MW_R, T_L, T_R, Y_L, Y_R) firstprivate(Re_size_loc1, Re_size_loc2)
872# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
873#elif defined(MFC_OpenMP)
874# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
875
876# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
877
878# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
879
880# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
881!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(i, j, k, l, &
882# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
883!$omp& alpha_rho_L, alpha_rho_R, vel_L, vel_R, alpha_L, alpha_R, Re_L, Re_R, rho_avg, h_avg, gamma_avg, s_L, s_R, s_S, Ys_L, Ys_R, Cp_iL, Cp_iR, Xs_L, Xs_R, Gamma_iL, Gamma_iR, Yi_avg, Phi_avg, h_iL, &
884# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
885!$omp& h_iR, h_avg_2, pcorr, zcoef, vel_grad_L, vel_grad_R, idx_right_phys, vel_L_rms, vel_R_rms, vel_avg_rms, vel_L_tmp, vel_R_tmp, Ms_L, Ms_R, pres_SL, pres_SR, alpha_L_sum, alpha_R_sum, c_avg, &
886# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
887!$omp& pres_L, pres_R, rho_L, rho_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, qv_L, qv_R, c_L, c_R, E_L, E_R, H_L, H_R, ptilde_L, ptilde_R, s_M, s_P, xi_M, xi_P, Cp_avg, Cv_avg, T_avg, eps, c_sum_Yi_Phi, &
888# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
889!$omp& Cp_L, Cp_R, Cv_L, Cv_R, R_gas_L, R_gas_R, MW_L, MW_R, T_L, T_R, Y_L, Y_R) firstprivate(Re_size_loc1, Re_size_loc2)
890# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
891#endif
892# 113 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
893 do l = is3%beg, is3%end
894 do k = is1%beg, is1%end
895 do j = is2%beg, is2%end
896
897# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
898#if defined(MFC_OpenACC)
899# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
900!$acc loop seq
901# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
902#elif defined(MFC_OpenMP)
903# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
904
905# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
906#endif
907 do i = 1, eqn_idx%cont%end
908 alpha_rho_l(i) = ql_prim_rsx_vf(j, k, l, i)
909 alpha_rho_r(i) = qr_prim_rsx_vf(j, k + 1, l, i)
910 end do
911
912 vel_l_rms = 0._wp; vel_r_rms = 0._wp
913
914
915# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
916#if defined(MFC_OpenACC)
917# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
918!$acc loop seq
919# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
920#elif defined(MFC_OpenMP)
921# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
922
923# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
924#endif
925 do i = 1, num_vels
926 vel_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%cont%end + i)
927 vel_r(i) = qr_prim_rsx_vf(j, k + 1, l, eqn_idx%cont%end + i)
928 vel_l_rms = vel_l_rms + vel_l(i)**2._wp
929 vel_r_rms = vel_r_rms + vel_r(i)**2._wp
930 end do
931
932
933# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
934#if defined(MFC_OpenACC)
935# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
936!$acc loop seq
937# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
938#elif defined(MFC_OpenMP)
939# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
940
941# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
942#endif
943 do i = 1, num_fluids
944 alpha_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%E + i)
945 alpha_r(i) = qr_prim_rsx_vf(j, k + 1, l, eqn_idx%E + i)
946 end do
947
948 pres_l = ql_prim_rsx_vf(j, k, l, eqn_idx%E)
949 pres_r = qr_prim_rsx_vf(j, k + 1, l, eqn_idx%E)
950
951 rho_l = 0._wp
952 gamma_l = 0._wp
953 pi_inf_l = 0._wp
954 qv_l = 0._wp
955
956 rho_r = 0._wp
957 gamma_r = 0._wp
958 pi_inf_r = 0._wp
959 qv_r = 0._wp
960
961 alpha_l_sum = 0._wp
962 alpha_r_sum = 0._wp
963
964 if (mpp_lim) then
965
966# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
967#if defined(MFC_OpenACC)
968# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
969!$acc loop seq
970# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
971#elif defined(MFC_OpenMP)
972# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
973
974# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
975#endif
976 do i = 1, num_fluids
977 alpha_rho_l(i) = max(0._wp, alpha_rho_l(i))
978 alpha_l(i) = min(max(0._wp, alpha_l(i)), 1._wp)
979 alpha_l_sum = alpha_l_sum + alpha_l(i)
980 alpha_rho_r(i) = max(0._wp, alpha_rho_r(i))
981 alpha_r(i) = min(max(0._wp, alpha_r(i)), 1._wp)
982 alpha_r_sum = alpha_r_sum + alpha_r(i)
983 end do
984
985 alpha_l = alpha_l/max(alpha_l_sum, sgm_eps)
986 alpha_r = alpha_r/max(alpha_r_sum, sgm_eps)
987 end if
988
989 call s_accumulate_mixture_properties(num_fluids, alpha_rho_l, alpha_l, rho_l, gamma_l, pi_inf_l, qv_l)
990 call s_accumulate_mixture_properties(num_fluids, alpha_rho_r, alpha_r, rho_r, gamma_r, pi_inf_r, qv_r)
991
992 if (viscous) then
993 call s_compute_interface_reynolds(alpha_l, re_l, re_size_loc1, re_size_loc2)
994 call s_compute_interface_reynolds(alpha_r, re_r, re_size_loc1, re_size_loc2)
995 end if
996
997 if (chemistry) then
998
999# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1000#if defined(MFC_OpenACC)
1001# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1002!$acc loop seq
1003# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1004#elif defined(MFC_OpenMP)
1005# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1006
1007# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1008#endif
1009 do i = eqn_idx%species%beg, eqn_idx%species%end
1010 ys_l(i - eqn_idx%species%beg + 1) = ql_prim_rsx_vf(j, k, l, i)
1011 ys_r(i - eqn_idx%species%beg + 1) = qr_prim_rsx_vf(j, k + 1, l, i)
1012 end do
1013
1014 call get_mixture_molecular_weight(ys_l, mw_l)
1015 call get_mixture_molecular_weight(ys_r, mw_r)
1016
1017 xs_l(:) = ys_l(:)*mw_l/molecular_weights(:)
1018 xs_r(:) = ys_r(:)*mw_r/molecular_weights(:)
1019
1020 r_gas_l = gas_constant/mw_l
1021 r_gas_r = gas_constant/mw_r
1022 t_l = pres_l/rho_l/r_gas_l
1023 t_r = pres_r/rho_r/r_gas_r
1024
1025 call get_species_specific_heats_r(t_l, cp_il)
1026 call get_species_specific_heats_r(t_r, cp_ir)
1027
1028 if (chem_params%gamma_method == 1) then
1029 ! gamma_method = 1: Ref. Section 2.3.1 Formulation of doi:10.7907/ZKW8-ES97.
1030 gamma_il = cp_il/(cp_il - 1.0_wp)
1031 gamma_ir = cp_ir/(cp_ir - 1.0_wp)
1032
1033 gamma_l = sum(xs_l(:)/(gamma_il(:) - 1.0_wp))
1034 gamma_r = sum(xs_r(:)/(gamma_ir(:) - 1.0_wp))
1035 else if (chem_params%gamma_method == 2) then
1036 ! gamma_method = 2: c_p / c_v where c_p, c_v are specific heats.
1037 call get_mixture_specific_heat_cp_mass(t_l, ys_l, cp_l)
1038 call get_mixture_specific_heat_cp_mass(t_r, ys_r, cp_r)
1039 call get_mixture_specific_heat_cv_mass(t_l, ys_l, cv_l)
1040 call get_mixture_specific_heat_cv_mass(t_r, ys_r, cv_r)
1041
1042 gamm_l = cp_l/cv_l
1043 gamma_l = 1.0_wp/(gamm_l - 1.0_wp)
1044 gamm_r = cp_r/cv_r
1045 gamma_r = 1.0_wp/(gamm_r - 1.0_wp)
1046 end if
1047
1048 call get_mixture_energy_mass(t_l, ys_l, e_l)
1049 call get_mixture_energy_mass(t_r, ys_r, e_r)
1050
1051 e_l = rho_l*e_l + 5.e-1*rho_l*vel_l_rms
1052 e_r = rho_r*e_r + 5.e-1*rho_r*vel_r_rms
1053 h_l = (e_l + pres_l)/rho_l
1054 h_r = (e_r + pres_r)/rho_r
1055 else
1056 e_l = gamma_l*pres_l + pi_inf_l + 5.e-1*rho_l*vel_l_rms + qv_l
1057 e_r = gamma_r*pres_r + pi_inf_r + 5.e-1*rho_r*vel_r_rms + qv_r
1058 h_l = (e_l + pres_l)/rho_l
1059 h_r = (e_r + pres_r)/rho_r
1060 end if
1061
1062 call s_compute_speed_of_sound(pres_l, rho_l, gamma_l, pi_inf_l, h_l, alpha_l, vel_l_rms, 0._wp, c_l, &
1063 & qv_l)
1064
1065 call s_compute_speed_of_sound(pres_r, rho_r, gamma_r, pi_inf_r, h_r, alpha_r, vel_r_rms, 0._wp, c_r, &
1066 & qv_r)
1067
1068 s_l = 0._wp; s_r = 0._wp
1069
1070
1071# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1072#if defined(MFC_OpenACC)
1073# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1074!$acc loop seq
1075# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1076#elif defined(MFC_OpenMP)
1077# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1078
1079# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1080#endif
1081 do i = 1, num_dims
1082 s_l = s_l + vel_l(i)**2._wp
1083 s_r = s_r + vel_r(i)**2._wp
1084 end do
1085
1086 s_l = sqrt(s_l)
1087 s_r = sqrt(s_r)
1088
1089 s_p = max(s_l, s_r) + max(c_l, c_r)
1090 s_m = -s_p
1091
1092 s_l = s_m
1093 s_r = s_p
1094
1095 ! Low Mach correction
1096 if (low_mach == 1) then
1097 if (riemann_solver == riemann_solver_hll .or. riemann_solver == riemann_solver_lax_friedrichs) then
1098# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1099 zcoef = min(1._wp, max(vel_l_rms**5.e-1_wp/c_l, vel_r_rms**5.e-1_wp/c_r))
1100# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1101 pcorr = 0._wp
1102# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1103
1104# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1105 if (low_mach == 1) then
1106# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1107 pcorr = -(s_p - s_m)*(rho_l + rho_r)/8._wp*(zcoef - 1._wp)
1108# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1109 end if
1110# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1111 else if (riemann_solver == riemann_solver_hllc) then
1112# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1113 zcoef = min(1._wp, max(vel_l_rms**5.e-1_wp/c_l, vel_r_rms**5.e-1_wp/c_r))
1114# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1115 pcorr = 0._wp
1116# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1117
1118# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1119 if (low_mach == 1) then
1120# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1121 pcorr = rho_l*rho_r*(s_l - vel_l(dir_idx(1)))*(s_r - vel_r(dir_idx(1)))*(vel_r(dir_idx(1)) - vel_l(dir_idx(1))) &
1122# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1123 & /(rho_r*(s_r - vel_r(dir_idx(1))) - rho_l*(s_l - vel_l(dir_idx(1))))*(zcoef - 1._wp)
1124# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1125 else if (low_mach == 2) then
1126# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1127 vel_l_tmp = 5.e-1_wp*((vel_l(dir_idx(1)) + vel_r(dir_idx(1))) + zcoef*(vel_l(dir_idx(1)) - vel_r(dir_idx(1))))
1128# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1129 vel_r_tmp = 5.e-1_wp*((vel_l(dir_idx(1)) + vel_r(dir_idx(1))) + zcoef*(vel_r(dir_idx(1)) - vel_l(dir_idx(1))))
1130# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1131 vel_l(dir_idx(1)) = vel_l_tmp
1132# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1133 vel_r(dir_idx(1)) = vel_r_tmp
1134# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1135 end if
1136# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1137 end if
1138 else
1139 pcorr = 0._wp
1140 end if
1141
1142 ! Mass
1143
1144# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1145#if defined(MFC_OpenACC)
1146# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1147!$acc loop seq
1148# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1149#elif defined(MFC_OpenMP)
1150# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1151
1152# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1153#endif
1154 do i = 1, eqn_idx%cont%end
1155 flux_rsx_vf(j, k, l, &
1156 & i) = (s_m*alpha_rho_r(i)*vel_r(norm_dir) - s_p*alpha_rho_l(i)*vel_l(norm_dir) &
1157 & + s_m*s_p*(alpha_rho_l(i) - alpha_rho_r(i)))/(s_m - s_p)
1158 end do
1159
1160 ! Momentum
1161 if (bubbles_euler) then
1162
1163# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1164#if defined(MFC_OpenACC)
1165# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1166!$acc loop seq
1167# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1168#elif defined(MFC_OpenMP)
1169# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1170
1171# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1172#endif
1173 do i = 1, num_vels
1174 flux_rsx_vf(j, k, l, &
1175 & eqn_idx%cont%end + dir_idx(i)) = (s_m*(rho_r*vel_r(dir_idx(1))*vel_r(dir_idx(i)) &
1176 & + dir_flg(dir_idx(i))*(pres_r - ptilde_r)) - s_p*(rho_l*vel_l(dir_idx(1)) &
1177 & *vel_l(dir_idx(i)) + dir_flg(dir_idx(i))*(pres_l - ptilde_l)) &
1178 & + s_m*s_p*(rho_l*vel_l(dir_idx(i)) - rho_r*vel_r(dir_idx(i))))/(s_m - s_p) &
1179 & + (s_m/s_l)*(s_p/s_r)*pcorr*(vel_r(dir_idx(i)) - vel_l(dir_idx(i)))
1180 end do
1181 else
1182
1183# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1184#if defined(MFC_OpenACC)
1185# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1186!$acc loop seq
1187# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1188#elif defined(MFC_OpenMP)
1189# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1190
1191# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1192#endif
1193 do i = 1, num_vels
1194 flux_rsx_vf(j, k, l, &
1195 & eqn_idx%cont%end + dir_idx(i)) = (s_m*(rho_r*vel_r(dir_idx(1))*vel_r(dir_idx(i)) &
1196 & + dir_flg(dir_idx(i))*pres_r) - s_p*(rho_l*vel_l(dir_idx(1))*vel_l(dir_idx(i)) &
1197 & + dir_flg(dir_idx(i))*pres_l) + s_m*s_p*(rho_l*vel_l(dir_idx(i)) &
1198 & - rho_r*vel_r(dir_idx(i))))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r) &
1199 & *pcorr*(vel_r(dir_idx(i)) - vel_l(dir_idx(i)))
1200 end do
1201 end if
1202
1203 ! Energy
1204 if (bubbles_euler) then
1205 flux_rsx_vf(j, k, l, &
1206 & eqn_idx%E) = (s_m*vel_r(dir_idx(1))*(e_r + pres_r - ptilde_r) - s_p*vel_l(dir_idx(1) &
1207 & )*(e_l + pres_l - ptilde_l) + s_m*s_p*(e_l - e_r))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r) &
1208 & *pcorr*(vel_r_rms - vel_l_rms)/2._wp
1209 else
1210 flux_rsx_vf(j, k, l, &
1211 & eqn_idx%E) = (s_m*vel_r(dir_idx(1))*(e_r + pres_r) - s_p*vel_l(dir_idx(1))*(e_l &
1212 & + pres_l) + s_m*s_p*(e_l - e_r))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r)*pcorr*(vel_r_rms &
1213 & - vel_l_rms)/2._wp
1214 end if
1215
1216 ! Advection flux and source: interface velocity for volume fraction transport
1217
1218# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1219#if defined(MFC_OpenACC)
1220# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1221!$acc loop seq
1222# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1223#elif defined(MFC_OpenMP)
1224# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1225
1226# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1227#endif
1228 do i = eqn_idx%adv%beg, eqn_idx%adv%end
1229 flux_rsx_vf(j, k, l, i) = (ql_prim_rsx_vf(j, k, l, i) - qr_prim_rsx_vf(j, k + 1, l, &
1230 & i))*s_m*s_p/(s_m - s_p)
1231 flux_src_rsx_vf(j, k, l, i) = (s_m*qr_prim_rsx_vf(j, k + 1, l, &
1232 & i) - s_p*ql_prim_rsx_vf(j, k, l, i))/(s_m - s_p)
1233 end do
1234
1235 if (bubbles_euler) then
1236 ! From HLLC: Kills mass transport @ bubble gas density
1237 if (num_fluids > 1) then
1238 flux_rsx_vf(j, k, l, eqn_idx%cont%end) = 0._wp
1239 end if
1240 end if
1241
1242 if (chemistry) then
1243
1244# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1245#if defined(MFC_OpenACC)
1246# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1247!$acc loop seq
1248# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1249#elif defined(MFC_OpenMP)
1250# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1251
1252# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1253#endif
1254 do i = eqn_idx%species%beg, eqn_idx%species%end
1255 y_l = ql_prim_rsx_vf(j, k, l, i)
1256 y_r = qr_prim_rsx_vf(j, k + 1, l, i)
1257
1258 flux_rsx_vf(j, k, l, &
1259 & i) = (s_m*y_r*rho_r*vel_r(dir_idx(1)) - s_p*y_l*rho_l*vel_l(dir_idx(1)) &
1260 & + s_m*s_p*(y_l*rho_l - y_r*rho_r))/(s_m - s_p)
1261 flux_src_rsx_vf(j, k, l, i) = 0._wp
1262 end do
1263 end if
1264
1265# 336 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1266 if (cyl_coord) then
1267 ! Substituting the advective flux into the inviscid geometrical source flux
1268
1269# 338 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1270#if defined(MFC_OpenACC)
1271# 338 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1272!$acc loop seq
1273# 338 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1274#elif defined(MFC_OpenMP)
1275# 338 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1276
1277# 338 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1278#endif
1279 do i = 1, eqn_idx%E
1280 flux_gsrc_rsx_vf(j, k, l, i) = flux_rsx_vf(j, k, l, i)
1281 end do
1282 ! Recalculating the radial momentum geometric source flux
1283 flux_gsrc_rsx_vf(j, k, l, eqn_idx%cont%end + 2) = flux_rsx_vf(j, k, l, &
1284 & eqn_idx%cont%end + 2) - (s_m*pres_r - s_p*pres_l)/(s_m - s_p)
1285 ! Geometrical source of the void fraction(s) is zero
1286
1287# 346 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1288#if defined(MFC_OpenACC)
1289# 346 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1290!$acc loop seq
1291# 346 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1292#elif defined(MFC_OpenMP)
1293# 346 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1294
1295# 346 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1296#endif
1297 do i = eqn_idx%adv%beg, eqn_idx%adv%end
1298 flux_gsrc_rsx_vf(j, k, l, i) = flux_rsx_vf(j, k, l, i)
1299 end do
1300 end if
1301# 352 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1302 end do
1303 end do
1304 end do
1305
1306# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1307#if defined(MFC_OpenACC)
1308# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1309!$acc end parallel loop
1310# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1311#elif defined(MFC_OpenMP)
1312# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1313
1314# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1315!$omp end target teams loop
1316# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1317#endif
1318 end if
1319# 102 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1320# 103 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1321# 104 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1322 if (norm_dir == 3) then
1323
1324# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1325
1326# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1327#if defined(MFC_OpenACC)
1328# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1329!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, l, alpha_rho_L, alpha_rho_R, vel_L, vel_R, alpha_L, alpha_R, Re_L, Re_R, rho_avg, h_avg, gamma_avg, s_L, s_R, s_S, Ys_L, &
1330# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1331!$acc& Ys_R, Cp_iL, Cp_iR, Xs_L, Xs_R, Gamma_iL, Gamma_iR, Yi_avg, Phi_avg, h_iL, h_iR, h_avg_2, pcorr, zcoef, vel_grad_L, vel_grad_R, idx_right_phys, vel_L_rms, vel_R_rms, vel_avg_rms, vel_L_tmp, &
1332# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1333!$acc& vel_R_tmp, Ms_L, Ms_R, pres_SL, pres_SR, alpha_L_sum, alpha_R_sum, c_avg, pres_L, pres_R, rho_L, rho_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, qv_L, qv_R, c_L, c_R, E_L, E_R, H_L, H_R, ptilde_L, &
1334# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1335!$acc& ptilde_R, s_M, s_P, xi_M, xi_P, Cp_avg, Cv_avg, T_avg, eps, c_sum_Yi_Phi, Cp_L, Cp_R, Cv_L, Cv_R, R_gas_L, R_gas_R, MW_L, MW_R, T_L, T_R, Y_L, Y_R) firstprivate(Re_size_loc1, Re_size_loc2)
1336# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1337#elif defined(MFC_OpenMP)
1338# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1339
1340# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1341
1342# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1343
1344# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1345!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(i, j, k, l, &
1346# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1347!$omp& alpha_rho_L, alpha_rho_R, vel_L, vel_R, alpha_L, alpha_R, Re_L, Re_R, rho_avg, h_avg, gamma_avg, s_L, s_R, s_S, Ys_L, Ys_R, Cp_iL, Cp_iR, Xs_L, Xs_R, Gamma_iL, Gamma_iR, Yi_avg, Phi_avg, h_iL, &
1348# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1349!$omp& h_iR, h_avg_2, pcorr, zcoef, vel_grad_L, vel_grad_R, idx_right_phys, vel_L_rms, vel_R_rms, vel_avg_rms, vel_L_tmp, vel_R_tmp, Ms_L, Ms_R, pres_SL, pres_SR, alpha_L_sum, alpha_R_sum, c_avg, &
1350# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1351!$omp& pres_L, pres_R, rho_L, rho_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, qv_L, qv_R, c_L, c_R, E_L, E_R, H_L, H_R, ptilde_L, ptilde_R, s_M, s_P, xi_M, xi_P, Cp_avg, Cv_avg, T_avg, eps, c_sum_Yi_Phi, &
1352# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1353!$omp& Cp_L, Cp_R, Cv_L, Cv_R, R_gas_L, R_gas_R, MW_L, MW_R, T_L, T_R, Y_L, Y_R) firstprivate(Re_size_loc1, Re_size_loc2)
1354# 105 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1355#endif
1356# 113 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1357 do l = is1%beg, is1%end
1358 do k = is2%beg, is2%end
1359 do j = is3%beg, is3%end
1360
1361# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1362#if defined(MFC_OpenACC)
1363# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1364!$acc loop seq
1365# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1366#elif defined(MFC_OpenMP)
1367# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1368
1369# 116 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1370#endif
1371 do i = 1, eqn_idx%cont%end
1372 alpha_rho_l(i) = ql_prim_rsx_vf(j, k, l, i)
1373 alpha_rho_r(i) = qr_prim_rsx_vf(j, k, l + 1, i)
1374 end do
1375
1376 vel_l_rms = 0._wp; vel_r_rms = 0._wp
1377
1378
1379# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1380#if defined(MFC_OpenACC)
1381# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1382!$acc loop seq
1383# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1384#elif defined(MFC_OpenMP)
1385# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1386
1387# 124 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1388#endif
1389 do i = 1, num_vels
1390 vel_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%cont%end + i)
1391 vel_r(i) = qr_prim_rsx_vf(j, k, l + 1, eqn_idx%cont%end + i)
1392 vel_l_rms = vel_l_rms + vel_l(i)**2._wp
1393 vel_r_rms = vel_r_rms + vel_r(i)**2._wp
1394 end do
1395
1396
1397# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1398#if defined(MFC_OpenACC)
1399# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1400!$acc loop seq
1401# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1402#elif defined(MFC_OpenMP)
1403# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1404
1405# 132 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1406#endif
1407 do i = 1, num_fluids
1408 alpha_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%E + i)
1409 alpha_r(i) = qr_prim_rsx_vf(j, k, l + 1, eqn_idx%E + i)
1410 end do
1411
1412 pres_l = ql_prim_rsx_vf(j, k, l, eqn_idx%E)
1413 pres_r = qr_prim_rsx_vf(j, k, l + 1, eqn_idx%E)
1414
1415 rho_l = 0._wp
1416 gamma_l = 0._wp
1417 pi_inf_l = 0._wp
1418 qv_l = 0._wp
1419
1420 rho_r = 0._wp
1421 gamma_r = 0._wp
1422 pi_inf_r = 0._wp
1423 qv_r = 0._wp
1424
1425 alpha_l_sum = 0._wp
1426 alpha_r_sum = 0._wp
1427
1428 if (mpp_lim) then
1429
1430# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1431#if defined(MFC_OpenACC)
1432# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1433!$acc loop seq
1434# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1435#elif defined(MFC_OpenMP)
1436# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1437
1438# 155 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1439#endif
1440 do i = 1, num_fluids
1441 alpha_rho_l(i) = max(0._wp, alpha_rho_l(i))
1442 alpha_l(i) = min(max(0._wp, alpha_l(i)), 1._wp)
1443 alpha_l_sum = alpha_l_sum + alpha_l(i)
1444 alpha_rho_r(i) = max(0._wp, alpha_rho_r(i))
1445 alpha_r(i) = min(max(0._wp, alpha_r(i)), 1._wp)
1446 alpha_r_sum = alpha_r_sum + alpha_r(i)
1447 end do
1448
1449 alpha_l = alpha_l/max(alpha_l_sum, sgm_eps)
1450 alpha_r = alpha_r/max(alpha_r_sum, sgm_eps)
1451 end if
1452
1453 call s_accumulate_mixture_properties(num_fluids, alpha_rho_l, alpha_l, rho_l, gamma_l, pi_inf_l, qv_l)
1454 call s_accumulate_mixture_properties(num_fluids, alpha_rho_r, alpha_r, rho_r, gamma_r, pi_inf_r, qv_r)
1455
1456 if (viscous) then
1457 call s_compute_interface_reynolds(alpha_l, re_l, re_size_loc1, re_size_loc2)
1458 call s_compute_interface_reynolds(alpha_r, re_r, re_size_loc1, re_size_loc2)
1459 end if
1460
1461 if (chemistry) then
1462
1463# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1464#if defined(MFC_OpenACC)
1465# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1466!$acc loop seq
1467# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1468#elif defined(MFC_OpenMP)
1469# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1470
1471# 178 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1472#endif
1473 do i = eqn_idx%species%beg, eqn_idx%species%end
1474 ys_l(i - eqn_idx%species%beg + 1) = ql_prim_rsx_vf(j, k, l, i)
1475 ys_r(i - eqn_idx%species%beg + 1) = qr_prim_rsx_vf(j, k, l + 1, i)
1476 end do
1477
1478 call get_mixture_molecular_weight(ys_l, mw_l)
1479 call get_mixture_molecular_weight(ys_r, mw_r)
1480
1481 xs_l(:) = ys_l(:)*mw_l/molecular_weights(:)
1482 xs_r(:) = ys_r(:)*mw_r/molecular_weights(:)
1483
1484 r_gas_l = gas_constant/mw_l
1485 r_gas_r = gas_constant/mw_r
1486 t_l = pres_l/rho_l/r_gas_l
1487 t_r = pres_r/rho_r/r_gas_r
1488
1489 call get_species_specific_heats_r(t_l, cp_il)
1490 call get_species_specific_heats_r(t_r, cp_ir)
1491
1492 if (chem_params%gamma_method == 1) then
1493 ! gamma_method = 1: Ref. Section 2.3.1 Formulation of doi:10.7907/ZKW8-ES97.
1494 gamma_il = cp_il/(cp_il - 1.0_wp)
1495 gamma_ir = cp_ir/(cp_ir - 1.0_wp)
1496
1497 gamma_l = sum(xs_l(:)/(gamma_il(:) - 1.0_wp))
1498 gamma_r = sum(xs_r(:)/(gamma_ir(:) - 1.0_wp))
1499 else if (chem_params%gamma_method == 2) then
1500 ! gamma_method = 2: c_p / c_v where c_p, c_v are specific heats.
1501 call get_mixture_specific_heat_cp_mass(t_l, ys_l, cp_l)
1502 call get_mixture_specific_heat_cp_mass(t_r, ys_r, cp_r)
1503 call get_mixture_specific_heat_cv_mass(t_l, ys_l, cv_l)
1504 call get_mixture_specific_heat_cv_mass(t_r, ys_r, cv_r)
1505
1506 gamm_l = cp_l/cv_l
1507 gamma_l = 1.0_wp/(gamm_l - 1.0_wp)
1508 gamm_r = cp_r/cv_r
1509 gamma_r = 1.0_wp/(gamm_r - 1.0_wp)
1510 end if
1511
1512 call get_mixture_energy_mass(t_l, ys_l, e_l)
1513 call get_mixture_energy_mass(t_r, ys_r, e_r)
1514
1515 e_l = rho_l*e_l + 5.e-1*rho_l*vel_l_rms
1516 e_r = rho_r*e_r + 5.e-1*rho_r*vel_r_rms
1517 h_l = (e_l + pres_l)/rho_l
1518 h_r = (e_r + pres_r)/rho_r
1519 else
1520 e_l = gamma_l*pres_l + pi_inf_l + 5.e-1*rho_l*vel_l_rms + qv_l
1521 e_r = gamma_r*pres_r + pi_inf_r + 5.e-1*rho_r*vel_r_rms + qv_r
1522 h_l = (e_l + pres_l)/rho_l
1523 h_r = (e_r + pres_r)/rho_r
1524 end if
1525
1526 call s_compute_speed_of_sound(pres_l, rho_l, gamma_l, pi_inf_l, h_l, alpha_l, vel_l_rms, 0._wp, c_l, &
1527 & qv_l)
1528
1529 call s_compute_speed_of_sound(pres_r, rho_r, gamma_r, pi_inf_r, h_r, alpha_r, vel_r_rms, 0._wp, c_r, &
1530 & qv_r)
1531
1532 s_l = 0._wp; s_r = 0._wp
1533
1534
1535# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1536#if defined(MFC_OpenACC)
1537# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1538!$acc loop seq
1539# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1540#elif defined(MFC_OpenMP)
1541# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1542
1543# 240 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1544#endif
1545 do i = 1, num_dims
1546 s_l = s_l + vel_l(i)**2._wp
1547 s_r = s_r + vel_r(i)**2._wp
1548 end do
1549
1550 s_l = sqrt(s_l)
1551 s_r = sqrt(s_r)
1552
1553 s_p = max(s_l, s_r) + max(c_l, c_r)
1554 s_m = -s_p
1555
1556 s_l = s_m
1557 s_r = s_p
1558
1559 ! Low Mach correction
1560 if (low_mach == 1) then
1561 if (riemann_solver == riemann_solver_hll .or. riemann_solver == riemann_solver_lax_friedrichs) then
1562# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1563 zcoef = min(1._wp, max(vel_l_rms**5.e-1_wp/c_l, vel_r_rms**5.e-1_wp/c_r))
1564# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1565 pcorr = 0._wp
1566# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1567
1568# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1569 if (low_mach == 1) then
1570# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1571 pcorr = -(s_p - s_m)*(rho_l + rho_r)/8._wp*(zcoef - 1._wp)
1572# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1573 end if
1574# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1575 else if (riemann_solver == riemann_solver_hllc) then
1576# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1577 zcoef = min(1._wp, max(vel_l_rms**5.e-1_wp/c_l, vel_r_rms**5.e-1_wp/c_r))
1578# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1579 pcorr = 0._wp
1580# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1581
1582# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1583 if (low_mach == 1) then
1584# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1585 pcorr = rho_l*rho_r*(s_l - vel_l(dir_idx(1)))*(s_r - vel_r(dir_idx(1)))*(vel_r(dir_idx(1)) - vel_l(dir_idx(1))) &
1586# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1587 & /(rho_r*(s_r - vel_r(dir_idx(1))) - rho_l*(s_l - vel_l(dir_idx(1))))*(zcoef - 1._wp)
1588# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1589 else if (low_mach == 2) then
1590# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1591 vel_l_tmp = 5.e-1_wp*((vel_l(dir_idx(1)) + vel_r(dir_idx(1))) + zcoef*(vel_l(dir_idx(1)) - vel_r(dir_idx(1))))
1592# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1593 vel_r_tmp = 5.e-1_wp*((vel_l(dir_idx(1)) + vel_r(dir_idx(1))) + zcoef*(vel_r(dir_idx(1)) - vel_l(dir_idx(1))))
1594# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1595 vel_l(dir_idx(1)) = vel_l_tmp
1596# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1597 vel_r(dir_idx(1)) = vel_r_tmp
1598# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1599 end if
1600# 257 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1601 end if
1602 else
1603 pcorr = 0._wp
1604 end if
1605
1606 ! Mass
1607
1608# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1609#if defined(MFC_OpenACC)
1610# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1611!$acc loop seq
1612# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1613#elif defined(MFC_OpenMP)
1614# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1615
1616# 263 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1617#endif
1618 do i = 1, eqn_idx%cont%end
1619 flux_rsx_vf(j, k, l, &
1620 & i) = (s_m*alpha_rho_r(i)*vel_r(norm_dir) - s_p*alpha_rho_l(i)*vel_l(norm_dir) &
1621 & + s_m*s_p*(alpha_rho_l(i) - alpha_rho_r(i)))/(s_m - s_p)
1622 end do
1623
1624 ! Momentum
1625 if (bubbles_euler) then
1626
1627# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1628#if defined(MFC_OpenACC)
1629# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1630!$acc loop seq
1631# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1632#elif defined(MFC_OpenMP)
1633# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1634
1635# 272 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1636#endif
1637 do i = 1, num_vels
1638 flux_rsx_vf(j, k, l, &
1639 & eqn_idx%cont%end + dir_idx(i)) = (s_m*(rho_r*vel_r(dir_idx(1))*vel_r(dir_idx(i)) &
1640 & + dir_flg(dir_idx(i))*(pres_r - ptilde_r)) - s_p*(rho_l*vel_l(dir_idx(1)) &
1641 & *vel_l(dir_idx(i)) + dir_flg(dir_idx(i))*(pres_l - ptilde_l)) &
1642 & + s_m*s_p*(rho_l*vel_l(dir_idx(i)) - rho_r*vel_r(dir_idx(i))))/(s_m - s_p) &
1643 & + (s_m/s_l)*(s_p/s_r)*pcorr*(vel_r(dir_idx(i)) - vel_l(dir_idx(i)))
1644 end do
1645 else
1646
1647# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1648#if defined(MFC_OpenACC)
1649# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1650!$acc loop seq
1651# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1652#elif defined(MFC_OpenMP)
1653# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1654
1655# 282 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1656#endif
1657 do i = 1, num_vels
1658 flux_rsx_vf(j, k, l, &
1659 & eqn_idx%cont%end + dir_idx(i)) = (s_m*(rho_r*vel_r(dir_idx(1))*vel_r(dir_idx(i)) &
1660 & + dir_flg(dir_idx(i))*pres_r) - s_p*(rho_l*vel_l(dir_idx(1))*vel_l(dir_idx(i)) &
1661 & + dir_flg(dir_idx(i))*pres_l) + s_m*s_p*(rho_l*vel_l(dir_idx(i)) &
1662 & - rho_r*vel_r(dir_idx(i))))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r) &
1663 & *pcorr*(vel_r(dir_idx(i)) - vel_l(dir_idx(i)))
1664 end do
1665 end if
1666
1667 ! Energy
1668 if (bubbles_euler) then
1669 flux_rsx_vf(j, k, l, &
1670 & eqn_idx%E) = (s_m*vel_r(dir_idx(1))*(e_r + pres_r - ptilde_r) - s_p*vel_l(dir_idx(1) &
1671 & )*(e_l + pres_l - ptilde_l) + s_m*s_p*(e_l - e_r))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r) &
1672 & *pcorr*(vel_r_rms - vel_l_rms)/2._wp
1673 else
1674 flux_rsx_vf(j, k, l, &
1675 & eqn_idx%E) = (s_m*vel_r(dir_idx(1))*(e_r + pres_r) - s_p*vel_l(dir_idx(1))*(e_l &
1676 & + pres_l) + s_m*s_p*(e_l - e_r))/(s_m - s_p) + (s_m/s_l)*(s_p/s_r)*pcorr*(vel_r_rms &
1677 & - vel_l_rms)/2._wp
1678 end if
1679
1680 ! Advection flux and source: interface velocity for volume fraction transport
1681
1682# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1683#if defined(MFC_OpenACC)
1684# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1685!$acc loop seq
1686# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1687#elif defined(MFC_OpenMP)
1688# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1689
1690# 307 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1691#endif
1692 do i = eqn_idx%adv%beg, eqn_idx%adv%end
1693 flux_rsx_vf(j, k, l, i) = (ql_prim_rsx_vf(j, k, l, i) - qr_prim_rsx_vf(j, k, l + 1, &
1694 & i))*s_m*s_p/(s_m - s_p)
1695 flux_src_rsx_vf(j, k, l, i) = (s_m*qr_prim_rsx_vf(j, k, l + 1, &
1696 & i) - s_p*ql_prim_rsx_vf(j, k, l, i))/(s_m - s_p)
1697 end do
1698
1699 if (bubbles_euler) then
1700 ! From HLLC: Kills mass transport @ bubble gas density
1701 if (num_fluids > 1) then
1702 flux_rsx_vf(j, k, l, eqn_idx%cont%end) = 0._wp
1703 end if
1704 end if
1705
1706 if (chemistry) then
1707
1708# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1709#if defined(MFC_OpenACC)
1710# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1711!$acc loop seq
1712# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1713#elif defined(MFC_OpenMP)
1714# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1715
1716# 323 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1717#endif
1718 do i = eqn_idx%species%beg, eqn_idx%species%end
1719 y_l = ql_prim_rsx_vf(j, k, l, i)
1720 y_r = qr_prim_rsx_vf(j, k, l + 1, i)
1721
1722 flux_rsx_vf(j, k, l, &
1723 & i) = (s_m*y_r*rho_r*vel_r(dir_idx(1)) - s_p*y_l*rho_l*vel_l(dir_idx(1)) &
1724 & + s_m*s_p*(y_l*rho_l - y_r*rho_r))/(s_m - s_p)
1725 flux_src_rsx_vf(j, k, l, i) = 0._wp
1726 end do
1727 end if
1728
1729# 352 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1730 end do
1731 end do
1732 end do
1733
1734# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1735#if defined(MFC_OpenACC)
1736# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1737!$acc end parallel loop
1738# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1739#elif defined(MFC_OpenMP)
1740# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1741
1742# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1743!$omp end target teams loop
1744# 355 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1745#endif
1746 end if
1747# 358 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1748
1749 if (viscous) then
1750
1751# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1752
1753# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1754#if defined(MFC_OpenACC)
1755# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1756!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, l, idx_right_phys, vel_grad_L, vel_grad_R, alpha_L, alpha_R, vel_L, vel_R, Re_L, Re_R) &
1757# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1758!$acc& firstprivate(Re_size_loc1, Re_size_loc2) copyin(norm_dir)
1759# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1760#elif defined(MFC_OpenMP)
1761# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1762
1763# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1764
1765# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1766
1767# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1768!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) &
1769# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1770!$omp& private(i, j, k, l, idx_right_phys, vel_grad_L, vel_grad_R, alpha_L, alpha_R, vel_L, vel_R, Re_L, Re_R) firstprivate(Re_size_loc1, Re_size_loc2) map(to:norm_dir)
1771# 360 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1772#endif
1773# 362 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1774 do l = isz%beg, isz%end
1775 do k = isy%beg, isy%end
1776 do j = isx%beg, isx%end
1777 idx_right_phys(1) = j
1778 idx_right_phys(2) = k
1779 idx_right_phys(3) = l
1780 idx_right_phys(norm_dir) = idx_right_phys(norm_dir) + 1
1781
1782 if (norm_dir == 1) then
1783
1784# 371 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1785#if defined(MFC_OpenACC)
1786# 371 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1787!$acc loop seq
1788# 371 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1789#elif defined(MFC_OpenMP)
1790# 371 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1791
1792# 371 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1793#endif
1794 do i = 1, num_fluids
1795 alpha_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%E + i)
1796 alpha_r(i) = qr_prim_rsx_vf(j + 1, k, l, eqn_idx%E + i)
1797 end do
1798
1799
1800# 377 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1801#if defined(MFC_OpenACC)
1802# 377 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1803!$acc loop seq
1804# 377 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1805#elif defined(MFC_OpenMP)
1806# 377 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1807
1808# 377 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1809#endif
1810 do i = 1, num_dims
1811 vel_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%mom%beg + i - 1)
1812 vel_r(i) = qr_prim_rsx_vf(j + 1, k, l, eqn_idx%mom%beg + i - 1)
1813 end do
1814 else if (norm_dir == 2) then
1815
1816# 383 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1817#if defined(MFC_OpenACC)
1818# 383 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1819!$acc loop seq
1820# 383 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1821#elif defined(MFC_OpenMP)
1822# 383 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1823
1824# 383 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1825#endif
1826 do i = 1, num_fluids
1827 alpha_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%E + i)
1828 alpha_r(i) = qr_prim_rsx_vf(j, k + 1, l, eqn_idx%E + i)
1829 end do
1830
1831# 388 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1832#if defined(MFC_OpenACC)
1833# 388 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1834!$acc loop seq
1835# 388 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1836#elif defined(MFC_OpenMP)
1837# 388 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1838
1839# 388 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1840#endif
1841 do i = 1, num_dims
1842 vel_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%mom%beg + i - 1)
1843 vel_r(i) = qr_prim_rsx_vf(j, k + 1, l, eqn_idx%mom%beg + i - 1)
1844 end do
1845 else
1846
1847# 394 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1848#if defined(MFC_OpenACC)
1849# 394 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1850!$acc loop seq
1851# 394 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1852#elif defined(MFC_OpenMP)
1853# 394 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1854
1855# 394 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1856#endif
1857 do i = 1, num_fluids
1858 alpha_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%E + i)
1859 alpha_r(i) = qr_prim_rsx_vf(j, k, l + 1, eqn_idx%E + i)
1860 end do
1861
1862
1863# 400 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1864#if defined(MFC_OpenACC)
1865# 400 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1866!$acc loop seq
1867# 400 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1868#elif defined(MFC_OpenMP)
1869# 400 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1870
1871# 400 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1872#endif
1873 do i = 1, num_dims
1874 vel_l(i) = ql_prim_rsx_vf(j, k, l, eqn_idx%mom%beg + i - 1)
1875 vel_r(i) = qr_prim_rsx_vf(j, k, l + 1, eqn_idx%mom%beg + i - 1)
1876 end do
1877 end if
1878
1879 call s_compute_interface_reynolds(alpha_l, re_l, re_size_loc1, re_size_loc2)
1880 call s_compute_interface_reynolds(alpha_r, re_r, re_size_loc1, re_size_loc2)
1881
1882 if (shear_stress) then
1883
1884# 411 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1885#if defined(MFC_OpenACC)
1886# 411 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1887!$acc loop seq
1888# 411 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1889#elif defined(MFC_OpenMP)
1890# 411 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1891
1892# 411 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1893#endif
1894 do i = 1, num_dims
1895 vel_grad_l(i, 1) = (dql_prim_dx_vf(eqn_idx%mom%beg + i - 1)%sf(j, k, l)/re_l(1))
1896 vel_grad_r(i, 1) = (dqr_prim_dx_vf(eqn_idx%mom%beg + i - 1)%sf(idx_right_phys(1), &
1897 & idx_right_phys(2), idx_right_phys(3))/re_r(1))
1898# 417 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1899 if (num_dims > 1) then
1900 vel_grad_l(i, 2) = (dql_prim_dy_vf(eqn_idx%mom%beg + i - 1)%sf(j, k, l)/re_l(1))
1901 vel_grad_r(i, 2) = (dqr_prim_dy_vf(eqn_idx%mom%beg + i - 1)%sf(idx_right_phys(1), &
1902 & idx_right_phys(2), idx_right_phys(3))/re_r(1))
1903 end if
1904# 423 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1905 if (num_dims > 2) then
1906 vel_grad_l(i, 3) = (dql_prim_dz_vf(eqn_idx%mom%beg + i - 1)%sf(j, k, l)/re_l(1))
1907 vel_grad_r(i, 3) = (dqr_prim_dz_vf(eqn_idx%mom%beg + i - 1)%sf(idx_right_phys(1), &
1908 & idx_right_phys(2), idx_right_phys(3))/re_r(1))
1909 end if
1910# 429 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1911# 430 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1912 end do
1913
1914 if (norm_dir == 1) then
1915 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
1916 & l) - (4._wp/3._wp)*0.5_wp*(vel_grad_l(1, 1) + vel_grad_r(1, 1))
1917 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1918 & l) - (4._wp/3._wp)*0.5_wp*(vel_grad_l(1, 1)*vel_l(1) + vel_grad_r(1, 1)*vel_r(1))
1919# 438 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1920 if (num_dims > 1) then
1921 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
1922 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(2, 2) + vel_grad_r(2, 2))
1923 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1924 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(2, 2)*vel_l(1) + vel_grad_r(2, &
1925 & 2)*vel_r(1))
1926
1927 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, &
1928 & l) - 0.5_wp*(vel_grad_l(1, 2) + vel_grad_r(1, 2)) - 0.5_wp*(vel_grad_l(2, &
1929 & 1) + vel_grad_r(2, 1))
1930 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1931 & l) - 0.5_wp*(vel_grad_l(1, 2)*vel_l(2) + vel_grad_r(1, &
1932 & 2)*vel_r(2)) - 0.5_wp*(vel_grad_l(2, 1)*vel_l(2) + vel_grad_r(2, 1)*vel_r(2))
1933# 452 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1934 if (num_dims > 2) then
1935 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
1936 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(3, 3) + vel_grad_r(3, 3))
1937 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1938 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(3, &
1939 & 3)*vel_l(1) + vel_grad_r(3, 3)*vel_r(1))
1940
1941 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
1942 & l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
1943 & l) - 0.5_wp*(vel_grad_l(1, 3) + vel_grad_r(1, &
1944 & 3)) - 0.5_wp*(vel_grad_l(3, 1) + vel_grad_r(3, 1))
1945 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1946 & l) - 0.5_wp*(vel_grad_l(1, 3)*vel_l(3) + vel_grad_r(1, &
1947 & 3)*vel_r(3)) - 0.5_wp*(vel_grad_l(3, 1)*vel_l(3) + vel_grad_r(3, &
1948 & 1)*vel_r(3))
1949 end if
1950# 469 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1951 end if
1952# 471 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1953 else if (norm_dir == 2) then
1954# 473 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1955 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, &
1956 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(1, 1) + vel_grad_r(1, 1))
1957 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1958 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(1, 1)*vel_l(2) + vel_grad_r(1, 1)*vel_r(2))
1959
1960 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, &
1961 & l) - (4._wp/3._wp)*0.5_wp*(vel_grad_l(2, 2) + vel_grad_r(2, 2))
1962 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1963 & l) - (4._wp/3._wp)*0.5_wp*(vel_grad_l(2, 2)*vel_l(2) + vel_grad_r(2, 2)*vel_r(2))
1964
1965 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
1966 & l) - 0.5_wp*(vel_grad_l(1, 2) + vel_grad_r(1, 2)) - 0.5_wp*(vel_grad_l(2, &
1967 & 1) + vel_grad_r(2, 1))
1968 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1969 & l) - 0.5_wp*(vel_grad_l(1, 2)*vel_l(1) + vel_grad_r(1, &
1970 & 2)*vel_r(1)) - 0.5_wp*(vel_grad_l(2, 1)*vel_l(1) + vel_grad_r(2, 1)*vel_r(1))
1971# 490 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1972 if (num_dims > 2) then
1973 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, &
1974 & k, l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(3, 3) + vel_grad_r(3, 3))
1975 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1976 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(3, 3)*vel_l(2) + vel_grad_r(3, &
1977 & 3)*vel_r(2))
1978
1979 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, &
1980 & k, l) - 0.5_wp*(vel_grad_l(2, 3) + vel_grad_r(2, &
1981 & 3)) - 0.5_wp*(vel_grad_l(3, 2) + vel_grad_r(3, 2))
1982 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1983 & l) - 0.5_wp*(vel_grad_l(2, 3)*vel_l(3) + vel_grad_r(2, &
1984 & 3)*vel_r(3)) - 0.5_wp*(vel_grad_l(3, 2)*vel_l(3) + vel_grad_r(3, &
1985 & 2)*vel_r(3))
1986 end if
1987# 506 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1988# 507 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1989 else
1990# 509 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
1991 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
1992 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(1, 1) + vel_grad_r(1, 1))
1993 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1994 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(1, 1)*vel_l(3) + vel_grad_r(1, 1)*vel_r(3))
1995
1996 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
1997 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(2, 2) + vel_grad_r(2, 2))
1998 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
1999 & l) - (-2._wp/3._wp)*0.5_wp*(vel_grad_l(2, 2)*vel_l(3) + vel_grad_r(2, 2)*vel_r(3))
2000
2001 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
2002 & l) - 0.5_wp*(vel_grad_l(1, 3) + vel_grad_r(1, 3)) - 0.5_wp*(vel_grad_l(3, &
2003 & 1) + vel_grad_r(3, 1))
2004 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2005 & l) - 0.5_wp*(vel_grad_l(1, 3)*vel_l(1) + vel_grad_r(1, &
2006 & 3)*vel_r(1)) - 0.5_wp*(vel_grad_l(3, 1)*vel_l(1) + vel_grad_r(3, 1)*vel_r(1))
2007
2008 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
2009 & l) - (4._wp/3._wp)*0.5_wp*(vel_grad_l(3, 3) + vel_grad_r(3, 3))
2010 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2011 & l) - (4._wp/3._wp)*0.5_wp*(vel_grad_l(3, 3)*vel_l(3) + vel_grad_r(3, 3)*vel_r(3))
2012
2013 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, &
2014 & l) - 0.5_wp*(vel_grad_l(2, 3) + vel_grad_r(2, 3)) - 0.5_wp*(vel_grad_l(3, &
2015 & 2) + vel_grad_r(3, 2))
2016 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2017 & l) - 0.5_wp*(vel_grad_l(2, 3)*vel_l(2) + vel_grad_r(2, &
2018 & 3)*vel_r(2)) - 0.5_wp*(vel_grad_l(3, 2)*vel_l(2) + vel_grad_r(3, 2)*vel_r(2))
2019# 538 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2020 end if
2021 end if
2022
2023 if (bulk_stress) then
2024
2025# 542 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2026#if defined(MFC_OpenACC)
2027# 542 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2028!$acc loop seq
2029# 542 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2030#elif defined(MFC_OpenMP)
2031# 542 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2032
2033# 542 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2034#endif
2035 do i = 1, num_dims
2036 vel_grad_l(i, 1) = (dql_prim_dx_vf(eqn_idx%mom%beg + i - 1)%sf(j, k, l)/re_l(2))
2037 vel_grad_r(i, 1) = (dqr_prim_dx_vf(eqn_idx%mom%beg + i - 1)%sf(idx_right_phys(1), &
2038 & idx_right_phys(2), idx_right_phys(3))/re_r(2))
2039# 548 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2040 if (num_dims > 1) then
2041 vel_grad_l(i, 2) = (dql_prim_dy_vf(eqn_idx%mom%beg + i - 1)%sf(j, k, l)/re_l(2))
2042 vel_grad_r(i, 2) = (dqr_prim_dy_vf(eqn_idx%mom%beg + i - 1)%sf(idx_right_phys(1), &
2043 & idx_right_phys(2), idx_right_phys(3))/re_r(2))
2044 end if
2045# 554 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2046# 555 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2047 if (num_dims > 2) then
2048 vel_grad_l(i, 3) = (dql_prim_dz_vf(eqn_idx%mom%beg + i - 1)%sf(j, k, l)/re_l(2))
2049 vel_grad_r(i, 3) = (dqr_prim_dz_vf(eqn_idx%mom%beg + i - 1)%sf(idx_right_phys(1), &
2050 & idx_right_phys(2), idx_right_phys(3))/re_r(2))
2051 end if
2052# 561 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2053 end do
2054
2055 if (norm_dir == 1) then
2056 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
2057 & l) - 0.5_wp*(vel_grad_l(1, 1) + vel_grad_r(1, 1))
2058 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, l) - 0.5_wp*(vel_grad_l(1, &
2059 & 1)*vel_l(1) + vel_grad_r(1, 1)*vel_r(1))
2060# 569 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2061 if (num_dims > 1) then
2062 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
2063 & l) - 0.5_wp*(vel_grad_l(2, 2) + vel_grad_r(2, 2))
2064 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2065 & l) - 0.5_wp*(vel_grad_l(2, 2)*vel_l(1) + vel_grad_r(2, 2)*vel_r(1))
2066
2067# 576 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2068 if (num_dims > 2) then
2069 flux_src_vf(eqn_idx%mom%beg)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg)%sf(j, k, &
2070 & l) - 0.5_wp*(vel_grad_l(3, 3) + vel_grad_r(3, 3))
2071 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2072 & l) - 0.5_wp*(vel_grad_l(3, 3)*vel_l(1) + vel_grad_r(3, 3)*vel_r(1))
2073 end if
2074# 583 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2075 end if
2076# 585 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2077 else if (norm_dir == 2) then
2078# 587 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2079 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, &
2080 & l) - 0.5_wp*(vel_grad_l(1, 1) + vel_grad_r(1, 1))
2081 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2082 & l) - 0.5_wp*(vel_grad_l(1, 1)*vel_l(2) + vel_grad_r(1, 1)*vel_r(2))
2083
2084 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, &
2085 & l) - 0.5_wp*(vel_grad_l(2, 2) + vel_grad_r(2, 2))
2086 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2087 & l) - 0.5_wp*(vel_grad_l(2, 2)*vel_l(2) + vel_grad_r(2, 2)*vel_r(2))
2088
2089# 598 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2090 if (num_dims > 2) then
2091 flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 1)%sf(j, &
2092 & k, l) - 0.5_wp*(vel_grad_l(3, 3) + vel_grad_r(3, 3))
2093 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2094 & l) - 0.5_wp*(vel_grad_l(3, 3)*vel_l(2) + vel_grad_r(3, 3)*vel_r(2))
2095 end if
2096# 605 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2097# 606 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2098 else
2099# 608 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2100 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
2101 & l) - 0.5_wp*(vel_grad_l(1, 1) + vel_grad_r(1, 1))
2102 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2103 & l) - 0.5_wp*(vel_grad_l(1, 1)*vel_l(3) + vel_grad_r(1, 1)*vel_r(3))
2104
2105 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
2106 & l) - 0.5_wp*(vel_grad_l(2, 2) + vel_grad_r(2, 2))
2107 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2108 & l) - 0.5_wp*(vel_grad_l(2, 2)*vel_l(3) + vel_grad_r(2, 2)*vel_r(3))
2109
2110 flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, l) = flux_src_vf(eqn_idx%mom%beg + 2)%sf(j, k, &
2111 & l) - 0.5_wp*(vel_grad_l(3, 3) + vel_grad_r(3, 3))
2112 flux_src_vf(eqn_idx%E)%sf(j, k, l) = flux_src_vf(eqn_idx%E)%sf(j, k, &
2113 & l) - 0.5_wp*(vel_grad_l(3, 3)*vel_l(3) + vel_grad_r(3, 3)*vel_r(3))
2114# 623 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2115 end if
2116 end if
2117 end do
2118 end do
2119 end do
2120
2121# 628 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2122#if defined(MFC_OpenACC)
2123# 628 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2124!$acc end parallel loop
2125# 628 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2126#elif defined(MFC_OpenMP)
2127# 628 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2128
2129# 628 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2130!$omp end target teams loop
2131# 628 "/home/runner/work/MFC/MFC/src/simulation/m_riemann_solver_lf.fpp"
2132#endif
2133 end if
2134
2135 call s_finalize_riemann_solver(flux_vf, flux_src_vf, flux_gsrc_vf, norm_dir)
2136
2137 end subroutine s_lf_riemann_solver
2138
2139end module m_riemann_solver_lf
Compile-time constant parameters: default values, tolerances, and physical constants.
integer, parameter riemann_solver_hll
real(wp), parameter sgm_eps
Segmentation tolerance.
integer, parameter riemann_solver_hllc
integer, parameter riemann_solver_lax_friedrichs
Shared derived types for field data, patch geometry, bubble dynamics, and MPI I/O structures.
Global parameters for the computational domain, fluid properties, and simulation algorithm configurat...
logical bulk_stress
Bulk stresses.
integer, dimension(3) dir_idx
real(wp), dimension(3) dir_flg
logical shear_stress
Shear stresses.
Lax-Friedrichs (Rusanov) approximate Riemann solver.
subroutine s_lf_riemann_solver(ql_prim_rsx_vf, dql_prim_dx_vf, dql_prim_dy_vf, dql_prim_dz_vf, ql_prim_vf, qr_prim_rsx_vf, dqr_prim_dx_vf, dqr_prim_dy_vf, dqr_prim_dz_vf, qr_prim_vf, q_prim_vf, flux_vf, flux_src_vf, flux_gsrc_vf, norm_dir, ix, iy, iz)
Lax-Friedrichs (Rusanov) approximate Riemann solver.
Shared Riemann-solver module state and the per-sweep setup, state-buffer population,...
real(wp), dimension(:,:,:,:), allocatable flux_src_rsx_vf
type(int_bounds_info) isz
type(int_bounds_info) isx
type(int_bounds_info) is3
real(wp), dimension(:,:,:,:), allocatable flux_rsx_vf
The cell-boundary values of the fluxes (src - source) that are computed through the chosen Riemann pr...
subroutine s_initialize_riemann_solver(flux_src_vf, norm_dir)
Set up the chosen Riemann solver algorithm for the current direction.
subroutine s_populate_riemann_states_variables_buffers(ql_prim_rsx_vf, dql_prim_dx_vf, dql_prim_dy_vf, dql_prim_dz_vf, qr_prim_rsx_vf, dqr_prim_dx_vf, dqr_prim_dy_vf, dqr_prim_dz_vf, norm_dir, ix, iy, iz)
Populate the left and right Riemann state variable buffers based on boundary conditions.
type(int_bounds_info) isy
type(int_bounds_info) is2
subroutine s_finalize_riemann_solver(flux_vf, flux_src_vf, flux_gsrc_vf, norm_dir)
Reshape and copy the Riemann-solver flux buffers back to the physical-space output arrays for the sel...
subroutine s_compute_interface_reynolds(alpha_k, re_k, re_size_loc1, re_size_loc2)
Compute the shear and volume Reynolds numbers of one Riemann state by inverse-weighting the fluid Rey...
real(wp), dimension(:,:,:,:), allocatable flux_gsrc_rsx_vf
The cell-boundary values of the geometrical source flux that are computed through the chosen Riemann ...
subroutine s_accumulate_mixture_properties(nf, alpha_rho_k, alpha_k, rho_k, gamma_k, pi_inf_k, qv_k)
Accumulate the mixture density, specific heat ratio function, liquid stiffness function,...
type(int_bounds_info) is1
Conservative-to-primitive variable conversion, mixture property evaluation, and pressure computation.
subroutine, public s_compute_speed_of_sound(pres, rho, gamma, pi_inf, h, adv, vel_sum, c_c, c, qv)
Compute the speed of sound from thermodynamic state variables, supporting multiple equation-of-state ...
Integer bounds for variables.
Derived type annexing a scalar field (SF).