MFC
Exascale flow solver
Loading...
Searching...
No Matches
m_acoustic_src.fpp.f90
Go to the documentation of this file.
1# 1 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
2!>
3!! @file
4!! @brief Contains module m_acoustic_src
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_acoustic_src.fpp" 2
333
334!> @brief One-way acoustic source injection, Maeda and Colonius JCP (2017)
336
339 use m_bubbles
342 use m_constants
343
344 implicit none
345
347
348 integer, allocatable, dimension(:) :: pulse, support
349
350# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
351#if defined(MFC_OpenACC)
352# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
353!$acc declare create(pulse, support)
354# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
355#elif defined(MFC_OpenMP)
356# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
357!$omp declare target (pulse, support)
358# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
359#endif
360
361 logical, allocatable, dimension(:) :: dipole
362
363# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
364#if defined(MFC_OpenACC)
365# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
366!$acc declare create(dipole)
367# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
368#elif defined(MFC_OpenMP)
369# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
370!$omp declare target (dipole)
371# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
372#endif
373
374 real(wp), allocatable, target, dimension(:,:) :: loc_acoustic
375
376# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
377#if defined(MFC_OpenACC)
378# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
379!$acc declare create(loc_acoustic)
380# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
381#elif defined(MFC_OpenMP)
382# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
383!$omp declare target (loc_acoustic)
384# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
385#endif
386
387 real(wp), allocatable, dimension(:) :: mag, length, height, wavelength, frequency
388 real(wp), allocatable, dimension(:) :: gauss_sigma_dist, gauss_sigma_time, npulse, dir, delay
389
390# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
391#if defined(MFC_OpenACC)
392# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
393!$acc declare create(mag, length, height, wavelength, frequency)
394# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
395#elif defined(MFC_OpenMP)
396# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
397!$omp declare target (mag, length, height, wavelength, frequency)
398# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
399#endif
400
401# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
402#if defined(MFC_OpenACC)
403# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
404!$acc declare create(gauss_sigma_dist, gauss_sigma_time, npulse, dir, delay)
405# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
406#elif defined(MFC_OpenMP)
407# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
408!$omp declare target (gauss_sigma_dist, gauss_sigma_time, npulse, dir, delay)
409# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
410#endif
411
412 real(wp), allocatable, dimension(:) :: foc_length, aperture
413
414# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
415#if defined(MFC_OpenACC)
416# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
417!$acc declare create(foc_length, aperture)
418# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
419#elif defined(MFC_OpenMP)
420# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
421!$omp declare target (foc_length, aperture)
422# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
423#endif
424
425 real(wp), allocatable, dimension(:) :: element_spacing_angle, element_polygon_ratio, rotate_angle
426
427# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
428#if defined(MFC_OpenACC)
429# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
430!$acc declare create(element_spacing_angle, element_polygon_ratio, rotate_angle)
431# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
432#elif defined(MFC_OpenMP)
433# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
434!$omp declare target (element_spacing_angle, element_polygon_ratio, rotate_angle)
435# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
436#endif
437
438 real(wp), allocatable, dimension(:) :: bb_bandwidth, bb_lowest_freq
439
440# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
441#if defined(MFC_OpenACC)
442# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
443!$acc declare create(bb_bandwidth, bb_lowest_freq)
444# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
445#elif defined(MFC_OpenMP)
446# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
447!$omp declare target (bb_bandwidth, bb_lowest_freq)
448# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
449#endif
450
451 integer, allocatable, dimension(:) :: num_elements, element_on, bb_num_freq
452
453# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
454#if defined(MFC_OpenACC)
455# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
456!$acc declare create(num_elements, element_on, bb_num_freq)
457# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
458#elif defined(MFC_OpenMP)
459# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
460!$omp declare target (num_elements, element_on, bb_num_freq)
461# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
462#endif
463
464 !> @name Acoustic source terms
465 !> @{
466 real(wp), allocatable, dimension(:,:,:) :: mass_src, e_src
467 real(wp), allocatable, dimension(:,:,:,:) :: mom_src
468 !> @}
469
470# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
471#if defined(MFC_OpenACC)
472# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
473!$acc declare create(mass_src, e_src, mom_src)
474# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
475#elif defined(MFC_OpenMP)
476# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
477!$omp declare target (mass_src, e_src, mom_src)
478# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
479#endif
480
481 integer, dimension(:), allocatable :: source_spatials_num_points !< Number of non-zero source grid points for each source
482
483# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
484#if defined(MFC_OpenACC)
485# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
486!$acc declare create(source_spatials_num_points)
487# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
488#elif defined(MFC_OpenMP)
489# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
490!$omp declare target (source_spatials_num_points)
491# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
492#endif
493
494 type(source_spatial_type), dimension(:), allocatable :: source_spatials !< Data of non-zero source grid points for each source
495
496# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
497#if defined(MFC_OpenACC)
498# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
499!$acc declare create(source_spatials)
500# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
501#elif defined(MFC_OpenMP)
502# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
503!$omp declare target (source_spatials)
504# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
505#endif
506
507contains
508
509 !> Initialize the acoustic source module
511
512 integer :: i, j !< generic loop variables
513
514#ifdef MFC_DEBUG
515# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
516 block
517# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
518 use iso_fortran_env, only: output_unit
519# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
520
521# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
522 print *, 'm_acoustic_src.fpp:67: ', '@:ALLOCATE(loc_acoustic(1:3, 1:num_source), mag(1:num_source), dipole(1:num_source), support(1:num_source), length(1:num_source), height(1:num_source), wavelength(1:num_source), frequency(1:num_source), gauss_sigma_dist(1:num_source), gauss_sigma_time(1:num_source), foc_length(1:num_source), aperture(1:num_source), npulse(1:num_source), pulse(1:num_source), dir(1:num_source), delay(1:num_source), element_polygon_ratio(1:num_source), rotate_angle(1:num_source), element_spacing_angle(1:num_source), num_elements(1:num_source), element_on(1:num_source), bb_num_freq(1:num_source), bb_bandwidth(1:num_source), bb_lowest_freq(1:num_source))'
523# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
524
525# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
526 call flush (output_unit)
527# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
528 end block
529# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
530#endif
531# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
532 allocate (loc_acoustic(1:3, 1:num_source), mag(1:num_source), dipole(1:num_source), support(1:num_source), length(1:num_source), height(1:num_source), wavelength(1:num_source), frequency(1:num_source), gauss_sigma_dist(1:num_source), gauss_sigma_time(1:num_source), foc_length(1:num_source), aperture(1:num_source), npulse(1:num_source), pulse(1:num_source), dir(1:num_source), delay(1:num_source), element_polygon_ratio(1:num_source), rotate_angle(1:num_source), element_spacing_angle(1:num_source), num_elements(1:num_source), element_on(1:num_source), bb_num_freq(1:num_source), bb_bandwidth(1:num_source), bb_lowest_freq(1:num_source))
533# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
534
535# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
536
537# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
538
539# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
540
541# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
542
543# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
544
545# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
546
547# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
548
549# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
550
551# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
552
553# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
554
555# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
556
557# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
558
559# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
560
561# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
562
563# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
564
565# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
566
567# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
568
569# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
570
571# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
572
573# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
574
575# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
576
577# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
578
579# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
580
581# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
582
583# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
584#if defined(MFC_OpenACC)
585# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
586!$acc enter data create(loc_acoustic, mag, dipole, support, length, height, wavelength, frequency, gauss_sigma_dist, gauss_sigma_time, foc_length, aperture, npulse, pulse, dir, delay, element_polygon_ratio, rotate_angle, element_spacing_angle, num_elements, element_on, bb_num_freq, bb_bandwidth, bb_lowest_freq)
587# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
588#elif defined(MFC_OpenMP)
589# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
590!$omp target enter data map(always,alloc:loc_acoustic, mag, dipole, support, length, height, wavelength, frequency, gauss_sigma_dist, gauss_sigma_time, foc_length, aperture, npulse, pulse, dir, delay, element_polygon_ratio, rotate_angle, element_spacing_angle, num_elements, element_on, bb_num_freq, bb_bandwidth, bb_lowest_freq)
591# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
592#endif
593# 73 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
594
595 do i = 1, num_source
596 do j = 1, 3
597 loc_acoustic(j, i) = acoustic(i)%loc(j)
598 end do
599 mag(i) = acoustic(i)%mag
600 dipole(i) = acoustic(i)%dipole
601 support(i) = acoustic(i)%support
602 length(i) = acoustic(i)%length
603 height(i) = acoustic(i)%height
604 wavelength(i) = acoustic(i)%wavelength
605 frequency(i) = acoustic(i)%frequency
606 gauss_sigma_dist(i) = acoustic(i)%gauss_sigma_dist
607 gauss_sigma_time(i) = acoustic(i)%gauss_sigma_time
608 foc_length(i) = acoustic(i)%foc_length
609 aperture(i) = acoustic(i)%aperture
610 npulse(i) = acoustic(i)%npulse
611 pulse(i) = acoustic(i)%pulse
612 dir(i) = acoustic(i)%dir
613 element_spacing_angle(i) = acoustic(i)%element_spacing_angle
614 element_polygon_ratio(i) = acoustic(i)%element_polygon_ratio
615 num_elements(i) = acoustic(i)%num_elements
616 bb_num_freq(i) = acoustic(i)%bb_num_freq
617 bb_bandwidth(i) = acoustic(i)%bb_bandwidth
618 bb_lowest_freq(i) = acoustic(i)%bb_lowest_freq
619
620 if (acoustic(i)%element_on == dflt_int) then
621 element_on(i) = 0
622 else
623 element_on(i) = acoustic(i)%element_on
624 end if
625 if (f_is_default(acoustic(i)%rotate_angle)) then
626 rotate_angle(i) = 0._wp
627 else
628 rotate_angle(i) = acoustic(i)%rotate_angle
629 end if
630 if (f_is_default(acoustic(i)%delay)) then ! m_checker guarantees acoustic(i)%delay is set for pulse = 2 (Gaussian)
631 delay(i) = 0._wp ! Defaults to zero for sine and square waves
632 else
633 delay(i) = acoustic(i)%delay
634 end if
635 end do
636
637# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
638#if defined(MFC_OpenACC)
639# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
640!$acc update device(loc_acoustic, mag, dipole, support, length, height, wavelength, frequency, gauss_sigma_dist, gauss_sigma_time, foc_length, aperture, npulse, pulse, dir, delay, element_polygon_ratio, rotate_angle, element_spacing_angle, num_elements, element_on, bb_num_freq, bb_bandwidth, bb_lowest_freq)
641# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
642#elif defined(MFC_OpenMP)
643# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
644!$omp target update to(loc_acoustic, mag, dipole, support, length, height, wavelength, frequency, gauss_sigma_dist, gauss_sigma_time, foc_length, aperture, npulse, pulse, dir, delay, element_polygon_ratio, rotate_angle, element_spacing_angle, num_elements, element_on, bb_num_freq, bb_bandwidth, bb_lowest_freq)
645# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
646#endif
647# 118 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
648
649#ifdef MFC_DEBUG
650# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
651 block
652# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
653 use iso_fortran_env, only: output_unit
654# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
655
656# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
657 print *, 'm_acoustic_src.fpp:119: ', '@:ALLOCATE(mass_src(0:m, 0:n, 0:p))'
658# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
659
660# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
661 call flush (output_unit)
662# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
663 end block
664# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
665#endif
666# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
667 allocate (mass_src(0:m, 0:n, 0:p))
668# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
669
670# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
671
672# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
673#if defined(MFC_OpenACC)
674# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
675!$acc enter data create(mass_src)
676# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
677#elif defined(MFC_OpenMP)
678# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
679!$omp target enter data map(always,alloc:mass_src)
680# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
681#endif
682#ifdef MFC_DEBUG
683# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
684 block
685# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
686 use iso_fortran_env, only: output_unit
687# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
688
689# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
690 print *, 'm_acoustic_src.fpp:120: ', '@:ALLOCATE(mom_src(1:num_vels, 0:m, 0:n, 0:p))'
691# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
692
693# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
694 call flush (output_unit)
695# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
696 end block
697# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
698#endif
699# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
700 allocate (mom_src(1:num_vels, 0:m, 0:n, 0:p))
701# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
702
703# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
704
705# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
706#if defined(MFC_OpenACC)
707# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
708!$acc enter data create(mom_src)
709# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
710#elif defined(MFC_OpenMP)
711# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
712!$omp target enter data map(always,alloc:mom_src)
713# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
714#endif
715#ifdef MFC_DEBUG
716# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
717 block
718# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
719 use iso_fortran_env, only: output_unit
720# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
721
722# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
723 print *, 'm_acoustic_src.fpp:121: ', '@:ALLOCATE(E_src(0:m, 0:n, 0:p))'
724# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
725
726# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
727 call flush (output_unit)
728# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
729 end block
730# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
731#endif
732# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
733 allocate (e_src(0:m, 0:n, 0:p))
734# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
735
736# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
737
738# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
739#if defined(MFC_OpenACC)
740# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
741!$acc enter data create(E_src)
742# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
743#elif defined(MFC_OpenMP)
744# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
745!$omp target enter data map(always,alloc:E_src)
746# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
747#endif
748
749 end subroutine s_initialize_acoustic_src
750
751 !> Compute mass, momentum, and energy acoustic source terms and add to the RHS
752 impure subroutine s_acoustic_src_calculations(q_cons_vf, q_prim_vf, rhs_vf)
753
754 type(scalar_field), dimension(sys_size), intent(inout) :: q_cons_vf !< Conservative variables
755 type(scalar_field), dimension(sys_size), intent(inout) :: q_prim_vf !< Primitive variables
756 type(scalar_field), dimension(sys_size), intent(inout) :: rhs_vf
757
758# 135 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
759 real(wp), dimension(num_fluids) :: myalpha, myalpha_rho
760# 137 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
761 real(wp) :: myrho, pi_inf_mix, qv_dummy
762 real(wp) :: sim_time, c, gamma_mix
763 real(wp) :: blkmod_q, pres_q, alpha_q, alpha_rho_q
764 real(wp) :: frequency_local, gauss_sigma_time_local
765 real(wp) :: mass_src_diff, mom_src_diff
766 real(wp) :: source_temporal
767 real(wp) :: period_bb !< period of each sine wave in broadband source
768 real(wp) :: sl_bb !< spectral level at each frequency
769 real(wp) :: ffre_bb !< source term corresponding to each frequency
770 real(wp) :: sum_bb !< total source term for the broadband wave
771 real(wp), allocatable, dimension(:) :: phi_rn !< random phase shift for each frequency
772 integer :: i, j, k, l, q !< generic loop variables
773 integer :: ai !< acoustic source index
774 integer :: num_points
775 logical :: freq_conv_flag, gauss_conv_flag
776 integer, parameter :: mass_label = 1, mom_label = 2
777
778 sim_time = mytime ! Accumulated time, correct under adaptive dt
779
780
781# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
782
783# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
784#if defined(MFC_OpenACC)
785# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
786!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l)
787# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
788#elif defined(MFC_OpenMP)
789# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
790
791# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
792
793# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
794
795# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
796!$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)
797# 156 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
798#endif
799 do l = 0, p
800 do k = 0, n
801 do j = 0, m
802 mass_src(j, k, l) = 0._wp
803 mom_src(1, j, k, l) = 0._wp
804 e_src(j, k, l) = 0._wp
805 if (n > 0) mom_src(2, j, k, l) = 0._wp
806 if (p > 0) mom_src(3, j, k, l) = 0._wp
807 end do
808 end do
809 end do
810
811# 168 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
812#if defined(MFC_OpenACC)
813# 168 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
814!$acc end parallel loop
815# 168 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
816#elif defined(MFC_OpenMP)
817# 168 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
818
819# 168 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
820!$omp end target teams loop
821# 168 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
822#endif
823
824 ! Keep outer loop sequential because different sources can have very different number of points
825 do ai = 1, num_source
826 ! Skip if the pulse has not started yet for sine and square waves
827 if (.not. (sim_time < delay(ai) .and. (pulse(ai) == 1 .or. pulse(ai) == 3))) then
828 ! Decide if frequency need to be converted from wavelength
829 freq_conv_flag = f_is_default(frequency(ai))
830 gauss_conv_flag = f_is_default(gauss_sigma_time(ai))
831
832 num_points = source_spatials_num_points(ai) ! Use scalar to force firstprivate to prevent GPU bug
833
834 ! Calculate the broadband source
835 period_bb = 0._wp
836 sl_bb = 0._wp
837 ffre_bb = 0._wp
838 sum_bb = 0._wp
839
840 ! Allocate buffers for random phase shift
841 allocate (phi_rn(1:bb_num_freq(ai)))
842 phi_rn(1:bb_num_freq(ai)) = 0._wp
843
844 if (pulse(ai) == 4) then
845 call random_number(phi_rn(1:bb_num_freq(ai)))
846 ! Ensure all the ranks have the same random phase shift
847 call s_mpi_send_random_number(phi_rn, bb_num_freq(ai))
848 end if
849
850 do k = 1, bb_num_freq(ai)
851 ! Acoustic period of the wave at each discrete frequency
852 period_bb = 1._wp/(bb_lowest_freq(ai) + k*bb_bandwidth(ai))
853 ! Spectral level at each frequency
854 sl_bb = broadband_spectral_level_constant*mag(ai) + k*mag(ai)/broadband_spectral_level_growth_rate
855 ! Source term corresponding to each frequencies
856 ffre_bb = sqrt((2._wp*sl_bb*bb_bandwidth(ai)))*cos((sim_time)*2._wp*pi/period_bb + 2._wp*pi*phi_rn(k))
857 ! Sum up the source term of each frequency to obtain the total source term for broadband wave
858 sum_bb = sum_bb + ffre_bb
859 end do
860
861 deallocate (phi_rn)
862
863
864# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
865
866# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
867#if defined(MFC_OpenACC)
868# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
869!$acc parallel loop gang vector default(present) private(myalpha, myalpha_rho, myRho, pi_inf_mix, qv_dummy, c, blkmod_q, pres_q, alpha_q, alpha_rho_q, gamma_mix, frequency_local, &
870# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
871!$acc& gauss_sigma_time_local, mass_src_diff, mom_src_diff, source_temporal, j, k, l, q) copyin(sum_BB, freq_conv_flag, gauss_conv_flag, sim_time)
872# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
873#elif defined(MFC_OpenMP)
874# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
875
876# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
877
878# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
879
880# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
881!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) private(myalpha, myalpha_rho, myRho, &
882# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
883!$omp& pi_inf_mix, qv_dummy, c, blkmod_q, pres_q, alpha_q, alpha_rho_q, gamma_mix, frequency_local, gauss_sigma_time_local, mass_src_diff, mom_src_diff, source_temporal, j, k, l, q) &
884# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
885!$omp& map(to:sum_BB, freq_conv_flag, gauss_conv_flag, sim_time)
886# 209 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
887#endif
888# 213 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
889 do i = 1, num_points
890 j = source_spatials(ai)%coord(1, i)
891 k = source_spatials(ai)%coord(2, i)
892 l = source_spatials(ai)%coord(3, i)
893
894 ! Compute speed of sound
895
896# 219 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
897#if defined(MFC_OpenACC)
898# 219 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
899!$acc loop seq
900# 219 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
901#elif defined(MFC_OpenMP)
902# 219 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
903
904# 219 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
905#endif
906 do q = 1, num_fluids
907 myalpha_rho(q) = q_cons_vf(q)%sf(j, k, l)
908 myalpha(q) = q_cons_vf(eqn_idx%adv%beg + q - 1)%sf(j, k, l)
909 end do
910
911 call s_compute_mixture_coefficients(myalpha_rho, myalpha, myrho, gamma_mix, pi_inf_mix, qv_dummy)
912
913 ! Mixture sound speed, as s_compute_speed_of_sound computes it. Inlined because that
914 ! call faults CCE OpenACC from this loop (#1794); fold it back once that is isolated.
915 if (any_state_dependent_eos) then ! frozen mixing of per-phase moduli, as in s_compute_speed_of_sound
916 c = 0._wp
917
918# 231 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
919#if defined(MFC_OpenACC)
920# 231 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
921!$acc loop seq
922# 231 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
923#elif defined(MFC_OpenMP)
924# 231 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
925
926# 231 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
927#endif
928 do q = 1, num_fluids
929 pres_q = q_prim_vf(eqn_idx%E)%sf(j, k, l)
930 alpha_q = myalpha(q)
931 alpha_rho_q = myalpha_rho(q)
932 call s_phase_bulk_modulus(pres_q, alpha_q, alpha_rho_q, q, blkmod_q)
933 c = c + myalpha(q)*blkmod_q
934 end do
935 c = c/myrho
936 else
937 c = f_bulk_modulus(q_prim_vf(eqn_idx%E)%sf(j, k, l), gamma_mix, pi_inf_mix)/myrho
938 end if
939 if (model_eqns == model_eqns_5eq .and. bubbles_euler .and. .not. (mpp_lim .and. num_fluids > 1)) then
940 c = c/(1._wp - myalpha(num_fluids))
941 end if
942 c = sqrt(c)
943
944 ! Wavelength to frequency conversion
945 if (pulse(ai) == 1 .or. pulse(ai) == 3) frequency_local = f_frequency_local(freq_conv_flag, ai, c)
946 if (pulse(ai) == 2) gauss_sigma_time_local = f_gauss_sigma_time_local(gauss_conv_flag, ai, c)
947
948 ! Update momentum source term
949 call s_source_temporal(sim_time, c, ai, mom_label, frequency_local, gauss_sigma_time_local, source_temporal, &
950 & sum_bb)
951 mom_src_diff = source_temporal*source_spatials(ai)%val(i)
952
953 if (dipole(ai)) then ! Double amplitude & No momentum source term (only works for Planar)
954 mass_src(j, k, l) = mass_src(j, k, l) + 2._wp*mom_src_diff/c
955 e_src(j, k, l) = e_src(j, k, l) + 2._wp*mom_src_diff*c*gamma_mix
956 cycle
957 end if
958
959 if (n == 0) then ! 1D
960 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*sign(1._wp, dir(ai)) ! Left or right-going wave
961 else if (p == 0) then ! 2D
962 if (support(ai) < 5) then ! Planar
963 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*cos(dir(ai))
964 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*sin(dir(ai))
965 else
966 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*cos(source_spatials(ai)%angle(i))
967 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*sin(source_spatials(ai)%angle(i))
968 end if
969 else ! 3D
970 if (support(ai) < 5) then ! Planar
971 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*cos(dir(ai))
972 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*sin(dir(ai))
973 else
974 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*source_spatials(ai)%xyz_to_r_ratios(1, i)
975 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*source_spatials(ai)%xyz_to_r_ratios(2, i)
976 mom_src(3, j, k, l) = mom_src(3, j, k, l) + mom_src_diff*source_spatials(ai)%xyz_to_r_ratios(3, i)
977 end if
978 end if
979
980 ! Update mass source term
981 if (support(ai) < 5) then ! Planar
982 mass_src_diff = mom_src_diff/c
983 else ! Spherical or cylindrical support
984 ! Mass source term must be calculated differently using a correction term for spherical and cylindrical
985 ! support
986 call s_source_temporal(sim_time, c, ai, mass_label, frequency_local, gauss_sigma_time_local, &
987 & source_temporal, sum_bb)
988 mass_src_diff = source_temporal*source_spatials(ai)%val(i)
989 end if
990 mass_src(j, k, l) = mass_src(j, k, l) + mass_src_diff
991
992 ! Update energy source term
993 e_src(j, k, l) = e_src(j, k, l) + mass_src_diff*c**2._wp*gamma_mix
994 end do
995
996# 299 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
997#if defined(MFC_OpenACC)
998# 299 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
999!$acc end parallel loop
1000# 299 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1001#elif defined(MFC_OpenMP)
1002# 299 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1003
1004# 299 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1005!$omp end target teams loop
1006# 299 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1007#endif
1008 end if
1009 end do
1010
1011 ! Update the rhs variables
1012
1013# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1014
1015# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1016#if defined(MFC_OpenACC)
1017# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1018!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l)
1019# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1020#elif defined(MFC_OpenMP)
1021# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1022
1023# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1024
1025# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1026
1027# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1028!$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)
1029# 304 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1030#endif
1031 do l = 0, p
1032 do k = 0, n
1033 do j = 0, m
1034
1035# 308 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1036#if defined(MFC_OpenACC)
1037# 308 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1038!$acc loop seq
1039# 308 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1040#elif defined(MFC_OpenMP)
1041# 308 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1042
1043# 308 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1044#endif
1045 do q = eqn_idx%cont%beg, eqn_idx%cont%end
1046 rhs_vf(q)%sf(j, k, l) = rhs_vf(q)%sf(j, k, l) + mass_src(j, k, l)
1047 end do
1048
1049# 312 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1050#if defined(MFC_OpenACC)
1051# 312 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1052!$acc loop seq
1053# 312 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1054#elif defined(MFC_OpenMP)
1055# 312 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1056
1057# 312 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1058#endif
1059 do q = eqn_idx%mom%beg, eqn_idx%mom%end
1060 rhs_vf(q)%sf(j, k, l) = rhs_vf(q)%sf(j, k, l) + mom_src(q - eqn_idx%cont%end, j, k, l)
1061 end do
1062 rhs_vf(eqn_idx%E)%sf(j, k, l) = rhs_vf(eqn_idx%E)%sf(j, k, l) + e_src(j, k, l)
1063 end do
1064 end do
1065 end do
1066
1067# 320 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1068#if defined(MFC_OpenACC)
1069# 320 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1070!$acc end parallel loop
1071# 320 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1072#elif defined(MFC_OpenMP)
1073# 320 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1074
1075# 320 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1076!$omp end target teams loop
1077# 320 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1078#endif
1079
1080 end subroutine s_acoustic_src_calculations
1081
1082 !> Compute the temporally varying amplitude of the pulse
1083 elemental subroutine s_source_temporal(sim_time, c, ai, term_index, frequency_local, gauss_sigma_time_local, source, sum_BB)
1084
1085
1086# 327 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1087#if MFC_OpenACC
1088# 327 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1089!$acc routine seq
1090# 327 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1091#elif MFC_OpenMP
1092# 327 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1093
1094# 327 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1095
1096# 327 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1097!$omp declare target device_type(any)
1098# 327 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1099#endif
1100 integer, intent(in) :: ai, term_index
1101 real(wp), intent(in) :: sim_time, c, sum_bb
1102 real(wp), intent(in) :: frequency_local, gauss_sigma_time_local
1103 real(wp), intent(out) :: source
1104 real(wp) :: omega !< angular frequency
1105 real(wp) :: sine_wave !< sine function for square wave
1106 real(wp) :: foc_length_factor !< Scale amplitude with radius for spherical support
1107 ! i.e. Spherical support -> 1/r scaling; Cylindrical support -> 1/sqrt(r) [empirical correction: ^-0.5 -> ^-0.85]
1108 integer, parameter :: mass_label = 1
1109
1110 ! An unfocused source leaves foc_length at its unset default (dflt_real < 0, or 0): apply no
1111 ! radial amplitude scaling. This guards every branch against the silent-corruption the default
1112 ! would otherwise cause - foc_length(ai)**(-0.85) = pow(negative, non-integer) = NaN on some
1113 ! compilers (e.g. nvhpc <= 24.3) in 2D, and 1/foc_length = wrong-sign or Inf in 3D/cylindrical.
1114 if (n == 0 .or. foc_length(ai) <= 0._wp) then
1115 foc_length_factor = 1._wp
1116 else if (p == 0 .and. (.not. cyl_coord)) then ! 2D Cartesian: cylindrical (line-source) spreading
1117 foc_length_factor = foc_length(ai)**(-0.85_wp) ! Empirical correction to 1/sqrt(r)
1118 else ! 3D or axisymmetric: spherical spreading
1119 foc_length_factor = 1._wp/foc_length(ai)
1120 end if
1121
1122 source = 0._wp
1123
1124 ! Temporal waveform: sine, Gaussian pulse, square wave, or broadband
1125 if (pulse(ai) == 1) then ! Sine wave
1126 if ((sim_time - delay(ai))*frequency_local > npulse(ai)) return
1127
1128 omega = 2._wp*pi*frequency_local
1129 source = mag(ai)*sin((sim_time - delay(ai))*omega)
1130
1131 if (term_index == mass_label) then
1132 source = source/c + foc_length_factor*mag(ai)*(cos((sim_time - delay(ai))*omega) - 1._wp)/omega
1133 end if
1134 else if (pulse(ai) == 2) then ! Gaussian pulse
1135 source = mag(ai)*exp(-0.5_wp*((sim_time - delay(ai))**2._wp)/(gauss_sigma_time_local**2._wp))
1136
1137 if (term_index == mass_label) then
1138 source = source/c - foc_length_factor*mag(ai)*sqrt(pi/2)*gauss_sigma_time_local*(erf((sim_time - delay(ai)) &
1139 & /(sqrt(2._wp)*gauss_sigma_time_local)) + 1)
1140 end if
1141 else if (pulse(ai) == 3) then ! Square wave
1142 if ((sim_time - delay(ai))*frequency_local > npulse(ai)) return
1143
1144 omega = 2._wp*pi*frequency_local
1145 sine_wave = sin((sim_time - delay(ai))*omega)
1146 source = mag(ai)*sign(1._wp, sine_wave)
1147
1148 ! Prevent max-norm differences due to compilers to pass CI
1149 if (abs(sine_wave) < 1.e-2_wp) then
1150 source = mag(ai)*sine_wave*1.e2_wp
1151 end if
1152 else if (pulse(ai) == 4) then ! Broadband wave
1153 source = sum_bb
1154 end if
1155
1156 end subroutine s_source_temporal
1157
1158 !> Pre-compute non-zero spatial source weights before time-stepping
1160
1161 integer :: j, k, l, ai
1162 integer :: count
1163 integer :: dim
1164 real(wp) :: source_spatial, angle, xyz_to_r_ratios(3)
1165 real(wp), parameter :: threshold = 1.e-10_wp
1166
1167 if (n == 0) then
1168 dim = 1
1169 else if (p == 0) then
1170 dim = 2
1171 else
1172 dim = 3
1173 end if
1174
1175#ifdef MFC_DEBUG
1176# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1177 block
1178# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1179 use iso_fortran_env, only: output_unit
1180# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1181
1182# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1183 print *, 'm_acoustic_src.fpp:403: ', '@:ALLOCATE(source_spatials_num_points(1:num_source))'
1184# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1185
1186# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1187 call flush (output_unit)
1188# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1189 end block
1190# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1191#endif
1192# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1193 allocate (source_spatials_num_points(1:num_source))
1194# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1195
1196# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1197
1198# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1199#if defined(MFC_OpenACC)
1200# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1201!$acc enter data create(source_spatials_num_points)
1202# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1203#elif defined(MFC_OpenMP)
1204# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1205!$omp target enter data map(always,alloc:source_spatials_num_points)
1206# 403 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1207#endif
1208#ifdef MFC_DEBUG
1209# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1210 block
1211# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1212 use iso_fortran_env, only: output_unit
1213# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1214
1215# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1216 print *, 'm_acoustic_src.fpp:404: ', '@:ALLOCATE(source_spatials(1:num_source))'
1217# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1218
1219# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1220 call flush (output_unit)
1221# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1222 end block
1223# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1224#endif
1225# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1226 allocate (source_spatials(1:num_source))
1227# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1228
1229# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1230
1231# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1232#if defined(MFC_OpenACC)
1233# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1234!$acc enter data create(source_spatials)
1235# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1236#elif defined(MFC_OpenMP)
1237# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1238!$omp target enter data map(always,alloc:source_spatials)
1239# 404 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1240#endif
1241
1242 do ai = 1, num_source
1243 ! First pass: Count the number of points for each source
1244 count = 0
1245 do l = 0, p
1246 do k = 0, n
1247 do j = 0, m
1248 call s_source_spatial(j, k, l, loc_acoustic(:,ai), ai, source_spatial, angle, xyz_to_r_ratios)
1249 if (abs(source_spatial) < threshold) cycle
1250 count = count + 1
1251 end do
1252 end do
1253 end do
1254 source_spatials_num_points(ai) = count
1255
1256 ! Allocate arrays with the correct size
1257
1258#ifdef MFC_DEBUG
1259# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1260 block
1261# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1262 use iso_fortran_env, only: output_unit
1263# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1264
1265# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1266 print *, 'm_acoustic_src.fpp:422: ', '@:ALLOCATE(source_spatials(ai)%coord(1:3, 1:count))'
1267# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1268
1269# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1270 call flush (output_unit)
1271# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1272 end block
1273# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1274#endif
1275# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1276 allocate (source_spatials(ai)%coord(1:3, 1:count))
1277# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1278
1279# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1280
1281# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1282#if defined(MFC_OpenACC)
1283# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1284!$acc enter data create(source_spatials(ai)%coord)
1285# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1286#elif defined(MFC_OpenMP)
1287# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1288!$omp target enter data map(always,alloc:source_spatials(ai)%coord)
1289# 422 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1290#endif
1291#ifdef MFC_DEBUG
1292# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1293 block
1294# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1295 use iso_fortran_env, only: output_unit
1296# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1297
1298# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1299 print *, 'm_acoustic_src.fpp:423: ', '@:ALLOCATE(source_spatials(ai)%val(1:count))'
1300# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1301
1302# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1303 call flush (output_unit)
1304# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1305 end block
1306# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1307#endif
1308# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1309 allocate (source_spatials(ai)%val(1:count))
1310# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1311
1312# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1313
1314# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1315#if defined(MFC_OpenACC)
1316# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1317!$acc enter data create(source_spatials(ai)%val)
1318# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1319#elif defined(MFC_OpenMP)
1320# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1321!$omp target enter data map(always,alloc:source_spatials(ai)%val)
1322# 423 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1323#endif
1324#ifdef MFC_DEBUG
1325# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1326 block
1327# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1328 use iso_fortran_env, only: output_unit
1329# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1330
1331# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1332 print *, 'm_acoustic_src.fpp:424: ', '@:ALLOCATE(source_spatials(ai)%angle(1:count))'
1333# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1334
1335# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1336 call flush (output_unit)
1337# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1338 end block
1339# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1340#endif
1341# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1342 allocate (source_spatials(ai)%angle(1:count))
1343# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1344
1345# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1346
1347# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1348#if defined(MFC_OpenACC)
1349# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1350!$acc enter data create(source_spatials(ai)%angle)
1351# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1352#elif defined(MFC_OpenMP)
1353# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1354!$omp target enter data map(always,alloc:source_spatials(ai)%angle)
1355# 424 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1356#endif
1357#ifdef MFC_DEBUG
1358# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1359 block
1360# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1361 use iso_fortran_env, only: output_unit
1362# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1363
1364# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1365 print *, 'm_acoustic_src.fpp:425: ', '@:ALLOCATE(source_spatials(ai)%xyz_to_r_ratios(1:3, 1:count))'
1366# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1367
1368# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1369 call flush (output_unit)
1370# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1371 end block
1372# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1373#endif
1374# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1375 allocate (source_spatials(ai)%xyz_to_r_ratios(1:3, 1:count))
1376# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1377
1378# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1379
1380# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1381#if defined(MFC_OpenACC)
1382# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1383!$acc enter data create(source_spatials(ai)%xyz_to_r_ratios)
1384# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1385#elif defined(MFC_OpenMP)
1386# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1387!$omp target enter data map(always,alloc:source_spatials(ai)%xyz_to_r_ratios)
1388# 425 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1389#endif
1390
1391#ifdef _CRAYFTN
1392# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1393 block
1394# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1395#ifdef MFC_DEBUG
1396# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1397 block
1398# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1399 use iso_fortran_env, only: output_unit
1400# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1401
1402# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1403 print *, 'm_acoustic_src.fpp:427: ', '@:ACC_SETUP_source_spatials(source_spatials(ai))'
1404# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1405
1406# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1407 call flush (output_unit)
1408# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1409 end block
1410# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1411#endif
1412# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1413
1414# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1415
1416# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1417#if defined(MFC_OpenACC)
1418# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1419!$acc enter data copyin(source_spatials(ai))
1420# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1421#elif defined(MFC_OpenMP)
1422# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1423!$omp target enter data map(to:source_spatials(ai))
1424# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1425#endif
1426# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1427 if (associated(source_spatials(ai)%coord)) then
1428# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1429
1430# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1431#if defined(MFC_OpenACC)
1432# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1433!$acc enter data copyin(source_spatials(ai)%coord)
1434# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1435#elif defined(MFC_OpenMP)
1436# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1437!$omp target enter data map(to:source_spatials(ai)%coord)
1438# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1439#endif
1440# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1441 end if
1442# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1443 if (associated(source_spatials(ai)%val)) then
1444# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1445
1446# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1447#if defined(MFC_OpenACC)
1448# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1449!$acc enter data copyin(source_spatials(ai)%val)
1450# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1451#elif defined(MFC_OpenMP)
1452# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1453!$omp target enter data map(to:source_spatials(ai)%val)
1454# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1455#endif
1456# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1457 end if
1458# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1459 if (associated(source_spatials(ai)%angle)) then
1460# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1461
1462# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1463#if defined(MFC_OpenACC)
1464# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1465!$acc enter data copyin(source_spatials(ai)%angle)
1466# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1467#elif defined(MFC_OpenMP)
1468# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1469!$omp target enter data map(to:source_spatials(ai)%angle)
1470# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1471#endif
1472# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1473 end if
1474# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1475 if (associated(source_spatials(ai)%xyz_to_r_ratios)) then
1476# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1477
1478# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1479#if defined(MFC_OpenACC)
1480# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1481!$acc enter data copyin(source_spatials(ai)%xyz_to_r_ratios)
1482# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1483#elif defined(MFC_OpenMP)
1484# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1485!$omp target enter data map(to:source_spatials(ai)%xyz_to_r_ratios)
1486# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1487#endif
1488# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1489 end if
1490# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1491 end block
1492# 427 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1493#endif
1494
1495 ! Second pass: Store the values
1496 count = 0 ! Reset counter
1497 do l = 0, p
1498 do k = 0, n
1499 do j = 0, m
1500 call s_source_spatial(j, k, l, loc_acoustic(:,ai), ai, source_spatial, angle, xyz_to_r_ratios)
1501 if (abs(source_spatial) < threshold) cycle
1502 count = count + 1
1503 source_spatials(ai)%coord(1, count) = j
1504 source_spatials(ai)%coord(2, count) = k
1505 source_spatials(ai)%coord(3, count) = l
1506 source_spatials(ai)%val(count) = source_spatial
1507 if (support(ai) >= 5) then
1508 if (dim == 2) source_spatials(ai)%angle(count) = angle
1509 if (dim == 3) source_spatials(ai)%xyz_to_r_ratios(1:3,count) = xyz_to_r_ratios
1510 end if
1511 end do
1512 end do
1513 end do
1514
1515 if (source_spatials_num_points(ai) /= count) then
1516 call s_mpi_abort('Fatal Error: Inconsistent allocation of source_spatials')
1517 end if
1518
1519 if (count > 0) then
1520
1521# 454 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1522#if defined(MFC_OpenACC)
1523# 454 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1524!$acc update device(source_spatials(ai)%coord)
1525# 454 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1526#elif defined(MFC_OpenMP)
1527# 454 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1528!$omp target update to(source_spatials(ai)%coord)
1529# 454 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1530#endif
1531
1532# 455 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1533#if defined(MFC_OpenACC)
1534# 455 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1535!$acc update device(source_spatials(ai)%val)
1536# 455 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1537#elif defined(MFC_OpenMP)
1538# 455 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1539!$omp target update to(source_spatials(ai)%val)
1540# 455 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1541#endif
1542 if (support(ai) >= 5) then
1543 if (dim == 2) then
1544
1545# 458 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1546#if defined(MFC_OpenACC)
1547# 458 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1548!$acc update device(source_spatials(ai)%angle)
1549# 458 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1550#elif defined(MFC_OpenMP)
1551# 458 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1552!$omp target update to(source_spatials(ai)%angle)
1553# 458 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1554#endif
1555 end if
1556 if (dim == 3) then
1557
1558# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1559#if defined(MFC_OpenACC)
1560# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1561!$acc update device(source_spatials(ai)%xyz_to_r_ratios)
1562# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1563#elif defined(MFC_OpenMP)
1564# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1565!$omp target update to(source_spatials(ai)%xyz_to_r_ratios)
1566# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1567#endif
1568 end if
1569 end if
1570 end if
1571 end do
1572
1573#ifdef MFC_DEBUG
1574 do ai = 1, num_source
1575 write (*, '(A,I2,A,I8,A)') 'Acoustic source ', ai, ' has ', source_spatials_num_points(ai), &
1576 & ' grid points with non-zero source term'
1577 end do
1578#endif
1579
1581
1582 !> Compute the spatial support of the acoustic source
1583 subroutine s_source_spatial(j, k, l, loc, ai, source, angle, xyz_to_r_ratios)
1584
1585 integer, intent(in) :: j, k, l, ai
1586 real(wp), dimension(3), intent(in) :: loc
1587 real(wp), intent(out) :: source, angle, xyz_to_r_ratios(3)
1588 real(wp) :: sig, r(3)
1589
1590 ! Calculate sig spatial support width
1591
1592 if (n == 0) then
1593 sig = dx(j)
1594 else if (p == 0) then
1595 sig = maxval((/dx(j), dy(k)/))
1596 else
1597 sig = maxval((/dx(j), dy(k), dz(l)/))
1598 end if
1599 sig = sig*acoustic_spatial_support_width
1600
1601 ! Calculate displacement from acoustic source location
1602 r(1) = x_cc(j) - loc(1)
1603 if (n /= 0) r(2) = y_cc(k) - loc(2)
1604 if (p /= 0) r(3) = z_cc(l) - loc(3)
1605
1606 if (any(support(ai) == (/1, 2, 3, 4/))) then
1607 call s_source_spatial_planar(ai, sig, r, source)
1608 else if (any(support(ai) == (/5, 6, 7/))) then
1609 call s_source_spatial_transducer(ai, sig, r, source, angle, xyz_to_r_ratios)
1610 else if (any(support(ai) == (/9, 10, 11/))) then
1611 call s_source_spatial_transducer_array(ai, sig, r, source, angle, xyz_to_r_ratios)
1612 end if
1613
1614 end subroutine s_source_spatial
1615
1616 !> Compute the spatial support for planar acoustic sources in 1D, 2D, and 3D
1617 subroutine s_source_spatial_planar(ai, sig, r, source)
1618
1619 integer, intent(in) :: ai
1620 real(wp), intent(in) :: sig, r(3)
1621 real(wp), intent(out) :: source
1622 real(wp) :: dist
1623
1624 source = 0._wp
1625
1626 ! Gaussian spatial pulse profile: exp(-0.5 * (d / sigma)^2) / (sqrt(2*pi) * sigma)
1627 if (support(ai) == 1) then ! 1D
1628 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(r(1)/(sig/2._wp))**2._wp)
1629 else if (support(ai) == 2 .or. support(ai) == 3) then ! 2D or 3D
1630 ! If we let unit vector e = (cos(dir), sin(dir)),
1631 dist = r(1)*cos(dir(ai)) + r(2)*sin(dir(ai)) ! dot(r,e)
1632 ! |r - dist*e| < length/2
1633 if ((r(1) - dist*cos(dir(ai)))**2._wp + (r(2) - dist*sin(dir(ai)))**2._wp < 0.25_wp*length(ai)**2._wp) then
1634 if (support(ai) /= 3 .or. abs(r(3)) < 0.25_wp*height(ai)) then ! additional height constraint for 3D
1635 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)
1636 end if
1637 end if
1638 end if
1639
1640 end subroutine s_source_spatial_planar
1641
1642 !> Compute the spatial support for a single transducer in 2D, 2D axisymmetric, and 3D
1643 subroutine s_source_spatial_transducer(ai, sig, r, source, angle, xyz_to_r_ratios)
1644
1645 integer, intent(in) :: ai
1646 real(wp), intent(in) :: sig, r(3)
1647 real(wp), intent(out) :: source, angle, xyz_to_r_ratios(3)
1648 real(wp) :: current_angle, angle_half_aperture, dist, norm
1649
1650 source = 0._wp ! If not affected by transducer
1651 angle = 0._wp
1652 xyz_to_r_ratios = 0._wp
1653
1654 if (support(ai) == 5 .or. support(ai) == 6) then ! 2D or 2D axisymmetric
1655 current_angle = -atan(r(2)/(foc_length(ai) - r(1)))
1656 angle_half_aperture = asin((aperture(ai)/2._wp)/(foc_length(ai)))
1657
1658 if (abs(current_angle) < angle_half_aperture .and. r(1) < foc_length(ai)) then
1659 dist = foc_length(ai) - sqrt(r(2)**2._wp + (foc_length(ai) - r(1))**2._wp)
1660 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)
1661 angle = -atan(r(2)/(foc_length(ai) - r(1)))
1662 end if
1663 else if (support(ai) == 7) then ! 3D
1664 current_angle = -atan(sqrt(r(2)**2 + r(3)**2)/(foc_length(ai) - r(1)))
1665 angle_half_aperture = asin((aperture(ai)/2._wp)/(foc_length(ai)))
1666
1667 if (abs(current_angle) < angle_half_aperture .and. r(1) < foc_length(ai)) then
1668 dist = foc_length(ai) - sqrt(r(2)**2._wp + r(3)**2._wp + (foc_length(ai) - r(1))**2._wp)
1669 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)
1670
1671 norm = sqrt(r(2)**2._wp + r(3)**2._wp + (foc_length(ai) - r(1))**2._wp)
1672 xyz_to_r_ratios(1) = -(r(1) - foc_length(ai))/norm
1673 xyz_to_r_ratios(2) = -r(2)/norm
1674 xyz_to_r_ratios(3) = -r(3)/norm
1675 end if
1676 end if
1677
1678 end subroutine s_source_spatial_transducer
1679
1680 !> Compute the spatial support for multiple transducers in 2D, 2D axisymmetric, and 3D
1681 subroutine s_source_spatial_transducer_array(ai, sig, r, source, angle, xyz_to_r_ratios)
1682
1683 integer, intent(in) :: ai
1684 real(wp), intent(in) :: sig, r(3)
1685 real(wp), intent(out) :: source, angle, xyz_to_r_ratios(3)
1686 integer :: elem, elem_min, elem_max
1687 real(wp) :: current_angle, angle_half_aperture, angle_per_elem, dist
1688 real(wp) :: angle_min, angle_max, norm
1689 real(wp) :: poly_side_length, aperture_element_3D, angle_elem
1690 real(wp) :: x2, y2, z2, x3, y3, z3, C, f, half_apert, dist_interp_to_elem_center
1691
1692 if (element_on(ai) == 0) then ! Full transducer
1693 elem_min = 1
1694 elem_max = num_elements(ai)
1695 else ! Transducer element specified
1696 elem_min = element_on(ai)
1697 elem_max = element_on(ai)
1698 end if
1699
1700 source = 0._wp ! If not affected by any transducer element
1701 angle = 0._wp
1702 xyz_to_r_ratios = 0._wp
1703
1704 if (support(ai) == 9 .or. support(ai) == 10) then ! 2D or 2D axisymmetric
1705 current_angle = -atan(r(2)/(foc_length(ai) - r(1)))
1706 angle_half_aperture = asin((aperture(ai)/2._wp)/(foc_length(ai)))
1707 angle_per_elem = (2._wp*angle_half_aperture - (num_elements(ai) - 1._wp)*element_spacing_angle(ai))/num_elements(ai)
1708 dist = foc_length(ai) - sqrt(r(2)**2._wp + (foc_length(ai) - r(1))**2._wp)
1709
1710 do elem = elem_min, elem_max
1711 angle_max = angle_half_aperture - (element_spacing_angle(ai) + angle_per_elem)*(elem - 1._wp)
1712 angle_min = angle_max - angle_per_elem
1713
1714 if (current_angle > angle_min .and. current_angle < angle_max .and. r(1) < foc_length(ai)) then
1715 source = exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)/(sqrt(2._wp*pi)*sig/2._wp)
1716 angle = current_angle
1717 exit ! Assume elements don't overlap
1718 end if
1719 end do
1720 else if (support(ai) == 11) then ! 3D
1721 poly_side_length = aperture(ai)*sin(pi/num_elements(ai))
1722 aperture_element_3d = poly_side_length*element_polygon_ratio(ai)
1723 f = foc_length(ai)
1724 half_apert = aperture(ai)/2._wp
1725
1726 do elem = elem_min, elem_max
1727 angle_elem = 2._wp*pi*real(elem, wp)/real(num_elements(ai), wp) + rotate_angle(ai)
1728
1729 ! Point 2 is the elem center
1730 x2 = f - sqrt(f**2 - half_apert**2)
1731 y2 = half_apert*cos(angle_elem)
1732 z2 = half_apert*sin(angle_elem)
1733
1734 ! Construct a plane normal to the line from the focal point to the elem center, Point 3 is the intercept of the
1735 ! plane and the line from the focal point to the current location
1736 c = f**2._wp/((r(1) - f)*(x2 - f) + r(2)*y2 + r(3)*z2) ! Constant for intermediate step
1737 x3 = c*(r(1) - f) + f
1738 y3 = c*r(2)
1739 z3 = c*r(3)
1740
1741 dist_interp_to_elem_center = sqrt((x2 - x3)**2._wp + (y2 - y3)**2._wp + (z2 - z3)**2._wp)
1742 if ((dist_interp_to_elem_center < aperture_element_3d/2._wp) .and. (r(1) < f)) then
1743 dist = sqrt((x3 - r(1))**2._wp + (y3 - r(2))**2._wp + (z3 - r(3))**2._wp)
1744 source = exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)/(sqrt(2._wp*pi)*sig/2._wp)
1745
1746 norm = sqrt(r(2)**2._wp + r(3)**2._wp + (f - r(1))**2._wp)
1747 xyz_to_r_ratios(1) = -(r(1) - f)/norm
1748 xyz_to_r_ratios(2) = -r(2)/norm
1749 xyz_to_r_ratios(3) = -r(3)/norm
1750 end if
1751 end do
1752 end if
1753
1755
1756 !> Convert wavelength to frequency
1757 elemental function f_frequency_local(freq_conv_flag, ai, c)
1758
1759
1760# 653 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1761#if MFC_OpenACC
1762# 653 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1763!$acc routine seq
1764# 653 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1765#elif MFC_OpenMP
1766# 653 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1767
1768# 653 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1769
1770# 653 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1771!$omp declare target device_type(any)
1772# 653 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1773#endif
1774 logical, intent(in) :: freq_conv_flag
1775 integer, intent(in) :: ai
1776 real(wp), intent(in) :: c
1777 real(wp) :: f_frequency_local
1778
1779 if (freq_conv_flag) then
1781 else
1783 end if
1784
1785 end function f_frequency_local
1786
1787 !> Convert Gaussian sigma from distance to time
1788 function f_gauss_sigma_time_local(gauss_conv_flag, ai, c)
1789
1790
1791# 670 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1792#if MFC_OpenACC
1793# 670 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1794!$acc routine seq
1795# 670 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1796#elif MFC_OpenMP
1797# 670 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1798
1799# 670 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1800
1801# 670 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1802!$omp declare target device_type(any)
1803# 670 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1804#endif
1805 logical, intent(in) :: gauss_conv_flag
1806 integer, intent(in) :: ai
1807 real(wp), intent(in) :: c
1808 real(wp) :: f_gauss_sigma_time_local
1809
1810 if (gauss_conv_flag) then
1812 else
1814 end if
1815
1816 end function f_gauss_sigma_time_local
1817
1818end module m_acoustic_src
type(scalar_field), dimension(sys_size), intent(inout) q_cons_vf
integer, intent(in) k
integer, intent(in) j
integer, intent(in) l
One-way acoustic source injection, Maeda and Colonius JCP (2017).
real(wp), dimension(:), allocatable gauss_sigma_time
elemental real(wp) function f_frequency_local(freq_conv_flag, ai, c)
Convert wavelength to frequency.
type(source_spatial_type), dimension(:), allocatable source_spatials
Data of non-zero source grid points for each source.
subroutine s_source_spatial_transducer(ai, sig, r, source, angle, xyz_to_r_ratios)
Compute the spatial support for a single transducer in 2D, 2D axisymmetric, and 3D.
real(wp), dimension(:), allocatable height
integer, dimension(:), allocatable pulse
integer, dimension(:), allocatable bb_num_freq
real(wp), dimension(:), allocatable dir
impure subroutine, public s_precalculate_acoustic_spatial_sources
Pre-compute non-zero spatial source weights before time-stepping.
real(wp), dimension(:), allocatable length
real(wp), dimension(:), allocatable rotate_angle
integer, dimension(:), allocatable support
real(wp), dimension(:), allocatable wavelength
impure subroutine, public s_acoustic_src_calculations(q_cons_vf, q_prim_vf, rhs_vf)
Compute mass, momentum, and energy acoustic source terms and add to the RHS.
real(wp), dimension(:), allocatable npulse
integer, dimension(:), allocatable element_on
real(wp), dimension(:), allocatable aperture
real(wp), dimension(:), allocatable foc_length
subroutine s_source_spatial_planar(ai, sig, r, source)
Compute the spatial support for planar acoustic sources in 1D, 2D, and 3D.
subroutine s_source_spatial_transducer_array(ai, sig, r, source, angle, xyz_to_r_ratios)
Compute the spatial support for multiple transducers in 2D, 2D axisymmetric, and 3D.
real(wp), dimension(:), allocatable bb_lowest_freq
integer, dimension(:), allocatable num_elements
impure subroutine, public s_initialize_acoustic_src
Initialize the acoustic source module.
real(wp), dimension(:,:), allocatable, target loc_acoustic
logical, dimension(:), allocatable dipole
real(wp), dimension(:,:,:), allocatable mass_src
real(wp), dimension(:), allocatable mag
real(wp), dimension(:), allocatable frequency
elemental subroutine s_source_temporal(sim_time, c, ai, term_index, frequency_local, gauss_sigma_time_local, source, sum_bb)
Compute the temporally varying amplitude of the pulse.
real(wp) function f_gauss_sigma_time_local(gauss_conv_flag, ai, c)
Convert Gaussian sigma from distance to time.
real(wp), dimension(:,:,:), allocatable e_src
subroutine s_source_spatial(j, k, l, loc, ai, source, angle, xyz_to_r_ratios)
Compute the spatial support of the acoustic source.
real(wp), dimension(:), allocatable element_spacing_angle
integer, dimension(:), allocatable source_spatials_num_points
Number of non-zero source grid points for each source.
real(wp), dimension(:), allocatable element_polygon_ratio
real(wp), dimension(:), allocatable bb_bandwidth
real(wp), dimension(:), allocatable gauss_sigma_dist
real(wp), dimension(:), allocatable delay
real(wp), dimension(:,:,:,:), allocatable mom_src
Bubble-dynamics procedures for ensemble- and volume-averaged model.
Compile-time constant parameters: default values, tolerances, and physical constants.
integer, parameter dflt_int
Default integer value.
Shared derived types for field data, patch geometry, bubble dynamics, and MPI I/O structures.
Global parameters for the computational domain, fluid properties, and simulation algorithm configurat...
Basic floating-point utilities: approximate equality, default detection, and coordinate bounds.
logical elemental function, public f_is_default(var)
Check if a real(wp) variable is of default value.
Conservative-to-primitive variable conversion, mixture property evaluation, and pressure computation.
Acoustic source source_spatial pre-calculated values.