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