MFC
Exascale flow solver
Loading...
Searching...
No Matches
m_bubbles_EL_kernels.fpp.f90
Go to the documentation of this file.
1# 1 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
2!>
3!! @file
4!! @brief Contains module @ref m_bubbles_el_kernels "m_bubbles_EL_kernels"
5
6# 1 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 1
7# 1 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 1
8# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
9# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
10# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
11# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
12# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
13# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
14
15# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
16# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
17# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
18
19# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
20# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
21# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
22
23# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
24
25# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
26
27# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
28
29# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
30
31# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
32
33# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
34
35# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
36
37# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
38! New line at end of file is required for FYPP
39# 2 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
40# 1 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 1
41# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
42# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
43# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
44# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
45# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
46# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
47
48# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
49# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
50# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
51
52# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
53# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
54# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
55
56# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
57
58# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
59
60# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
61
62# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
63
64# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
65
66# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
67
68# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
69
70# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
71! New line at end of file is required for FYPP
72# 2 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 2
73
74# 4 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
75# 5 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
76# 6 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
77# 7 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
78# 8 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
79
80# 20 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
81
82# 43 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
83
84# 48 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
85
86# 53 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
87
88# 58 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
89
90# 63 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
91
92# 68 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
93
94# 81 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
95
96# 86 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
97
98# 91 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
99
100# 96 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
101
102# 101 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
103
104# 106 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
105
106# 111 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
107
108# 116 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
109
110# 121 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
111
112# 126 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
113
114# 156 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
115
116# 197 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
117
118# 211 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
119
120# 236 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
121
122# 247 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
123
124# 249 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
125# 260 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
126
127# 310 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
128
129# 320 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
130
131# 330 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
132
133# 339 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
134
135# 356 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
136
137# 366 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
138
139# 373 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
140
141# 379 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
142
143# 385 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
144
145# 391 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
146
147# 397 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
148
149# 403 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
150! New line at end of file is required for FYPP
151# 3 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
152# 1 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 1
153# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
154# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
155# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
156# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
157# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
158# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
159
160# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
161# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
162# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
163
164# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
165# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
166# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
167
168# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
169
170# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
171
172# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
173
174# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
175
176# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
177
178# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
179
180# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
181
182# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
183! New line at end of file is required for FYPP
184# 2 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 2
185
186# 7 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
187
188# 17 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
189
190# 22 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
191
192# 27 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
193
194# 32 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
195
196# 37 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
197
198# 42 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
199
200# 47 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
201
202# 52 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
203
204# 57 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
205
206# 62 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
207
208# 73 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
209
210# 78 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
211
212# 83 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
213
214# 88 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
215
216# 103 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
217
218# 131 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
219
220# 160 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
221
222# 175 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
223
224# 193 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
225
226# 215 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
227
228# 244 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
229
230# 259 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
231
232# 269 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
233
234# 278 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
235
236# 294 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
237
238# 304 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
239
240# 311 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
241! New line at end of file is required for FYPP
242# 4 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
243
244! GPU parallel region (scalar reductions, maxval/minval)
245# 23 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
246
247! GPU parallel loop over threads (most common GPU macro)
248# 43 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
249
250! Required closing for GPU_PARALLEL_LOOP
251# 55 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
252
253! Mark routine for device compilation
254# 112 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
255
256! Declare device-resident data
257# 130 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
258
259! Inner loop within a GPU parallel region
260# 145 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
261
262! Scoped GPU data region
263# 164 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
264
265! Host code with device pointers (for MPI with GPU buffers)
266# 193 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
267
268! Allocate device memory (unscoped)
269# 207 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
270
271! Free device memory
272# 219 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
273
274! Atomic operation on device
275# 231 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
276
277! End atomic capture block
278# 242 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
279
280! Copy data between host and device
281# 254 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
282
283! Synchronization barrier
284# 266 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
285
286! Import GPU library module (openacc or omp_lib)
287# 275 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
288
289! Emit code only for AMD compiler
290# 282 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
291
292! Emit code for non-Cray compilers
293# 289 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
294
295! Emit code only for Cray compiler
296# 296 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
297
298! Emit code for non-NVIDIA compilers
299# 303 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
300
301# 305 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
302# 306 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
303! New line at end of file is required for FYPP
304# 2 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 2
305
306# 14 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
307
308! Caution: This macro requires the use of a binding script to set CUDA_VISIBLE_DEVICES, such that we have one GPU device per MPI
309! rank. That's because for both cudaMemAdvise (preferred location) and cudaMemPrefetchAsync we use location = device_id = 0. For an
310! example see misc/nvidia_uvm/bind.sh.
311# 52 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
312
313! Allocate and create GPU device memory
314# 72 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
315
316! Free GPU device memory and deallocate
317# 80 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
318
319! Cray-specific GPU pointer setup for vector fields
320# 104 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
321
322! Cray-specific GPU pointer setup for scalar fields
323# 120 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
324
325! Cray-specific GPU pointer setup for acoustic source spatials
326# 145 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
327
328# 151 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
329
330# 158 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
331! New line at end of file is required for FYPP
332# 6 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp" 2
333
334!> @brief Kernel functions (Gaussian, delta) that smear Lagrangian bubble effects onto the Eulerian grid
336
337 use m_mpi_proxy
338
339 implicit none
340
341 ! Cell-centered pressure gradients (precomputed for translational motion)
342 real(wp), allocatable, dimension(:,:,:) :: grad_p_x, grad_p_y, grad_p_z
343
344# 16 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
345#if defined(MFC_OpenACC)
346# 16 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
347!$acc declare create(grad_p_x, grad_p_y, grad_p_z)
348# 16 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
349#elif defined(MFC_OpenMP)
350# 16 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
351!$omp declare target (grad_p_x, grad_p_y, grad_p_z)
352# 16 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
353#endif
354
355 ! Finite-difference coefficients for pressure gradient computation
356 real(wp), allocatable, dimension(:,:) :: fd_coeff_x_pgrad
357 real(wp), allocatable, dimension(:,:) :: fd_coeff_y_pgrad
358 real(wp), allocatable, dimension(:,:) :: fd_coeff_z_pgrad
359
360# 22 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
361#if defined(MFC_OpenACC)
362# 22 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
363!$acc declare create(fd_coeff_x_pgrad, fd_coeff_y_pgrad, fd_coeff_z_pgrad)
364# 22 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
365#elif defined(MFC_OpenMP)
366# 22 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
367!$omp declare target (fd_coeff_x_pgrad, fd_coeff_y_pgrad, fd_coeff_z_pgrad)
368# 22 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
369#endif
370
371 ! Cell list for bubble-to-cell mapping (rebuilt each RK stage before smearing)
372 integer, allocatable, dimension(:,:,:) :: cell_list_start ! (0:m, 0:n, 0:p)
373 integer, allocatable, dimension(:,:,:) :: cell_list_count ! (0:m, 0:n, 0:p)
374 integer, allocatable, dimension(:) :: cell_list_idx ! (1:nBubs_glb) sorted bubble indices
375
376# 28 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
377#if defined(MFC_OpenACC)
378# 28 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
379!$acc declare create(cell_list_start, cell_list_count, cell_list_idx)
380# 28 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
381#elif defined(MFC_OpenMP)
382# 28 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
383!$omp declare target (cell_list_start, cell_list_count, cell_list_idx)
384# 28 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
385#endif
386
387contains
388
389 !> Smear the Lagrangian bubble effects onto the Eulerian grid using the selected kernel
390 subroutine s_smoothfunction(nBubs, lbk_rad, lbk_vel, lbk_s, lbk_pos, updatedvar, kcomp)
391
392 integer, intent(in) :: nBubs
393 real(wp), dimension(1:lag_params%nBubs_glb,1:3,1:2), intent(in) :: lbk_s, lbk_pos
394 real(wp), dimension(1:lag_params%nBubs_glb,1:2), intent(in) :: lbk_rad, lbk_vel
395 type(scalar_field), dimension(:), intent(inout) :: updatedvar
396 type(scalar_field), dimension(:), intent(inout) :: kcomp
397
398 smoothfunc:select case(lag_params%smooth_type)
399 case (1)
400 call s_gaussian(nbubs, lbk_rad, lbk_vel, lbk_s, lbk_pos, updatedvar, kcomp)
401 case (2)
402 call s_deltafunc(nbubs, lbk_rad, lbk_vel, lbk_s, updatedvar, kcomp)
403 end select smoothfunc
404
405 end subroutine s_smoothfunction
406
407 !> Builds a sorted cell list mapping each interior cell (0:m,0:n,0:p) to its resident bubbles. Uses a counting-sort on the host
408 !! (O(nBubs + N_cells)). Must be called before s_gaussian each RK stage.
409 subroutine s_build_cell_list(nBubs, lbk_s)
410
411 integer, intent(in) :: nBubs
412 real(wp), dimension(1:lag_params%nBubs_glb,1:3,1:2), intent(in) :: lbk_s
413 integer :: l, ci, cj, ck, idx
414 real(wp), dimension(3) :: s_coord
415
416 ! Bring current bubble positions to host
417
418
419# 61 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
420#if defined(MFC_OpenACC)
421# 61 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
422!$acc update host(lbk_s)
423# 61 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
424#elif defined(MFC_OpenMP)
425# 61 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
426!$omp target update from(lbk_s)
427# 61 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
428#endif
429
430 ! Pass 1: zero counts and count bubbles per cell
432 do l = 1, nbubs
433 s_coord(1:3) = lbk_s(l,1:3,2)
434 ci = int(s_coord(1))
435 cj = int(s_coord(2))
436 ck = int(s_coord(3))
437 ! Clamp to interior (bubbles should already be in [0:m,0:n,0:p])
438 ci = max(0, min(ci, m))
439 cj = max(0, min(cj, n))
440 ck = max(0, min(ck, p))
441 cell_list_count(ci, cj, ck) = cell_list_count(ci, cj, ck) + 1
442 end do
443
444 ! Prefix sum to compute start indices (1-based into cell_list_idx)
445 idx = 1
446 do ck = 0, p
447 do cj = 0, n
448 do ci = 0, m
449 cell_list_start(ci, cj, ck) = idx
450 idx = idx + cell_list_count(ci, cj, ck)
451 end do
452 end do
453 end do
454
455 ! Pass 2: place bubble indices into cell_list_idx Temporarily reuse cell_list_count as a running offset
457 do l = 1, nbubs
458 s_coord(1:3) = lbk_s(l,1:3,2)
459 ci = int(s_coord(1))
460 cj = int(s_coord(2))
461 ck = int(s_coord(3))
462 ci = max(0, min(ci, m))
463 cj = max(0, min(cj, n))
464 ck = max(0, min(ck, p))
465 cell_list_idx(cell_list_start(ci, cj, ck) + cell_list_count(ci, cj, ck)) = l
466 cell_list_count(ci, cj, ck) = cell_list_count(ci, cj, ck) + 1
467 end do
468
469 ! Send cell list arrays to GPU
470
471# 103 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
472#if defined(MFC_OpenACC)
473# 103 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
474!$acc update device(cell_list_start, cell_list_count, cell_list_idx)
475# 103 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
476#elif defined(MFC_OpenMP)
477# 103 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
478!$omp target update to(cell_list_start, cell_list_count, cell_list_idx)
479# 103 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
480#endif
481
482 end subroutine s_build_cell_list
483
484 !> Cell-centric delta-function smearing using the cell list (no GPU atomics). Each bubble only affects the cell it resides in.
485 !! The outer GPU loop iterates over interior cells and sums contributions from resident bubbles.
486 subroutine s_deltafunc(nBubs, lbk_rad, lbk_vel, lbk_s, updatedvar, kcomp)
487
488 integer, intent(in) :: nBubs
489 real(wp), dimension(1:lag_params%nBubs_glb,1:3,1:2), intent(in) :: lbk_s
490 real(wp), dimension(1:lag_params%nBubs_glb,1:2), intent(in) :: lbk_rad, lbk_vel
491 type(scalar_field), dimension(:), intent(inout) :: updatedvar
492 type(scalar_field), dimension(:), intent(inout) :: kcomp
493 real(wp) :: strength_vel, strength_vol
494 real(wp) :: volpart, Vol
495 real(wp) :: y_kahan, t_kahan
496 integer :: i, j, k, lb, bub_idx
497
498
499# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
500
501# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
502#if defined(MFC_OpenACC)
503# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
504!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, lb, bub_idx, volpart, Vol, strength_vel, strength_vol, y_kahan, t_kahan)
505# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
506#elif defined(MFC_OpenMP)
507# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
508
509# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
510
511# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
512
513# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
514!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) &
515# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
516!$omp& private(i, j, k, lb, bub_idx, volpart, Vol, strength_vel, strength_vol, y_kahan, t_kahan)
517# 121 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
518#endif
519# 123 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
520 do k = 0, p
521 do j = 0, n
522 do i = 0, m
523 ! Cell volume
524 if (num_dims == 2) then
525 vol = dx(i)*dy(j)*lag_params%charwidth
526 if (cyl_coord) vol = dx(i)*dy(j)*y_cc(j)*2._wp*pi
527 else
528 vol = dx(i)*dy(j)*dz(k)
529 end if
530
531 ! Loop over bubbles in this cell
532
533# 135 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
534#if defined(MFC_OpenACC)
535# 135 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
536!$acc loop seq
537# 135 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
538#elif defined(MFC_OpenMP)
539# 135 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
540
541# 135 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
542#endif
543 do lb = cell_list_start(i, j, k), cell_list_start(i, j, k) + cell_list_count(i, j, k) - 1
544 bub_idx = cell_list_idx(lb)
545
546 volpart = 4._wp/3._wp*pi*lbk_rad(bub_idx, 2)**3._wp
547 strength_vol = volpart
548 strength_vel = 4._wp*pi*lbk_rad(bub_idx, 2)**2._wp*lbk_vel(bub_idx, 2)
549
550 if (lag_params%kahan_summation) then
551 ! Kahan summation for void fraction
552 y_kahan = real(strength_vol/vol, kind=wp) - kcomp(1)%sf(i, j, k)
553 t_kahan = updatedvar(1)%sf(i, j, k) + y_kahan
554 kcomp(1)%sf(i, j, k) = (t_kahan - updatedvar(1)%sf(i, j, k)) - y_kahan
555 updatedvar(1)%sf(i, j, k) = t_kahan
556
557 ! Kahan summation for time derivative of void fraction
558 y_kahan = real(strength_vel/vol, kind=wp) - kcomp(2)%sf(i, j, k)
559 t_kahan = updatedvar(2)%sf(i, j, k) + y_kahan
560 kcomp(2)%sf(i, j, k) = (t_kahan - updatedvar(2)%sf(i, j, k)) - y_kahan
561 updatedvar(2)%sf(i, j, k) = t_kahan
562 else
563 updatedvar(1)%sf(i, j, k) = updatedvar(1)%sf(i, j, k) + real(strength_vol/vol, kind=wp)
564 updatedvar(2)%sf(i, j, k) = updatedvar(2)%sf(i, j, k) + real(strength_vel/vol, kind=wp)
565 end if
566
567 ! Product of two smeared functions
568 if (lag_params%kahan_summation .and. lag_params%cluster_type >= 4) then
569 y_kahan = real((strength_vol*strength_vel)/vol, kind=wp) - kcomp(5)%sf(i, j, k)
570 t_kahan = updatedvar(5)%sf(i, j, k) + y_kahan
571 kcomp(5)%sf(i, j, k) = (t_kahan - updatedvar(5)%sf(i, j, k)) - y_kahan
572 updatedvar(5)%sf(i, j, k) = t_kahan
573 else if (lag_params%cluster_type >= 4) then
574 updatedvar(5)%sf(i, j, k) = updatedvar(5)%sf(i, j, k) + real((strength_vol*strength_vel)/vol, kind=wp)
575 end if
576 end do
577 end do
578 end do
579 end do
580
581# 173 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
582#if defined(MFC_OpenACC)
583# 173 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
584!$acc end parallel loop
585# 173 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
586#elif defined(MFC_OpenMP)
587# 173 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
588
589# 173 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
590!$omp end target teams loop
591# 173 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
592#endif
593
594 end subroutine s_deltafunc
595
596 !> Cell-centric gaussian smearing using the cell list (no GPU atomics). Each grid cell accumulates contributions from nearby
597 !! bubbles looked up via cell_list_start/count/idx.
598 subroutine s_gaussian(nBubs, lbk_rad, lbk_vel, lbk_s, lbk_pos, updatedvar, kcomp)
599
600 integer, intent(in) :: nBubs
601 real(wp), dimension(1:lag_params%nBubs_glb,1:3,1:2), intent(in) :: lbk_s, lbk_pos
602 real(wp), dimension(1:lag_params%nBubs_glb,1:2), intent(in) :: lbk_rad, lbk_vel
603 type(scalar_field), dimension(:), intent(inout) :: updatedvar
604 type(scalar_field), dimension(:), intent(inout) :: kcomp
605 real(wp), dimension(3) :: center, nodecoord, s_coord
606 integer, dimension(3) :: cell, cellijk
607 real(wp) :: stddsv, volpart
608 real(wp) :: strength_vel, strength_vol
609 real(wp) :: func, func2
610 real(wp) :: y_kahan, t_kahan
611 integer :: i, j, k, di, dj, dk, lb, bub_idx
612 integer :: di_beg, di_end, dj_beg, dj_end, dk_beg, dk_end
613 integer :: smear_x_beg, smear_x_end
614 integer :: smear_y_beg, smear_y_end
615 integer :: smear_z_beg, smear_z_end
616
617 ! Extended grid range for smearing (includes buffer cells for MPI communication)
618
619 smear_x_beg = -mapcells - 1
620 smear_x_end = m + mapcells + 1
621 smear_y_beg = merge(-mapcells - 1, 0, n > 0)
622 smear_y_end = merge(n + mapcells + 1, n, n > 0)
623 smear_z_beg = merge(-mapcells - 1, 0, p > 0)
624 smear_z_end = merge(p + mapcells + 1, p, p > 0)
625
626
627# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
628
629# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
630#if defined(MFC_OpenACC)
631# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
632!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, di, dj, dk, lb, bub_idx, center, nodecoord, s_coord, cell, cellijk, stddsv, volpart, strength_vel, strength_vol, func, &
633# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
634!$acc& func2, y_kahan, t_kahan, di_beg, di_end, dj_beg, dj_end, dk_beg, dk_end) copyin(smear_x_beg, smear_x_end, smear_y_beg, smear_y_end, smear_z_beg, smear_z_end)
635# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
636#elif defined(MFC_OpenMP)
637# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
638
639# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
640
641# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
642
643# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
644!$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, di, dj, &
645# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
646!$omp& dk, lb, bub_idx, center, nodecoord, s_coord, cell, cellijk, stddsv, volpart, strength_vel, strength_vol, func, func2, y_kahan, t_kahan, di_beg, di_end, dj_beg, dj_end, dk_beg, dk_end) &
647# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
648!$omp& map(to:smear_x_beg, smear_x_end, smear_y_beg, smear_y_end, smear_z_beg, smear_z_end)
649# 207 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
650#endif
651# 210 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
652 do k = smear_z_beg, smear_z_end
653 do j = smear_y_beg, smear_y_end
654 do i = smear_x_beg, smear_x_end
655 cellijk(1) = i
656 cellijk(2) = j
657 cellijk(3) = k
658
659 nodecoord(1) = x_cc(i)
660 nodecoord(2) = y_cc(j)
661 nodecoord(3) = 0._wp
662 if (p > 0) nodecoord(3) = z_cc(k)
663
664 ! Neighbor cell range clamped to interior [0:m, 0:n, 0:p]
665 di_beg = max(i - mapcells, 0)
666 di_end = min(i + mapcells, m)
667 dj_beg = max(j - mapcells, 0)
668 dj_end = min(j + mapcells, n)
669 dk_beg = max(k - mapcells, 0)
670 dk_end = min(k + mapcells, p)
671
672
673# 230 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
674#if defined(MFC_OpenACC)
675# 230 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
676!$acc loop seq
677# 230 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
678#elif defined(MFC_OpenMP)
679# 230 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
680
681# 230 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
682#endif
683 do dk = dk_beg, dk_end
684
685# 232 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
686#if defined(MFC_OpenACC)
687# 232 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
688!$acc loop seq
689# 232 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
690#elif defined(MFC_OpenMP)
691# 232 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
692
693# 232 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
694#endif
695 do dj = dj_beg, dj_end
696
697# 234 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
698#if defined(MFC_OpenACC)
699# 234 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
700!$acc loop seq
701# 234 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
702#elif defined(MFC_OpenMP)
703# 234 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
704
705# 234 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
706#endif
707 do di = di_beg, di_end
708
709# 236 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
710#if defined(MFC_OpenACC)
711# 236 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
712!$acc loop seq
713# 236 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
714#elif defined(MFC_OpenMP)
715# 236 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
716
717# 236 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
718#endif
719 do lb = cell_list_start(di, dj, dk), cell_list_start(di, dj, dk) + cell_list_count(di, dj, dk) - 1
720 bub_idx = cell_list_idx(lb)
721
722 ! Bubble properties
723 volpart = 4._wp/3._wp*pi*lbk_rad(bub_idx, 2)**3._wp
724 s_coord(1:3) = lbk_s(bub_idx,1:3,2)
725 call s_get_cell(s_coord, cell)
726 call s_compute_stddsv(cell, volpart, stddsv)
727
728 strength_vol = volpart
729 strength_vel = 4._wp*pi*lbk_rad(bub_idx, 2)**2._wp*lbk_vel(bub_idx, 2)
730
731 center(1:2) = lbk_pos(bub_idx,1:2,2)
732 center(3) = 0._wp
733 if (p > 0) center(3) = lbk_pos(bub_idx, 3, 2)
734
735 call s_applygaussian(center, cellijk, nodecoord, stddsv, 0._wp, func)
736
737 ! Kahan summation for void fraction
738 y_kahan = real(func*strength_vol, kind=wp) - kcomp(1)%sf(i, j, k)
739 t_kahan = updatedvar(1)%sf(i, j, k) + y_kahan
740 kcomp(1)%sf(i, j, k) = (t_kahan - updatedvar(1)%sf(i, j, k)) - y_kahan
741 updatedvar(1)%sf(i, j, k) = t_kahan
742
743 ! Kahan summation for time derivative of void fraction
744 y_kahan = real(func*strength_vel, kind=wp) - kcomp(2)%sf(i, j, k)
745 t_kahan = updatedvar(2)%sf(i, j, k) + y_kahan
746 kcomp(2)%sf(i, j, k) = (t_kahan - updatedvar(2)%sf(i, j, k)) - y_kahan
747 updatedvar(2)%sf(i, j, k) = t_kahan
748
749 if (lag_params%cluster_type >= 4) then
750 call s_applygaussian(center, cellijk, nodecoord, stddsv, 1._wp, func2)
751 y_kahan = real(func2*strength_vol*strength_vel, kind=wp) - kcomp(5)%sf(i, j, k)
752 t_kahan = updatedvar(5)%sf(i, j, k) + y_kahan
753 kcomp(5)%sf(i, j, k) = (t_kahan - updatedvar(5)%sf(i, j, k)) - y_kahan
754 updatedvar(5)%sf(i, j, k) = t_kahan
755 end if
756 end do
757 end do
758 end do
759 end do
760 end do
761 end do
762 end do
763
764# 281 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
765#if defined(MFC_OpenACC)
766# 281 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
767!$acc end parallel loop
768# 281 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
769#elif defined(MFC_OpenMP)
770# 281 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
771
772# 281 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
773!$omp end target teams loop
774# 281 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
775#endif
776
777 end subroutine s_gaussian
778
779 !> Evaluate the Gaussian kernel at a grid node for a given bubble center
780 subroutine s_applygaussian(center, cellaux, nodecoord, stddsv, strength_idx, func)
781
782
783# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
784#ifdef _CRAYFTN
785# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
786#if MFC_OpenACC
787# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
788!$acc routine seq
789# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
790#elif MFC_OpenMP
791# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
792
793# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
794
795# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
796!$omp declare target device_type(any)
797# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
798#else
799# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
800!DIR$ INLINEALWAYS s_applygaussian
801# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
802#endif
803# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
804#elif MFC_OpenACC
805# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
806!$acc routine seq
807# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
808#elif MFC_OpenMP
809# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
810
811# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
812
813# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
814!$omp declare target device_type(any)
815# 288 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
816#endif
817
818 real(wp), dimension(3), intent(in) :: center
819 integer, dimension(3), intent(in) :: cellaux
820 real(wp), dimension(3), intent(in) :: nodecoord
821 real(wp), intent(in) :: stddsv
822 real(wp), intent(in) :: strength_idx
823 real(wp), intent(out) :: func
824 integer :: i
825 real(wp) :: distance
826 real(wp) :: theta, dtheta, L2, dzp, Lz2, zc
827 real(wp) :: Nr, Nr_count
828
829 distance = sqrt((center(1) - nodecoord(1))**2._wp + (center(2) - nodecoord(2))**2._wp + (center(3) - nodecoord(3))**2._wp)
830
831 if (num_dims == 3) then
832 !> 3D gaussian function
833 func = exp(-0.5_wp*(distance/stddsv)**2._wp)/(sqrt(2._wp*pi)*stddsv)**3._wp
834 else
835 if (cyl_coord) then
836 !> 2D cylindrical function:
837 ! We smear particles in the azimuthal direction for given r
838 theta = 0._wp
839 nr = ceiling(2._wp*pi*nodecoord(2)/(y_cb(cellaux(2)) - y_cb(cellaux(2) - 1)))
840 dtheta = 2._wp*pi/nr
841 l2 = center(2)**2._wp + nodecoord(2)**2._wp - 2._wp*center(2)*nodecoord(2)*cos(theta)
842 distance = sqrt((center(1) - nodecoord(1))**2._wp + l2)
843 ! Factor 2._wp is for symmetry (upper half of the 2D field (+r) is considered)
844 func = dtheta/2._wp/pi*exp(-0.5_wp*(distance/stddsv)**2._wp)/(sqrt(2._wp*pi)*stddsv)**3._wp
845 nr_count = 0._wp
846 do while (nr_count < nr - 1._wp)
847 nr_count = nr_count + 1._wp
848 theta = nr_count*dtheta
849 ! trigonometric relation
850 l2 = center(2)**2._wp + nodecoord(2)**2._wp - 2._wp*center(2)*nodecoord(2)*cos(theta)
851 distance = sqrt((center(1) - nodecoord(1))**2._wp + l2)
852 ! nodecoord(2)*dtheta is the azimuthal width of the cell
853 func = func + dtheta/2._wp/pi*exp(-0.5_wp*(distance/stddsv)**2._wp)/(sqrt(2._wp*pi)*stddsv) &
854 & **(3._wp*(strength_idx + 1._wp))
855 end do
856 else
857 !> 2D cartesian function: Equation (48) from Maeda and Colonius 2018
858 ! We smear particles considering a virtual depth (lag_params%charwidth) with lag_params%charNz cells
859 dzp = (lag_params%charwidth/(lag_params%charNz + 1._wp))
860
861 func = 0._wp
862 do i = 0, lag_params%charNz
863 zc = (-lag_params%charwidth/2._wp + dzp*(0.5_wp + i)) ! Center of virtual cell i in z-direction
864 lz2 = (center(3) - zc)**2._wp
865 distance = sqrt((center(1) - nodecoord(1))**2._wp + (center(2) - nodecoord(2))**2._wp + lz2)
866 func = func + dzp/lag_params%charwidth*exp(-0.5_wp*(distance/stddsv)**2._wp)/(sqrt(2._wp*pi)*stddsv)**3._wp
867 end do
868 end if
869 end if
870
871 end subroutine s_applygaussian
872
873 !> Check if the current cell is outside the computational domain including ghost cells
874 subroutine s_check_celloutside(cellaux, celloutside)
875
876
877# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
878#ifdef _CRAYFTN
879# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
880#if MFC_OpenACC
881# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
882!$acc routine seq
883# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
884#elif MFC_OpenMP
885# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
886
887# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
888
889# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
890!$omp declare target device_type(any)
891# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
892#else
893# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
894!DIR$ INLINEALWAYS s_check_celloutside
895# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
896#endif
897# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
898#elif MFC_OpenACC
899# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
900!$acc routine seq
901# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
902#elif MFC_OpenMP
903# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
904
905# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
906
907# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
908!$omp declare target device_type(any)
909# 348 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
910#endif
911
912 integer, dimension(3), intent(inout) :: cellaux
913 logical, intent(out) :: celloutside
914
915 celloutside = .false.
916
917 if (num_dims == 2) then
918 if ((cellaux(1) < -buff_size) .or. (cellaux(2) < -buff_size)) then
919 celloutside = .true.
920 end if
921 if (cyl_coord .and. y_cc(cellaux(2)) < 0._wp) then
922 celloutside = .true.
923 end if
924 if ((cellaux(2) > n + buff_size) .or. (cellaux(1) > m + buff_size)) then
925 celloutside = .true.
926 end if
927 else
928 if ((cellaux(3) < -buff_size) .or. (cellaux(1) < -buff_size) .or. (cellaux(2) < -buff_size)) then
929 celloutside = .true.
930 end if
931
932 if ((cellaux(3) > p + buff_size) .or. (cellaux(2) > n + buff_size) .or. (cellaux(1) > m + buff_size)) then
933 celloutside = .true.
934 end if
935 end if
936
937 end subroutine s_check_celloutside
938
939 !> Relocate cells that intersect a symmetric boundary
940 subroutine s_shift_cell_symmetric_bc(cellaux, cell)
941
942
943# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
944#ifdef _CRAYFTN
945# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
946#if MFC_OpenACC
947# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
948!$acc routine seq
949# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
950#elif MFC_OpenMP
951# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
952
953# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
954
955# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
956!$omp declare target device_type(any)
957# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
958#else
959# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
960!DIR$ INLINEALWAYS s_shift_cell_symmetric_bc
961# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
962#endif
963# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
964#elif MFC_OpenACC
965# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
966!$acc routine seq
967# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
968#elif MFC_OpenMP
969# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
970
971# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
972
973# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
974!$omp declare target device_type(any)
975# 380 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
976#endif
977
978 integer, dimension(3), intent(inout) :: cellaux
979 integer, dimension(3), intent(in) :: cell
980
981 ! x-dir
982 if (bc_x%beg == bc_reflective .and. (cell(1) <= mapcells - 1)) then
983 cellaux(1) = abs(cellaux(1)) - 1
984 end if
985 if (bc_x%end == bc_reflective .and. (cell(1) >= m + 1 - mapcells)) then
986 cellaux(1) = cellaux(1) - (2*(cellaux(1) - m) - 1)
987 end if
988
989 ! y-dir
990 if (bc_y%beg == bc_reflective .and. (cell(2) <= mapcells - 1)) then
991 cellaux(2) = abs(cellaux(2)) - 1
992 end if
993 if (bc_y%end == bc_reflective .and. (cell(2) >= n + 1 - mapcells)) then
994 cellaux(2) = cellaux(2) - (2*(cellaux(2) - n) - 1)
995 end if
996
997 if (p > 0) then
998 ! z-dir
999 if (bc_z%beg == bc_reflective .and. (cell(3) <= mapcells - 1)) then
1000 cellaux(3) = abs(cellaux(3)) - 1
1001 end if
1002 if (bc_z%end == bc_reflective .and. (cell(3) >= p + 1 - mapcells)) then
1003 cellaux(3) = cellaux(3) - (2*(cellaux(3) - p) - 1)
1004 end if
1005 end if
1006
1007 end subroutine s_shift_cell_symmetric_bc
1008
1009 !> Calculates the standard deviation of the bubble being smeared in the Eulerian framework.
1010 subroutine s_compute_stddsv(cell, volpart, stddsv)
1011
1012
1013# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1014#ifdef _CRAYFTN
1015# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1016#if MFC_OpenACC
1017# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1018!$acc routine seq
1019# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1020#elif MFC_OpenMP
1021# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1022
1023# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1024
1025# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1026!$omp declare target device_type(any)
1027# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1028#else
1029# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1030!DIR$ INLINEALWAYS s_compute_stddsv
1031# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1032#endif
1033# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1034#elif MFC_OpenACC
1035# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1036!$acc routine seq
1037# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1038#elif MFC_OpenMP
1039# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1040
1041# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1042
1043# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1044!$omp declare target device_type(any)
1045# 416 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1046#endif
1047
1048 integer, dimension(3), intent(in) :: cell
1049 real(wp), intent(in) :: volpart
1050 real(wp), intent(out) :: stddsv
1051 real(wp) :: chardist, charvol
1052 real(wp) :: rad
1053
1054 !> Compute characteristic distance
1055 chardist = sqrt(dx(cell(1))*dy(cell(2)))
1056 if (p > 0) chardist = (dx(cell(1))*dy(cell(2))*dz(cell(3)))**(1._wp/3._wp)
1057
1058 !> Compute characteristic volume
1059 if (p > 0) then
1060 charvol = dx(cell(1))*dy(cell(2))*dz(cell(3))
1061 else
1062 if (cyl_coord) then
1063 charvol = dx(cell(1))*dy(cell(2))*y_cc(cell(2))*2._wp*pi
1064 else
1065 charvol = dx(cell(1))*dy(cell(2))*lag_params%charwidth
1066 end if
1067 end if
1068
1069 !> Compute Standard deviaton
1070 if ((volpart/charvol) > 0.5_wp*lag_params%valmaxvoid .or. (lag_params%smooth_type == 1)) then
1071 rad = (3._wp*volpart/(4._wp*pi))**(1._wp/3._wp)
1072 stddsv = 1._wp*lag_params%epsilonb*max(chardist, rad)
1073 else
1074 stddsv = 0._wp
1075 end if
1076
1077 end subroutine s_compute_stddsv
1078
1079 !> Compute the characteristic cell volume
1080 subroutine s_get_char_vol(cellx, celly, cellz, Charvol)
1081
1082
1083# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1084#ifdef _CRAYFTN
1085# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1086#if MFC_OpenACC
1087# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1088!$acc routine seq
1089# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1090#elif MFC_OpenMP
1091# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1092
1093# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1094
1095# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1096!$omp declare target device_type(any)
1097# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1098#else
1099# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1100!DIR$ INLINEALWAYS s_get_char_vol
1101# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1102#endif
1103# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1104#elif MFC_OpenACC
1105# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1106!$acc routine seq
1107# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1108#elif MFC_OpenMP
1109# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1110
1111# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1112
1113# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1114!$omp declare target device_type(any)
1115# 452 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1116#endif
1117
1118 integer, intent(in) :: cellx, celly, cellz
1119 real(wp), intent(out) :: Charvol
1120
1121 if (p > 0) then
1122 charvol = dx(cellx)*dy(celly)*dz(cellz)
1123 else
1124 if (cyl_coord) then
1125 charvol = dx(cellx)*dy(celly)*y_cc(celly)*2._wp*pi
1126 else
1127 charvol = dx(cellx)*dy(celly)*lag_params%charwidth
1128 end if
1129 end if
1130
1131 end subroutine s_get_char_vol
1132
1133 !> Convert bubble computational coordinates from real to integer cell indices
1134 subroutine s_get_cell(s_cell, get_cell)
1135
1136
1137# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1138#ifdef _CRAYFTN
1139# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1140#if MFC_OpenACC
1141# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1142!$acc routine seq
1143# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1144#elif MFC_OpenMP
1145# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1146
1147# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1148
1149# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1150!$omp declare target device_type(any)
1151# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1152#else
1153# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1154!DIR$ INLINEALWAYS s_get_cell
1155# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1156#endif
1157# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1158#elif MFC_OpenACC
1159# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1160!$acc routine seq
1161# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1162#elif MFC_OpenMP
1163# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1164
1165# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1166
1167# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1168!$omp declare target device_type(any)
1169# 472 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1170#endif
1171
1172 real(wp), dimension(3), intent(in) :: s_cell
1173 integer, dimension(3), intent(out) :: get_cell
1174 integer :: i
1175
1176 get_cell(:) = int(s_cell(:))
1177 do i = 1, num_dims
1178 if (s_cell(i) < 0._wp) get_cell(i) = get_cell(i) - 1
1179 end do
1180
1181 end subroutine s_get_cell
1182
1183 !> Precompute cell-centered pressure gradients (dp/dx, dp/dy, dp/dz)
1184 subroutine s_compute_pressure_gradients(q_prim_vf)
1185
1186 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf
1187 integer :: i, j, k, r
1188
1189 ! dp/dx at all cell centers
1190
1191
1192# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1193
1194# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1195#if defined(MFC_OpenACC)
1196# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1197!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, r)
1198# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1199#elif defined(MFC_OpenMP)
1200# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1201
1202# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1203
1204# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1205
1206# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1207!$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, r)
1208# 493 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1209#endif
1210 do k = 0, p
1211 do j = 0, n
1212 do i = 0, m
1213 grad_p_x(i, j, k) = 0._wp
1214
1215# 498 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1216#if defined(MFC_OpenACC)
1217# 498 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1218!$acc loop seq
1219# 498 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1220#elif defined(MFC_OpenMP)
1221# 498 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1222
1223# 498 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1224#endif
1225 do r = -fd_number, fd_number
1226 grad_p_x(i, j, k) = grad_p_x(i, j, k) + q_prim_vf(eqn_idx%E)%sf(i + r, j, k)*fd_coeff_x_pgrad(r, i)
1227 end do
1228 end do
1229 end do
1230 end do
1231
1232# 505 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1233#if defined(MFC_OpenACC)
1234# 505 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1235!$acc end parallel loop
1236# 505 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1237#elif defined(MFC_OpenMP)
1238# 505 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1239
1240# 505 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1241!$omp end target teams loop
1242# 505 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1243#endif
1244
1245 ! dp/dy at all cell centers
1246 if (n > 0) then
1247
1248# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1249
1250# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1251#if defined(MFC_OpenACC)
1252# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1253!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, r)
1254# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1255#elif defined(MFC_OpenMP)
1256# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1257
1258# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1259
1260# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1261
1262# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1263!$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, r)
1264# 509 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1265#endif
1266 do k = 0, p
1267 do j = 0, n
1268 do i = 0, m
1269 grad_p_y(i, j, k) = 0._wp
1270
1271# 514 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1272#if defined(MFC_OpenACC)
1273# 514 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1274!$acc loop seq
1275# 514 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1276#elif defined(MFC_OpenMP)
1277# 514 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1278
1279# 514 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1280#endif
1281 do r = -fd_number, fd_number
1282 grad_p_y(i, j, k) = grad_p_y(i, j, k) + q_prim_vf(eqn_idx%E)%sf(i, j + r, k)*fd_coeff_y_pgrad(r, j)
1283 end do
1284 end do
1285 end do
1286 end do
1287
1288# 521 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1289#if defined(MFC_OpenACC)
1290# 521 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1291!$acc end parallel loop
1292# 521 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1293#elif defined(MFC_OpenMP)
1294# 521 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1295
1296# 521 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1297!$omp end target teams loop
1298# 521 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1299#endif
1300 end if
1301
1302 ! dp/dz at all cell centers
1303 if (p > 0) then
1304
1305# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1306
1307# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1308#if defined(MFC_OpenACC)
1309# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1310!$acc parallel loop collapse(3) gang vector default(present) private(i, j, k, r)
1311# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1312#elif defined(MFC_OpenMP)
1313# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1314
1315# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1316
1317# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1318
1319# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1320!$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, r)
1321# 526 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1322#endif
1323 do k = 0, p
1324 do j = 0, n
1325 do i = 0, m
1326 grad_p_z(i, j, k) = 0._wp
1327
1328# 531 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1329#if defined(MFC_OpenACC)
1330# 531 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1331!$acc loop seq
1332# 531 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1333#elif defined(MFC_OpenMP)
1334# 531 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1335
1336# 531 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1337#endif
1338 do r = -fd_number, fd_number
1339 grad_p_z(i, j, k) = grad_p_z(i, j, k) + q_prim_vf(eqn_idx%E)%sf(i, j, k + r)*fd_coeff_z_pgrad(r, k)
1340 end do
1341 end do
1342 end do
1343 end do
1344
1345# 538 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1346#if defined(MFC_OpenACC)
1347# 538 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1348!$acc end parallel loop
1349# 538 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1350#elif defined(MFC_OpenMP)
1351# 538 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1352
1353# 538 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1354!$omp end target teams loop
1355# 538 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1356#endif
1357 end if
1358
1359 end subroutine s_compute_pressure_gradients
1360
1361 !! Interpolate the velocity of Eulerian field at the position of the bubble.
1362 function f_interpolate_velocity(pos, cell, i, q_prim_vf) result(v)
1363
1364
1365# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1366#ifdef _CRAYFTN
1367# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1368#if MFC_OpenACC
1369# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1370!$acc routine seq
1371# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1372#elif MFC_OpenMP
1373# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1374
1375# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1376
1377# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1378!$omp declare target device_type(any)
1379# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1380#else
1381# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1382!DIR$ NOINLINE f_interpolate_velocity
1383# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1384#endif
1385# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1386#elif MFC_OpenACC
1387# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1388!$acc routine seq
1389# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1390#elif MFC_OpenMP
1391# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1392
1393# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1394
1395# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1396!$omp declare target device_type(any)
1397# 546 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1398#endif
1399
1400 real(wp), intent(in) :: pos
1401 integer, dimension(3), intent(in) :: cell
1402 integer, intent(in) :: i
1403 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf
1404 real(wp) :: v
1405 real(wp), dimension(5) :: xi, eta, l
1406
1407 if (fd_order == 2) then
1408 if (i == 1) then
1409 xi(1) = x_cc(cell(1) - 1)
1410 eta(1) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1) - 1, cell(2), cell(3))
1411 xi(2) = x_cc(cell(1))
1412 eta(2) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1), cell(2), cell(3))
1413 xi(3) = x_cc(cell(1) + 1)
1414 eta(3) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1) + 1, cell(2), cell(3))
1415 else if (i == 2) then
1416 xi(1) = y_cc(cell(2) - 1)
1417 eta(1) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2) - 1, cell(3))
1418 xi(2) = y_cc(cell(2))
1419 eta(2) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2), cell(3))
1420 xi(3) = y_cc(cell(2) + 1)
1421 eta(3) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2) + 1, cell(3))
1422 else if (i == 3) then
1423 xi(1) = z_cc(cell(3) - 1)
1424 eta(1) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3) - 1)
1425 xi(2) = z_cc(cell(3))
1426 eta(2) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3))
1427 xi(3) = z_cc(cell(3) + 1)
1428 eta(3) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3) + 1)
1429 end if
1430
1431 l(1) = ((pos - xi(2))*(pos - xi(3)))/((xi(1) - xi(2))*(xi(1) - xi(3)))
1432 l(2) = ((pos - xi(1))*(pos - xi(3)))/((xi(2) - xi(1))*(xi(2) - xi(3)))
1433 l(3) = ((pos - xi(1))*(pos - xi(2)))/((xi(3) - xi(1))*(xi(3) - xi(2)))
1434
1435 v = l(1)*eta(1) + l(2)*eta(2) + l(3)*eta(3)
1436 else if (fd_order == 4) then
1437 if (i == 1) then
1438 xi(1) = x_cc(cell(1) - 2)
1439 eta(1) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1) - 2, cell(2), cell(3))
1440 xi(2) = x_cc(cell(1) - 1)
1441 eta(2) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1) - 1, cell(2), cell(3))
1442 xi(3) = x_cc(cell(1))
1443 eta(3) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1), cell(2), cell(3))
1444 xi(4) = x_cc(cell(1) + 1)
1445 eta(4) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1) + 1, cell(2), cell(3))
1446 xi(5) = x_cc(cell(1) + 2)
1447 eta(5) = q_prim_vf(eqn_idx%mom%beg)%sf(cell(1) + 2, cell(2), cell(3))
1448 else if (i == 2) then
1449 xi(1) = y_cc(cell(2) - 2)
1450 eta(1) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2) - 2, cell(3))
1451 xi(2) = y_cc(cell(2) - 1)
1452 eta(2) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2) - 1, cell(3))
1453 xi(3) = y_cc(cell(2))
1454 eta(3) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2), cell(3))
1455 xi(4) = y_cc(cell(2) + 1)
1456 eta(4) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2) + 1, cell(3))
1457 xi(5) = y_cc(cell(2) + 2)
1458 eta(5) = q_prim_vf(eqn_idx%mom%beg + 1)%sf(cell(1), cell(2) + 2, cell(3))
1459 else if (i == 3) then
1460 xi(1) = z_cc(cell(3) - 2)
1461 eta(1) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3) - 2)
1462 xi(2) = z_cc(cell(3) - 1)
1463 eta(2) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3) - 1)
1464 xi(3) = z_cc(cell(3))
1465 eta(3) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3))
1466 xi(4) = z_cc(cell(3) + 1)
1467 eta(4) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3) + 1)
1468 xi(5) = z_cc(cell(3) + 2)
1469 eta(5) = q_prim_vf(eqn_idx%mom%end)%sf(cell(1), cell(2), cell(3) + 2)
1470 end if
1471
1472 l(1) = ((pos - xi(2))*(pos - xi(3))*(pos - xi(4))*(pos - xi(5)))/((xi(1) - xi(2))*(xi(1) - xi(3))*(xi(1) - xi(4)) &
1473 & *(xi(1) - xi(5)))
1474 l(2) = ((pos - xi(1))*(pos - xi(3))*(pos - xi(4))*(pos - xi(5)))/((xi(2) - xi(1))*(xi(2) - xi(3))*(xi(2) - xi(4)) &
1475 & *(xi(2) - xi(5)))
1476 l(3) = ((pos - xi(1))*(pos - xi(2))*(pos - xi(4))*(pos - xi(5)))/((xi(3) - xi(1))*(xi(3) - xi(2))*(xi(3) - xi(4)) &
1477 & *(xi(3) - xi(5)))
1478 l(4) = ((pos - xi(1))*(pos - xi(2))*(pos - xi(3))*(pos - xi(5)))/((xi(4) - xi(1))*(xi(4) - xi(2))*(xi(4) - xi(3)) &
1479 & *(xi(4) - xi(5)))
1480 l(5) = ((pos - xi(1))*(pos - xi(2))*(pos - xi(3))*(pos - xi(4)))/((xi(5) - xi(1))*(xi(5) - xi(2))*(xi(5) - xi(3)) &
1481 & *(xi(5) - xi(4)))
1482
1483 v = l(1)*eta(1) + l(2)*eta(2) + l(3)*eta(3) + l(4)*eta(4) + l(5)*eta(5)
1484 end if
1485
1486 end function f_interpolate_velocity
1487
1488 !! Calculate the force on a bubble based on the pressure gradient, velocity, and drag model.
1489 !! @param mg Mass of the gas in the bubble
1490 !! @param mv Mass of the liquid in the bubble
1491 !! @param Re Reynolds number
1492 !! @param i Direction of the velocity (1: x, 2: y, 3: z)
1493 function f_get_bubble_force(pos, rad, rdot, vel, mg, mv, Re, rho, cell, i, q_prim_vf) result(force)
1494
1495
1496# 643 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1497#if MFC_OpenACC
1498# 643 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1499!$acc routine seq
1500# 643 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1501#elif MFC_OpenMP
1502# 643 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1503
1504# 643 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1505
1506# 643 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1507!$omp declare target device_type(any)
1508# 643 "/home/runner/work/MFC/MFC/src/simulation/m_bubbles_EL_kernels.fpp"
1509#endif
1510 real(wp), intent(in) :: pos, rad, rdot, mg, mv, re, rho, vel
1511 integer, dimension(3), intent(in) :: cell
1512 integer, intent(in) :: i
1513 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf
1514 real(wp) :: dp, vol, force
1515 real(wp) :: v_rel
1516
1517 if (fd_order > 1) then
1518 v_rel = vel - f_interpolate_velocity(pos, cell, i, q_prim_vf)
1519 else
1520 v_rel = vel - q_prim_vf(eqn_idx%mom%beg + i - 1)%sf(cell(1), cell(2), cell(3))
1521 end if
1522
1523 force = 0._wp
1524
1525 if (lag_params%drag_model == 1) then ! Free slip Stokes drag
1526 force = force - (4._wp*pi*rad*v_rel)/re
1527 else if (lag_params%drag_model == 2) then ! No slip Stokes drag
1528 force = force - (6._wp*pi*rad*v_rel)/re
1529 else if (lag_params%drag_model == 3) then ! Levich drag
1530 force = force - (12._wp*pi*rad*v_rel)/re
1531 end if
1532
1533 if (lag_pressure_force) then
1534 ! Use precomputed cell-centered pressure gradients
1535 if (i == 1) then
1536 dp = grad_p_x(cell(1), cell(2), cell(3))
1537 else if (i == 2) then
1538 dp = grad_p_y(cell(1), cell(2), cell(3))
1539 else if (i == 3) then
1540 dp = grad_p_z(cell(1), cell(2), cell(3))
1541 end if
1542
1543 vol = (4._wp/3._wp)*pi*(rad**3._wp)
1544 force = force - vol*dp
1545 end if
1546
1547 if (lag_params%gravity_force) then
1548 force = force + (mg + mv)*accel_bf(i)
1549 end if
1550
1551 end function f_get_bubble_force
1552
1553end module m_bubbles_el_kernels
integer, intent(in) l
Kernel functions (Gaussian, delta) that smear Lagrangian bubble effects onto the Eulerian grid.
subroutine s_build_cell_list(nbubs, lbk_s)
Builds a sorted cell list mapping each interior cell (0:m,0:n,0:p) to its resident bubbles....
real(wp), dimension(:,:), allocatable fd_coeff_z_pgrad
subroutine s_deltafunc(nbubs, lbk_rad, lbk_vel, lbk_s, updatedvar, kcomp)
Cell-centric delta-function smearing using the cell list (no GPU atomics). Each bubble only affects t...
subroutine s_applygaussian(center, cellaux, nodecoord, stddsv, strength_idx, func)
Evaluate the Gaussian kernel at a grid node for a given bubble center.
subroutine s_check_celloutside(cellaux, celloutside)
Check if the current cell is outside the computational domain including ghost cells.
real(wp), dimension(:,:,:), allocatable grad_p_y
subroutine s_compute_pressure_gradients(q_prim_vf)
Precompute cell-centered pressure gradients (dp/dx, dp/dy, dp/dz).
real(wp) function f_get_bubble_force(pos, rad, rdot, vel, mg, mv, re, rho, cell, i, q_prim_vf)
integer, dimension(:), allocatable cell_list_idx
real(wp) function f_interpolate_velocity(pos, cell, i, q_prim_vf)
subroutine s_shift_cell_symmetric_bc(cellaux, cell)
Relocate cells that intersect a symmetric boundary.
integer, dimension(:,:,:), allocatable cell_list_start
subroutine s_compute_stddsv(cell, volpart, stddsv)
Calculates the standard deviation of the bubble being smeared in the Eulerian framework.
real(wp), dimension(:,:,:), allocatable grad_p_x
real(wp), dimension(:,:), allocatable fd_coeff_y_pgrad
real(wp), dimension(:,:,:), allocatable grad_p_z
subroutine s_get_char_vol(cellx, celly, cellz, charvol)
Compute the characteristic cell volume.
subroutine s_smoothfunction(nbubs, lbk_rad, lbk_vel, lbk_s, lbk_pos, updatedvar, kcomp)
Smear the Lagrangian bubble effects onto the Eulerian grid using the selected kernel.
subroutine s_gaussian(nbubs, lbk_rad, lbk_vel, lbk_s, lbk_pos, updatedvar, kcomp)
Cell-centric gaussian smearing using the cell list (no GPU atomics). Each grid cell accumulates contr...
subroutine s_get_cell(s_cell, get_cell)
Convert bubble computational coordinates from real to integer cell indices.
integer, dimension(:,:,:), allocatable cell_list_count
real(wp), dimension(:,:), allocatable fd_coeff_x_pgrad
MPI halo exchange, domain decomposition, and buffer packing/unpacking for the simulation solver.