MFC
Exascale flow solver
Loading...
Searching...
No Matches
m_thinc.fpp.f90
Go to the documentation of this file.
1# 1 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
2!>
3!! @file
4!! @brief Contains module m_thinc
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# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
20
21# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
22
23# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
24
25# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
26
27# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
28
29# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
30
31# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
32! New line at end of file is required for FYPP
33# 2 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
34# 1 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 1
35# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
36# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
37# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
38# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
39# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
40# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
41
42# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
43# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
44# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
45
46# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
47
48# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
49
50# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
51
52# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
53
54# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
55
56# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
57
58# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
59! New line at end of file is required for FYPP
60# 2 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 2
61
62# 4 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
63# 5 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
64# 6 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
65# 7 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
66# 8 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
67
68# 20 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
69
70# 43 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
71
72# 48 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
73
74# 53 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
75
76# 58 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
77
78# 63 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
79
80# 68 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
81
82# 76 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
83
84# 81 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
85
86# 86 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
87
88# 91 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
89
90# 96 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
91
92# 101 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
93
94# 106 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
95
96# 111 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
97
98# 116 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
99
100# 121 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
101
102# 151 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
103
104# 192 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
105
106# 206 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
107
108# 231 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
109
110# 242 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
111
112# 244 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
113# 255 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
114
115# 284 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
116
117# 294 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
118
119# 304 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
120
121# 313 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
122
123# 330 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
124
125# 340 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
126
127# 347 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
128
129# 353 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
130
131# 359 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
132
133# 365 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
134
135# 371 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
136
137# 377 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
138! New line at end of file is required for FYPP
139# 3 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
140# 1 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 1
141# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
142# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
143# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
144# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
145# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
146# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
147
148# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
149# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
150# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
151
152# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
153
154# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
155
156# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
157
158# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
159
160# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
161
162# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
163
164# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
165! New line at end of file is required for FYPP
166# 2 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 2
167
168# 7 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
169
170# 17 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
171
172# 22 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
173
174# 27 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
175
176# 32 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
177
178# 37 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
179
180# 42 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
181
182# 47 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
183
184# 52 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
185
186# 57 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
187
188# 62 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
189
190# 73 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
191
192# 78 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
193
194# 83 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
195
196# 88 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
197
198# 103 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
199
200# 131 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
201
202# 160 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
203
204# 175 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
205
206# 193 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
207
208# 215 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
209
210# 244 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
211
212# 259 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
213
214# 269 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
215
216# 278 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
217
218# 294 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
219
220# 304 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
221
222# 311 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
223! New line at end of file is required for FYPP
224# 4 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
225
226! GPU parallel region (scalar reductions, maxval/minval)
227# 23 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
228
229! GPU parallel loop over threads (most common GPU macro)
230# 43 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
231
232! Required closing for GPU_PARALLEL_LOOP
233# 55 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
234
235! Mark routine for device compilation
236# 112 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
237
238! Declare device-resident data
239# 130 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
240
241! Inner loop within a GPU parallel region
242# 145 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
243
244! Scoped GPU data region
245# 164 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
246
247! Host code with device pointers (for MPI with GPU buffers)
248# 193 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
249
250! Allocate device memory (unscoped)
251# 207 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
252
253! Free device memory
254# 219 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
255
256! Atomic operation on device
257# 231 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
258
259! End atomic capture block
260# 242 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
261
262! Copy data between host and device
263# 254 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
264
265! Synchronization barrier
266# 266 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
267
268! Import GPU library module (openacc or omp_lib)
269# 275 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
270
271! Emit code only for AMD compiler
272# 282 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
273
274! Emit code for non-Cray compilers
275# 289 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
276
277! Emit code only for Cray compiler
278# 296 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
279
280! Emit code for non-NVIDIA compilers
281# 303 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
282
283# 305 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
284# 306 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
285! New line at end of file is required for FYPP
286# 2 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 2
287
288# 14 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
289
290! Caution: This macro requires the use of a binding script to set CUDA_VISIBLE_DEVICES, such that we have one GPU device per MPI
291! rank. That's because for both cudaMemAdvise (preferred location) and cudaMemPrefetchAsync we use location = device_id = 0. For an
292! example see misc/nvidia_uvm/bind.sh. NVIDIA unified memory page placement hint
293# 57 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
294
295! Allocate and create GPU device memory
296# 77 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
297
298! Free GPU device memory and deallocate
299# 85 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
300
301! Cray-specific GPU pointer setup for vector fields
302# 109 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
303
304! Cray-specific GPU pointer setup for scalar fields
305# 125 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
306
307! Cray-specific GPU pointer setup for acoustic source spatials
308# 150 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
309
310# 156 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
311
312# 163 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
313! New line at end of file is required for FYPP
314# 6 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp" 2
315
316!> @brief THINC and MTHINC interface compression for volume fraction sharpening. THINC (int_comp=1): 1D directional interface
317!! compression applied after MUSCL/WENO reconstruction. Uses hyperbolic tangent profile per dimension. MTHINC (int_comp=2):
318!! Multi-dimensional THINC that reconstructs a tanh profile oriented along the interface normal (computed from the gradient of
319!! alpha), then integrates that profile over each cell face using Gaussian quadrature. A Newton iteration enforces the conservation
320!! constraint (cell-averaged alpha is preserved). Reference: B. Xie and F. Xiao, "Toward efficient and accurate interface capturing
321!! on arbitrary hybrid unstructured grids: The THINC method with quadratic surface representation and Gaussian quadrature," Journal
322!! of Computational Physics, vol. 349, pp. 415-440, 2017.
323module m_thinc
324
325 use m_derived_types
326 use m_global_parameters
327 use m_helper
328
329#ifdef MFC_OpenACC
330 use openacc
331#endif
332
333 implicit none
334
335 private; public :: s_initialize_thinc_module, s_thinc_compression, s_compute_mthinc_normals, s_finalize_thinc_module
336
337 !> 3-point Gauss-Legendre quadrature on [-1/2, 1/2]
338 !> Node locations: +-sqrt(3/5)/2, 0
339 real(wp) :: gq3_pts(3) = [-5e-1_wp*0.7745966692414834_wp, 0._wp, 5e-1_wp*0.7745966692414834_wp]
340 !> Weights: 5/18, 8/18, 5/18
341 real(wp) :: gq3_wts(3) = [5._wp/18._wp, 8._wp/18._wp, 5._wp/18._wp]
342 !> ln(2)
343 real(wp) :: ln2 = 0.6931471805599453_wp
344
345# 35 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
346#if defined(MFC_OpenACC)
347# 35 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
348!$acc declare create(gq3_pts, gq3_wts, ln2)
349# 35 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
350#elif defined(MFC_OpenMP)
351# 35 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
352!$omp declare target (gq3_pts, gq3_wts, ln2)
353# 35 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
354#endif
355
356 real(wp), allocatable, dimension(:,:,:,:) :: mthinc_nhat !> Unit normal vector
357 real(wp), allocatable, dimension(:,:,:) :: mthinc_d !> Interface position parameter
358
359# 39 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
360#if defined(MFC_OpenACC)
361# 39 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
362!$acc declare create(mthinc_nhat, mthinc_d)
363# 39 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
364#elif defined(MFC_OpenMP)
365# 39 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
366!$omp declare target (mthinc_nhat, mthinc_d)
367# 39 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
368#endif
369
370contains
371
372 !> @brief Stable computation of ln(cosh(x))
373 function f_log_cosh(x) result(res)
374
375
376# 46 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
377#if MFC_OpenACC
378# 46 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
379!$acc routine seq
380# 46 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
381#elif MFC_OpenMP
382# 46 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
383
384# 46 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
385
386# 46 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
387!$omp declare target device_type(any)
388# 46 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
389#endif
390 real(wp), intent(in) :: x
391 real(wp) :: res, ax
392
393 ax = abs(x)
394 if (ax > 20._wp) then
395 res = ax - ln2
396 else
397 res = ax + log(1._wp + exp(-2._wp*ax)) - ln2
398 end if
399
400 end function f_log_cosh
401
402 !> @brief Stable difference: log_cosh(a+h) - log_cosh(a-h) = 2*atanh(tanh(a)*tanh(h)). Avoids catastrophic cancellation when h
403 !! is small relative to a.
404 function f_log_cosh_diff(a, h) result(res)
405
406
407# 63 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
408#if MFC_OpenACC
409# 63 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
410!$acc routine seq
411# 63 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
412#elif MFC_OpenMP
413# 63 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
414
415# 63 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
416
417# 63 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
418!$omp declare target device_type(any)
419# 63 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
420#endif
421 real(wp), intent(in) :: a, h
422 real(wp) :: res, t
423
424 t = tanh(a)*tanh(h)
425 res = 2._wp*atanh(sign(min(abs(t), 1._wp - epsilon(1._wp)), t))
426
427 end function f_log_cosh_diff
428
429 !> @brief Analytical 1-D integral of the THINC function
430 function f_thinc_integral_1d(a, b) result(res)
431
432
433# 75 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
434#if MFC_OpenACC
435# 75 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
436!$acc routine seq
437# 75 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
438#elif MFC_OpenMP
439# 75 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
440
441# 75 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
442
443# 75 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
444!$omp declare target device_type(any)
445# 75 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
446#endif
447 real(wp), intent(in) :: a, b
448 real(wp) :: res
449
450 if (abs(b) < verysmall) then
451 res = 5e-1_wp*(1._wp + tanh(a))
452 else
453 res = 5e-1_wp + f_log_cosh_diff(a, 5e-1_wp*b)/(2._wp*b)
454 end if
455
456 end function f_thinc_integral_1d
457
458 !> @brief Volume integral of H(xi) = 0.5*(1 + tanh(beta*(n.xi + d))) over the cell [-1/2, 1/2]^ndim
459 function f_mthinc_volume_integral(n1, n2, n3, d, beta, ndim) result(res)
460
461
462# 90 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
463#if MFC_OpenACC
464# 90 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
465!$acc routine seq
466# 90 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
467#elif MFC_OpenMP
468# 90 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
469
470# 90 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
471
472# 90 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
473!$omp declare target device_type(any)
474# 90 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
475#endif
476 real(wp), intent(in) :: n1, n2, n3, d, beta
477 integer, intent(in) :: ndim
478 real(wp) :: res, a
479 integer :: q1, q2
480
481 if (ndim == 1) then
482 res = f_thinc_integral_1d(beta*d, beta*n1)
483 else if (ndim == 2) then
484 res = 0._wp
485 do q1 = 1, 3
486 a = beta*(n2*gq3_pts(q1) + d)
487 res = res + gq3_wts(q1)*f_thinc_integral_1d(a, beta*n1)
488 end do
489 else
490 res = 0._wp
491 do q1 = 1, 3
492 do q2 = 1, 3
493 a = beta*(n2*gq3_pts(q1) + n3*gq3_pts(q2) + d)
494 res = res + gq3_wts(q1)*gq3_wts(q2)*f_thinc_integral_1d(a, beta*n1)
495 end do
496 end do
497 end if
498
499 end function f_mthinc_volume_integral
500
501 !> @brief Derivative dV/dd of the volume integral (for Newton iteration)
502 function f_mthinc_volume_integral_dd(n1, n2, n3, d, beta, ndim) result(res)
503
504
505# 119 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
506#if MFC_OpenACC
507# 119 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
508!$acc routine seq
509# 119 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
510#elif MFC_OpenMP
511# 119 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
512
513# 119 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
514
515# 119 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
516!$omp declare target device_type(any)
517# 119 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
518#endif
519 real(wp), intent(in) :: n1, n2, n3, d, beta
520 integer, intent(in) :: ndim
521 real(wp) :: res, th
522 integer :: q1, q2, q3
523
524 res = 0._wp
525 if (ndim == 1) then
526 do q1 = 1, 3
527 th = tanh(beta*(n1*gq3_pts(q1) + d))
528 res = res + gq3_wts(q1)*(1._wp - th*th)
529 end do
530 else if (ndim == 2) then
531 do q1 = 1, 3
532 do q2 = 1, 3
533 th = tanh(beta*(n1*gq3_pts(q1) + n2*gq3_pts(q2) + d))
534 res = res + gq3_wts(q1)*gq3_wts(q2)*(1._wp - th*th)
535 end do
536 end do
537 else
538 do q1 = 1, 3
539 do q2 = 1, 3
540 do q3 = 1, 3
541 th = tanh(beta*(n1*gq3_pts(q1) + n2*gq3_pts(q2) + n3*gq3_pts(q3) + d))
542 res = res + gq3_wts(q1)*gq3_wts(q2)*gq3_wts(q3)*(1._wp - th*th)
543 end do
544 end do
545 end do
546 end if
547 res = 5e-1_wp*beta*res
548
549 end function f_mthinc_volume_integral_dd
550
551 !> @brief Solve for the interface-position parameter d
552 function f_mthinc_solve_d(n1, n2, n3, beta, alpha_cell, ndim) result(d)
553
554
555# 155 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
556#if MFC_OpenACC
557# 155 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
558!$acc routine seq
559# 155 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
560#elif MFC_OpenMP
561# 155 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
562
563# 155 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
564
565# 155 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
566!$omp declare target device_type(any)
567# 155 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
568#endif
569 real(wp), intent(in) :: n1, n2, n3, beta, alpha_cell
570 integer, intent(in) :: ndim
571 real(wp) :: d, v, residual, dv
572 integer :: iter
573
574 d = 0._wp
575 do iter = 1, 30
576 v = f_mthinc_volume_integral(n1, n2, n3, d, beta, ndim)
577 residual = v - alpha_cell
578 if (abs(residual) < verysmall) exit
579 dv = f_mthinc_volume_integral_dd(n1, n2, n3, d, beta, ndim)
580 if (abs(dv) < verysmall) exit
581 d = d - residual/dv
582 end do
583
584 end function f_mthinc_solve_d
585
586 !> @brief Face-averaged THINC function at a cell face. face_dir: 1=x, 2=y, 3=z. face_pos: -0.5 (low) or +0.5 (high). The face
587 !! coordinate in face_dir is fixed; remaining directions are integrated over [-1/2, 1/2] (analytically along one, Gauss for
588 !! others).
589 function f_mthinc_face_average(n1, n2, n3, d, beta, face_dir, face_pos, ndim) result(res)
590
591
592# 178 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
593#if MFC_OpenACC
594# 178 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
595!$acc routine seq
596# 178 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
597#elif MFC_OpenMP
598# 178 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
599
600# 178 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
601
602# 178 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
603!$omp declare target device_type(any)
604# 178 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
605#endif
606 real(wp), intent(in) :: n1, n2, n3, d, beta, face_pos
607 integer, intent(in) :: face_dir, ndim
608 real(wp) :: res, a, n_face, n_t0, n_t1
609 integer :: q
610
611 if (ndim == 1) then
612 res = 5e-1_wp*(1._wp + tanh(beta*(n1*face_pos + d)))
613 else if (ndim == 2) then
614 ! One transverse direction
615 if (face_dir == 1) then
616 n_face = n1; n_t0 = n2
617 else
618 n_face = n2; n_t0 = n1
619 end if
620 a = beta*(n_face*face_pos + d)
621 res = f_thinc_integral_1d(a, beta*n_t0)
622 else
623 ! Two transverse directions
624 if (face_dir == 1) then
625 n_face = n1; n_t0 = n2; n_t1 = n3
626 else if (face_dir == 2) then
627 n_face = n2; n_t0 = n1; n_t1 = n3
628 else
629 n_face = n3; n_t0 = n1; n_t1 = n2
630 end if
631 res = 0._wp
632 do q = 1, 3
633 a = beta*(n_face*face_pos + n_t0*gq3_pts(q) + d)
634 res = res + gq3_wts(q)*f_thinc_integral_1d(a, beta*n_t1)
635 end do
636 end if
637
638 end function f_mthinc_face_average
639
640 subroutine s_initialize_thinc_module()
641
642 if (int_comp == 2) then
643#ifdef MFC_DEBUG
644# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
645 block
646# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
647 use iso_fortran_env, only: output_unit
648# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
649
650# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
651 print *, 'm_thinc.fpp:216: ', '@:ALLOCATE(mthinc_nhat(1:3, idwbuff(1)%beg:idwbuff(1)%end, idwbuff(2)%beg:idwbuff(2)%end, idwbuff(3)%beg:idwbuff(3)%end))'
652# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
653
654# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
655 call flush (output_unit)
656# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
657 end block
658# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
659#endif
660# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
661 allocate (mthinc_nhat(1:3, idwbuff(1)%beg:idwbuff(1)%end, idwbuff(2)%beg:idwbuff(2)%end, idwbuff(3)%beg:idwbuff(3)%end))
662# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
663
664# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
665
666# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
667#if defined(MFC_OpenACC)
668# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
669!$acc enter data create(mthinc_nhat)
670# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
671#elif defined(MFC_OpenMP)
672# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
673!$omp target enter data map(always,alloc:mthinc_nhat)
674# 216 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
675#endif
676# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
677#ifdef MFC_DEBUG
678# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
679 block
680# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
681 use iso_fortran_env, only: output_unit
682# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
683
684# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
685 print *, 'm_thinc.fpp:218: ', '@:ALLOCATE(mthinc_d(idwbuff(1)%beg:idwbuff(1)%end, idwbuff(2)%beg:idwbuff(2)%end, idwbuff(3)%beg:idwbuff(3)%end))'
686# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
687
688# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
689 call flush (output_unit)
690# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
691 end block
692# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
693#endif
694# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
695 allocate (mthinc_d(idwbuff(1)%beg:idwbuff(1)%end, idwbuff(2)%beg:idwbuff(2)%end, idwbuff(3)%beg:idwbuff(3)%end))
696# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
697
698# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
699
700# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
701#if defined(MFC_OpenACC)
702# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
703!$acc enter data create(mthinc_d)
704# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
705#elif defined(MFC_OpenMP)
706# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
707!$omp target enter data map(always,alloc:mthinc_d)
708# 218 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
709#endif
710 end if
711
712 end subroutine s_initialize_thinc_module
713
714 !> @brief Compute the unit normal and interface-position parameter d at each interface cell
715 subroutine s_compute_mthinc_normals(v_vf)
716
717 type(scalar_field), dimension(:), intent(in) :: v_vf
718 integer :: j, k, l
719 real(wp) :: nr_x, nr_y, nr_z, nmag, nmax, ac
720 type(int_bounds_info), dimension(3) :: id_norm
721
722
723# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
724
725# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
726#if defined(MFC_OpenACC)
727# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
728!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l)
729# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
730#elif defined(MFC_OpenMP)
731# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
732
733# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
734
735# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
736
737# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
738!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(j, k, l)
739# 231 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
740#endif
741 do l = idwbuff(3)%beg, idwbuff(3)%end
742 do k = idwbuff(2)%beg, idwbuff(2)%end
743 do j = idwbuff(1)%beg, idwbuff(1)%end
744 mthinc_nhat(1, j, k, l) = 0._wp
745 mthinc_nhat(2, j, k, l) = 0._wp
746 mthinc_nhat(3, j, k, l) = 0._wp
747 mthinc_d(j, k, l) = 0._wp
748 end do
749 end do
750 end do
751
752# 242 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
753#if defined(MFC_OpenACC)
754# 242 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
755!$acc end parallel loop
756# 242 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
757#elif defined(MFC_OpenMP)
758# 242 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
759
760# 242 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
761!$omp end target teams loop
762# 242 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
763#endif
764
765 id_norm(1)%beg = idwbuff(1)%beg + 1; id_norm(1)%end = idwbuff(1)%end - 1
766 id_norm(2)%beg = 0; id_norm(2)%end = 0
767 id_norm(3)%beg = 0; id_norm(3)%end = 0
768 if (n > 0) then
769 id_norm(2)%beg = idwbuff(2)%beg + 1; id_norm(2)%end = idwbuff(2)%end - 1
770 end if
771 if (p > 0) then
772 id_norm(3)%beg = idwbuff(3)%beg + 1; id_norm(3)%end = idwbuff(3)%end - 1
773 end if
774
775 ! Compute unit normal and solve for d at interior
776
777# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
778
779# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
780#if defined(MFC_OpenACC)
781# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
782!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l, nr_x, nr_y, nr_z, nmag, nmax, ac) copyin(id_norm)
783# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
784#elif defined(MFC_OpenMP)
785# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
786
787# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
788
789# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
790
791# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
792!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(j, k, l, nr_x, nr_y, nr_z, nmag, nmax, ac) map(to:id_norm)
793# 255 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
794#endif
795 do l = id_norm(3)%beg, id_norm(3)%end
796 do k = id_norm(2)%beg, id_norm(2)%end
797 do j = id_norm(1)%beg, id_norm(1)%end
798 ac = v_vf(eqn_idx%adv%beg)%sf(j, k, l)
799
800 if (ac >= ic_eps .and. ac <= 1._wp - ic_eps) then
801 nr_x = (v_vf(eqn_idx%adv%beg)%sf(j + 1, k, l) - v_vf(eqn_idx%adv%beg)%sf(j - 1, k, l))*5e-1_wp
802
803 nr_y = 0._wp
804 if (n > 0) then
805 nr_y = (v_vf(eqn_idx%adv%beg)%sf(j, k + 1, l) - v_vf(eqn_idx%adv%beg)%sf(j, k - 1, l))*5e-1_wp
806 end if
807
808 nr_z = 0._wp
809 if (p > 0) then
810 nr_z = (v_vf(eqn_idx%adv%beg)%sf(j, k, l + 1) - v_vf(eqn_idx%adv%beg)%sf(j, k, l - 1))*5e-1_wp
811 end if
812
813 nmag = sqrt(nr_x*nr_x + nr_y*nr_y + nr_z*nr_z)
814
815 if (nmag > verysmall) then
816 nr_x = nr_x/nmag
817 nr_y = nr_y/nmag
818 nr_z = nr_z/nmag
819
820 ! Snap near-grid-aligned normals to exact alignment
821 nmax = max(abs(nr_x), abs(nr_y), abs(nr_z))
822 if (abs(nr_x) < mthinc_align_tol*nmax) nr_x = 0._wp
823 if (abs(nr_y) < mthinc_align_tol*nmax) nr_y = 0._wp
824 if (abs(nr_z) < mthinc_align_tol*nmax) nr_z = 0._wp
825 nmag = sqrt(nr_x*nr_x + nr_y*nr_y + nr_z*nr_z)
826 if (nmag > verysmall) then
827 nr_x = nr_x/nmag
828 nr_y = nr_y/nmag
829 nr_z = nr_z/nmag
830
831 mthinc_nhat(1, j, k, l) = nr_x
832 mthinc_nhat(2, j, k, l) = nr_y
833 mthinc_nhat(3, j, k, l) = nr_z
834
835 mthinc_d(j, k, l) = f_mthinc_solve_d(nr_x, nr_y, nr_z, ic_beta, ac, num_dims)
836 end if
837 end if
838 end if
839 end do
840 end do
841 end do
842
843# 303 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
844#if defined(MFC_OpenACC)
845# 303 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
846!$acc end parallel loop
847# 303 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
848#elif defined(MFC_OpenMP)
849# 303 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
850
851# 303 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
852!$omp end target teams loop
853# 303 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
854#endif
855
856 end subroutine s_compute_mthinc_normals
857
858 !> @brief Applies THINC (int_comp=1) or MTHINC (int_comp=2) interface compression to sharpen volume-fraction and density
859 !! reconstructions at material interfaces. Called after WENO/MUSCL reconstruction per direction.
860 subroutine s_thinc_compression(v_rs_ws, vL_rs_vf_x, vL_rs_vf_y, vL_rs_vf_z, vR_rs_vf_x, vR_rs_vf_y, vR_rs_vf_z, recon_dir, &
861 & is1_d, is2_d, is3_d)
862
863 real(wp), dimension(idwbuff(1)%beg:,idwbuff(2)%beg:,idwbuff(3)%beg:,1:), intent(in) :: v_rs_ws
864 real(wp), dimension(idwbuff(1)%beg:,idwbuff(2)%beg:,idwbuff(3)%beg:,1:), intent(inout) :: vl_rs_vf_x, vl_rs_vf_y, &
865 & vL_rs_vf_z, vR_rs_vf_x, vR_rs_vf_y, vR_rs_vf_z
866 integer, intent(in) :: recon_dir
867 type(int_bounds_info), intent(in) :: is1_d, is2_d, is3_d
868 integer :: j, k, l, ix, iy, iz
869 real(wp) :: acl, acr, ac, athinc, qmin, qmax, a, b, c
870 real(wp) :: sgn, moncon, beta_eff
871 real(wp) :: nh1, nh2, nh3, d_local, rho1, rho2
872 real(wp) :: rho_b, rho_e
873
874# 324 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
875 if (recon_dir == 1) then
876
877# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
878
879# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
880#if defined(MFC_OpenACC)
881# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
882!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l, ix, iy, iz, aCL, aC, aCR, aTHINC, moncon, sgn, qmin, qmax, A, B, C, beta_eff, nh1, nh2, nh3, d_local, rho1, rho2, rho_b, rho_e)
883# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
884#elif defined(MFC_OpenMP)
885# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
886
887# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
888
889# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
890
891# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
892!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(j, k, l, ix, iy, iz, aCL, aC, aCR, aTHINC, moncon, sgn, qmin, qmax, A, B, C, beta_eff, nh1, nh2, nh3, d_local, rho1, rho2, rho_b, rho_e)
893# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
894#endif
895# 327 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
896 do l = is3_d%beg, is3_d%end
897 do k = is2_d%beg, is2_d%end
898 do j = is1_d%beg, is1_d%end
899 acl = v_rs_ws(j - 1, k, l, eqn_idx%adv%beg)
900 ac = v_rs_ws(j, k, l, eqn_idx%adv%beg)
901 acr = v_rs_ws(j + 1, k, l, eqn_idx%adv%beg)
902
903 if (ac >= ic_eps .and. ac <= 1._wp - ic_eps) then
904 if (int_comp == 2 .and. n > 0) then ! MTHINC
905 ! Map reshaped (j,k,l) to physical (ix,iy,iz)
906# 338 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
907 ix = j; iy = k; iz = l
908# 344 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
909
910 nh1 = mthinc_nhat(1, ix, iy, iz)
911 nh2 = mthinc_nhat(2, ix, iy, iz)
912 nh3 = mthinc_nhat(3, ix, iy, iz)
913 d_local = mthinc_d(ix, iy, iz)
914
915 ! Skip if no valid normal was computed
916 if (nh1*nh1 + nh2*nh2 + nh3*nh3 > 5e-1_wp) then
917 rho1 = v_rs_ws(j, k, l, eqn_idx%cont%beg)/ac
918 rho2 = v_rs_ws(j, k, l, eqn_idx%cont%end)/(1._wp - ac)
919
920 ! Left face
921 athinc = f_mthinc_face_average(nh1, nh2, nh3, d_local, ic_beta, 1, -5e-1_wp, &
922 & num_dims)
923 if (athinc < ic_eps) athinc = ic_eps
924 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
925 vl_rs_vf_x(j, k, l, eqn_idx%cont%beg) = rho1*athinc
926 vl_rs_vf_x(j, k, l, eqn_idx%cont%end) = rho2*(1._wp - athinc)
927 vl_rs_vf_x(j, k, l, eqn_idx%adv%beg) = athinc
928 vl_rs_vf_x(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
929
930 ! Right face
931 athinc = f_mthinc_face_average(nh1, nh2, nh3, d_local, ic_beta, 1, 5e-1_wp, &
932 & num_dims)
933 if (athinc < ic_eps) athinc = ic_eps
934 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
935 vr_rs_vf_x(j, k, l, eqn_idx%cont%beg) = rho1*athinc
936 vr_rs_vf_x(j, k, l, eqn_idx%cont%end) = rho2*(1._wp - athinc)
937 vr_rs_vf_x(j, k, l, eqn_idx%adv%beg) = athinc
938 vr_rs_vf_x(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
939 end if
940 else ! THINC
941 moncon = (acr - ac)*(ac - acl)
942
943 if (moncon > moncon_cutoff) then
944 if (acr - acl > 0._wp) then
945 sgn = 1._wp
946 else
947 sgn = -1._wp
948 end if
949
950 beta_eff = ic_beta
951
952 qmin = min(acr, acl)
953 qmax = max(acr, acl) - qmin
954
955 c = (ac - qmin + sgm_eps)/(qmax + sgm_eps)
956 b = exp(sgn*beta_eff*(2._wp*c - 1._wp))
957 a = (b/cosh(beta_eff) - 1._wp)/tanh(beta_eff)
958
959 rho_b = v_rs_ws(j, k, l, eqn_idx%cont%beg)/ac
960 rho_e = v_rs_ws(j, k, l, eqn_idx%cont%end)/(1._wp - ac)
961
962 ! Left face
963 athinc = qmin + 5e-1_wp*qmax*(1._wp + sgn*a)
964 if (athinc < ic_eps) athinc = ic_eps
965 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
966 vl_rs_vf_x(j, k, l, eqn_idx%cont%beg) = rho_b*athinc
967 vl_rs_vf_x(j, k, l, eqn_idx%cont%end) = rho_e*(1._wp - athinc)
968 vl_rs_vf_x(j, k, l, eqn_idx%adv%beg) = athinc
969 vl_rs_vf_x(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
970
971 ! Right face
972 athinc = qmin + 5e-1_wp*qmax*(1._wp + sgn*(tanh(beta_eff) + a)/(1._wp + a*tanh(beta_eff)))
973 if (athinc < ic_eps) athinc = ic_eps
974 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
975 vr_rs_vf_x(j, k, l, eqn_idx%cont%beg) = rho_b*athinc
976 vr_rs_vf_x(j, k, l, eqn_idx%cont%end) = rho_e*(1._wp - athinc)
977 vr_rs_vf_x(j, k, l, eqn_idx%adv%beg) = athinc
978 vr_rs_vf_x(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
979 end if
980 end if
981 end if
982 end do
983 end do
984 end do
985
986# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
987#if defined(MFC_OpenACC)
988# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
989!$acc end parallel loop
990# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
991#elif defined(MFC_OpenMP)
992# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
993
994# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
995!$omp end target teams loop
996# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
997#endif
998 end if
999# 324 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1000 if (recon_dir == 2) then
1001
1002# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1003
1004# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1005#if defined(MFC_OpenACC)
1006# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1007!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l, ix, iy, iz, aCL, aC, aCR, aTHINC, moncon, sgn, qmin, qmax, A, B, C, beta_eff, nh1, nh2, nh3, d_local, rho1, rho2, rho_b, rho_e)
1008# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1009#elif defined(MFC_OpenMP)
1010# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1011
1012# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1013
1014# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1015
1016# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1017!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(j, k, l, ix, iy, iz, aCL, aC, aCR, aTHINC, moncon, sgn, qmin, qmax, A, B, C, beta_eff, nh1, nh2, nh3, d_local, rho1, rho2, rho_b, rho_e)
1018# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1019#endif
1020# 327 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1021 do l = is3_d%beg, is3_d%end
1022 do k = is2_d%beg, is2_d%end
1023 do j = is1_d%beg, is1_d%end
1024 acl = v_rs_ws(j - 1, k, l, eqn_idx%adv%beg)
1025 ac = v_rs_ws(j, k, l, eqn_idx%adv%beg)
1026 acr = v_rs_ws(j + 1, k, l, eqn_idx%adv%beg)
1027
1028 if (ac >= ic_eps .and. ac <= 1._wp - ic_eps) then
1029 if (int_comp == 2 .and. n > 0) then ! MTHINC
1030 ! Map reshaped (j,k,l) to physical (ix,iy,iz)
1031# 340 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1032 ix = k; iy = j; iz = l
1033# 344 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1034
1035 nh1 = mthinc_nhat(1, ix, iy, iz)
1036 nh2 = mthinc_nhat(2, ix, iy, iz)
1037 nh3 = mthinc_nhat(3, ix, iy, iz)
1038 d_local = mthinc_d(ix, iy, iz)
1039
1040 ! Skip if no valid normal was computed
1041 if (nh1*nh1 + nh2*nh2 + nh3*nh3 > 5e-1_wp) then
1042 rho1 = v_rs_ws(j, k, l, eqn_idx%cont%beg)/ac
1043 rho2 = v_rs_ws(j, k, l, eqn_idx%cont%end)/(1._wp - ac)
1044
1045 ! Left face
1046 athinc = f_mthinc_face_average(nh1, nh2, nh3, d_local, ic_beta, 2, -5e-1_wp, &
1047 & num_dims)
1048 if (athinc < ic_eps) athinc = ic_eps
1049 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1050 vl_rs_vf_y(j, k, l, eqn_idx%cont%beg) = rho1*athinc
1051 vl_rs_vf_y(j, k, l, eqn_idx%cont%end) = rho2*(1._wp - athinc)
1052 vl_rs_vf_y(j, k, l, eqn_idx%adv%beg) = athinc
1053 vl_rs_vf_y(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1054
1055 ! Right face
1056 athinc = f_mthinc_face_average(nh1, nh2, nh3, d_local, ic_beta, 2, 5e-1_wp, &
1057 & num_dims)
1058 if (athinc < ic_eps) athinc = ic_eps
1059 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1060 vr_rs_vf_y(j, k, l, eqn_idx%cont%beg) = rho1*athinc
1061 vr_rs_vf_y(j, k, l, eqn_idx%cont%end) = rho2*(1._wp - athinc)
1062 vr_rs_vf_y(j, k, l, eqn_idx%adv%beg) = athinc
1063 vr_rs_vf_y(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1064 end if
1065 else ! THINC
1066 moncon = (acr - ac)*(ac - acl)
1067
1068 if (moncon > moncon_cutoff) then
1069 if (acr - acl > 0._wp) then
1070 sgn = 1._wp
1071 else
1072 sgn = -1._wp
1073 end if
1074
1075 beta_eff = ic_beta
1076
1077 qmin = min(acr, acl)
1078 qmax = max(acr, acl) - qmin
1079
1080 c = (ac - qmin + sgm_eps)/(qmax + sgm_eps)
1081 b = exp(sgn*beta_eff*(2._wp*c - 1._wp))
1082 a = (b/cosh(beta_eff) - 1._wp)/tanh(beta_eff)
1083
1084 rho_b = v_rs_ws(j, k, l, eqn_idx%cont%beg)/ac
1085 rho_e = v_rs_ws(j, k, l, eqn_idx%cont%end)/(1._wp - ac)
1086
1087 ! Left face
1088 athinc = qmin + 5e-1_wp*qmax*(1._wp + sgn*a)
1089 if (athinc < ic_eps) athinc = ic_eps
1090 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1091 vl_rs_vf_y(j, k, l, eqn_idx%cont%beg) = rho_b*athinc
1092 vl_rs_vf_y(j, k, l, eqn_idx%cont%end) = rho_e*(1._wp - athinc)
1093 vl_rs_vf_y(j, k, l, eqn_idx%adv%beg) = athinc
1094 vl_rs_vf_y(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1095
1096 ! Right face
1097 athinc = qmin + 5e-1_wp*qmax*(1._wp + sgn*(tanh(beta_eff) + a)/(1._wp + a*tanh(beta_eff)))
1098 if (athinc < ic_eps) athinc = ic_eps
1099 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1100 vr_rs_vf_y(j, k, l, eqn_idx%cont%beg) = rho_b*athinc
1101 vr_rs_vf_y(j, k, l, eqn_idx%cont%end) = rho_e*(1._wp - athinc)
1102 vr_rs_vf_y(j, k, l, eqn_idx%adv%beg) = athinc
1103 vr_rs_vf_y(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1104 end if
1105 end if
1106 end if
1107 end do
1108 end do
1109 end do
1110
1111# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1112#if defined(MFC_OpenACC)
1113# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1114!$acc end parallel loop
1115# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1116#elif defined(MFC_OpenMP)
1117# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1118
1119# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1120!$omp end target teams loop
1121# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1122#endif
1123 end if
1124# 324 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1125 if (recon_dir == 3) then
1126
1127# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1128
1129# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1130#if defined(MFC_OpenACC)
1131# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1132!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l, ix, iy, iz, aCL, aC, aCR, aTHINC, moncon, sgn, qmin, qmax, A, B, C, beta_eff, nh1, nh2, nh3, d_local, rho1, rho2, rho_b, rho_e)
1133# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1134#elif defined(MFC_OpenMP)
1135# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1136
1137# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1138
1139# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1140
1141# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1142!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(j, k, l, ix, iy, iz, aCL, aC, aCR, aTHINC, moncon, sgn, qmin, qmax, A, B, C, beta_eff, nh1, nh2, nh3, d_local, rho1, rho2, rho_b, rho_e)
1143# 325 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1144#endif
1145# 327 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1146 do l = is3_d%beg, is3_d%end
1147 do k = is2_d%beg, is2_d%end
1148 do j = is1_d%beg, is1_d%end
1149 acl = v_rs_ws(j - 1, k, l, eqn_idx%adv%beg)
1150 ac = v_rs_ws(j, k, l, eqn_idx%adv%beg)
1151 acr = v_rs_ws(j + 1, k, l, eqn_idx%adv%beg)
1152
1153 if (ac >= ic_eps .and. ac <= 1._wp - ic_eps) then
1154 if (int_comp == 2 .and. n > 0) then ! MTHINC
1155 ! Map reshaped (j,k,l) to physical (ix,iy,iz)
1156# 342 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1157 ix = l; iy = k; iz = j
1158# 344 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1159
1160 nh1 = mthinc_nhat(1, ix, iy, iz)
1161 nh2 = mthinc_nhat(2, ix, iy, iz)
1162 nh3 = mthinc_nhat(3, ix, iy, iz)
1163 d_local = mthinc_d(ix, iy, iz)
1164
1165 ! Skip if no valid normal was computed
1166 if (nh1*nh1 + nh2*nh2 + nh3*nh3 > 5e-1_wp) then
1167 rho1 = v_rs_ws(j, k, l, eqn_idx%cont%beg)/ac
1168 rho2 = v_rs_ws(j, k, l, eqn_idx%cont%end)/(1._wp - ac)
1169
1170 ! Left face
1171 athinc = f_mthinc_face_average(nh1, nh2, nh3, d_local, ic_beta, 3, -5e-1_wp, &
1172 & num_dims)
1173 if (athinc < ic_eps) athinc = ic_eps
1174 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1175 vl_rs_vf_z(j, k, l, eqn_idx%cont%beg) = rho1*athinc
1176 vl_rs_vf_z(j, k, l, eqn_idx%cont%end) = rho2*(1._wp - athinc)
1177 vl_rs_vf_z(j, k, l, eqn_idx%adv%beg) = athinc
1178 vl_rs_vf_z(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1179
1180 ! Right face
1181 athinc = f_mthinc_face_average(nh1, nh2, nh3, d_local, ic_beta, 3, 5e-1_wp, &
1182 & num_dims)
1183 if (athinc < ic_eps) athinc = ic_eps
1184 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1185 vr_rs_vf_z(j, k, l, eqn_idx%cont%beg) = rho1*athinc
1186 vr_rs_vf_z(j, k, l, eqn_idx%cont%end) = rho2*(1._wp - athinc)
1187 vr_rs_vf_z(j, k, l, eqn_idx%adv%beg) = athinc
1188 vr_rs_vf_z(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1189 end if
1190 else ! THINC
1191 moncon = (acr - ac)*(ac - acl)
1192
1193 if (moncon > moncon_cutoff) then
1194 if (acr - acl > 0._wp) then
1195 sgn = 1._wp
1196 else
1197 sgn = -1._wp
1198 end if
1199
1200 beta_eff = ic_beta
1201
1202 qmin = min(acr, acl)
1203 qmax = max(acr, acl) - qmin
1204
1205 c = (ac - qmin + sgm_eps)/(qmax + sgm_eps)
1206 b = exp(sgn*beta_eff*(2._wp*c - 1._wp))
1207 a = (b/cosh(beta_eff) - 1._wp)/tanh(beta_eff)
1208
1209 rho_b = v_rs_ws(j, k, l, eqn_idx%cont%beg)/ac
1210 rho_e = v_rs_ws(j, k, l, eqn_idx%cont%end)/(1._wp - ac)
1211
1212 ! Left face
1213 athinc = qmin + 5e-1_wp*qmax*(1._wp + sgn*a)
1214 if (athinc < ic_eps) athinc = ic_eps
1215 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1216 vl_rs_vf_z(j, k, l, eqn_idx%cont%beg) = rho_b*athinc
1217 vl_rs_vf_z(j, k, l, eqn_idx%cont%end) = rho_e*(1._wp - athinc)
1218 vl_rs_vf_z(j, k, l, eqn_idx%adv%beg) = athinc
1219 vl_rs_vf_z(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1220
1221 ! Right face
1222 athinc = qmin + 5e-1_wp*qmax*(1._wp + sgn*(tanh(beta_eff) + a)/(1._wp + a*tanh(beta_eff)))
1223 if (athinc < ic_eps) athinc = ic_eps
1224 if (athinc > 1._wp - ic_eps) athinc = 1._wp - ic_eps
1225 vr_rs_vf_z(j, k, l, eqn_idx%cont%beg) = rho_b*athinc
1226 vr_rs_vf_z(j, k, l, eqn_idx%cont%end) = rho_e*(1._wp - athinc)
1227 vr_rs_vf_z(j, k, l, eqn_idx%adv%beg) = athinc
1228 vr_rs_vf_z(j, k, l, eqn_idx%adv%end) = 1._wp - athinc
1229 end if
1230 end if
1231 end if
1232 end do
1233 end do
1234 end do
1235
1236# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1237#if defined(MFC_OpenACC)
1238# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1239!$acc end parallel loop
1240# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1241#elif defined(MFC_OpenMP)
1242# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1243
1244# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1245!$omp end target teams loop
1246# 420 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1247#endif
1248 end if
1249# 423 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1250
1251 end subroutine s_thinc_compression
1252
1253 subroutine s_finalize_thinc_module()
1254
1255 if (int_comp == 2) then
1256#ifdef MFC_DEBUG
1257# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1258 block
1259# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1260 use iso_fortran_env, only: output_unit
1261# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1262
1263# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1264 print *, 'm_thinc.fpp:429: ', '@:DEALLOCATE(mthinc_nhat)'
1265# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1266
1267# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1268 call flush (output_unit)
1269# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1270 end block
1271# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1272#endif
1273# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1274
1275# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1276#if defined(MFC_OpenACC)
1277# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1278!$acc exit data delete(mthinc_nhat)
1279# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1280#elif defined(MFC_OpenMP)
1281# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1282!$omp target exit data map(release:mthinc_nhat)
1283# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1284#endif
1285# 429 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1286 deallocate (mthinc_nhat)
1287#ifdef MFC_DEBUG
1288# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1289 block
1290# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1291 use iso_fortran_env, only: output_unit
1292# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1293
1294# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1295 print *, 'm_thinc.fpp:430: ', '@:DEALLOCATE(mthinc_d)'
1296# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1297
1298# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1299 call flush (output_unit)
1300# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1301 end block
1302# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1303#endif
1304# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1305
1306# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1307#if defined(MFC_OpenACC)
1308# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1309!$acc exit data delete(mthinc_d)
1310# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1311#elif defined(MFC_OpenMP)
1312# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1313!$omp target exit data map(release:mthinc_d)
1314# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1315#endif
1316# 430 "/home/runner/work/MFC/MFC/src/simulation/m_thinc.fpp"
1317 deallocate (mthinc_d)
1318 end if
1319
1320 end subroutine s_finalize_thinc_module
1321
1322end module m_thinc
k
integer, intent(in) k
Definition m_phase_change.fpp.f90:1067
j
integer, intent(in) j
Definition m_phase_change.fpp.f90:1067
l
integer, intent(in) l
Definition m_phase_change.fpp.f90:1067
m_derived_types
Shared derived types for field data, patch geometry, bubble dynamics, and MPI I/O structures.
Definition m_derived_types.fpp.f90:317
m_global_parameters
Global parameters for the computational domain, fluid properties, and simulation algorithm configurat...
Definition m_global_parameters.fpp.f90:328
m_global_parameters::int_comp
integer int_comp
Interface compression: 0=off, 1=THINC, 2=MTHINC.
Definition m_global_parameters.fpp.f90:490
m_global_parameters::idwbuff
type(int_bounds_info), dimension(1:3) idwbuff
Definition m_global_parameters.fpp.f90:856
m_helper
Utility routines for bubble model setup, coordinate transforms, array sampling, and special functions...
Definition m_helper.fpp.f90:328
m_thinc
THINC and MTHINC interface compression for volume fraction sharpening. THINC (int_comp=1): 1D directi...
Definition m_thinc.fpp.f90:323
m_thinc::f_thinc_integral_1d
real(wp) function f_thinc_integral_1d(a, b)
Analytical 1-D integral of the THINC function.
Definition m_thinc.fpp.f90:431
m_thinc::gq3_pts
real(wp), dimension(3) gq3_pts
3-point Gauss-Legendre quadrature on [-1/2, 1/2] Node locations: +-sqrt(3/5)/2, 0
Definition m_thinc.fpp.f90:339
m_thinc::ln2
real(wp) ln2
ln(2)
Definition m_thinc.fpp.f90:343
m_thinc::s_finalize_thinc_module
subroutine, public s_finalize_thinc_module()
Definition m_thinc.fpp.f90:1254
m_thinc::mthinc_d
real(wp), dimension(:,:,:), allocatable mthinc_d
Definition m_thinc.fpp.f90:357
m_thinc::s_initialize_thinc_module
subroutine, public s_initialize_thinc_module()
Definition m_thinc.fpp.f90:641
m_thinc::s_thinc_compression
subroutine, public s_thinc_compression(v_rs_ws, vl_rs_vf_x, vl_rs_vf_y, vl_rs_vf_z, vr_rs_vf_x, vr_rs_vf_y, vr_rs_vf_z, recon_dir, is1_d, is2_d, is3_d)
Applies THINC (int_comp=1) or MTHINC (int_comp=2) interface compression to sharpen volume-fraction an...
Definition m_thinc.fpp.f90:862
m_thinc::f_log_cosh
real(wp) function f_log_cosh(x)
Stable computation of ln(cosh(x)).
Definition m_thinc.fpp.f90:374
m_thinc::f_mthinc_solve_d
real(wp) function f_mthinc_solve_d(n1, n2, n3, beta, alpha_cell, ndim)
Solve for the interface-position parameter d.
Definition m_thinc.fpp.f90:553
m_thinc::mthinc_nhat
real(wp), dimension(:,:,:,:), allocatable mthinc_nhat
Definition m_thinc.fpp.f90:356
m_thinc::f_mthinc_volume_integral
real(wp) function f_mthinc_volume_integral(n1, n2, n3, d, beta, ndim)
Volume integral of H(xi) = 0.5*(1 + tanh(beta*(n.xi + d))) over the cell [-1/2, 1/2]^ndim.
Definition m_thinc.fpp.f90:460
m_thinc::s_compute_mthinc_normals
subroutine, public s_compute_mthinc_normals(v_vf)
Compute the unit normal and interface-position parameter d at each interface cell.
Definition m_thinc.fpp.f90:716
m_thinc::gq3_wts
real(wp), dimension(3) gq3_wts
Weights: 5/18, 8/18, 5/18.
Definition m_thinc.fpp.f90:341
m_thinc::f_log_cosh_diff
real(wp) function f_log_cosh_diff(a, h)
Stable difference: log_cosh(a+h) - log_cosh(a-h) = 2*atanh(tanh(a)*tanh(h)). Avoids catastrophic canc...
Definition m_thinc.fpp.f90:405
m_thinc::f_mthinc_face_average
real(wp) function f_mthinc_face_average(n1, n2, n3, d, beta, face_dir, face_pos, ndim)
Face-averaged THINC function at a cell face. face_dir: 1=x, 2=y, 3=z. face_pos: -0....
Definition m_thinc.fpp.f90:590
m_thinc::f_mthinc_volume_integral_dd
real(wp) function f_mthinc_volume_integral_dd(n1, n2, n3, d, beta, ndim)
Derivative dV/dd of the volume integral (for Newton iteration).
Definition m_thinc.fpp.f90:503