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