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# 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
33# 167 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
34! New line at end of file is required for FYPP
35# 2 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
36# 1 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 1
37# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
38# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
39# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
40# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
41# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
42# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
43
44# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
45# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
46# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
47
48# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
49
50# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
51
52# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
53
54# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
55
56# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
57
58# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
59
60# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
61
62# 167 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
63! New line at end of file is required for FYPP
64# 2 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 2
65
66# 4 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
67# 5 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
68# 6 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
69# 7 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
70# 8 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
71
72# 20 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
73
74# 43 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
75
76# 48 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
77
78# 53 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
79
80# 58 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
81
82# 63 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
83
84# 68 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
85
86# 76 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
87
88# 81 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
89
90# 86 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
91
92# 91 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
93
94# 96 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
95
96# 101 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
97
98# 106 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
99
100# 111 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
101
102# 116 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
103
104# 121 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
105
106# 151 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
107
108# 192 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
109
110# 206 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
111
112# 231 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
113
114# 242 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
115
116# 244 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
117# 255 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
118
119# 284 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
120
121# 294 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
122
123# 304 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
124
125# 313 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
126
127# 330 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
128
129# 340 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
130
131# 347 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
132
133# 353 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
134
135# 359 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
136
137# 365 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
138
139# 371 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
140
141# 377 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
142! New line at end of file is required for FYPP
143# 3 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
144# 1 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 1
145# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
146# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
147# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
148# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
149# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
150# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
151
152# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
153# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
154# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
155
156# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
157
158# 46 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
159
160# 58 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
161
162# 68 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
163
164# 98 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
165
166# 110 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
167
168# 120 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
169
170# 167 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
171! New line at end of file is required for FYPP
172# 2 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 2
173
174# 7 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
175
176# 17 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
177
178# 22 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
179
180# 27 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
181
182# 32 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
183
184# 37 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
185
186# 42 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
187
188# 47 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
189
190# 52 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
191
192# 57 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
193
194# 62 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
195
196# 73 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
197
198# 78 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
199
200# 83 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
201
202# 88 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
203
204# 103 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
205
206# 131 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
207
208# 160 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
209
210# 175 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
211
212# 193 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
213
214# 215 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
215
216# 244 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
217
218# 259 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
219
220# 269 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
221
222# 278 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
223
224# 294 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
225
226# 304 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
227
228# 311 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
229! New line at end of file is required for FYPP
230# 4 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
231
232! GPU parallel region (scalar reductions, maxval/minval)
233# 23 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
234
235! GPU parallel loop over threads (most common GPU macro)
236# 43 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
237
238! Required closing for GPU_PARALLEL_LOOP
239# 55 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
240
241! Mark routine for device compilation
242# 112 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
243
244! Declare device-resident data
245# 130 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
246
247! Inner loop within a GPU parallel region
248# 145 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
249
250! Scoped GPU data region
251# 164 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
252
253! Host code with device pointers (for MPI with GPU buffers)
254# 193 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
255
256! Allocate device memory (unscoped)
257# 207 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
258
259! Free device memory
260# 219 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
261
262! Atomic operation on device
263# 231 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
264
265! End atomic capture block
266# 242 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
267
268! Copy data between host and device
269# 254 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
270
271! Synchronization barrier
272# 266 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
273
274! Import GPU library module (openacc or omp_lib)
275# 275 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
276
277! Emit code only for AMD compiler
278# 282 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
279
280! Emit code for non-Cray compilers
281# 289 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
282
283! Emit code only for Cray compiler
284# 296 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
285
286! Emit code for non-NVIDIA compilers
287# 303 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
288
289# 305 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
290# 306 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
291! New line at end of file is required for FYPP
292# 2 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 2
293
294# 14 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
295
296! Caution: This macro requires the use of a binding script to set CUDA_VISIBLE_DEVICES, such that we have one GPU device per MPI
297! rank. That's because for both cudaMemAdvise (preferred location) and cudaMemPrefetchAsync we use location = device_id = 0. For an
298! example see misc/nvidia_uvm/bind.sh.
299# 55 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
300
301! Allocate and create GPU device memory
302# 75 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
303
304! Free GPU device memory and deallocate
305# 83 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
306
307! Cray-specific GPU pointer setup for vector fields
308# 107 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
309
310! Cray-specific GPU pointer setup for scalar fields
311# 123 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
312
313! Cray-specific GPU pointer setup for acoustic source spatials
314# 148 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
315
316# 154 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
317
318# 161 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
319! New line at end of file is required for FYPP
320# 6 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp" 2
321
322!> @brief One-way acoustic source injection, Maeda and Colonius JCP (2017)
324
327 use m_bubbles
330 use m_constants
331
332 implicit none
333
335
336 integer, allocatable, dimension(:) :: pulse, support
337
338# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
339#if defined(MFC_OpenACC)
340# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
341!$acc declare create(pulse, support)
342# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
343#elif defined(MFC_OpenMP)
344# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
345!$omp declare target (pulse, support)
346# 22 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
347#endif
348
349 logical, allocatable, dimension(:) :: dipole
350
351# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
352#if defined(MFC_OpenACC)
353# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
354!$acc declare create(dipole)
355# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
356#elif defined(MFC_OpenMP)
357# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
358!$omp declare target (dipole)
359# 25 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
360#endif
361
362 real(wp), allocatable, target, dimension(:,:) :: loc_acoustic
363
364# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
365#if defined(MFC_OpenACC)
366# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
367!$acc declare create(loc_acoustic)
368# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
369#elif defined(MFC_OpenMP)
370# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
371!$omp declare target (loc_acoustic)
372# 28 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
373#endif
374
375 real(wp), allocatable, dimension(:) :: mag, length, height, wavelength, frequency
376 real(wp), allocatable, dimension(:) :: gauss_sigma_dist, gauss_sigma_time, npulse, dir, delay
377
378# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
379#if defined(MFC_OpenACC)
380# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
381!$acc declare create(mag, length, height, wavelength, frequency)
382# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
383#elif defined(MFC_OpenMP)
384# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
385!$omp declare target (mag, length, height, wavelength, frequency)
386# 32 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
387#endif
388
389# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
390#if defined(MFC_OpenACC)
391# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
392!$acc declare create(gauss_sigma_dist, gauss_sigma_time, npulse, dir, delay)
393# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
394#elif defined(MFC_OpenMP)
395# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
396!$omp declare target (gauss_sigma_dist, gauss_sigma_time, npulse, dir, delay)
397# 33 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
398#endif
399
400 real(wp), allocatable, dimension(:) :: foc_length, aperture
401
402# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
403#if defined(MFC_OpenACC)
404# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
405!$acc declare create(foc_length, aperture)
406# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
407#elif defined(MFC_OpenMP)
408# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
409!$omp declare target (foc_length, aperture)
410# 36 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
411#endif
412
413 real(wp), allocatable, dimension(:) :: element_spacing_angle, element_polygon_ratio, rotate_angle
414
415# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
416#if defined(MFC_OpenACC)
417# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
418!$acc declare create(element_spacing_angle, element_polygon_ratio, rotate_angle)
419# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
420#elif defined(MFC_OpenMP)
421# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
422!$omp declare target (element_spacing_angle, element_polygon_ratio, rotate_angle)
423# 39 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
424#endif
425
426 real(wp), allocatable, dimension(:) :: bb_bandwidth, bb_lowest_freq
427
428# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
429#if defined(MFC_OpenACC)
430# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
431!$acc declare create(bb_bandwidth, bb_lowest_freq)
432# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
433#elif defined(MFC_OpenMP)
434# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
435!$omp declare target (bb_bandwidth, bb_lowest_freq)
436# 42 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
437#endif
438
439 integer, allocatable, dimension(:) :: num_elements, element_on, bb_num_freq
440
441# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
442#if defined(MFC_OpenACC)
443# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
444!$acc declare create(num_elements, element_on, bb_num_freq)
445# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
446#elif defined(MFC_OpenMP)
447# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
448!$omp declare target (num_elements, element_on, bb_num_freq)
449# 45 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
450#endif
451
452 !> @name Acoustic source terms
453 !> @{
454 real(wp), allocatable, dimension(:,:,:) :: mass_src, e_src
455 real(wp), allocatable, dimension(:,:,:,:) :: mom_src
456 !> @}
457
458# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
459#if defined(MFC_OpenACC)
460# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
461!$acc declare create(mass_src, e_src, mom_src)
462# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
463#elif defined(MFC_OpenMP)
464# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
465!$omp declare target (mass_src, e_src, mom_src)
466# 52 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
467#endif
468
469 integer, dimension(:), allocatable :: source_spatials_num_points !< Number of non-zero source grid points for each source
470
471# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
472#if defined(MFC_OpenACC)
473# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
474!$acc declare create(source_spatials_num_points)
475# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
476#elif defined(MFC_OpenMP)
477# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
478!$omp declare target (source_spatials_num_points)
479# 55 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
480#endif
481
482 type(source_spatial_type), dimension(:), allocatable :: source_spatials !< Data of non-zero source grid points for each source
483
484# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
485#if defined(MFC_OpenACC)
486# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
487!$acc declare create(source_spatials)
488# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
489#elif defined(MFC_OpenMP)
490# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
491!$omp declare target (source_spatials)
492# 58 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
493#endif
494
495contains
496
497 !> Initialize the acoustic source module
499
500 integer :: i, j !< generic loop variables
501
502#ifdef MFC_DEBUG
503# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
504 block
505# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
506 use iso_fortran_env, only: output_unit
507# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
508
509# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
510 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))'
511# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
512
513# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
514 call flush (output_unit)
515# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
516 end block
517# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
518#endif
519# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
520 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))
521# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
522
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
527# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
528
529# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
530
531# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
532
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#if defined(MFC_OpenACC)
573# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
574!$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)
575# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
576#elif defined(MFC_OpenMP)
577# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
578!$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)
579# 67 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
580#endif
581# 73 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
582
583 do i = 1, num_source
584 do j = 1, 3
585 loc_acoustic(j, i) = acoustic(i)%loc(j)
586 end do
587 mag(i) = acoustic(i)%mag
588 dipole(i) = acoustic(i)%dipole
589 support(i) = acoustic(i)%support
590 length(i) = acoustic(i)%length
591 height(i) = acoustic(i)%height
592 wavelength(i) = acoustic(i)%wavelength
593 frequency(i) = acoustic(i)%frequency
594 gauss_sigma_dist(i) = acoustic(i)%gauss_sigma_dist
595 gauss_sigma_time(i) = acoustic(i)%gauss_sigma_time
596 foc_length(i) = acoustic(i)%foc_length
597 aperture(i) = acoustic(i)%aperture
598 npulse(i) = acoustic(i)%npulse
599 pulse(i) = acoustic(i)%pulse
600 dir(i) = acoustic(i)%dir
601 element_spacing_angle(i) = acoustic(i)%element_spacing_angle
602 element_polygon_ratio(i) = acoustic(i)%element_polygon_ratio
603 num_elements(i) = acoustic(i)%num_elements
604 bb_num_freq(i) = acoustic(i)%bb_num_freq
605 bb_bandwidth(i) = acoustic(i)%bb_bandwidth
606 bb_lowest_freq(i) = acoustic(i)%bb_lowest_freq
607
608 if (acoustic(i)%element_on == dflt_int) then
609 element_on(i) = 0
610 else
611 element_on(i) = acoustic(i)%element_on
612 end if
613 if (f_is_default(acoustic(i)%rotate_angle)) then
614 rotate_angle(i) = 0._wp
615 else
616 rotate_angle(i) = acoustic(i)%rotate_angle
617 end if
618 if (f_is_default(acoustic(i)%delay)) then ! m_checker guarantees acoustic(i)%delay is set for pulse = 2 (Gaussian)
619 delay(i) = 0._wp ! Defaults to zero for sine and square waves
620 else
621 delay(i) = acoustic(i)%delay
622 end if
623 end do
624
625# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
626#if defined(MFC_OpenACC)
627# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
628!$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)
629# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
630#elif defined(MFC_OpenMP)
631# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
632!$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)
633# 115 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
634#endif
635# 118 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
636
637#ifdef MFC_DEBUG
638# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
639 block
640# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
641 use iso_fortran_env, only: output_unit
642# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
643
644# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
645 print *, 'm_acoustic_src.fpp:119: ', '@:ALLOCATE(mass_src(0:m, 0:n, 0:p))'
646# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
647
648# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
649 call flush (output_unit)
650# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
651 end block
652# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
653#endif
654# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
655 allocate (mass_src(0:m, 0:n, 0:p))
656# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
657
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#if defined(MFC_OpenACC)
662# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
663!$acc enter data create(mass_src)
664# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
665#elif defined(MFC_OpenMP)
666# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
667!$omp target enter data map(always,alloc:mass_src)
668# 119 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
669#endif
670#ifdef MFC_DEBUG
671# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
672 block
673# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
674 use iso_fortran_env, only: output_unit
675# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
676
677# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
678 print *, 'm_acoustic_src.fpp:120: ', '@:ALLOCATE(mom_src(1:num_vels, 0:m, 0:n, 0:p))'
679# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
680
681# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
682 call flush (output_unit)
683# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
684 end block
685# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
686#endif
687# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
688 allocate (mom_src(1:num_vels, 0:m, 0:n, 0:p))
689# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
690
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#if defined(MFC_OpenACC)
695# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
696!$acc enter data create(mom_src)
697# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
698#elif defined(MFC_OpenMP)
699# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
700!$omp target enter data map(always,alloc:mom_src)
701# 120 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
702#endif
703#ifdef MFC_DEBUG
704# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
705 block
706# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
707 use iso_fortran_env, only: output_unit
708# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
709
710# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
711 print *, 'm_acoustic_src.fpp:121: ', '@:ALLOCATE(E_src(0:m, 0:n, 0:p))'
712# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
713
714# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
715 call flush (output_unit)
716# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
717 end block
718# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
719#endif
720# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
721 allocate (e_src(0:m, 0:n, 0:p))
722# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
723
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#if defined(MFC_OpenACC)
728# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
729!$acc enter data create(E_src)
730# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
731#elif defined(MFC_OpenMP)
732# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
733!$omp target enter data map(always,alloc:E_src)
734# 121 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
735#endif
736
737 end subroutine s_initialize_acoustic_src
738
739 !> Compute mass, momentum, and energy acoustic source terms and add to the RHS
740 impure subroutine s_acoustic_src_calculations(q_cons_vf, q_prim_vf, rhs_vf)
741
742 type(scalar_field), dimension(sys_size), intent(inout) :: q_cons_vf !< Conservative variables
743 type(scalar_field), dimension(sys_size), intent(inout) :: q_prim_vf !< Primitive variables
744 type(scalar_field), dimension(sys_size), intent(inout) :: rhs_vf
745
746# 135 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
747 real(wp), dimension(num_fluids) :: myalpha, myalpha_rho
748# 137 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
749 real(wp) :: myrho, b_tait
750 real(wp) :: sim_time, c, small_gamma
751 real(wp) :: frequency_local, gauss_sigma_time_local
752 real(wp) :: mass_src_diff, mom_src_diff
753 real(wp) :: source_temporal
754 real(wp) :: period_bb !< period of each sine wave in broadband source
755 real(wp) :: sl_bb !< spectral level at each frequency
756 real(wp) :: ffre_bb !< source term corresponding to each frequency
757 real(wp) :: sum_bb !< total source term for the broadband wave
758 real(wp), allocatable, dimension(:) :: phi_rn !< random phase shift for each frequency
759 integer :: i, j, k, l, q !< generic loop variables
760 integer :: ai !< acoustic source index
761 integer :: num_points
762 logical :: freq_conv_flag, gauss_conv_flag
763 integer, parameter :: mass_label = 1, mom_label = 2
764
765 sim_time = mytime ! Accumulated time, correct under adaptive dt
766
767
768# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
769
770# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
771#if defined(MFC_OpenACC)
772# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
773!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l)
774# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
775#elif defined(MFC_OpenMP)
776# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
777
778# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
779
780# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
781
782# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
783!$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)
784# 155 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
785#endif
786 do l = 0, p
787 do k = 0, n
788 do j = 0, m
789 mass_src(j, k, l) = 0._wp
790 mom_src(1, j, k, l) = 0._wp
791 e_src(j, k, l) = 0._wp
792 if (n > 0) mom_src(2, j, k, l) = 0._wp
793 if (p > 0) mom_src(3, j, k, l) = 0._wp
794 end do
795 end do
796 end do
797
798# 167 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
799#if defined(MFC_OpenACC)
800# 167 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
801!$acc end parallel loop
802# 167 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
803#elif defined(MFC_OpenMP)
804# 167 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
805
806# 167 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
807!$omp end target teams loop
808# 167 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
809#endif
810
811 ! Keep outer loop sequential because different sources can have very different number of points
812 do ai = 1, num_source
813 ! Skip if the pulse has not started yet for sine and square waves
814 if (.not. (sim_time < delay(ai) .and. (pulse(ai) == 1 .or. pulse(ai) == 3))) then
815 ! Decide if frequency need to be converted from wavelength
816 freq_conv_flag = f_is_default(frequency(ai))
817 gauss_conv_flag = f_is_default(gauss_sigma_time(ai))
818
819 num_points = source_spatials_num_points(ai) ! Use scalar to force firstprivate to prevent GPU bug
820
821 ! Calculate the broadband source
822 period_bb = 0._wp
823 sl_bb = 0._wp
824 ffre_bb = 0._wp
825 sum_bb = 0._wp
826
827 ! Allocate buffers for random phase shift
828 allocate (phi_rn(1:bb_num_freq(ai)))
829 phi_rn(1:bb_num_freq(ai)) = 0._wp
830
831 if (pulse(ai) == 4) then
832 call random_number(phi_rn(1:bb_num_freq(ai)))
833 ! Ensure all the ranks have the same random phase shift
834 call s_mpi_send_random_number(phi_rn, bb_num_freq(ai))
835 end if
836
837 do k = 1, bb_num_freq(ai)
838 ! Acoustic period of the wave at each discrete frequency
839 period_bb = 1._wp/(bb_lowest_freq(ai) + k*bb_bandwidth(ai))
840 ! Spectral level at each frequency
841 sl_bb = broadband_spectral_level_constant*mag(ai) + k*mag(ai)/broadband_spectral_level_growth_rate
842 ! Source term corresponding to each frequencies
843 ffre_bb = sqrt((2._wp*sl_bb*bb_bandwidth(ai)))*cos((sim_time)*2._wp*pi/period_bb + 2._wp*pi*phi_rn(k))
844 ! Sum up the source term of each frequency to obtain the total source term for broadband wave
845 sum_bb = sum_bb + ffre_bb
846 end do
847
848 deallocate (phi_rn)
849
850
851# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
852
853# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
854#if defined(MFC_OpenACC)
855# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
856!$acc parallel loop gang vector default(present) &
857# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
858!$acc& private(myalpha, myalpha_rho, myRho, B_tait, c, small_gamma, frequency_local, gauss_sigma_time_local, mass_src_diff, mom_src_diff, source_temporal, j, k, l, q) &
859# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
860!$acc& copyin(sum_BB, freq_conv_flag, gauss_conv_flag, sim_time)
861# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
862#elif defined(MFC_OpenMP)
863# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
864
865# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
866
867# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
868
869# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
870!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) &
871# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
872!$omp& private(myalpha, myalpha_rho, myRho, B_tait, c, small_gamma, frequency_local, gauss_sigma_time_local, mass_src_diff, mom_src_diff, source_temporal, j, k, l, q) &
873# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
874!$omp& map(to:sum_BB, freq_conv_flag, gauss_conv_flag, sim_time)
875# 208 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
876#endif
877# 211 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
878 do i = 1, num_points
879 j = source_spatials(ai)%coord(1, i)
880 k = source_spatials(ai)%coord(2, i)
881 l = source_spatials(ai)%coord(3, i)
882
883 ! Compute speed of sound
884 myrho = 0._wp
885 b_tait = 0._wp
886 small_gamma = 0._wp
887
888
889# 221 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
890#if defined(MFC_OpenACC)
891# 221 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
892!$acc loop seq
893# 221 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
894#elif defined(MFC_OpenMP)
895# 221 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
896
897# 221 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
898#endif
899 do q = 1, num_fluids
900 myalpha_rho(q) = q_cons_vf(q)%sf(j, k, l)
901 myalpha(q) = q_cons_vf(eqn_idx%adv%beg + q - 1)%sf(j, k, l)
902 end do
903
904 if (bubbles_euler) then
905 if (num_fluids > 2) then
906
907# 229 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
908#if defined(MFC_OpenACC)
909# 229 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
910!$acc loop seq
911# 229 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
912#elif defined(MFC_OpenMP)
913# 229 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
914
915# 229 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
916#endif
917 do q = 1, num_fluids - 1
918 myrho = myrho + myalpha_rho(q)
919 b_tait = b_tait + myalpha(q)*pi_infs(q)
920 small_gamma = small_gamma + myalpha(q)*gammas(q)
921 end do
922 else
923 myrho = myalpha_rho(1)
924 b_tait = pi_infs(1)
925 small_gamma = gammas(1)
926 end if
927 end if
928
929 if ((.not. bubbles_euler) .or. (mpp_lim .and. (num_fluids > 2))) then
930
931# 243 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
932#if defined(MFC_OpenACC)
933# 243 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
934!$acc loop seq
935# 243 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
936#elif defined(MFC_OpenMP)
937# 243 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
938
939# 243 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
940#endif
941 do q = 1, num_fluids
942 myrho = myrho + myalpha_rho(q)
943 b_tait = b_tait + myalpha(q)*pi_infs(q)
944 small_gamma = small_gamma + myalpha(q)*gammas(q)
945 end do
946 end if
947
948 small_gamma = 1._wp/small_gamma + 1._wp
949 c = sqrt(small_gamma*(q_prim_vf(eqn_idx%E)%sf(j, k, l) + ((small_gamma - 1._wp)/small_gamma)*b_tait)/myrho)
950
951 ! Wavelength to frequency conversion
952 if (pulse(ai) == 1 .or. pulse(ai) == 3) frequency_local = f_frequency_local(freq_conv_flag, ai, c)
953 if (pulse(ai) == 2) gauss_sigma_time_local = f_gauss_sigma_time_local(gauss_conv_flag, ai, c)
954
955 ! Update momentum source term
956 call s_source_temporal(sim_time, c, ai, mom_label, frequency_local, gauss_sigma_time_local, source_temporal, &
957 & sum_bb)
958 mom_src_diff = source_temporal*source_spatials(ai)%val(i)
959
960 if (dipole(ai)) then ! Double amplitude & No momentum source term (only works for Planar)
961 mass_src(j, k, l) = mass_src(j, k, l) + 2._wp*mom_src_diff/c
962 e_src(j, k, l) = e_src(j, k, l) + 2._wp*mom_src_diff*c/(small_gamma - 1._wp)
963 cycle
964 end if
965
966 if (n == 0) then ! 1D
967 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
968 else if (p == 0) then ! 2D
969 if (support(ai) < 5) then ! Planar
970 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*cos(dir(ai))
971 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*sin(dir(ai))
972 else
973 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*cos(source_spatials(ai)%angle(i))
974 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*sin(source_spatials(ai)%angle(i))
975 end if
976 else ! 3D
977 if (support(ai) < 5) then ! Planar
978 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*cos(dir(ai))
979 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*sin(dir(ai))
980 else
981 mom_src(1, j, k, l) = mom_src(1, j, k, l) + mom_src_diff*source_spatials(ai)%xyz_to_r_ratios(1, i)
982 mom_src(2, j, k, l) = mom_src(2, j, k, l) + mom_src_diff*source_spatials(ai)%xyz_to_r_ratios(2, i)
983 mom_src(3, j, k, l) = mom_src(3, j, k, l) + mom_src_diff*source_spatials(ai)%xyz_to_r_ratios(3, i)
984 end if
985 end if
986
987 ! Update mass source term
988 if (support(ai) < 5) then ! Planar
989 mass_src_diff = mom_src_diff/c
990 else ! Spherical or cylindrical support
991 ! Mass source term must be calculated differently using a correction term for spherical and cylindrical
992 ! support
993 call s_source_temporal(sim_time, c, ai, mass_label, frequency_local, gauss_sigma_time_local, &
994 & source_temporal, sum_bb)
995 mass_src_diff = source_temporal*source_spatials(ai)%val(i)
996 end if
997 mass_src(j, k, l) = mass_src(j, k, l) + mass_src_diff
998
999 ! Update energy source term
1000 e_src(j, k, l) = e_src(j, k, l) + mass_src_diff*c**2._wp/(small_gamma - 1._wp)
1001 end do
1002
1003# 305 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1004#if defined(MFC_OpenACC)
1005# 305 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1006!$acc end parallel loop
1007# 305 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1008#elif defined(MFC_OpenMP)
1009# 305 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1010
1011# 305 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1012!$omp end target teams loop
1013# 305 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1014#endif
1015 end if
1016 end do
1017
1018 ! Update the rhs variables
1019
1020# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1021
1022# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1023#if defined(MFC_OpenACC)
1024# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1025!$acc parallel loop collapse(3) gang vector default(present) private(j, k, l)
1026# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1027#elif defined(MFC_OpenMP)
1028# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1029
1030# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1031
1032# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1033
1034# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1035!$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)
1036# 310 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1037#endif
1038 do l = 0, p
1039 do k = 0, n
1040 do j = 0, m
1041
1042# 314 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1043#if defined(MFC_OpenACC)
1044# 314 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1045!$acc loop seq
1046# 314 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1047#elif defined(MFC_OpenMP)
1048# 314 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1049
1050# 314 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1051#endif
1052 do q = eqn_idx%cont%beg, eqn_idx%cont%end
1053 rhs_vf(q)%sf(j, k, l) = rhs_vf(q)%sf(j, k, l) + mass_src(j, k, l)
1054 end do
1055
1056# 318 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1057#if defined(MFC_OpenACC)
1058# 318 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1059!$acc loop seq
1060# 318 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1061#elif defined(MFC_OpenMP)
1062# 318 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1063
1064# 318 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1065#endif
1066 do q = eqn_idx%mom%beg, eqn_idx%mom%end
1067 rhs_vf(q)%sf(j, k, l) = rhs_vf(q)%sf(j, k, l) + mom_src(q - eqn_idx%cont%end, j, k, l)
1068 end do
1069 rhs_vf(eqn_idx%E)%sf(j, k, l) = rhs_vf(eqn_idx%E)%sf(j, k, l) + e_src(j, k, l)
1070 end do
1071 end do
1072 end do
1073
1074# 326 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1075#if defined(MFC_OpenACC)
1076# 326 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1077!$acc end parallel loop
1078# 326 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1079#elif defined(MFC_OpenMP)
1080# 326 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1081
1082# 326 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1083!$omp end target teams loop
1084# 326 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1085#endif
1086
1087 end subroutine s_acoustic_src_calculations
1088
1089 !> Compute the temporally varying amplitude of the pulse
1090 elemental subroutine s_source_temporal(sim_time, c, ai, term_index, frequency_local, gauss_sigma_time_local, source, sum_BB)
1091
1092
1093# 333 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1094#if MFC_OpenACC
1095# 333 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1096!$acc routine seq
1097# 333 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1098#elif MFC_OpenMP
1099# 333 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1100
1101# 333 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1102
1103# 333 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1104!$omp declare target device_type(any)
1105# 333 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1106#endif
1107 integer, intent(in) :: ai, term_index
1108 real(wp), intent(in) :: sim_time, c, sum_bb
1109 real(wp), intent(in) :: frequency_local, gauss_sigma_time_local
1110 real(wp), intent(out) :: source
1111 real(wp) :: omega !< angular frequency
1112 real(wp) :: sine_wave !< sine function for square wave
1113 real(wp) :: foc_length_factor !< Scale amplitude with radius for spherical support
1114 ! i.e. Spherical support -> 1/r scaling; Cylindrical support -> 1/sqrt(r) [empirical correction: ^-0.5 -> ^-0.85]
1115 integer, parameter :: mass_label = 1
1116
1117 ! An unfocused source leaves foc_length at its unset default (dflt_real < 0, or 0): apply no
1118 ! radial amplitude scaling. This guards every branch against the silent-corruption the default
1119 ! would otherwise cause - foc_length(ai)**(-0.85) = pow(negative, non-integer) = NaN on some
1120 ! compilers (e.g. nvhpc <= 24.3) in 2D, and 1/foc_length = wrong-sign or Inf in 3D/cylindrical.
1121 if (n == 0 .or. foc_length(ai) <= 0._wp) then
1122 foc_length_factor = 1._wp
1123 else if (p == 0 .and. (.not. cyl_coord)) then ! 2D Cartesian: cylindrical (line-source) spreading
1124 foc_length_factor = foc_length(ai)**(-0.85_wp) ! Empirical correction to 1/sqrt(r)
1125 else ! 3D or axisymmetric: spherical spreading
1126 foc_length_factor = 1._wp/foc_length(ai)
1127 end if
1128
1129 source = 0._wp
1130
1131 ! Temporal waveform: sine, Gaussian pulse, square wave, or broadband
1132 if (pulse(ai) == 1) then ! Sine wave
1133 if ((sim_time - delay(ai))*frequency_local > npulse(ai)) return
1134
1135 omega = 2._wp*pi*frequency_local
1136 source = mag(ai)*sin((sim_time - delay(ai))*omega)
1137
1138 if (term_index == mass_label) then
1139 source = source/c + foc_length_factor*mag(ai)*(cos((sim_time - delay(ai))*omega) - 1._wp)/omega
1140 end if
1141 else if (pulse(ai) == 2) then ! Gaussian pulse
1142 source = mag(ai)*exp(-0.5_wp*((sim_time - delay(ai))**2._wp)/(gauss_sigma_time_local**2._wp))
1143
1144 if (term_index == mass_label) then
1145 source = source/c - foc_length_factor*mag(ai)*sqrt(pi/2)*gauss_sigma_time_local*(erf((sim_time - delay(ai)) &
1146 & /(sqrt(2._wp)*gauss_sigma_time_local)) + 1)
1147 end if
1148 else if (pulse(ai) == 3) then ! Square wave
1149 if ((sim_time - delay(ai))*frequency_local > npulse(ai)) return
1150
1151 omega = 2._wp*pi*frequency_local
1152 sine_wave = sin((sim_time - delay(ai))*omega)
1153 source = mag(ai)*sign(1._wp, sine_wave)
1154
1155 ! Prevent max-norm differences due to compilers to pass CI
1156 if (abs(sine_wave) < 1.e-2_wp) then
1157 source = mag(ai)*sine_wave*1.e2_wp
1158 end if
1159 else if (pulse(ai) == 4) then ! Broadband wave
1160 source = sum_bb
1161 end if
1162
1163 end subroutine s_source_temporal
1164
1165 !> Pre-compute non-zero spatial source weights before time-stepping
1167
1168 integer :: j, k, l, ai
1169 integer :: count
1170 integer :: dim
1171 real(wp) :: source_spatial, angle, xyz_to_r_ratios(3)
1172 real(wp), parameter :: threshold = 1.e-10_wp
1173
1174 if (n == 0) then
1175 dim = 1
1176 else if (p == 0) then
1177 dim = 2
1178 else
1179 dim = 3
1180 end if
1181
1182#ifdef MFC_DEBUG
1183# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1184 block
1185# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1186 use iso_fortran_env, only: output_unit
1187# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1188
1189# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1190 print *, 'm_acoustic_src.fpp:409: ', '@:ALLOCATE(source_spatials_num_points(1:num_source))'
1191# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1192
1193# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1194 call flush (output_unit)
1195# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1196 end block
1197# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1198#endif
1199# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1200 allocate (source_spatials_num_points(1:num_source))
1201# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1202
1203# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1204
1205# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1206#if defined(MFC_OpenACC)
1207# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1208!$acc enter data create(source_spatials_num_points)
1209# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1210#elif defined(MFC_OpenMP)
1211# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1212!$omp target enter data map(always,alloc:source_spatials_num_points)
1213# 409 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1214#endif
1215#ifdef MFC_DEBUG
1216# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1217 block
1218# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1219 use iso_fortran_env, only: output_unit
1220# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1221
1222# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1223 print *, 'm_acoustic_src.fpp:410: ', '@:ALLOCATE(source_spatials(1:num_source))'
1224# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1225
1226# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1227 call flush (output_unit)
1228# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1229 end block
1230# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1231#endif
1232# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1233 allocate (source_spatials(1:num_source))
1234# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1235
1236# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1237
1238# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1239#if defined(MFC_OpenACC)
1240# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1241!$acc enter data create(source_spatials)
1242# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1243#elif defined(MFC_OpenMP)
1244# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1245!$omp target enter data map(always,alloc:source_spatials)
1246# 410 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1247#endif
1248
1249 do ai = 1, num_source
1250 ! First pass: Count the number of points for each source
1251 count = 0
1252 do l = 0, p
1253 do k = 0, n
1254 do j = 0, m
1255 call s_source_spatial(j, k, l, loc_acoustic(:,ai), ai, source_spatial, angle, xyz_to_r_ratios)
1256 if (abs(source_spatial) < threshold) cycle
1257 count = count + 1
1258 end do
1259 end do
1260 end do
1261 source_spatials_num_points(ai) = count
1262
1263 ! Allocate arrays with the correct size
1264
1265#ifdef MFC_DEBUG
1266# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1267 block
1268# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1269 use iso_fortran_env, only: output_unit
1270# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1271
1272# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1273 print *, 'm_acoustic_src.fpp:428: ', '@:ALLOCATE(source_spatials(ai)%coord(1:3, 1:count))'
1274# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1275
1276# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1277 call flush (output_unit)
1278# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1279 end block
1280# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1281#endif
1282# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1283 allocate (source_spatials(ai)%coord(1:3, 1:count))
1284# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1285
1286# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1287
1288# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1289#if defined(MFC_OpenACC)
1290# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1291!$acc enter data create(source_spatials(ai)%coord)
1292# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1293#elif defined(MFC_OpenMP)
1294# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1295!$omp target enter data map(always,alloc:source_spatials(ai)%coord)
1296# 428 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1297#endif
1298#ifdef MFC_DEBUG
1299# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1300 block
1301# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1302 use iso_fortran_env, only: output_unit
1303# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1304
1305# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1306 print *, 'm_acoustic_src.fpp:429: ', '@:ALLOCATE(source_spatials(ai)%val(1:count))'
1307# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1308
1309# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1310 call flush (output_unit)
1311# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1312 end block
1313# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1314#endif
1315# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1316 allocate (source_spatials(ai)%val(1:count))
1317# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1318
1319# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1320
1321# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1322#if defined(MFC_OpenACC)
1323# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1324!$acc enter data create(source_spatials(ai)%val)
1325# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1326#elif defined(MFC_OpenMP)
1327# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1328!$omp target enter data map(always,alloc:source_spatials(ai)%val)
1329# 429 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1330#endif
1331#ifdef MFC_DEBUG
1332# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1333 block
1334# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1335 use iso_fortran_env, only: output_unit
1336# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1337
1338# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1339 print *, 'm_acoustic_src.fpp:430: ', '@:ALLOCATE(source_spatials(ai)%angle(1:count))'
1340# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1341
1342# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1343 call flush (output_unit)
1344# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1345 end block
1346# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1347#endif
1348# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1349 allocate (source_spatials(ai)%angle(1:count))
1350# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1351
1352# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1353
1354# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1355#if defined(MFC_OpenACC)
1356# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1357!$acc enter data create(source_spatials(ai)%angle)
1358# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1359#elif defined(MFC_OpenMP)
1360# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1361!$omp target enter data map(always,alloc:source_spatials(ai)%angle)
1362# 430 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1363#endif
1364#ifdef MFC_DEBUG
1365# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1366 block
1367# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1368 use iso_fortran_env, only: output_unit
1369# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1370
1371# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1372 print *, 'm_acoustic_src.fpp:431: ', '@:ALLOCATE(source_spatials(ai)%xyz_to_r_ratios(1:3, 1:count))'
1373# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1374
1375# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1376 call flush (output_unit)
1377# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1378 end block
1379# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1380#endif
1381# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1382 allocate (source_spatials(ai)%xyz_to_r_ratios(1:3, 1:count))
1383# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1384
1385# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1386
1387# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1388#if defined(MFC_OpenACC)
1389# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1390!$acc enter data create(source_spatials(ai)%xyz_to_r_ratios)
1391# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1392#elif defined(MFC_OpenMP)
1393# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1394!$omp target enter data map(always,alloc:source_spatials(ai)%xyz_to_r_ratios)
1395# 431 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1396#endif
1397
1398#ifdef _CRAYFTN
1399# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1400 block
1401# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1402#ifdef MFC_DEBUG
1403# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1404 block
1405# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1406 use iso_fortran_env, only: output_unit
1407# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1408
1409# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1410 print *, 'm_acoustic_src.fpp:433: ', '@:ACC_SETUP_source_spatials(source_spatials(ai))'
1411# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1412
1413# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1414 call flush (output_unit)
1415# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1416 end block
1417# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1418#endif
1419# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1420
1421# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1422
1423# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1424#if defined(MFC_OpenACC)
1425# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1426!$acc enter data copyin(source_spatials(ai))
1427# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1428#elif defined(MFC_OpenMP)
1429# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1430!$omp target enter data map(to:source_spatials(ai))
1431# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1432#endif
1433# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1434 if (associated(source_spatials(ai)%coord)) then
1435# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1436
1437# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1438#if defined(MFC_OpenACC)
1439# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1440!$acc enter data copyin(source_spatials(ai)%coord)
1441# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1442#elif defined(MFC_OpenMP)
1443# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1444!$omp target enter data map(to:source_spatials(ai)%coord)
1445# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1446#endif
1447# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1448 end if
1449# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1450 if (associated(source_spatials(ai)%val)) then
1451# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1452
1453# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1454#if defined(MFC_OpenACC)
1455# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1456!$acc enter data copyin(source_spatials(ai)%val)
1457# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1458#elif defined(MFC_OpenMP)
1459# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1460!$omp target enter data map(to:source_spatials(ai)%val)
1461# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1462#endif
1463# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1464 end if
1465# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1466 if (associated(source_spatials(ai)%angle)) then
1467# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1468
1469# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1470#if defined(MFC_OpenACC)
1471# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1472!$acc enter data copyin(source_spatials(ai)%angle)
1473# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1474#elif defined(MFC_OpenMP)
1475# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1476!$omp target enter data map(to:source_spatials(ai)%angle)
1477# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1478#endif
1479# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1480 end if
1481# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1482 if (associated(source_spatials(ai)%xyz_to_r_ratios)) then
1483# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1484
1485# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1486#if defined(MFC_OpenACC)
1487# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1488!$acc enter data copyin(source_spatials(ai)%xyz_to_r_ratios)
1489# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1490#elif defined(MFC_OpenMP)
1491# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1492!$omp target enter data map(to:source_spatials(ai)%xyz_to_r_ratios)
1493# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1494#endif
1495# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1496 end if
1497# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1498 end block
1499# 433 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1500#endif
1501
1502 ! Second pass: Store the values
1503 count = 0 ! Reset counter
1504 do l = 0, p
1505 do k = 0, n
1506 do j = 0, m
1507 call s_source_spatial(j, k, l, loc_acoustic(:,ai), ai, source_spatial, angle, xyz_to_r_ratios)
1508 if (abs(source_spatial) < threshold) cycle
1509 count = count + 1
1510 source_spatials(ai)%coord(1, count) = j
1511 source_spatials(ai)%coord(2, count) = k
1512 source_spatials(ai)%coord(3, count) = l
1513 source_spatials(ai)%val(count) = source_spatial
1514 if (support(ai) >= 5) then
1515 if (dim == 2) source_spatials(ai)%angle(count) = angle
1516 if (dim == 3) source_spatials(ai)%xyz_to_r_ratios(1:3,count) = xyz_to_r_ratios
1517 end if
1518 end do
1519 end do
1520 end do
1521
1522 if (source_spatials_num_points(ai) /= count) then
1523 call s_mpi_abort('Fatal Error: Inconsistent allocation of source_spatials')
1524 end if
1525
1526 if (count > 0) then
1527
1528# 460 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1529#if defined(MFC_OpenACC)
1530# 460 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1531!$acc update device(source_spatials(ai)%coord)
1532# 460 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1533#elif defined(MFC_OpenMP)
1534# 460 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1535!$omp target update to(source_spatials(ai)%coord)
1536# 460 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1537#endif
1538
1539# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1540#if defined(MFC_OpenACC)
1541# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1542!$acc update device(source_spatials(ai)%val)
1543# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1544#elif defined(MFC_OpenMP)
1545# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1546!$omp target update to(source_spatials(ai)%val)
1547# 461 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1548#endif
1549 if (support(ai) >= 5) then
1550 if (dim == 2) then
1551
1552# 464 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1553#if defined(MFC_OpenACC)
1554# 464 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1555!$acc update device(source_spatials(ai)%angle)
1556# 464 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1557#elif defined(MFC_OpenMP)
1558# 464 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1559!$omp target update to(source_spatials(ai)%angle)
1560# 464 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1561#endif
1562 end if
1563 if (dim == 3) then
1564
1565# 467 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1566#if defined(MFC_OpenACC)
1567# 467 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1568!$acc update device(source_spatials(ai)%xyz_to_r_ratios)
1569# 467 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1570#elif defined(MFC_OpenMP)
1571# 467 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1572!$omp target update to(source_spatials(ai)%xyz_to_r_ratios)
1573# 467 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1574#endif
1575 end if
1576 end if
1577 end if
1578 end do
1579
1580#ifdef MFC_DEBUG
1581 do ai = 1, num_source
1582 write (*, '(A,I2,A,I8,A)') 'Acoustic source ', ai, ' has ', source_spatials_num_points(ai), &
1583 & ' grid points with non-zero source term'
1584 end do
1585#endif
1586
1588
1589 !> Compute the spatial support of the acoustic source
1590 subroutine s_source_spatial(j, k, l, loc, ai, source, angle, xyz_to_r_ratios)
1591
1592 integer, intent(in) :: j, k, l, ai
1593 real(wp), dimension(3), intent(in) :: loc
1594 real(wp), intent(out) :: source, angle, xyz_to_r_ratios(3)
1595 real(wp) :: sig, r(3)
1596
1597 ! Calculate sig spatial support width
1598
1599 if (n == 0) then
1600 sig = dx(j)
1601 else if (p == 0) then
1602 sig = maxval((/dx(j), dy(k)/))
1603 else
1604 sig = maxval((/dx(j), dy(k), dz(l)/))
1605 end if
1606 sig = sig*acoustic_spatial_support_width
1607
1608 ! Calculate displacement from acoustic source location
1609 r(1) = x_cc(j) - loc(1)
1610 if (n /= 0) r(2) = y_cc(k) - loc(2)
1611 if (p /= 0) r(3) = z_cc(l) - loc(3)
1612
1613 if (any(support(ai) == (/1, 2, 3, 4/))) then
1614 call s_source_spatial_planar(ai, sig, r, source)
1615 else if (any(support(ai) == (/5, 6, 7/))) then
1616 call s_source_spatial_transducer(ai, sig, r, source, angle, xyz_to_r_ratios)
1617 else if (any(support(ai) == (/9, 10, 11/))) then
1618 call s_source_spatial_transducer_array(ai, sig, r, source, angle, xyz_to_r_ratios)
1619 end if
1620
1621 end subroutine s_source_spatial
1622
1623 !> Compute the spatial support for planar acoustic sources in 1D, 2D, and 3D
1624 subroutine s_source_spatial_planar(ai, sig, r, source)
1625
1626 integer, intent(in) :: ai
1627 real(wp), intent(in) :: sig, r(3)
1628 real(wp), intent(out) :: source
1629 real(wp) :: dist
1630
1631 source = 0._wp
1632
1633 ! Gaussian spatial pulse profile: exp(-0.5 * (d / sigma)^2) / (sqrt(2*pi) * sigma)
1634 if (support(ai) == 1) then ! 1D
1635 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(r(1)/(sig/2._wp))**2._wp)
1636 else if (support(ai) == 2 .or. support(ai) == 3) then ! 2D or 3D
1637 ! If we let unit vector e = (cos(dir), sin(dir)),
1638 dist = r(1)*cos(dir(ai)) + r(2)*sin(dir(ai)) ! dot(r,e)
1639 ! |r - dist*e| < length/2
1640 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
1641 if (support(ai) /= 3 .or. abs(r(3)) < 0.25_wp*height(ai)) then ! additional height constraint for 3D
1642 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)
1643 end if
1644 end if
1645 end if
1646
1647 end subroutine s_source_spatial_planar
1648
1649 !> Compute the spatial support for a single transducer in 2D, 2D axisymmetric, and 3D
1650 subroutine s_source_spatial_transducer(ai, sig, r, source, angle, xyz_to_r_ratios)
1651
1652 integer, intent(in) :: ai
1653 real(wp), intent(in) :: sig, r(3)
1654 real(wp), intent(out) :: source, angle, xyz_to_r_ratios(3)
1655 real(wp) :: current_angle, angle_half_aperture, dist, norm
1656
1657 source = 0._wp ! If not affected by transducer
1658 angle = 0._wp
1659 xyz_to_r_ratios = 0._wp
1660
1661 if (support(ai) == 5 .or. support(ai) == 6) then ! 2D or 2D axisymmetric
1662 current_angle = -atan(r(2)/(foc_length(ai) - r(1)))
1663 angle_half_aperture = asin((aperture(ai)/2._wp)/(foc_length(ai)))
1664
1665 if (abs(current_angle) < angle_half_aperture .and. r(1) < foc_length(ai)) then
1666 dist = foc_length(ai) - sqrt(r(2)**2._wp + (foc_length(ai) - r(1))**2._wp)
1667 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)
1668 angle = -atan(r(2)/(foc_length(ai) - r(1)))
1669 end if
1670 else if (support(ai) == 7) then ! 3D
1671 current_angle = -atan(sqrt(r(2)**2 + r(3)**2)/(foc_length(ai) - r(1)))
1672 angle_half_aperture = asin((aperture(ai)/2._wp)/(foc_length(ai)))
1673
1674 if (abs(current_angle) < angle_half_aperture .and. r(1) < foc_length(ai)) then
1675 dist = foc_length(ai) - sqrt(r(2)**2._wp + r(3)**2._wp + (foc_length(ai) - r(1))**2._wp)
1676 source = 1._wp/(sqrt(2._wp*pi)*sig/2._wp)*exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)
1677
1678 norm = sqrt(r(2)**2._wp + r(3)**2._wp + (foc_length(ai) - r(1))**2._wp)
1679 xyz_to_r_ratios(1) = -(r(1) - foc_length(ai))/norm
1680 xyz_to_r_ratios(2) = -r(2)/norm
1681 xyz_to_r_ratios(3) = -r(3)/norm
1682 end if
1683 end if
1684
1685 end subroutine s_source_spatial_transducer
1686
1687 !> Compute the spatial support for multiple transducers in 2D, 2D axisymmetric, and 3D
1688 subroutine s_source_spatial_transducer_array(ai, sig, r, source, angle, xyz_to_r_ratios)
1689
1690 integer, intent(in) :: ai
1691 real(wp), intent(in) :: sig, r(3)
1692 real(wp), intent(out) :: source, angle, xyz_to_r_ratios(3)
1693 integer :: elem, elem_min, elem_max
1694 real(wp) :: current_angle, angle_half_aperture, angle_per_elem, dist
1695 real(wp) :: angle_min, angle_max, norm
1696 real(wp) :: poly_side_length, aperture_element_3D, angle_elem
1697 real(wp) :: x2, y2, z2, x3, y3, z3, C, f, half_apert, dist_interp_to_elem_center
1698
1699 if (element_on(ai) == 0) then ! Full transducer
1700 elem_min = 1
1701 elem_max = num_elements(ai)
1702 else ! Transducer element specified
1703 elem_min = element_on(ai)
1704 elem_max = element_on(ai)
1705 end if
1706
1707 source = 0._wp ! If not affected by any transducer element
1708 angle = 0._wp
1709 xyz_to_r_ratios = 0._wp
1710
1711 if (support(ai) == 9 .or. support(ai) == 10) then ! 2D or 2D axisymmetric
1712 current_angle = -atan(r(2)/(foc_length(ai) - r(1)))
1713 angle_half_aperture = asin((aperture(ai)/2._wp)/(foc_length(ai)))
1714 angle_per_elem = (2._wp*angle_half_aperture - (num_elements(ai) - 1._wp)*element_spacing_angle(ai))/num_elements(ai)
1715 dist = foc_length(ai) - sqrt(r(2)**2._wp + (foc_length(ai) - r(1))**2._wp)
1716
1717 do elem = elem_min, elem_max
1718 angle_max = angle_half_aperture - (element_spacing_angle(ai) + angle_per_elem)*(elem - 1._wp)
1719 angle_min = angle_max - angle_per_elem
1720
1721 if (current_angle > angle_min .and. current_angle < angle_max .and. r(1) < foc_length(ai)) then
1722 source = exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)/(sqrt(2._wp*pi)*sig/2._wp)
1723 angle = current_angle
1724 exit ! Assume elements don't overlap
1725 end if
1726 end do
1727 else if (support(ai) == 11) then ! 3D
1728 poly_side_length = aperture(ai)*sin(pi/num_elements(ai))
1729 aperture_element_3d = poly_side_length*element_polygon_ratio(ai)
1730 f = foc_length(ai)
1731 half_apert = aperture(ai)/2._wp
1732
1733 do elem = elem_min, elem_max
1734 angle_elem = 2._wp*pi*real(elem, wp)/real(num_elements(ai), wp) + rotate_angle(ai)
1735
1736 ! Point 2 is the elem center
1737 x2 = f - sqrt(f**2 - half_apert**2)
1738 y2 = half_apert*cos(angle_elem)
1739 z2 = half_apert*sin(angle_elem)
1740
1741 ! Construct a plane normal to the line from the focal point to the elem center, Point 3 is the intercept of the
1742 ! plane and the line from the focal point to the current location
1743 c = f**2._wp/((r(1) - f)*(x2 - f) + r(2)*y2 + r(3)*z2) ! Constant for intermediate step
1744 x3 = c*(r(1) - f) + f
1745 y3 = c*r(2)
1746 z3 = c*r(3)
1747
1748 dist_interp_to_elem_center = sqrt((x2 - x3)**2._wp + (y2 - y3)**2._wp + (z2 - z3)**2._wp)
1749 if ((dist_interp_to_elem_center < aperture_element_3d/2._wp) .and. (r(1) < f)) then
1750 dist = sqrt((x3 - r(1))**2._wp + (y3 - r(2))**2._wp + (z3 - r(3))**2._wp)
1751 source = exp(-0.5_wp*(dist/(sig/2._wp))**2._wp)/(sqrt(2._wp*pi)*sig/2._wp)
1752
1753 norm = sqrt(r(2)**2._wp + r(3)**2._wp + (f - r(1))**2._wp)
1754 xyz_to_r_ratios(1) = -(r(1) - f)/norm
1755 xyz_to_r_ratios(2) = -r(2)/norm
1756 xyz_to_r_ratios(3) = -r(3)/norm
1757 end if
1758 end do
1759 end if
1760
1762
1763 !> Convert wavelength to frequency
1764 elemental function f_frequency_local(freq_conv_flag, ai, c)
1765
1766
1767# 659 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1768#if MFC_OpenACC
1769# 659 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1770!$acc routine seq
1771# 659 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1772#elif MFC_OpenMP
1773# 659 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1774
1775# 659 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1776
1777# 659 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1778!$omp declare target device_type(any)
1779# 659 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1780#endif
1781 logical, intent(in) :: freq_conv_flag
1782 integer, intent(in) :: ai
1783 real(wp), intent(in) :: c
1784 real(wp) :: f_frequency_local
1785
1786 if (freq_conv_flag) then
1788 else
1790 end if
1791
1792 end function f_frequency_local
1793
1794 !> Convert Gaussian sigma from distance to time
1795 function f_gauss_sigma_time_local(gauss_conv_flag, ai, c)
1796
1797
1798# 676 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1799#if MFC_OpenACC
1800# 676 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1801!$acc routine seq
1802# 676 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1803#elif MFC_OpenMP
1804# 676 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1805
1806# 676 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1807
1808# 676 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1809!$omp declare target device_type(any)
1810# 676 "/home/runner/work/MFC/MFC/src/simulation/m_acoustic_src.fpp"
1811#endif
1812 logical, intent(in) :: gauss_conv_flag
1813 integer, intent(in) :: ai
1814 real(wp), intent(in) :: c
1815 real(wp) :: f_gauss_sigma_time_local
1816
1817 if (gauss_conv_flag) then
1819 else
1821 end if
1822
1823 end function f_gauss_sigma_time_local
1824
1825end 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.