MFC
Exascale flow solver
Loading...
Searching...
No Matches
m_data_output.fpp.f90
Go to the documentation of this file.
1# 1 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2!>
3!! @file
4!! @brief Contains module m_data_output
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# 12 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
20# 1 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp" 1
21! Default for Fypp runs outside the toolchain (e.g. documentation). Searched last, so the
22! build/include/<target>/thermochem.fpp the toolchain writes from m_thermochem takes precedence.
23# 4 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp"
24# 5 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp"
25# 13 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 2
26
27# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
28# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
29
30# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
31
32# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
33
34# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
35
36# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
37
38# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
39
40# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
41
42# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
43
44# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
45! New line at end of file is required for FYPP
46# 2 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
47# 1 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 1
48# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
49# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
50# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
51# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
52# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
53# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
54
55# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
56# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
57# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
58
59# 12 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
60# 1 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp" 1
61! Default for Fypp runs outside the toolchain (e.g. documentation). Searched last, so the
62! build/include/<target>/thermochem.fpp the toolchain writes from m_thermochem takes precedence.
63# 4 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp"
64# 5 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp"
65# 13 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 2
66
67# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
68# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
69
70# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
71
72# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
73
74# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
75
76# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
77
78# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
79
80# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
81
82# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
83
84# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
85! New line at end of file is required for FYPP
86# 2 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 2
87
88# 4 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
89# 5 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
90# 6 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
91# 7 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
92# 8 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
93
94# 20 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
95
96# 43 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
97
98# 48 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
99
100# 53 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
101
102# 58 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
103
104# 63 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
105
106# 68 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
107
108# 81 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
109
110# 95 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
111
112# 100 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
113
114# 105 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
115
116# 110 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
117
118# 115 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
119
120# 120 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
121
122# 125 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
123
124# 130 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
125
126# 135 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
127
128# 165 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
129
130# 219 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
131
132# 233 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
133
134# 258 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
135
136# 269 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
137
138# 271 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
139# 282 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
140
141# 332 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
142
143# 342 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
144
145# 352 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
146
147# 361 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
148
149# 378 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
150
151# 388 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
152
153# 395 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
154
155# 401 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
156
157# 407 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
158
159# 413 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
160
161# 419 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
162
163# 425 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
164! New line at end of file is required for FYPP
165# 3 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
166# 1 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 1
167# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
168# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
169# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
170# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
171# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
172# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
173
174# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
175# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
176# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
177
178# 12 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
179# 1 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp" 1
180! Default for Fypp runs outside the toolchain (e.g. documentation). Searched last, so the
181! build/include/<target>/thermochem.fpp the toolchain writes from m_thermochem takes precedence.
182# 4 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp"
183# 5 "/home/runner/work/MFC/MFC/src/common/include/defaults/thermochem.fpp"
184# 13 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 2
185
186# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
187# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
188
189# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
190
191# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
192
193# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
194
195# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
196
197# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
198
199# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
200
201# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
202
203# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
204! New line at end of file is required for FYPP
205# 2 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 2
206
207# 7 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
208
209# 17 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
210
211# 22 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
212
213# 27 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
214
215# 32 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
216
217# 37 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
218
219# 42 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
220
221# 47 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
222
223# 52 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
224
225# 57 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
226
227# 62 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
228
229# 73 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
230
231# 78 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
232
233# 83 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
234
235# 88 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
236
237# 103 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
238
239# 131 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
240
241# 160 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
242
243# 175 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
244
245# 193 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
246
247# 215 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
248
249# 244 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
250
251# 259 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
252
253# 269 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
254
255# 278 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
256
257# 294 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
258
259# 304 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
260
261# 311 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
262! New line at end of file is required for FYPP
263# 4 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
264
265! GPU parallel region (scalar reductions, maxval/minval)
266# 23 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
267
268! GPU parallel loop over threads (most common GPU macro)
269# 43 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
270
271! Required closing for GPU_PARALLEL_LOOP
272# 55 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
273
274! Mark routine for device compilation
275# 112 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
276
277! Declare device-resident data
278# 130 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
279
280! Inner loop within a GPU parallel region
281# 145 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
282
283! Scoped GPU data region
284# 164 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
285
286! Host code with device pointers (for MPI with GPU buffers)
287# 193 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
288
289! Allocate device memory (unscoped)
290# 207 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
291
292! Free device memory
293# 219 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
294
295! Atomic operation on device
296# 231 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
297
298! End atomic capture block
299# 242 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
300
301! Copy data between host and device
302# 254 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
303
304! Synchronization barrier
305# 266 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
306
307! Import GPU library module (openacc or omp_lib)
308# 275 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
309
310! Emit code only for AMD compiler
311# 282 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
312
313! Emit code for non-Cray compilers
314# 289 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
315
316! Emit code only for Cray compiler
317# 296 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
318
319! Emit code for non-NVIDIA compilers
320# 303 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
321
322# 305 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
323# 306 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
324! New line at end of file is required for FYPP
325# 2 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 2
326
327# 14 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
328
329! Caution: This macro requires the use of a binding script to set CUDA_VISIBLE_DEVICES, such that we have one GPU device per MPI
330! rank. That's because for both cudaMemAdvise (preferred location) and cudaMemPrefetchAsync we use location = device_id = 0. For an
331! example see misc/nvidia_uvm/bind.sh.
332# 52 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
333
334! Allocate and create GPU device memory
335# 72 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
336
337! Free GPU device memory and deallocate
338# 80 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
339
340! Cray-specific GPU pointer setup for vector fields
341# 104 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
342
343! Cray-specific GPU pointer setup for scalar fields
344# 120 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
345
346! Cray-specific GPU pointer setup for acoustic source spatials
347# 145 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
348
349# 151 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
350
351# 158 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
352! New line at end of file is required for FYPP
353# 6 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp" 2
354# 1 "/home/runner/work/MFC/MFC/src/common/include/case.fpp" 1
355! This file exists so that Fypp can be run without generating case.fpp files for
356! each target. This is useful when generating documentation, for example. This
357! should also let MFC be built with CMake directly, without invoking mfc.sh.
358
359! For pre-process.
360# 8 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
361
362! For moving immersed boundaries in simulation
363# 12 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
364# 7 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp" 2
365
366!> @brief Writes solution data, run-time stability diagnostics (ICFL, VCFL, CCFL, Rc), and probe/center-of-mass files
367module m_data_output
368
369 use m_derived_types
370 use m_global_parameters
371 use m_mpi_proxy
372 use m_variables_conversion
373 use m_compile_specific
374 use m_helper
375 use m_helper_basic
376 use m_sim_helpers
377 use m_delay_file_access
378 use m_ibm
379 use m_boundary_common
380 use m_constants, only: model_eqns_5eq, precision_single
381
382 implicit none
383
384 private
385 public :: s_initialize_data_output_module, s_open_run_time_information_file, s_open_probe_files, &
386 & s_write_run_time_information, s_write_data_files, s_write_serial_data_files, s_write_parallel_data_files, &
387 & s_write_ib_data_file, s_write_probe_files, s_write_ib_state_file, s_write_ib_force_history, s_close_ib_force_history, &
388 & s_close_run_time_information_file, s_close_probe_files, s_finalize_data_output_module
389
390 !> @name ICFL, VCFL, CCFL, TCFL, and Rc stability criteria extrema over all the time-steps
391 !> @{
392 real(wp) :: icfl_max !< ICFL criterion maximum
393 real(wp) :: vcfl_max !< VCFL criterion maximum
394 real(wp) :: ccfl_max !< CCFL criterion maximum
395 real(wp) :: tcfl_max !< TCFL criterion maximum
396 real(wp) :: rc_min !< Rc criterion maximum
397 !> @}
398
399 type(scalar_field), allocatable, dimension(:) :: q_cons_temp_ds
400
401 !> One fixed-width text record per immersed body per recorded step, in D/ib_forces.dat.
402 !!
403 !! IB_REC_FMT is fixed width by construction: every ES descriptor right-justifies in its field,
404 !! including for negatives, three-digit exponents, NaN and Inf, so a record is always
405 !! IB_REC_BODY characters. That is what lets a rank compute a byte offset for (step, body) and
406 !! write there directly, giving one shared text file with no gather and no per-rank shards.
407 !! The two must be edited together: widening the format without IB_REC_LEN shears the file.
408 character(len=*), parameter :: ib_rec_fmt = '(I10,19(1X,ES17.9E3))'
409 integer, parameter :: ib_rec_body = 10 + 19*18 !< characters the format emits
410 integer, parameter :: ib_rec_len = ib_rec_body + 1 !< plus the newline
411 integer :: ib_hist_file = -1 !< held open for the run; -1 until first write
412
413contains
414
415 !> Write data files. Dispatch subroutine that replaces procedure pointer.
416 impure subroutine s_write_data_files(q_cons_vf, q_T_sf, q_prim_vf, t_step, bc_type, beta)
417
418 type(scalar_field), dimension(sys_size), intent(inout) :: q_cons_vf
419 type(scalar_field), intent(inout) :: q_t_sf
420 type(scalar_field), dimension(sys_size), intent(inout) :: q_prim_vf
421 integer, intent(in) :: t_step
422 type(scalar_field), intent(inout), optional :: beta
423 type(integer_field), dimension(1:num_dims,-1:1), intent(in) :: bc_type
424
425 if (.not. parallel_io) then
426 call s_write_serial_data_files(q_cons_vf, q_t_sf, q_prim_vf, t_step, bc_type, beta)
427 else
428 call s_write_parallel_data_files(q_cons_vf, t_step, bc_type, beta)
429 end if
430
431 end subroutine s_write_data_files
432
433 !> Open the run-time information file and write the stability criteria table header
434 impure subroutine s_open_run_time_information_file
435
436 character(LEN=name_len), parameter :: file_name = 'run_time.inf' !< Name of the run-time information file
437 character(LEN=path_len + name_len) :: file_path !< Relative path to a file in the case directory
438 character(LEN=8) :: file_date !< Creation date of the run-time information file
439
440 file_path = trim(case_dir) // '/' // trim(file_name)
441
442 open (3, file=trim(file_path), form='formatted', status='replace')
443
444 write (3, '(A)') 'Description: Stability information at ' // 'each time-step of the simulation. This'
445 write (3, '(13X,A)') 'data is composed of the inviscid ' // 'Courant-Friedrichs-Lewy (ICFL)'
446 write (3, '(13X,A)') 'number, the viscous CFL (VCFL) number, ' // 'the capillary CFL (CCFL)'
447 write (3, '(13X,A)') 'number, the thermal diffusion CFL (TCFL) ' // 'number, and the cell Reynolds (Rc)'
448 write (3, '(13X,A)') 'number. Please note that only'
449 write (3, '(13X,A)') 'those stability conditions pertinent ' // 'to the physics included in'
450 write (3, '(13X,A)') 'the current computation are displayed.'
451 if (hypoelasticity) then
452 write (3, '(13X,A)') 'NOTE: the reported ICFL uses the acoustic ' // 'sound speed only; it may'
453 write (3, '(13X,A)') 'underestimate the elastic characteristic ' // 'speeds.'
454 end if
455
456 call date_and_time(date=file_date)
457
458 write (3, '(A)') 'Date: ' // file_date(5:6) // '/' // file_date(7:8) // '/' // file_date(3:4)
459
460 write (3, '(A)') ''; write (3, '(A)') ''
461
462 write (3, '(13X,A9,13X,A10,13X,A10,13X,A10)', advance="no") trim('Time-step'), trim('dt'), trim('Time'), trim('ICFL Max')
463
464 if (surface_tension) then
465 write (3, '(13X,A10)', advance="no") trim('CCFL Max')
466 end if
467
468 if (heat_conduction) then
469 write (3, '(13X,A10)', advance="no") trim('TCFL Max')
470 end if
471
472 if (viscous) then
473 write (3, '(13X,A10,13X,A16)', advance="no") trim('VCFL Max'), trim('Rc Min')
474 end if
475
476 if (bubbles_lagrange) then
477 write (3, '(13X,A10)', advance="no") trim('N Bubbles')
478 end if
479
480 write (3, *) ! new line
481
482 end subroutine s_open_run_time_information_file
483
484 !> Open flow probe data files for writing
485 impure subroutine s_open_probe_files
486
487 character(LEN=path_len + 3*name_len) :: file_path !< Relative path to the probe data file in the case directory
488 integer :: i !< Generic loop iterator
489 logical :: file_exist
490 logical :: fresh_start
491
492 ! A run continues from a checkpoint when t_step_start > 0, or under cfl_dt when n_start > 0 (t_step_start stays
493 ! at its default there); pick the criterion by mode -- OR-ing them made every cfl_dt run look fresh.
494
495 fresh_start = merge(n_start == 0, t_step_start == 0, cfl_dt)
496
497 do i = 1, num_probes
498 write (file_path, '(A,I0,A)') '/D/probe', i, '_prim.dat'
499 file_path = trim(case_dir) // trim(file_path)
500
501 inquire (file=trim(file_path), exist=file_exist)
502
503 ! Append only when continuing a run. A fresh start that appends splices the previous run's rows
504 ! onto this one's, and nothing in the file marks the join: the time column simply resets partway
505 ! down, and the two runs need not even share a grid. Readers see one monotonic series and are
506 ! silently wrong.
507 if (file_exist .and. .not. fresh_start) then
508 open (i + 30, file=trim(file_path), form='formatted', status='old', position='append')
509 else
510 open (i + 30, file=trim(file_path), form='formatted', status='replace')
511 end if
512 end do
513
514 end subroutine s_open_probe_files
515
516 !> Write stability criteria extrema to the run-time information file at the given time step
517 impure subroutine s_write_run_time_information(q_prim_vf, t_step)
518
519 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf
520 integer, intent(in) :: t_step
521 real(wp) :: rho !< Cell-avg. density
522
523# 169 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
524 real(wp), dimension(num_fluids) :: alpha, alpha_rho !< Cell-avg. volume fraction, partial density
525 real(wp), dimension(num_vels) :: vel !< Cell-avg. velocity
526# 172 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
527 real(wp) :: vel_sum !< Cell-avg. velocity sum
528 real(wp) :: pres !< Cell-avg. pressure
529 real(wp) :: gamma !< Cell-avg. sp. heat ratio
530 real(wp) :: pi_inf !< Cell-avg. liquid stiffness function
531 real(wp) :: qv !< Cell-avg. internal energy reference value
532 real(wp) :: c !< Cell-avg. sound speed
533 real(wp), dimension(2) :: re !< Cell-avg. Reynolds numbers
534 integer :: j, k, l
535 real(wp) :: icfl_max_loc, icfl_max_glb !< ICFL stability extrema on local and global grids
536 real(wp) :: vcfl_max_loc, vcfl_max_glb !< VCFL stability extrema on local and global grids
537 real(wp) :: ccfl_max_loc, ccfl_max_glb !< CCFL stability extrema on local and global grids
538 real(wp) :: tcfl_max_loc, tcfl_max_glb !< TCFL stability extrema on local and global grids
539 real(wp) :: rc_min_loc, rc_min_glb !< Rc stability extrema on local and global grids
540 real(wp) :: icfl, vcfl, ccfl, tcfl, rc
541 real(wp) :: mu_frac, mu_frac_max_loc, mu_frac_max_glb !< Compression as a fraction of the EOS limit
542 integer :: fl !< Fluid loop iterator
543 logical :: include_cell !< Cell is fluid, not ghost/inside an IB
544 real(wp), dimension(4) :: stab_max_loc, stab_max_glb !< Max-reduced criteria (ICFL, VCFL, CCFL, TCFL), packed
545 real(wp), dimension(1) :: stab_min_loc, stab_min_glb !< Min-reduced criteria (Rc), packed
546
547 icfl_max_loc = 0._wp
548 vcfl_max_loc = 0._wp
549 ccfl_max_loc = 0._wp
550 tcfl_max_loc = 0._wp
551 rc_min_loc = huge(1.0_wp)
552 mu_frac_max_loc = 0._wp
553 ! Computing Stability Criteria at Current Time-step
554
555# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
556
557# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
558#if defined(MFC_OpenACC)
559# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
560!$acc parallel loop collapse(3) gang vector default(present) &
561# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
562!$acc& private(j, k, l, vel, alpha, alpha_rho, Re, rho, vel_sum, pres, gamma, pi_inf, c, qv, icfl, vcfl, Rc, ccfl, tcfl, fl, mu_frac, include_cell) &
563# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
564!$acc& reduction(max:icfl_max_loc, vcfl_max_loc, ccfl_max_loc, tcfl_max_loc, mu_frac_max_loc) reduction(min:Rc_min_loc)
565# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
566#elif defined(MFC_OpenMP)
567# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
568
569# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
570
571# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
572
573# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
574!$omp target teams loop defaultmap(firstprivate:scalar) bind(teams,parallel) collapse(3) defaultmap(tofrom:aggregate) defaultmap(tofrom:allocatable) defaultmap(tofrom:pointer) &
575# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
576!$omp& private(j, k, l, vel, alpha, alpha_rho, Re, rho, vel_sum, pres, gamma, pi_inf, c, qv, icfl, vcfl, Rc, ccfl, tcfl, fl, mu_frac, include_cell) &
577# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
578!$omp& reduction(max:icfl_max_loc, vcfl_max_loc, ccfl_max_loc, tcfl_max_loc, mu_frac_max_loc) reduction(min:Rc_min_loc)
579# 199 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
580#endif
581# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
582 do l = 0, p
583 do k = 0, n
584 do j = 0, m
585 ! exclude cells inside of immersed boundaries
586 include_cell = .true.
587 if (ib) include_cell = (ib_markers%sf(j, k, l) == 0)
588 if (include_cell) then
589 call s_compute_cell_state(q_prim_vf, pres, rho, gamma, pi_inf, re, alpha, alpha_rho, vel, vel_sum, qv, j, &
590 & k, l)
591
592 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
593
594 ! How close each Mie-Gruneisen phase is to the compression its Hugoniot fit can represent.
595 ! Past 1 there is no shock state to find and the reference curve is fiction, so it is reduced
596 ! out of the kernel and turned into an abort on the host -- s_mpi_abort cannot be called here.
597 if (any_state_dependent_eos) then
598
599# 218 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
600#if defined(MFC_OpenACC)
601# 218 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
602!$acc loop seq
603# 218 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
604#elif defined(MFC_OpenMP)
605# 218 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
606
607# 218 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
608#endif
609# 219 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
610 do fl = 1, num_fluids
611 if (eoss(fl) == eos_mie_gruneisen) then
612 mu_frac = (alpha_rho(fl)/max(alpha(fl), &
613 & sgm_eps)/eos_coeffs(fl)%rho0 - 1._wp)/eos_coeffs(fl)%mu_max
614 mu_frac_max_loc = max(mu_frac_max_loc, mu_frac)
615 end if
616 end do
617 end if
618
619 if (any_non_newtonian) then
620 re(1) = 0._wp
621 do fl = 1, num_fluids
622 if (is_non_newtonian(fl)) then
623 re(1) = re(1) + alpha(fl)*hb_mu_max(fl)
624 else
625 re(1) = re(1) + alpha(fl)*fluid_inv_re(fl)
626 end if
627 end do
628 re(1) = 1._wp/max(re(1), sgm_eps)
629 end if
630
631 call s_compute_stability_from_dt(vel, c, rho, re, alpha, alpha_rho, j, k, l, icfl, vcfl, rc, ccfl, tcfl)
632
633 icfl_max_loc = max(icfl_max_loc, icfl)
634 vcfl_max_loc = max(vcfl_max_loc, merge(vcfl, 0.0_wp, viscous))
635 ccfl_max_loc = max(ccfl_max_loc, merge(ccfl, 0.0_wp, surface_tension))
636 tcfl_max_loc = max(tcfl_max_loc, merge(tcfl, 0.0_wp, heat_conduction))
637 rc_min_loc = min(rc_min_loc, merge(rc, huge(1.0_wp), viscous))
638 end if
639 end do
640 end do
641 end do
642
643# 251 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
644#if defined(MFC_OpenACC)
645# 251 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
646!$acc end parallel loop
647# 251 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
648#elif defined(MFC_OpenMP)
649# 251 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
650
651# 251 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
652!$omp end target teams loop
653# 251 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
654#endif
655# 252 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
656 ! end: Computing Stability Criteria at Current Time-step
657
658 if (num_procs > 1) then
659 stab_max_loc = (/icfl_max_loc, vcfl_max_loc, ccfl_max_loc, tcfl_max_loc/)
660 stab_min_loc = (/rc_min_loc/)
661 call s_mpi_reduce_stability_criteria_extrema(stab_max_loc, stab_min_loc, n_el_bubs_loc, stab_max_glb, stab_min_glb, &
662 & n_el_bubs_glb)
663 icfl_max_glb = stab_max_glb(1)
664 vcfl_max_glb = stab_max_glb(2)
665 ccfl_max_glb = stab_max_glb(3)
666 tcfl_max_glb = stab_max_glb(4)
667 rc_min_glb = stab_min_glb(1)
668 else
669 icfl_max_glb = icfl_max_loc
670 if (viscous) vcfl_max_glb = vcfl_max_loc
671 if (viscous) rc_min_glb = rc_min_loc
672 if (surface_tension) ccfl_max_glb = ccfl_max_loc
673 if (heat_conduction) tcfl_max_glb = tcfl_max_loc
674 if (bubbles_lagrange) n_el_bubs_glb = n_el_bubs_loc
675 end if
676
677 mu_frac_max_glb = mu_frac_max_loc
678 if (num_procs > 1) call s_mpi_allreduce_max(mu_frac_max_loc, mu_frac_max_glb)
679
680 if (icfl_max_glb > icfl_max) icfl_max = icfl_max_glb
681
682 if (surface_tension) then
683 if (ccfl_max_glb > ccfl_max) ccfl_max = ccfl_max_glb
684 end if
685
686 if (heat_conduction) then
687 if (tcfl_max_glb > tcfl_max) tcfl_max = tcfl_max_glb
688 end if
689
690 if (viscous) then
691 if (vcfl_max_glb > vcfl_max) vcfl_max = vcfl_max_glb
692 if (rc_min_glb < rc_min) rc_min = rc_min_glb
693 end if
694
695 ! Any rank whose own local extremum violates the limit is, by construction of the
696 ! max-reduction above, a rank that actually contains the offending cell(s).
697 if ((.not. f_approx_equal(icfl_max_loc, icfl_max_loc)) .or. icfl_max_loc > 1._wp) then
698 call s_report_icfl_violation(q_prim_vf)
699 end if
700 call s_mpi_barrier() ! ensure diagnostic output above is flushed before any rank aborts below
701
702 if (proc_rank == 0) then
703 write (3, '(13X,I9,13X,F10.6,13X,F10.6,13X,F10.6)', advance="no") t_step, dt, mytime, icfl_max_glb
704
705 if (surface_tension) then
706 write (3, '(13X,F10.6)', advance="no") ccfl_max_glb
707 end if
708
709 if (heat_conduction) then
710 write (3, '(13X,F10.6)', advance="no") tcfl_max_glb
711 end if
712
713 if (viscous) then
714 write (3, '(13X,F10.6,13X,ES16.6)', advance="no") vcfl_max_glb, rc_min_glb
715 end if
716
717 if (bubbles_lagrange) then
718 write (3, '(13X,I10)', advance="no") n_el_bubs_glb
719 end if
720
721 write (3, *) ! new line
722
723 if (mu_frac_max_glb > 1._wp) then
724 print *, 'compression as a fraction of the Hugoniot fit limit', mu_frac_max_glb
725 call s_mpi_abort('A Mie-Gruneisen phase is compressed past what its Hugoniot fit represents. Exiting.')
726 end if
727
728 if (.not. f_approx_equal(icfl_max_glb, icfl_max_glb)) then
729 call s_mpi_abort('ICFL is NaN. Exiting.')
730 else if (icfl_max_glb > 1._wp) then
731 print *, 'icfl', icfl_max_glb
732 call s_mpi_abort('ICFL is greater than 1.0. Exiting.')
733 end if
734
735 if (viscous) then
736 if (.not. f_approx_equal(vcfl_max_glb, vcfl_max_glb)) then
737 call s_mpi_abort('VCFL is NaN. Exiting.')
738 else if (vcfl_max_glb > 1._wp) then
739 print *, 'vcfl', vcfl_max_glb
740 call s_mpi_abort('VCFL is greater than 1.0. Exiting.')
741 end if
742 end if
743
744 if (bubbles_lagrange) then
745 if (n_el_bubs_glb == 0) then
746 call s_mpi_abort('No Lagrangian bubbles remain in the domain. Exiting.')
747 end if
748 end if
749 end if
750
751 call s_mpi_barrier()
752
753 end subroutine s_write_run_time_information
754
755 !> Locate the grid cell responsible for an ICFL violation on this rank and report its state plus the nearest immersed-boundary
756 !! particles, to aid debugging stability failures in particle-laden high-Mach cases.
757 impure subroutine s_report_icfl_violation(q_prim_vf)
758
759 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf
760 real(wp), dimension(num_fluids) :: alpha, alpha_rho
761 real(wp), dimension(num_vels) :: vel, vel_hit
762 real(wp), dimension(2) :: re
763 real(wp) :: rho, vel_sum, pres, gamma, pi_inf, qv, c
764 real(wp) :: rho_hit, pres_hit, c_hit
765 real(wp) :: icfl, vcfl, rc, ccfl, tcfl, icfl_hit
766 integer :: i, j, k, l, fl, j_hit, k_hit, l_hit
767 real(wp) :: x_hit, y_hit, z_hit, dist
768 logical :: nan_hit
769 integer :: near1_id, near2_id
770 real(wp) :: near1_dist, near2_dist
771
772 do i = 1, sys_size
773
774# 369 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
775#if defined(MFC_OpenACC)
776# 369 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
777!$acc update host(q_prim_vf(i)%sf(:, :, :))
778# 369 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
779#elif defined(MFC_OpenMP)
780# 369 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
781!$omp target update from(q_prim_vf(i)%sf(:, :, :))
782# 369 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
783#endif
784# 370 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
785 end do
786 if (ib) then
787
788# 372 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
789#if defined(MFC_OpenACC)
790# 372 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
791!$acc update host(ib_markers%sf)
792# 372 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
793#elif defined(MFC_OpenMP)
794# 372 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
795!$omp target update from(ib_markers%sf)
796# 372 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
797#endif
798# 373 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
799 end if
800
801 icfl_hit = -huge(1._wp)
802 nan_hit = .false.
803 j_hit = 0; k_hit = 0; l_hit = 0
804
805 scan: do l = 0, p
806 do k = 0, n
807 do j = 0, m
808 if (ib) then
809 if (ib_markers%sf(j, k, l) /= 0) cycle
810 end if
811
812 call s_compute_cell_state(q_prim_vf, pres, rho, gamma, pi_inf, re, alpha, alpha_rho, vel, vel_sum, qv, j, k, l)
813 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
814
815 if (any_non_newtonian) then
816 re(1) = 0._wp
817 do fl = 1, num_fluids
818 if (is_non_newtonian(fl)) then
819 re(1) = re(1) + alpha(fl)*hb_mu_max(fl)
820 else
821 re(1) = re(1) + alpha(fl)*fluid_inv_re(fl)
822 end if
823 end do
824 re(1) = 1._wp/max(re(1), sgm_eps)
825 end if
826
827 call s_compute_stability_from_dt(vel, c, rho, re, alpha, alpha_rho, j, k, l, icfl, vcfl, rc, ccfl, tcfl)
828
829 if (.not. f_approx_equal(icfl, icfl)) then
830 nan_hit = .true.
831 j_hit = j; k_hit = k; l_hit = l
832 rho_hit = rho; pres_hit = pres; c_hit = c; vel_hit = vel
833 exit scan
834 else if (icfl > icfl_hit) then
835 icfl_hit = icfl
836 j_hit = j; k_hit = k; l_hit = l
837 rho_hit = rho; pres_hit = pres; c_hit = c; vel_hit = vel
838 end if
839 end do
840 end do
841 end do scan
842
843 x_hit = x_cc(j_hit)
844 y_hit = 0._wp; if (n > 0) y_hit = y_cc(k_hit)
845 z_hit = 0._wp; if (p > 0) z_hit = z_cc(l_hit)
846
847 print '(A,I0,A,I0,A,I0,A,I0,A)', 'ICFL violation on rank ', proc_rank, ': cell (j,k,l) = (', j_hit, ',', k_hit, ',', &
848 & l_hit, ')'
849 if (nan_hit) then
850 print '(A)', ' icfl = NaN'
851 else
852 print '(A,ES16.6)', ' icfl = ', icfl_hit
853 end if
854 print '(A,3(ES16.6,1X))', ' position = ', x_hit, y_hit, z_hit
855 print '(A,ES16.6,A,ES16.6,A,ES16.6)', ' rho, pres, c = ', rho_hit, ', ', pres_hit, ', ', c_hit
856 print '(A,3(ES16.6,1X))', ' velocity = ', vel_hit
857 if (ib) print '(A,I0)', ' ib_markers = ', ib_markers%sf(j_hit, k_hit, l_hit)
858
859 if (ib .and. num_ibs > 0) then
860 near1_id = 0; near1_dist = huge(1._wp)
861 near2_id = 0; near2_dist = huge(1._wp)
862 do i = 1, num_ibs
863 dist = sqrt((x_hit - patch_ib(i)%x_centroid)**2 + (y_hit - patch_ib(i)%y_centroid)**2 + (z_hit &
864 & - patch_ib(i)%z_centroid)**2)
865 if (dist < near1_dist) then
866 near2_dist = near1_dist; near2_id = near1_id
867 near1_dist = dist; near1_id = i
868 else if (dist < near2_dist) then
869 near2_dist = dist; near2_id = i
870 end if
871 end do
872 if (near1_id > 0) then
873 print '(A,I0,A,ES16.6,A,ES16.6,A,3(ES16.6,1X))', ' nearest particle id=', near1_id, ' dist=', near1_dist, &
874 & ' gap=', near1_dist - patch_ib(near1_id)%radius, ' vel=', patch_ib(near1_id)%vel
875 print '(A,3(ES16.6,1X))', ' centroid = ', patch_ib(near1_id)%x_centroid, patch_ib(near1_id)%y_centroid, &
876 & patch_ib(near1_id)%z_centroid
877 print '(A,3(ES16.6,1X))', ' angular_vel = ', patch_ib(near1_id)%angular_vel
878 print '(A,3(ES16.6,1X))', ' force = ', patch_ib(near1_id)%force
879 print '(A,3(ES16.6,1X))', ' torque = ', patch_ib(near1_id)%torque
880 print '(A,I0,A,ES16.6,A,ES16.6)', ' moving_ibm = ', patch_ib(near1_id)%moving_ibm, ' mass=', &
881 & patch_ib(near1_id)%mass, ' moment=', patch_ib(near1_id)%moment
882 end if
883 if (near2_id > 0) then
884 print '(A,I0,A,ES16.6,A,ES16.6,A,3(ES16.6,1X))', ' 2nd nearest particle id=', near2_id, ' dist=', near2_dist, &
885 & ' gap=', near2_dist - patch_ib(near2_id)%radius, ' vel=', patch_ib(near2_id)%vel
886 print '(A,3(ES16.6,1X))', ' centroid = ', patch_ib(near2_id)%x_centroid, patch_ib(near2_id)%y_centroid, &
887 & patch_ib(near2_id)%z_centroid
888 print '(A,3(ES16.6,1X))', ' angular_vel = ', patch_ib(near2_id)%angular_vel
889 print '(A,3(ES16.6,1X))', ' force = ', patch_ib(near2_id)%force
890 end if
891 end if
892
893 ! Dump a small x-neighborhood around the violating cell (reaching into the ghost/halo region on either side) to
894 ! distinguish a sharp discontinuity at a processor boundary - the signature of stale or corrupted halo/IB state -
895 ! from a smoothly diverging field, which indicates a genuine physical/numerical instability.
896 print '(A)', ' x-neighborhood (dj, rho, pres, vel) around violating cell:'
897 do j = max(-buff_size, j_hit - 3), min(m + buff_size, j_hit + 3)
898 call s_compute_cell_state(q_prim_vf, pres, rho, gamma, pi_inf, re, alpha, alpha_rho, vel, vel_sum, qv, j, k_hit, l_hit)
899 print '(A,I0,A,ES16.6,A,ES16.6,A,3(ES16.6,1X))', ' dj=', j - j_hit, ' rho=', rho, ' pres=', pres, ' vel=', vel
900 end do
901
902 call flush (6)
903
904 end subroutine s_report_icfl_violation
905
906 !> Write grid and conservative variable data files in serial format
907 impure subroutine s_write_serial_data_files(q_cons_vf, q_T_sf, q_prim_vf, t_step, bc_type, beta)
908
909 type(scalar_field), dimension(sys_size), intent(inout) :: q_cons_vf
910 type(scalar_field), intent(inout) :: q_t_sf
911 type(scalar_field), dimension(sys_size), intent(inout) :: q_prim_vf
912 integer, intent(in) :: t_step
913 type(scalar_field), intent(inout), optional :: beta
914 type(integer_field), dimension(1:num_dims,-1:1), intent(in) :: bc_type
915 character(LEN=path_len + 2*name_len) :: t_step_dir !< Relative path to the current time-step directory
916 character(LEN=path_len + 3*name_len) :: file_path !< Relative path to the grid and conservative variables data files
917 logical :: file_exist !< Logical used to check existence of current time-step directory
918 character(LEN=15) :: fmt
919 integer :: i, j, k, l, r
920
921 write (t_step_dir, '(A,I0,A,I0)') trim(case_dir) // '/p_all'
922 write (t_step_dir, '(a,i0,a,i0)') trim(case_dir) // '/p_all/p', proc_rank, '/', t_step
923
924 file_path = trim(t_step_dir) // '/.'
925 call my_inquire(file_path, file_exist)
926 if (file_exist) call s_delete_directory(trim(t_step_dir))
927 call s_create_directory(trim(t_step_dir))
928
929 file_path = trim(t_step_dir) // '/x_cb.dat'
930
931 open (2, file=trim(file_path), form='unformatted', status='new')
932 write (2) x_cb(-1:m); close (2)
933
934 if (n > 0) then
935 file_path = trim(t_step_dir) // '/y_cb.dat'
936
937 open (2, file=trim(file_path), form='unformatted', status='new')
938 write (2) y_cb(-1:n); close (2)
939
940 if (p > 0) then
941 file_path = trim(t_step_dir) // '/z_cb.dat'
942
943 open (2, file=trim(file_path), form='unformatted', status='new')
944 write (2) z_cb(-1:p); close (2)
945 end if
946 end if
947
948 do i = 1, sys_size
949 write (file_path, '(A,I0,A)') trim(t_step_dir) // '/q_cons_vf', i, '.dat'
950
951 open (2, file=trim(file_path), form='unformatted', status='new')
952
953 write (2) q_cons_vf(i)%sf(0:m,0:n,0:p); close (2)
954 end do
955
956 ! Lagrangian beta (void fraction) written as q_cons_vf(sys_size+1) to match the parallel I/O path and allow post_process to
957 ! read it.
958 if (bubbles_lagrange) then
959 write (file_path, '(A,I0,A)') trim(t_step_dir) // '/q_cons_vf', sys_size + 1, '.dat'
960
961 open (2, file=trim(file_path), form='unformatted', status='new')
962
963 write (2) beta%sf(0:m,0:n,0:p); close (2)
964 end if
965
966 if (qbmm .and. .not. polytropic) then
967 do i = 1, nb
968 do r = 1, nnode
969 write (file_path, '(A,I0,A)') trim(t_step_dir) // '/pb', sys_size + (i - 1)*nnode + r, '.dat'
970
971 open (2, file=trim(file_path), form='unformatted', status='new')
972
973 write (2) pb_ts(1)%sf(0:m,0:n,0:p,r, i); close (2)
974 end do
975 end do
976
977 do i = 1, nb
978 do r = 1, nnode
979 write (file_path, '(A,I0,A)') trim(t_step_dir) // '/mv', sys_size + (i - 1)*nnode + r, '.dat'
980
981 open (2, file=trim(file_path), form='unformatted', status='new')
982
983 write (2) mv_ts(1)%sf(0:m,0:n,0:p,r, i); close (2)
984 end do
985 end do
986 end if
987
988 ! Writing the IB markers
989 if (ib) then
990 call s_write_serial_ib_data(t_step)
991 end if
992
993 if (precision == precision_single) then
994 fmt = "(2F30.3)"
995 else
996 fmt = "(2F40.14)"
997 end if
998
999 write (t_step_dir, '(A,I0,A,I0)') trim(case_dir) // '/D'
1000 file_path = trim(t_step_dir) // '/.'
1001
1002 inquire (file=trim(file_path), exist=file_exist)
1003
1004 if (.not. file_exist) call s_create_directory(trim(t_step_dir))
1005
1006 if ((prim_vars_wrt .or. (n == 0 .and. p == 0)) .and. (.not. igr)) then
1007 call s_convert_conservative_to_primitive_variables(q_cons_vf, q_t_sf, q_prim_vf, idwint)
1008 do i = 1, sys_size
1009
1010# 583 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1011#if defined(MFC_OpenACC)
1012# 583 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1013!$acc update host(q_prim_vf(i)%sf(:, :, :))
1014# 583 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1015#elif defined(MFC_OpenMP)
1016# 583 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1017!$omp target update from(q_prim_vf(i)%sf(:, :, :))
1018# 583 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1019#endif
1020# 584 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1021 end do
1022 ! q_prim_vf(eqn_idx%bub%beg) stores the value of nb needed in riemann solvers, so replace with true primitive value
1023 ! (=1._wp)
1024 if (qbmm) then
1025 q_prim_vf(eqn_idx%bub%beg)%sf = 1._wp
1026 end if
1027 end if
1028
1029 if (n == 0 .and. p == 0) then
1030 if (model_eqns == model_eqns_5eq .and. (.not. igr)) then
1031 do i = 1, sys_size
1032 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/prim.', i, '.', proc_rank, '.', t_step, '.dat'
1033
1034 open (2, file=trim(file_path))
1035 do j = 0, m
1036 ! todo: revisit change here
1037 if (((i >= eqn_idx%adv%beg) .and. (i <= eqn_idx%adv%end))) then
1038 write (2, fmt) x_cb(j), q_cons_vf(i)%sf(j, 0, 0)
1039 else
1040 write (2, fmt) x_cb(j), q_prim_vf(i)%sf(j, 0, 0)
1041 end if
1042 end do
1043 close (2)
1044 end do
1045 end if
1046
1047 do i = 1, sys_size
1048 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/cons.', i, '.', proc_rank, '.', t_step, '.dat'
1049
1050 open (2, file=trim(file_path))
1051 do j = 0, m
1052 write (2, fmt) x_cb(j), q_cons_vf(i)%sf(j, 0, 0)
1053 end do
1054 close (2)
1055 end do
1056
1057 if (qbmm .and. .not. polytropic) then
1058 do i = 1, nb
1059 do r = 1, nnode
1060 write (file_path, '(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/pres.', i, '.', r, '.', proc_rank, &
1061 & '.', t_step, '.dat'
1062
1063 open (2, file=trim(file_path))
1064 do j = 0, m
1065 write (2, fmt) x_cb(j), pb_ts(1)%sf(j, 0, 0, r, i)
1066 end do
1067 close (2)
1068 end do
1069 end do
1070 do i = 1, nb
1071 do r = 1, nnode
1072 write (file_path, '(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/mv.', i, '.', r, '.', proc_rank, &
1073 & '.', t_step, '.dat'
1074
1075 open (2, file=trim(file_path))
1076 do j = 0, m
1077 write (2, fmt) x_cb(j), mv_ts(1)%sf(j, 0, 0, r, i)
1078 end do
1079 close (2)
1080 end do
1081 end do
1082 end if
1083 end if
1084
1085 if (precision == precision_single) then
1086 fmt = "(3F30.7)"
1087 else
1088 fmt = "(3F40.14)"
1089 end if
1090
1091 if ((n > 0) .and. (p == 0)) then
1092 do i = 1, sys_size
1093 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/cons.', i, '.', proc_rank, '.', t_step, '.dat'
1094 open (2, file=trim(file_path))
1095 do j = 0, m
1096 do k = 0, n
1097 write (2, fmt) x_cb(j), y_cb(k), q_cons_vf(i)%sf(j, k, 0)
1098 end do
1099 write (2, *)
1100 end do
1101 close (2)
1102 end do
1103
1104 if (present(beta)) then
1105 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/beta.', i, '.', proc_rank, '.', t_step, '.dat'
1106 open (2, file=trim(file_path))
1107 do j = 0, m
1108 do k = 0, n
1109 write (2, fmt) x_cb(j), y_cb(k), beta%sf(j, k, 0)
1110 end do
1111 write (2, *)
1112 end do
1113 close (2)
1114 end if
1115
1116 if (qbmm .and. .not. polytropic) then
1117 do i = 1, nb
1118 do r = 1, nnode
1119 write (file_path, '(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/pres.', i, '.', r, '.', proc_rank, &
1120 & '.', t_step, '.dat'
1121
1122 open (2, file=trim(file_path))
1123 do j = 0, m
1124 do k = 0, n
1125 write (2, fmt) x_cb(j), y_cb(k), pb_ts(1)%sf(j, k, 0, r, i)
1126 end do
1127 end do
1128 close (2)
1129 end do
1130 end do
1131 do i = 1, nb
1132 do r = 1, nnode
1133 write (file_path, '(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/mv.', i, '.', r, '.', proc_rank, &
1134 & '.', t_step, '.dat'
1135
1136 open (2, file=trim(file_path))
1137 do j = 0, m
1138 do k = 0, n
1139 write (2, fmt) x_cb(j), y_cb(k), mv_ts(1)%sf(j, k, 0, r, i)
1140 end do
1141 end do
1142 close (2)
1143 end do
1144 end do
1145 end if
1146
1147 if (prim_vars_wrt .and. (.not. igr)) then
1148 do i = 1, sys_size
1149 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/prim.', i, '.', proc_rank, '.', t_step, '.dat'
1150
1151 open (2, file=trim(file_path))
1152
1153 do j = 0, m
1154 do k = 0, n
1155 if (((i >= eqn_idx%cont%beg) .and. (i <= eqn_idx%cont%end)) .or. ((i >= eqn_idx%adv%beg) &
1156 & .and. (i <= eqn_idx%adv%end))) then
1157 write (2, fmt) x_cb(j), y_cb(k), q_cons_vf(i)%sf(j, k, 0)
1158 else
1159 write (2, fmt) x_cb(j), y_cb(k), q_prim_vf(i)%sf(j, k, 0)
1160 end if
1161 end do
1162 write (2, *)
1163 end do
1164 close (2)
1165 end do
1166 end if
1167 end if
1168
1169 if (precision == precision_single) then
1170 fmt = "(4F30.7)"
1171 else
1172 fmt = "(4F40.14)"
1173 end if
1174
1175 if (p > 0) then
1176 do i = 1, sys_size
1177 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/cons.', i, '.', proc_rank, '.', t_step, '.dat'
1178 open (2, file=trim(file_path))
1179 do j = 0, m
1180 do k = 0, n
1181 do l = 0, p
1182 write (2, fmt) x_cb(j), y_cb(k), z_cb(l), q_cons_vf(i)%sf(j, k, l)
1183 end do
1184 write (2, *)
1185 end do
1186 write (2, *)
1187 end do
1188 close (2)
1189 end do
1190
1191 if (present(beta)) then
1192 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/beta.', i, '.', proc_rank, '.', t_step, '.dat'
1193 open (2, file=trim(file_path))
1194 do j = 0, m
1195 do k = 0, n
1196 do l = 0, p
1197 write (2, fmt) x_cb(j), y_cb(k), z_cb(l), beta%sf(j, k, l)
1198 end do
1199 write (2, *)
1200 end do
1201 write (2, *)
1202 end do
1203 close (2)
1204 end if
1205
1206 if (qbmm .and. .not. polytropic) then
1207 do i = 1, nb
1208 do r = 1, nnode
1209 write (file_path, '(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/pres.', i, '.', r, '.', proc_rank, &
1210 & '.', t_step, '.dat'
1211
1212 open (2, file=trim(file_path))
1213 do j = 0, m
1214 do k = 0, n
1215 do l = 0, p
1216 write (2, fmt) x_cb(j), y_cb(k), z_cb(l), pb_ts(1)%sf(j, k, l, r, i)
1217 end do
1218 end do
1219 end do
1220 close (2)
1221 end do
1222 end do
1223 do i = 1, nb
1224 do r = 1, nnode
1225 write (file_path, '(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/mv.', i, '.', r, '.', proc_rank, &
1226 & '.', t_step, '.dat'
1227
1228 open (2, file=trim(file_path))
1229 do j = 0, m
1230 do k = 0, n
1231 do l = 0, p
1232 write (2, fmt) x_cb(j), y_cb(k), z_cb(l), mv_ts(1)%sf(j, k, l, r, i)
1233 end do
1234 end do
1235 end do
1236 close (2)
1237 end do
1238 end do
1239 end if
1240
1241 if (prim_vars_wrt .and. (.not. igr)) then
1242 do i = 1, sys_size
1243 write (file_path, '(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/prim.', i, '.', proc_rank, '.', t_step, '.dat'
1244
1245 open (2, file=trim(file_path))
1246
1247 do j = 0, m
1248 do k = 0, n
1249 do l = 0, p
1250 if (((i >= eqn_idx%cont%beg) .and. (i <= eqn_idx%cont%end)) .or. ((i >= eqn_idx%adv%beg) &
1251 & .and. (i <= eqn_idx%adv%end)) .or. ((i >= eqn_idx%species%beg) &
1252 & .and. (i <= eqn_idx%species%end))) then
1253 write (2, fmt) x_cb(j), y_cb(k), z_cb(l), q_cons_vf(i)%sf(j, k, l)
1254 else
1255 write (2, fmt) x_cb(j), y_cb(k), z_cb(l), q_prim_vf(i)%sf(j, k, l)
1256 end if
1257 end do
1258 write (2, *)
1259 end do
1260 write (2, *)
1261 end do
1262 close (2)
1263 end do
1264 end if
1265 end if
1266
1267 end subroutine s_write_serial_data_files
1268
1269 !> Write grid and conservative variable data files in parallel via MPI I/O
1270 impure subroutine s_write_parallel_data_files(q_cons_vf, t_step, bc_type, beta, q_T_sf)
1271
1272 type(scalar_field), dimension(sys_size), intent(inout) :: q_cons_vf
1273 integer, intent(in) :: t_step
1274 type(scalar_field), intent(inout), optional :: beta
1275 type(integer_field), dimension(1:num_dims,-1:1), intent(in) :: bc_type
1276 type(scalar_field), intent(inout), optional :: q_t_sf
1277
1278#ifdef MFC_MPI
1279 integer :: ifile, ierr, data_size
1280 integer, dimension(MPI_STATUS_SIZE) :: status
1281 integer(kind=MPI_OFFSET_kind) :: disp
1282 integer(kind=MPI_OFFSET_kind) :: m_mok, n_mok, p_mok
1283 integer(kind=MPI_OFFSET_kind) :: wp_mok, var_mok, str_mok
1284 integer(kind=MPI_OFFSET_kind) :: nvars_mok
1285 integer(kind=MPI_OFFSET_kind) :: mok
1286 character(LEN=path_len + 2*name_len) :: file_loc
1287 logical :: file_exist, dir_check
1288 character(len=10) :: t_step_string
1289 integer :: i !< Generic loop iterator
1290 integer :: alt_sys !< Altered system size for the lagrangian subgrid bubble model
1291 ! Down sampling variables
1292 integer :: m_ds, n_ds, p_ds
1293 integer :: m_glb_ds, n_glb_ds, p_glb_ds
1294 integer :: m_glb_save, n_glb_save, p_glb_save !< Global save size
1295
1296 if (down_sample) then
1297 call s_downsample_data(q_cons_vf, q_cons_temp_ds, m_ds, n_ds, p_ds, m_glb_ds, n_glb_ds, p_glb_ds)
1298 end if
1299
1300 if (present(beta)) then
1301 alt_sys = sys_size + 1
1302 else
1303 alt_sys = sys_size
1304 end if
1305
1306 if (file_per_process) then
1307 call s_int_to_str(t_step, t_step_string)
1308
1309 if (down_sample) then
1310 call s_initialize_mpi_data_ds(m_ds, n_ds, p_ds)
1311 else
1312 if (ib) then
1313 call s_initialize_mpi_data(q_cons_vf, ib_markers=ib_markers, ib_mpi_data=mpi_io_ib_data, qbmm_pb=pb_ts(1), &
1314 & qbmm_mv=mv_ts(1))
1315 else
1316 call s_initialize_mpi_data(q_cons_vf, qbmm_pb=pb_ts(1), qbmm_mv=mv_ts(1))
1317 end if
1318 end if
1319
1320 if (proc_rank == 0) then
1321 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string)
1322 call my_inquire(file_loc, dir_check)
1323 if (dir_check .neqv. .true.) then
1324 call s_create_directory(trim(file_loc))
1325 end if
1326 call s_create_directory(trim(file_loc))
1327 end if
1328 call s_mpi_barrier()
1329 call s_delay_file_access(proc_rank)
1330
1331 call s_initialize_mpi_data(q_cons_vf, qbmm_pb=pb_ts(1), qbmm_mv=mv_ts(1))
1332
1333 write (file_loc, '(I0,A,i7.7,A)') t_step, '_', proc_rank, '.dat'
1334 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string) // trim(mpiiofs) // trim(file_loc)
1335 inquire (file=trim(file_loc), exist=file_exist)
1336 if (file_exist .and. proc_rank == 0) then
1337 call mpi_file_delete(file_loc, mpi_info_int, ierr)
1338 end if
1339 call mpi_file_open(mpi_comm_self, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1340 call s_check_mpi_file_open(ierr, file_loc)
1341
1342 if (down_sample) then
1343 data_size = (m_ds + 3)*(n_ds + 3)*(p_ds + 3)
1344 m_glb_save = m_glb_ds + 1
1345 n_glb_save = n_glb_ds + 1
1346 p_glb_save = p_glb_ds + 1
1347 else
1348 data_size = (m + 1)*(n + 1)*(p + 1)
1349 m_glb_save = m_glb + 1
1350 n_glb_save = n_glb + 1
1351 p_glb_save = p_glb + 1
1352 end if
1353
1354 m_mok = int(m_glb_save + 1, mpi_offset_kind)
1355 n_mok = int(n_glb_save + 1, mpi_offset_kind)
1356 p_mok = int(p_glb_save + 1, mpi_offset_kind)
1357 wp_mok = int(storage_size(0._stp)/8, mpi_offset_kind)
1358 mok = int(1._wp, mpi_offset_kind)
1359 str_mok = int(name_len, mpi_offset_kind)
1360 nvars_mok = int(sys_size, mpi_offset_kind)
1361
1362 if (bubbles_euler) then
1363 do i = 1, sys_size
1364 var_mok = int(i, mpi_offset_kind)
1365
1366 call mpi_file_write_all(ifile, mpi_io_data%var(i)%sf(:,:,:), data_size*mpi_io_type, mpi_io_p, status, ierr)
1367 end do
1368 if (qbmm .and. .not. polytropic) then
1369 do i = sys_size + 1, sys_size + 2*nb*nnode
1370 var_mok = int(i, mpi_offset_kind)
1371
1372 call mpi_file_write_all(ifile, mpi_io_data%var(i)%sf(:,:,:), data_size*mpi_io_type, mpi_io_p, status, ierr)
1373 end do
1374 end if
1375 else
1376 if (down_sample) then
1377 do i = 1, sys_size ! TODO: check if sys_size is correct
1378 var_mok = int(i, mpi_offset_kind)
1379
1380 call mpi_file_write_all(ifile, q_cons_temp_ds(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1381 end do
1382 else
1383 do i = 1, sys_size ! TODO: check if sys_size is correct
1384 var_mok = int(i, mpi_offset_kind)
1385
1386 call mpi_file_write_all(ifile, mpi_io_data%var(i)%sf(:,:,:), data_size*mpi_io_type, mpi_io_p, status, ierr)
1387 end do
1388 end if
1389 end if
1390
1391 call mpi_file_close(ifile, ierr)
1392
1393 if (ib) then
1394 call s_write_parallel_ib_data(t_step)
1395 end if
1396 else
1397 if (ib) then
1398 call s_initialize_mpi_data(q_cons_vf, ib_markers=ib_markers, ib_mpi_data=mpi_io_ib_data, qbmm_pb=pb_ts(1), &
1399 & qbmm_mv=mv_ts(1))
1400 else if (present(beta)) then
1401 call s_initialize_mpi_data(q_cons_vf, beta=beta, qbmm_pb=pb_ts(1), qbmm_mv=mv_ts(1))
1402 else
1403 call s_initialize_mpi_data(q_cons_vf, qbmm_pb=pb_ts(1), qbmm_mv=mv_ts(1))
1404 end if
1405
1406 write (file_loc, '(I0,A)') t_step, '.dat'
1407 file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // trim(file_loc)
1408 inquire (file=trim(file_loc), exist=file_exist)
1409 if (file_exist .and. proc_rank == 0) then
1410 call mpi_file_delete(file_loc, mpi_info_int, ierr)
1411 end if
1412 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1413 call s_check_mpi_file_open(ierr, file_loc)
1414
1415 data_size = (m + 1)*(n + 1)*(p + 1)
1416
1417 m_mok = int(m_glb + 1, mpi_offset_kind)
1418 n_mok = int(n_glb + 1, mpi_offset_kind)
1419 p_mok = int(p_glb + 1, mpi_offset_kind)
1420 wp_mok = int(storage_size(0._stp)/8, mpi_offset_kind)
1421 mok = int(1._wp, mpi_offset_kind)
1422 str_mok = int(name_len, mpi_offset_kind)
1423 nvars_mok = int(alt_sys, mpi_offset_kind)
1424
1425 if (bubbles_euler) then
1426 do i = 1, sys_size
1427 var_mok = int(i, mpi_offset_kind)
1428
1429 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1430
1431 call mpi_file_set_view(ifile, disp, mpi_p, mpi_io_data%view(i), 'native', mpi_info_int, ierr)
1432 call mpi_file_write_all(ifile, mpi_io_data%var(i)%sf(:,:,:), data_size*mpi_io_type, mpi_io_p, status, ierr)
1433 end do
1434 if (qbmm .and. .not. polytropic) then
1435 do i = sys_size + 1, sys_size + 2*nb*nnode
1436 var_mok = int(i, mpi_offset_kind)
1437
1438 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1439
1440 call mpi_file_set_view(ifile, disp, mpi_p, mpi_io_data%view(i), 'native', mpi_info_int, ierr)
1441 call mpi_file_write_all(ifile, mpi_io_data%var(i)%sf(:,:,:), data_size*mpi_io_type, mpi_io_p, status, ierr)
1442 end do
1443 end if
1444 else
1445 do i = 1, sys_size ! TODO: check if sys_size is correct
1446 var_mok = int(i, mpi_offset_kind)
1447
1448 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1449
1450 call mpi_file_set_view(ifile, disp, mpi_p, mpi_io_data%view(i), 'native', mpi_info_int, ierr)
1451 call mpi_file_write_all(ifile, mpi_io_data%var(i)%sf(:,:,:), data_size*mpi_io_type, mpi_io_p, status, ierr)
1452 end do
1453 end if
1454
1455 if (present(beta)) then
1456 var_mok = int(sys_size + 1, mpi_offset_kind)
1457
1458 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1459
1460 call mpi_file_set_view(ifile, disp, mpi_p, mpi_io_data%view(sys_size + 1), 'native', mpi_info_int, ierr)
1461 call mpi_file_write_all(ifile, mpi_io_data%var(sys_size + 1)%sf(:,:,:), data_size*mpi_io_type, mpi_io_p, status, &
1462 & ierr)
1463 end if
1464
1465 call mpi_file_close(ifile, ierr)
1466
1467 if (ib) then
1468 call s_write_parallel_ib_data(t_step)
1469 end if
1470 end if
1471#endif
1472
1473 end subroutine s_write_parallel_data_files
1474
1475 !> Write immersed boundary marker data to a serial (per-processor) unformatted file
1476 subroutine s_write_serial_ib_data(time_step)
1477
1478 integer, intent(in) :: time_step
1479 character(LEN=path_len + 2*name_len) :: file_path
1480 character(LEN=path_len + 2*name_len) :: t_step_dir
1481
1482 write (t_step_dir, '(A,I0,A,I0)') trim(case_dir) // '/p_all'
1483 write (t_step_dir, '(a,i0,a,i0)') trim(case_dir) // '/p_all/p', proc_rank, '/', time_step
1484 write (file_path, '(A,I0,A)') trim(t_step_dir) // '/ib_data.dat'
1485
1486 open (2, file=trim(file_path), form='unformatted', status='new')
1487
1488
1489# 1051 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1490#if defined(MFC_OpenACC)
1491# 1051 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1492!$acc update host(ib_markers%sf)
1493# 1051 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1494#elif defined(MFC_OpenMP)
1495# 1051 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1496!$omp target update from(ib_markers%sf)
1497# 1051 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1498#endif
1499# 1052 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1500 write (2) ib_markers%sf(0:m,0:n,0:p); close (2)
1501
1502 end subroutine s_write_serial_ib_data
1503
1504 !> Write immersed boundary marker data in parallel using MPI I/O
1505 subroutine s_write_parallel_ib_data(time_step)
1506
1507 integer, intent(in) :: time_step
1508
1509#ifdef MFC_MPI
1510 character(LEN=path_len + 2*name_len) :: file_loc
1511 integer(kind=MPI_OFFSET_kind) :: disp
1512 integer(kind=MPI_OFFSET_kind) :: m_MOK, n_MOK, p_MOK
1513 integer(kind=MPI_OFFSET_kind) :: WP_MOK, var_MOK, MOK
1514 integer :: ifile, ierr, data_size
1515 integer, dimension(MPI_STATUS_SIZE) :: status
1516 character(len=10) :: t_step_string
1517
1518
1519# 1070 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1520#if defined(MFC_OpenACC)
1521# 1070 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1522!$acc update host(ib_markers%sf)
1523# 1070 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1524#elif defined(MFC_OpenMP)
1525# 1070 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1526!$omp target update from(ib_markers%sf)
1527# 1070 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1528#endif
1529# 1071 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1530
1531 data_size = (m + 1)*(n + 1)*(p + 1)
1532 m_mok = int(m_glb + 1, mpi_offset_kind)
1533 n_mok = int(n_glb + 1, mpi_offset_kind)
1534 p_mok = int(p_glb + 1, mpi_offset_kind)
1535 wp_mok = int(storage_size(0._stp)/8, mpi_offset_kind)
1536 mok = int(1._wp, mpi_offset_kind)
1537
1538 if (file_per_process) then
1539 call s_int_to_str(time_step, t_step_string)
1540
1541 if (proc_rank == 0) then
1542 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string)
1543 call s_create_directory(trim(file_loc))
1544 end if
1545 call s_mpi_barrier()
1546 call s_delay_file_access(proc_rank)
1547
1548 write (file_loc, '(A,I0,A,i7.7,A)') 'ib_markers_', time_step, '_', proc_rank, '.dat'
1549 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string) // '/' // trim(file_loc)
1550
1551 call mpi_file_open(mpi_comm_self, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1552 call s_check_mpi_file_open(ierr, file_loc)
1553 call mpi_file_write_all(ifile, mpi_io_ib_data%var%sf, data_size, mpi_integer, status, ierr)
1554 call mpi_file_close(ifile, ierr)
1555 else
1556 write (file_loc, '(A)') 'ib.dat'
1557 file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // trim(file_loc)
1558
1559 call s_mpi_barrier()
1560 call s_delay_file_access(proc_rank)
1561
1562 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1563 call s_check_mpi_file_open(ierr, file_loc)
1564
1565 var_mok = int(sys_size + 1, mpi_offset_kind)
1566 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1 + int(time_step/t_step_save))
1567 if (time_step == 0) disp = 0
1568
1569 call mpi_file_set_view(ifile, disp, mpi_integer, mpi_io_ib_data%view, 'native', mpi_info_int, ierr)
1570 call mpi_file_write_all(ifile, mpi_io_ib_data%var%sf, data_size, mpi_integer, status, ierr)
1571 call mpi_file_close(ifile, ierr)
1572 end if
1573#endif
1574
1575 end subroutine s_write_parallel_ib_data
1576
1577 !> Dispatch immersed boundary data output to the serial or parallel writer
1578 subroutine s_write_ib_data_file(time_step)
1579
1580 integer, intent(in) :: time_step
1581
1582
1583# 1123 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1584#if defined(MFC_OpenACC)
1585# 1123 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1586!$acc update host(patch_ib(1:num_ibs))
1587# 1123 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1588#elif defined(MFC_OpenMP)
1589# 1123 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1590!$omp target update from(patch_ib(1:num_ibs))
1591# 1123 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1592#endif
1593# 1124 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1594
1595 if (parallel_io) then
1596 call s_write_parallel_ib_data(time_step)
1597 else
1598 call s_write_serial_ib_data(time_step)
1599 end if
1600
1601 end subroutine s_write_ib_data_file
1602
1603 !> Writes the IB state information out to file
1604 subroutine s_write_parallel_ib_state(t_step)
1605
1606 integer, intent(in) :: t_step
1607
1608#ifdef MFC_MPI
1609 character(LEN=path_len + 2*name_len) :: file_loc
1610 integer(kind=MPI_OFFSET_KIND) :: disp
1611 integer(kind=MPI_OFFSET_KIND) :: WP_MOK
1612 integer :: ifile, ierr
1613 integer, dimension(MPI_STATUS_SIZE) :: status
1614 logical :: file_exist, dir_check
1615 integer :: i, ib_idx
1616 integer, parameter :: NFIELDS_PER_IB = 20
1617 real(wp) :: ib_buf(NFIELDS_PER_IB)
1618 integer :: file_unit
1619 character(len=10) :: t_step_string
1620
1621 ! Partition IBs across ranks round-robin style
1622 integer :: ib_start, ib_end, nibs_per_rank, remainder
1623
1624 wp_mok = int(storage_size(0._wp)/8, mpi_offset_kind)
1625
1626 if (file_per_process) then
1627 call s_int_to_str(t_step, t_step_string)
1628
1629 if (proc_rank == 0) then
1630 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string)
1631 call s_create_directory(trim(file_loc))
1632 end if
1633 call s_mpi_barrier()
1634 call s_delay_file_access(proc_rank)
1635
1636 write (file_loc, '(A,I0,A,i7.7,A)') 'ib_state_', t_step, '_', proc_rank, '.dat'
1637 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string) // '/' // trim(file_loc)
1638
1639 inquire (file=trim(file_loc), exist=file_exist)
1640 if (file_exist) then
1641 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='replace')
1642 else
1643 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='new')
1644 end if
1645
1646 write (file_unit) num_local_ibs
1647 do i = 1, num_local_ibs
1648 ib_idx = local_ib_patch_ids(i)
1649 ib_buf(1) = mytime
1650 ib_buf(2:4) = patch_ib(ib_idx)%force(1:3)
1651 ib_buf(5:7) = patch_ib(ib_idx)%torque(1:3)
1652 ib_buf(8:10) = patch_ib(ib_idx)%vel(1:3)
1653 ib_buf(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1654 ib_buf(14:16) = patch_ib(ib_idx)%angles(1:3)
1655 ib_buf(17) = patch_ib(ib_idx)%x_centroid
1656 ib_buf(18) = patch_ib(ib_idx)%y_centroid
1657 ib_buf(19) = patch_ib(ib_idx)%z_centroid
1658 ib_buf(20) = patch_ib(ib_idx)%radius
1659
1660 write (file_unit) patch_ib(ib_idx)%gbl_patch_id
1661 write (file_unit) ib_buf
1662 end do
1663
1664 close (file_unit)
1665 else
1666 if (proc_rank == 0) then
1667 call s_create_directory(trim(case_dir) // '/restart_data')
1668 end if
1669 call s_mpi_barrier()
1670
1671 write (file_loc, '(A,I0,A)') '/restart_data/ib_state_', t_step, '.dat'
1672 file_loc = trim(case_dir) // trim(file_loc)
1673
1674 inquire (file=trim(file_loc), exist=file_exist)
1675 if (file_exist .and. proc_rank == 0) then
1676 call mpi_file_delete(file_loc, mpi_info_int, ierr)
1677 end if
1678 call s_mpi_barrier()
1679
1680 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1681 call s_check_mpi_file_open(ierr, file_loc)
1682
1683 do i = 1, num_local_ibs
1684 ib_idx = local_ib_patch_ids(i)
1685 ib_buf(1) = mytime
1686 ib_buf(2:4) = patch_ib(ib_idx)%force(1:3)
1687 ib_buf(5:7) = patch_ib(ib_idx)%torque(1:3)
1688 ib_buf(8:10) = patch_ib(ib_idx)%vel(1:3)
1689 ib_buf(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1690 ib_buf(14:16) = patch_ib(ib_idx)%angles(1:3)
1691 ib_buf(17) = patch_ib(ib_idx)%x_centroid
1692 ib_buf(18) = patch_ib(ib_idx)%y_centroid
1693 ib_buf(19) = patch_ib(ib_idx)%z_centroid
1694 ib_buf(20) = patch_ib(ib_idx)%radius
1695
1696 ! Global IB index determines position in file
1697 disp = int(patch_ib(ib_idx)%gbl_patch_id - 1, mpi_offset_kind)*int(nfields_per_ib, mpi_offset_kind)*wp_mok
1698
1699 call mpi_file_write_at(ifile, disp, ib_buf, nfields_per_ib, mpi_p, status, ierr)
1700 end do
1701
1702 call mpi_file_close(ifile, ierr)
1703 end if
1704#endif
1705
1706 end subroutine s_write_parallel_ib_state
1707
1708 !> Write IB state data to a per-timestep serial (unformatted) file
1709 subroutine s_write_serial_ib_state(t_step)
1710
1711 integer, intent(in) :: t_step
1712 character(LEN=path_len + 2*name_len) :: file_loc
1713 integer :: i, ios, file_unit
1714 integer, parameter :: NFIELDS_PER_IB = 20
1715 real(wp) :: ib_buf(NFIELDS_PER_IB)
1716
1717 call s_create_directory(trim(case_dir) // '/restart_data')
1718
1719 write (file_loc, '(A,I0,A)') '/restart_data/ib_state_', t_step, '.dat'
1720 file_loc = trim(case_dir) // trim(file_loc)
1721
1722 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='replace', iostat=ios)
1723 if (ios /= 0) call s_mpi_abort('Cannot open IB state output file: ' // trim(file_loc))
1724
1725 do i = 1, num_ibs
1726 ib_buf(1) = mytime
1727 ib_buf(2:4) = patch_ib(i)%force(1:3)
1728 ib_buf(5:7) = patch_ib(i)%torque(1:3)
1729 ib_buf(8:10) = patch_ib(i)%vel(1:3)
1730 ib_buf(11:13) = patch_ib(i)%angular_vel(1:3)
1731 ib_buf(14:16) = patch_ib(i)%angles(1:3)
1732 ib_buf(17) = patch_ib(i)%x_centroid
1733 ib_buf(18) = patch_ib(i)%y_centroid
1734 ib_buf(19) = patch_ib(i)%z_centroid
1735 ib_buf(20) = patch_ib(i)%radius
1736
1737 write (file_unit) ib_buf
1738 end do
1739
1740 close (file_unit)
1741
1742 end subroutine s_write_serial_ib_state
1743
1744 !> Record every immersed body's force, torque and kinematics for this step.
1745 !!
1746 !! One shared text file, D/ib_forces.dat, opened once for the run. Each rank writes only the
1747 !! bodies it owns, at a byte offset computed from the step and the global body id, so the file
1748 !! is byte-identical however the domain is decomposed and needs no merge step. Writing a file
1749 !! per body instead costs an inquire, open and close per body per rank per step, which is
1750 !! 3e5 filesystem metadata operations per step at 1000 ranks holding 100 bodies each.
1751 !!
1752 !! Layout: row r holds all num_gbl_ibs bodies in global id order, so body g occupies bytes
1753 !! ((r*num_gbl_ibs) + g - 1)*IB_REC_LEN. The file carries no header line, which would shift
1754 !! every offset after it; the columns are listed in docs/documentation/case.md.
1755 !!
1756 !! Rows are numbered from the first step this run records, not from t_step, so that row 0 is
1757 !! always written. Nothing pre-fills the file -- both open paths create it empty and every
1758 !! write lands at a computed offset -- so a row no rank ever writes stays a hole, and a hole
1759 !! reads back as NUL bytes rather than blanks. Counting from t_step would leave exactly such a
1760 !! hole wherever the run starts: the skip below means step t_step_start is never written, and
1761 !! on a restart every row beneath it would be missing as well.
1762 impure subroutine s_write_ib_force_history(t_step)
1763
1764 integer, intent(in) :: t_step
1765 character(LEN=IB_REC_LEN) :: rec
1766 character(LEN=path_len + 2*name_len) :: file_loc
1767 real(wp) :: fields(19)
1768 integer :: i, ib_idx, n_write, row
1769
1770#ifdef MFC_MPI
1771 integer(kind=MPI_OFFSET_KIND) :: disp
1772 integer :: ierr, status(mpi_status_size)
1773#endif
1774
1775 if (.not. ib_force_wrt) return
1776 if (mod(t_step, max(ib_force_stride, 1)) /= 0) return
1777 ! This runs at RK stage 1, before the step's force has been computed, so the row for step N carries the
1778 ! force from the end of step N-1. The first step of a run has no N-1: patch_ib%force is still zero and
1779 ! the row would record identically zero force. That is not a measurement, and on a run chained across a
1780 ! queue's walltime limit it lands once per restart -- in a six-wingbeat case, zeros at steps 20649,
1781 ! 34649 and 48649 sitting among neighbours of -0.134, +0.474 and -0.475, corrupting every per-beat
1782 ! trough taken over the joined trace.
1783 if (t_step == t_step_start) return
1784
1785 n_write = num_local_ibs
1786 if (num_procs == 1) n_write = num_ibs
1787 ! Relative to the first recorded step, so the first one written is row 0 and the file is
1788 ! dense. See the hole discussion above.
1789 row = t_step/max(ib_force_stride, 1) - t_step_start/max(ib_force_stride, 1) - 1
1790
1791
1792# 1321 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1793#if defined(MFC_OpenACC)
1794# 1321 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1795!$acc update host(patch_ib(1:num_ibs))
1796# 1321 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1797#elif defined(MFC_OpenMP)
1798# 1321 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1799!$omp target update from(patch_ib(1:num_ibs))
1800# 1321 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1801#endif
1802# 1322 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1803
1804 call s_open_ib_force_history()
1805
1806 do i = 1, n_write
1807 ib_idx = i
1808 if (num_procs > 1) ib_idx = local_ib_patch_ids(i)
1809
1810 fields(1) = mytime
1811 fields(2:4) = patch_ib(ib_idx)%force(1:3)
1812 fields(5:7) = patch_ib(ib_idx)%torque(1:3)
1813 fields(8:10) = patch_ib(ib_idx)%vel(1:3)
1814 fields(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1815 fields(14:16) = patch_ib(ib_idx)%angles(1:3)
1816 fields(17) = patch_ib(ib_idx)%x_centroid
1817 fields(18) = patch_ib(ib_idx)%y_centroid
1818 fields(19) = patch_ib(ib_idx)%z_centroid
1819
1820 write (rec, ib_rec_fmt) patch_ib(ib_idx)%gbl_patch_id, fields
1821 rec(ib_rec_len:ib_rec_len) = new_line('a')
1822
1823#ifdef MFC_MPI
1824 disp = (int(row, mpi_offset_kind)*int(num_gbl_ibs, mpi_offset_kind) + int(patch_ib(ib_idx)%gbl_patch_id - 1, &
1825 & mpi_offset_kind))*int(ib_rec_len, mpi_offset_kind)
1826 call mpi_file_write_at(ib_hist_file, disp, rec, ib_rec_len, mpi_character, status, ierr)
1827#else
1828 write (ib_hist_file, rec=row*num_gbl_ibs + patch_ib(ib_idx)%gbl_patch_id) rec
1829#endif
1830 end do
1831
1832 end subroutine s_write_ib_force_history
1833
1834 !> Open the shared history file. Done once for the run.
1835 impure subroutine s_open_ib_force_history
1836
1837 character(LEN=path_len + 2*name_len) :: file_loc
1838 character(LEN=IB_REC_LEN) :: probe
1839 integer :: i
1840
1841#ifdef MFC_MPI
1842 integer :: ierr
1843 logical :: file_exist
1844#endif
1845
1846 if (ib_hist_file /= -1) return
1847
1848 ! Every offset below assumes the format emits exactly IB_REC_BODY characters. Measure it once
1849 ! rather than trusting that the format and the constant were edited together: a format one
1850 ! character wider would shear every record past the first without any other symptom.
1851 write (probe, ib_rec_fmt) 0, [(0._wp, i=1, 19)]
1852 if (len_trim(probe) /= ib_rec_body) then
1853# 1371 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1854 call s_prohibit_abort("len_trim(probe) /= IB_REC_BODY", "IB force record width disagrees with IB_REC_BODY; IB_REC_FMT and IB_REC_BODY must be changed together")
1855# 1371 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1856 end if
1857# 1373 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1858
1859 file_loc = trim(case_dir) // '/D/ib_forces.dat'
1860#ifdef MFC_MPI
1861 ! MPI_MODE_CREATE does not truncate, so a shorter run following a longer one in the same
1862 ! directory would keep the old tail past its last record. Delete first, as the ib_state
1863 ! writer does, then barrier so no rank opens before the delete lands.
1864 inquire (file=trim(file_loc), exist=file_exist)
1865 if (file_exist .and. proc_rank == 0) call mpi_file_delete(file_loc, mpi_info_null, ierr)
1866
1867 ! MPI_INFO_NULL, not mpi_info_int: the latter is only created when parallel_io is on
1868 ! (m_global_parameters_common.fpp returns before MPI_INFO_CREATE otherwise), and the
1869 ! force history is written whatever parallel_io is set to. Passing the uninitialised
1870 ! handle aborted every IBM case that runs the solver with parallel_io = F, in
1871 ! MPI_Info_dup, at the first recorded step. The hint it carries only disables ROMIO
1872 ! write data sieving, which this writer does not depend on.
1873 ! Collective: every rank opens, including one holding no body this step.
1874 call s_mpi_barrier()
1875 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_null, ib_hist_file, ierr)
1876 call s_check_mpi_file_open(ierr, file_loc)
1877#else
1878 ! Unformatted: the record is already a formatted string, so this writes its bytes verbatim and
1879 ! produces the same file the MPI branch does. A formatted direct-access write would need a
1880 ! format and would pad rather than emit the string as-is.
1881 open (newunit=ib_hist_file, file=trim(file_loc), form='unformatted', access='direct', recl=ib_rec_len, status='replace')
1882#endif
1883
1884 end subroutine s_open_ib_force_history
1885
1886 !> Close the history file. Nothing is buffered, so there is nothing to flush first.
1887 impure subroutine s_close_ib_force_history
1888
1889#ifdef MFC_MPI
1890 integer :: ierr
1891#endif
1892
1893 if (ib_hist_file == -1) return
1894#ifdef MFC_MPI
1895 call mpi_file_close(ib_hist_file, ierr)
1896#else
1897 close (ib_hist_file)
1898#endif
1899 ib_hist_file = -1
1900
1901 end subroutine s_close_ib_force_history
1902
1903 !> @brief Writes IB state records to restart_data/ib_state.dat. Called on every rank.
1904 impure subroutine s_write_ib_state_file(time_step)
1905
1906 integer, intent(in) :: time_step
1907
1908
1909# 1423 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1910#if defined(MFC_OpenACC)
1911# 1423 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1912!$acc update host(patch_ib(1:num_ibs))
1913# 1423 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1914#elif defined(MFC_OpenMP)
1915# 1423 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1916!$omp target update from(patch_ib(1:num_ibs))
1917# 1423 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1918#endif
1919# 1424 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1920
1921 if (parallel_io) then
1922 call s_write_parallel_ib_state(time_step)
1923 else
1924 call s_write_serial_ib_state(time_step)
1925 end if
1926 call s_write_centroid_offsets(time_step)
1927
1928 end subroutine s_write_ib_state_file
1929
1930 !> Write each global patch's centroid_offset to restart_data/ib_offset_<step>.dat, so a restart continues about the same centre
1931 !! of mass instead of re-measuring it from the body voxelised at the restart attitude. Collective.
1932 impure subroutine s_write_centroid_offsets(step)
1933
1934 integer, intent(in) :: step
1935 character(len=path_len + 2*name_len) :: file_loc
1936 real(wp), dimension(:,:), allocatable :: off_loc, off_glb
1937 integer :: gid, i, ib_idx, n_own, file_unit
1938
1939 if (.not. centroid_offsets_active) return
1940 ! column 4 flags a patch that has an offset; only the owning rank contributes, so the sum is the value
1941 allocate (off_loc(num_gbl_ibs, 4), off_glb(num_gbl_ibs, 4))
1942 off_loc = 0._wp
1943 n_own = num_local_ibs
1944 if (num_procs == 1) n_own = num_ibs
1945 do i = 1, n_own
1946 ib_idx = i
1947 if (num_procs > 1) ib_idx = local_ib_patch_ids(i)
1948 if (.not. f_needs_centroid_offset(patch_ib(ib_idx))) cycle
1949 off_loc(patch_ib(ib_idx)%gbl_patch_id,:) = [patch_ib(ib_idx)%centroid_offset, 1._wp]
1950 end do
1951 call s_mpi_allreduce_vectors_sum(off_loc, off_glb, num_gbl_ibs, 4)
1952 if (proc_rank == 0) then
1953 write (file_loc, '(A,I0,A)') trim(case_dir) // '/restart_data/ib_offset_', step, '.dat'
1954 open (newunit=file_unit, file=trim(file_loc), status='replace', action='write')
1955 do gid = 1, num_gbl_ibs
1956 if (off_glb(gid, 4) > 0.5_wp) write (file_unit, '(I0,3(1X,ES24.16))') gid, off_glb(gid,1:3)
1957 end do
1958 close (file_unit)
1959 end if
1960 deallocate (off_loc, off_glb)
1961
1962 end subroutine s_write_centroid_offsets
1963
1964 !> Write flow probe data at the current time step
1965 impure subroutine s_write_probe_files(t_step, q_cons_vf, accel_mag)
1966
1967 integer, intent(in) :: t_step
1968 type(scalar_field), dimension(sys_size), intent(in) :: q_cons_vf
1969 real(wp), dimension(0:m,0:n,0:p), intent(in) :: accel_mag
1970 real(wp), dimension(-1:m) :: distx
1971 real(wp), dimension(-1:n) :: disty
1972 real(wp), dimension(-1:p) :: distz
1973
1974 ! The cell-averaged partial densities, density, velocity, pressure, volume fractions, specific heat ratio function, liquid
1975 ! stiffness function, and sound speed.
1976 real(wp) :: lit_gamma, nbub
1977 real(wp) :: rho
1978 real(wp), dimension(num_vels) :: vel
1979 real(wp) :: pres
1980 real(wp) :: ptilde
1981 real(wp) :: ptot
1982 real(wp) :: alf
1983 real(wp) :: alfgr
1984 real(wp), dimension(num_fluids) :: alpha, alpha_rho
1985 real(wp) :: gamma
1986 real(wp) :: pi_inf
1987 real(wp) :: qv
1988 real(wp) :: c
1989 real(wp) :: m00, m10, m01, m20, m02
1990 real(wp) :: varr, varv
1991 real(wp), dimension(Nb) :: nr, r, nrdot, rdot
1992 real(wp) :: nr3
1993 real(wp) :: accel
1994 real(wp) :: int_pres
1995 real(wp) :: max_pres
1996 real(wp), dimension(2) :: re
1997 real(wp), dimension(6) :: tau_e
1998 real(wp) :: g_undamaged, g_damaged
1999 real(wp) :: dyn_p, t
2000 real(wp) :: damage_state
2001 real(wp) :: solid_partial_density !< damageable-solid partial density at the probe cell
2002 integer :: i, j, k, l, s, d !< Generic loop iterator
2003 real(wp) :: nondim_time !< Non-dimensional time
2004 real(wp) :: tmp !< Temporary variable to store quantity for mpi_allreduce
2005 real(wp) :: rhoyks(1:num_species)
2006
2007 t = dflt_t_guess
2008
2009 if (time_stepper == 23) then
2010 nondim_time = mytime
2011 else
2012 if (t_step_old /= dflt_int) then
2013 nondim_time = real(t_step + t_step_old, wp)*dt
2014 else
2015 nondim_time = real(t_step, wp)*dt
2016 end if
2017 end if
2018
2019 do i = 1, num_probes
2020 rho = 0._wp
2021 do s = 1, num_vels
2022 vel(s) = 0._wp
2023 end do
2024 pres = 0._wp
2025 gamma = 0._wp
2026 pi_inf = 0._wp
2027 qv = 0._wp
2028 c = 0._wp
2029 accel = 0._wp
2030 nr = 0._wp; r = 0._wp
2031 nrdot = 0._wp; rdot = 0._wp
2032 nbub = 0._wp
2033 m00 = 0._wp
2034 m10 = 0._wp
2035 m01 = 0._wp
2036 m20 = 0._wp
2037 m02 = 0._wp
2038 varr = 0._wp; varv = 0._wp
2039 alf = 0._wp
2040 do s = 1, (num_dims*(num_dims + 1))/2
2041 tau_e(s) = 0._wp
2042 end do
2043 damage_state = 0._wp
2044 g_damaged = 0._wp
2045
2046 if (n == 0) then
2047 if (f_probe_owned(probe(i)%x, x_cb(-1), x_cb(m), 1)) then
2048 do s = -1, m
2049 distx(s) = x_cb(s) - probe(i)%x
2050 if (distx(s) < 0._wp) distx(s) = 1000._wp
2051 end do
2052 j = minloc(distx, 1)
2053 if (j == 1) j = 2 ! Pick first point if probe is at edge
2054 k = 0
2055 l = 0
2056
2057 if (chemistry) then
2058 do d = 1, num_species
2059 rhoyks(d) = q_cons_vf(eqn_idx%species%beg + d - 1)%sf(j - 2, k, l)
2060 end do
2061 end if
2062
2063 ! Computing/Sharing necessary state variables
2064 if (hypoelasticity) then
2065 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k, l, rho, gamma, pi_inf, qv, re, g_undamaged, &
2066 & fluid_pp(:)%G)
2067 else
2068 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k, l, rho, gamma, pi_inf, qv)
2069 end if
2070 do s = 1, num_vels
2071 vel(s) = q_cons_vf(eqn_idx%cont%end + s)%sf(j - 2, k, l)/rho
2072 end do
2073 do s = 1, num_fluids
2074 alpha(s) = q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(j - 2, k, l)
2075 alpha_rho(s) = q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, k, l)
2076 end do
2077
2078 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
2079
2080 if (hypoelasticity) then
2081 if (cont_damage) then
2082 ! Recover D = U_D/m_s, clamped to [0, 1]
2083 solid_partial_density = 0._wp
2084 do s = 1, num_fluids
2085 if (fluid_pp(s)%G > verysmall) then
2086 solid_partial_density = solid_partial_density + q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, &
2087 & k, l)
2088 end if
2089 end do
2090 damage_state = min(max(q_cons_vf(eqn_idx%damage)%sf(j - 2, k, l)/max(solid_partial_density, &
2091 & verysmall), 0._wp), 1._wp)
2092 end if
2093 g_damaged = g_undamaged*max(1._wp - damage_state, 0._wp)
2094
2095 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k, l), q_cons_vf(eqn_idx%alf)%sf(j - 2, k, l), &
2096 & dyn_p, pi_inf, gamma, rho, qv, rhoyks(:), pres, t, &
2097 & f_hypoelastic_energy(q_cons_vf, j - 2, k, l, rho, g_undamaged))
2098 else
2099 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k, l), q_cons_vf(eqn_idx%alf)%sf(j - 2, k, l), &
2100 & dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t)
2101 end if
2102
2103 if (hypoelasticity) then
2104 tau_e(1) = q_cons_vf(eqn_idx%stress%end)%sf(j - 2, k, l)/rho
2105 end if
2106
2107 if (bubbles_euler) then
2108 alf = q_cons_vf(eqn_idx%alf)%sf(j - 2, k, l)
2109 if (num_fluids == 3) then
2110 alfgr = q_cons_vf(eqn_idx%alf - 1)%sf(j - 2, k, l)
2111 end if
2112 do s = 1, nb
2113 nr(s) = q_cons_vf(qbmm_idx%rs(s))%sf(j - 2, k, l)
2114 nrdot(s) = q_cons_vf(qbmm_idx%vs(s))%sf(j - 2, k, l)
2115 end do
2116
2117 if (adv_n) then
2118 nbub = q_cons_vf(eqn_idx%n)%sf(j - 2, k, l)
2119 else
2120 nr3 = 0._wp
2121 do s = 1, nb
2122 nr3 = nr3 + weight(s)*(nr(s)**3._wp)
2123 end do
2124
2125 nbub = sqrt((4._wp*pi/3._wp)*nr3/alf)
2126 end if
2127#ifdef MFC_DEBUG
2128 print *, 'In probe, nbub: ', nbub
2129#endif
2130 if (qbmm) then
2131 m00 = q_cons_vf(qbmm_idx%moms(1, 1))%sf(j - 2, k, l)/nbub
2132 m10 = q_cons_vf(qbmm_idx%moms(1, 2))%sf(j - 2, k, l)/nbub
2133 m01 = q_cons_vf(qbmm_idx%moms(1, 3))%sf(j - 2, k, l)/nbub
2134 m20 = q_cons_vf(qbmm_idx%moms(1, 4))%sf(j - 2, k, l)/nbub
2135 m02 = q_cons_vf(qbmm_idx%moms(1, 6))%sf(j - 2, k, l)/nbub
2136
2137 m10 = m10/m00
2138 m01 = m01/m00
2139 m20 = m20/m00
2140 m02 = m02/m00
2141
2142 varr = m20 - m10**2._wp
2143 varv = m02 - m01**2._wp
2144 end if
2145 r(:) = nr(:)/nbub
2146 rdot(:) = nrdot(:)/nbub
2147
2148 ptilde = ptil(j - 2, k, l)
2149 ptot = pres - ptilde
2150 end if
2151
2152 ! Compute mixture sound Speed
2153 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
2154 if (hypoelasticity) c = sqrt(c*c + (4._wp/3._wp)*g_damaged/rho)
2155
2156 accel = accel_mag(j - 2, k, l)
2157 end if
2158 else if (p == 0) then
2159 if (chemistry) then
2160 do d = 1, num_species
2161 rhoyks(d) = q_cons_vf(eqn_idx%species%beg + d - 1)%sf(j - 2, k - 2, l)
2162 end do
2163 end if
2164
2165 if (f_probe_owned(probe(i)%x, x_cb(-1), x_cb(m), 1)) then
2166 if (f_probe_owned(probe(i)%y, y_cb(-1), y_cb(n), 2)) then
2167 do s = -1, m
2168 distx(s) = x_cb(s) - probe(i)%x
2169 if (distx(s) < 0._wp) distx(s) = 1000._wp
2170 end do
2171 do s = -1, n
2172 disty(s) = y_cb(s) - probe(i)%y
2173 if (disty(s) < 0._wp) disty(s) = 1000._wp
2174 end do
2175 j = minloc(distx, 1)
2176 k = minloc(disty, 1)
2177 if (j == 1) j = 2 ! Pick first point if probe is at edge
2178 if (k == 1) k = 2 ! Pick first point if probe is at edge
2179 l = 0
2180
2181 ! Computing/Sharing necessary state variables
2182 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k - 2, l, rho, gamma, pi_inf, qv, re, g_undamaged, &
2183 & fluid_pp(:)%G)
2184 do s = 1, num_vels
2185 vel(s) = q_cons_vf(eqn_idx%cont%end + s)%sf(j - 2, k - 2, l)/rho
2186 end do
2187 do s = 1, num_fluids
2188 alpha(s) = q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(j - 2, k - 2, l)
2189 alpha_rho(s) = q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, k - 2, l)
2190 end do
2191
2192 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
2193
2194 if (hypoelasticity) then
2195 if (cont_damage) then
2196 ! Recover D = U_D/m_s, clamped to [0, 1]
2197 solid_partial_density = 0._wp
2198 do s = 1, num_fluids
2199 if (fluid_pp(s)%G > verysmall) then
2200 solid_partial_density = solid_partial_density + q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j &
2201 & - 2, k - 2, l)
2202 end if
2203 end do
2204 damage_state = min(max(q_cons_vf(eqn_idx%damage)%sf(j - 2, k - 2, l)/max(solid_partial_density, &
2205 & verysmall), 0._wp), 1._wp)
2206 end if
2207 g_damaged = g_undamaged*max(1._wp - damage_state, 0._wp)
2208
2209 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k - 2, l), q_cons_vf(eqn_idx%alf)%sf(j - 2, &
2210 & k - 2, l), dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t, &
2211 & f_hypoelastic_energy(q_cons_vf, j - 2, k - 2, l, rho, g_undamaged))
2212 else
2213 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k - 2, l), q_cons_vf(eqn_idx%alf)%sf(j - 2, &
2214 & k - 2, l), dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t)
2215 end if
2216
2217 if (hypoelasticity) then
2218 do s = 1, 3
2219 tau_e(s) = q_cons_vf(eqn_idx%stress%beg + s - 1)%sf(j - 2, k - 2, l)/rho
2220 end do
2221 end if
2222
2223 if (bubbles_euler) then
2224 alf = q_cons_vf(eqn_idx%alf)%sf(j - 2, k - 2, l)
2225 do s = 1, nb
2226 nr(s) = q_cons_vf(qbmm_idx%rs(s))%sf(j - 2, k - 2, l)
2227 nrdot(s) = q_cons_vf(qbmm_idx%vs(s))%sf(j - 2, k - 2, l)
2228 end do
2229
2230 if (adv_n) then
2231 nbub = q_cons_vf(eqn_idx%n)%sf(j - 2, k - 2, l)
2232 else
2233 nr3 = 0._wp
2234 do s = 1, nb
2235 nr3 = nr3 + weight(s)*(nr(s)**3._wp)
2236 end do
2237
2238 nbub = sqrt((4._wp*pi/3._wp)*nr3/alf)
2239 end if
2240
2241 r(:) = nr(:)/nbub
2242 rdot(:) = nrdot(:)/nbub
2243 end if
2244 ! Compute mixture sound speed
2245 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
2246 if (hypoelasticity) c = sqrt(c*c + (4._wp/3._wp)*g_damaged/rho)
2247 end if
2248 end if
2249 else
2250 if (f_probe_owned(probe(i)%x, x_cb(-1), x_cb(m), 1)) then
2251 if (f_probe_owned(probe(i)%y, y_cb(-1), y_cb(n), 2)) then
2252 if (f_probe_owned(probe(i)%z, z_cb(-1), z_cb(p), 3)) then
2253 do s = -1, m
2254 distx(s) = x_cb(s) - probe(i)%x
2255 if (distx(s) < 0._wp) distx(s) = 1000._wp
2256 end do
2257 do s = -1, n
2258 disty(s) = y_cb(s) - probe(i)%y
2259 if (disty(s) < 0._wp) disty(s) = 1000._wp
2260 end do
2261 do s = -1, p
2262 distz(s) = z_cb(s) - probe(i)%z
2263 if (distz(s) < 0._wp) distz(s) = 1000._wp
2264 end do
2265 j = minloc(distx, 1)
2266 k = minloc(disty, 1)
2267 l = minloc(distz, 1)
2268 if (j == 1) j = 2 ! Pick first point if probe is at edge
2269 if (k == 1) k = 2 ! Pick first point if probe is at edge
2270 if (l == 1) l = 2 ! Pick first point if probe is at edge
2271
2272 ! Computing/Sharing necessary state variables
2273 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k - 2, l - 2, rho, gamma, pi_inf, qv, re, &
2274 & g_undamaged, fluid_pp(:)%G)
2275 do s = 1, num_vels
2276 vel(s) = q_cons_vf(eqn_idx%cont%end + s)%sf(j - 2, k - 2, l - 2)/rho
2277 end do
2278 do s = 1, num_fluids
2279 alpha(s) = q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(j - 2, k - 2, l - 2)
2280 alpha_rho(s) = q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, k - 2, l - 2)
2281 end do
2282
2283 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
2284
2285 if (chemistry) then
2286 do d = 1, num_species
2287 rhoyks(d) = q_cons_vf(eqn_idx%species%beg + d - 1)%sf(j - 2, k - 2, l - 2)
2288 end do
2289 end if
2290
2291 if (hypoelasticity) then
2292 if (cont_damage) then
2293 ! Recover D = U_D/m_s, clamped to [0, 1]
2294 solid_partial_density = 0._wp
2295 do s = 1, num_fluids
2296 if (fluid_pp(s)%G > verysmall) then
2297 solid_partial_density = solid_partial_density + q_cons_vf(eqn_idx%cont%beg + s &
2298 & - 1)%sf(j - 2, k - 2, l - 2)
2299 end if
2300 end do
2301 damage_state = min(max(q_cons_vf(eqn_idx%damage)%sf(j - 2, k - 2, &
2302 & l - 2)/max(solid_partial_density, verysmall), 0._wp), 1._wp)
2303 end if
2304 g_damaged = g_undamaged*max(1._wp - damage_state, 0._wp)
2305
2306 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k - 2, l - 2), &
2307 & q_cons_vf(eqn_idx%alf)%sf(j - 2, k - 2, l - 2), dyn_p, pi_inf, gamma, &
2308 & rho, qv, rhoyks, pres, t, f_hypoelastic_energy(q_cons_vf, j - 2, k - 2, &
2309 & l - 2, rho, g_undamaged))
2310 else
2311 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k - 2, l - 2), &
2312 & q_cons_vf(eqn_idx%alf)%sf(j - 2, k - 2, l - 2), dyn_p, pi_inf, gamma, &
2313 & rho, qv, rhoyks, pres, t)
2314 end if
2315
2316 if (hypoelasticity) then
2317 do s = 1, 6
2318 tau_e(s) = q_cons_vf(eqn_idx%stress%beg + s - 1)%sf(j - 2, k - 2, l - 2)/rho
2319 end do
2320 end if
2321
2322 ! Compute mixture sound speed
2323 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
2324 if (hypoelasticity) c = sqrt(c*c + (4._wp/3._wp)*g_damaged/rho)
2325
2326 accel = accel_mag(j - 2, k - 2, l - 2)
2327 end if
2328 end if
2329 end if
2330 end if
2331 if (num_procs > 1) then
2332# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2333 tmp = rho
2334# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2335 call s_mpi_allreduce_sum(tmp, rho)
2336# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2337# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2338 tmp = pres
2339# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2340 call s_mpi_allreduce_sum(tmp, pres)
2341# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2342# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2343 tmp = gamma
2344# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2345 call s_mpi_allreduce_sum(tmp, gamma)
2346# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2347# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2348 tmp = pi_inf
2349# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2350 call s_mpi_allreduce_sum(tmp, pi_inf)
2351# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2352# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2353 tmp = qv
2354# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2355 call s_mpi_allreduce_sum(tmp, qv)
2356# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2357# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2358 tmp = c
2359# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2360 call s_mpi_allreduce_sum(tmp, c)
2361# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2362# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2363 tmp = accel
2364# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2365 call s_mpi_allreduce_sum(tmp, accel)
2366# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2367# 1840 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2368
2369 do s = 1, num_vels
2370 tmp = vel(s)
2371 call s_mpi_allreduce_sum(tmp, vel(s))
2372 end do
2373
2374 if (bubbles_euler) then
2375# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2376 tmp = alf
2377# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2378 call s_mpi_allreduce_sum(tmp, alf)
2379# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2380# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2381 tmp = alfgr
2382# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2383 call s_mpi_allreduce_sum(tmp, alfgr)
2384# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2385# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2386 tmp = nbub
2387# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2388 call s_mpi_allreduce_sum(tmp, nbub)
2389# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2390# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2391 tmp = nr(1)
2392# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2393 call s_mpi_allreduce_sum(tmp, nr(1))
2394# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2395# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2396 tmp = nrdot(1)
2397# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2398 call s_mpi_allreduce_sum(tmp, nrdot(1))
2399# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2400# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2401 tmp = m00
2402# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2403 call s_mpi_allreduce_sum(tmp, m00)
2404# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2405# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2406 tmp = r(1)
2407# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2408 call s_mpi_allreduce_sum(tmp, r(1))
2409# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2410# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2411 tmp = rdot(1)
2412# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2413 call s_mpi_allreduce_sum(tmp, rdot(1))
2414# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2415# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2416 tmp = ptilde
2417# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2418 call s_mpi_allreduce_sum(tmp, ptilde)
2419# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2420# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2421 tmp = ptot
2422# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2423 call s_mpi_allreduce_sum(tmp, ptot)
2424# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2425# 1851 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2426
2427 if (qbmm) then
2428# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2429 tmp = varr
2430# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2431 call s_mpi_allreduce_sum(tmp, varr)
2432# 1856 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2433# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2434 tmp = varv
2435# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2436 call s_mpi_allreduce_sum(tmp, varv)
2437# 1856 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2438# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2439 tmp = m10
2440# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2441 call s_mpi_allreduce_sum(tmp, m10)
2442# 1856 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2443# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2444 tmp = m01
2445# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2446 call s_mpi_allreduce_sum(tmp, m01)
2447# 1856 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2448# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2449 tmp = m20
2450# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2451 call s_mpi_allreduce_sum(tmp, m20)
2452# 1856 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2453# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2454 tmp = m02
2455# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2456 call s_mpi_allreduce_sum(tmp, m02)
2457# 1856 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2458# 1857 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2459 end if
2460 end if
2461
2462 if (hypoelasticity) then
2463 do s = 1, (num_dims*(num_dims + 1))/2
2464 tmp = tau_e(s)
2465 call s_mpi_allreduce_sum(tmp, tau_e(s))
2466 end do
2467 end if
2468
2469 if (cont_damage) then
2470 tmp = damage_state
2471 call s_mpi_allreduce_sum(tmp, damage_state)
2472 end if
2473 end if
2474 if (proc_rank == 0) then
2475 if (n == 0) then
2476 if (bubbles_euler .and. (num_fluids <= 2)) then
2477 if (qbmm) then
2478 write (i + 30, '(6x,f12.6,14f28.16)') nondim_time, rho, vel(1), pres, alf, r(1), rdot(1), nr(1), &
2479 & nrdot(1), varr, varv, m10, m01, m20, m02
2480 else
2481 write (i + 30, '(6x,f12.6,8f24.8)') nondim_time, rho, vel(1), pres, alf, r(1), rdot(1), nr(1), nrdot(1)
2482 ! ptilde, & ptot
2483 end if
2484 else if (bubbles_euler .and. (num_fluids == 3)) then
2485 write (i + 30, &
2486 & '(6x,f12.6,f24.8,f24.8,f24.8,f24.8,f24.8,' // 'f24.8,f24.8,f24.8,f24.8,f24.8, f24.8)') &
2487 & nondim_time, rho, vel(1), pres, alf, alfgr, nr(1), nrdot(1), r(1), rdot(1), ptilde, ptot
2488 else if (bubbles_euler .and. num_fluids == 4) then
2489 write (i + 30, &
2490 & '(6x,f12.6,f24.8,f24.8,f24.8,f24.8,' // 'f24.8,f24.8,f24.8,f24.8,f24.8,f24.8,f24.8,f24.8,f24.8)') &
2491 & nondim_time, q_cons_vf(1)%sf(j - 2, 0, 0), q_cons_vf(2)%sf(j - 2, 0, 0), q_cons_vf(3)%sf(j - 2, &
2492 & 0, 0), q_cons_vf(4)%sf(j - 2, 0, 0), q_cons_vf(5)%sf(j - 2, 0, 0), q_cons_vf(6)%sf(j - 2, 0, 0), &
2493 & q_cons_vf(7)%sf(j - 2, 0, 0), q_cons_vf(8)%sf(j - 2, 0, 0), q_cons_vf(9)%sf(j - 2, 0, 0), &
2494 & q_cons_vf(10)%sf(j - 2, 0, 0), nbub, r(1), rdot(1)
2495 else if (hypoelasticity) then
2496 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), pres, tau_e(1)
2497 else
2498 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), pres
2499 end if
2500 else if (p == 0) then
2501 if (bubbles_euler) then
2502# 1901 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2503 write (i + 30, '(6X,10F24.8)') nondim_time, rho, vel(1), vel(2), pres, alf, nr(1), nrdot(1), r(1), &
2504 & rdot(1)
2505# 1904 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2506 else if (hypoelasticity) then
2507# 1906 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2508 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8,F24.8,' // 'F24.8,F24.8,F24.8)') nondim_time, rho, &
2509 & vel(1), vel(2), pres, tau_e(1), tau_e(2), tau_e(3)
2510# 1909 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2511 else
2512# 1911 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2513 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), vel(2), pres
2514# 1913 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2515 end if
2516 else
2517# 1916 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2518 if (hypoelasticity) then
2519 write (i + 30, '(6X,F12.6,16F24.8)') nondim_time, rho, vel(1), vel(2), vel(3), pres, gamma, pi_inf, &
2520 & qv, c, accel, tau_e(1), tau_e(2), tau_e(3), tau_e(4), tau_e(5), tau_e(6)
2521 else
2522 write (i + 30, &
2523 & '(6X,F12.6,F24.8,F24.8,F24.8,F24.8,' // 'F24.8,F24.8,F24.8,F24.8,F24.8,' // 'F24.8)') &
2524 & nondim_time, rho, vel(1), vel(2), vel(3), pres, gamma, pi_inf, qv, c, accel
2525 end if
2526# 1925 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2527 end if
2528 end if
2529 end do
2530
2531 end subroutine s_write_probe_files
2532
2533 !> Write footer with stability criteria extrema and run-time to the information file, then close it
2534 impure subroutine s_close_run_time_information_file
2535
2536 real(wp) :: run_time !< Run-time of the simulation
2537
2538 write (3, '(A)') ' '
2539 write (3, '(A)') ''
2540
2541 write (3, '(A,F9.6)') 'ICFL Max: ', icfl_max
2542 if (surface_tension) write (3, '(A,F9.6)') 'CCFL Max: ', ccfl_max
2543 if (heat_conduction) write (3, '(A,F9.6)') 'TCFL Max: ', tcfl_max
2544 if (viscous) write (3, '(A,F9.6)') 'VCFL Max: ', vcfl_max
2545 if (viscous) write (3, '(A,ES16.6)') 'Rc Min: ', rc_min
2546
2547 call cpu_time(run_time)
2548
2549 write (3, '(A)') ''
2550 write (3, '(A,I0,A)') 'Run-time: ', int(anint(run_time)), 's'
2551 write (3, '(A)') ' '
2552 close (3)
2553
2554 end subroutine s_close_run_time_information_file
2555
2556 !> Closes probe files
2557 impure subroutine s_close_probe_files
2558
2559 integer :: i !< Generic loop iterator
2560
2561 do i = 1, num_probes
2562 close (i + 30)
2563 end do
2564
2565 end subroutine s_close_probe_files
2566
2567 !> Initialize the data output module
2568 impure subroutine s_initialize_data_output_module
2569
2570 integer :: i, m_ds, n_ds, p_ds
2571
2572 if (run_time_info) then
2573 icfl_max = 0._wp
2574 if (surface_tension) then
2575 ccfl_max = 0._wp
2576 end if
2577 if (heat_conduction) then
2578 tcfl_max = 0._wp
2579 end if
2580 if (viscous) then
2581 vcfl_max = 0._wp
2582 rc_min = 1.e12_wp
2583 end if
2584 end if
2585
2586 if (down_sample) then
2587 m_ds = int((m + 1)/3) - 1
2588 n_ds = int((n + 1)/3) - 1
2589 p_ds = int((p + 1)/3) - 1
2590
2591 allocate (q_cons_temp_ds(1:sys_size))
2592 do i = 1, sys_size
2593 allocate (q_cons_temp_ds(i)%sf(-1:m_ds + 1,-1:n_ds + 1,-1:p_ds + 1))
2594 end do
2595 end if
2596
2597 end subroutine s_initialize_data_output_module
2598
2599 !> Module deallocation and/or disassociation procedures
2600 impure subroutine s_finalize_data_output_module
2601
2602 integer :: i
2603
2604 if (down_sample) then
2605 do i = 1, sys_size
2606 deallocate (q_cons_temp_ds(i)%sf)
2607 end do
2608 deallocate (q_cons_temp_ds)
2609 end if
2610
2611 end subroutine s_finalize_data_output_module
2612
2613 !> Half-open ownership lo < v <= hi, so a probe on a rank face is sampled once, by the rank holding the cell the serial code
2614 !! samples (left of the face); the first rank in a direction also owns its lower face.
2615 logical function f_probe_owned(v, lo, hi, dir)
2616
2617 real(wp), intent(in) :: v, lo, hi
2618 integer, intent(in) :: dir
2619 logical :: first
2620
2621 first = num_procs == 1
2622 if (.not. first) first = proc_coords(dir) == 0
2623 f_probe_owned = v <= hi .and. (v > lo .or. (first .and. v >= lo))
2624
2625 end function f_probe_owned
2626
2627end module m_data_output
q_cons_vf
type(scalar_field), dimension(sys_size), intent(inout) q_cons_vf
Definition m_phase_change.fpp.f90:1047
k
integer, intent(in) k
Definition m_phase_change.fpp.f90:1049
j
integer, intent(in) j
Definition m_phase_change.fpp.f90:1049
l
integer, intent(in) l
Definition m_phase_change.fpp.f90:1049
m_boundary_common
Noncharacteristic and processor boundary condition application for ghost cells and buffer regions.
Definition m_boundary_common.fpp.f90:367
m_compile_specific
Platform-specific file and directory operations: create, delete, inquire, getcwd, and basename.
Definition m_compile_specific.f90:6
m_compile_specific::s_delete_directory
impure subroutine s_delete_directory(dir_name)
Recursively delete a directory using a platform-specific system command.
Definition m_compile_specific.f90:43
m_compile_specific::my_inquire
impure subroutine my_inquire(fileloc, dircheck)
Inquire on the existence of a directory or file.
Definition m_compile_specific.f90:56
m_compile_specific::s_create_directory
impure subroutine s_create_directory(dir_name)
Create a directory and all its parents if it does not exist.
Definition m_compile_specific.f90:17
m_constants
Compile-time constant parameters: default values, tolerances, and physical constants.
Definition m_constants.fpp.f90:7
m_constants::model_eqns_5eq
integer, parameter model_eqns_5eq
Definition m_constants.fpp.f90:155
m_constants::name_len
integer, parameter name_len
Maximum name length.
Definition m_constants.fpp.f90:29
m_constants::dflt_t_guess
real(wp), parameter dflt_t_guess
Default guess for temperature (when a previous value is not available).
Definition m_constants.fpp.f90:66
m_constants::dflt_int
integer, parameter dflt_int
Default integer value.
Definition m_constants.fpp.f90:30
m_constants::eos_mie_gruneisen
integer, parameter eos_mie_gruneisen
Definition m_constants.fpp.f90:146
m_constants::sgm_eps
real(wp), parameter sgm_eps
Segmentation tolerance.
Definition m_constants.fpp.f90:13
m_constants::nnode
integer, parameter nnode
Number of QBMM nodes.
Definition m_constants.fpp.f90:46
m_constants::precision_single
integer, parameter precision_single
Definition m_constants.fpp.f90:165
m_constants::pi
real(wp), parameter pi
Pi.
Definition m_constants.fpp.f90:16
m_constants::verysmall
real(wp), parameter verysmall
Very small number.
Definition m_constants.fpp.f90:17
m_data_output
Writes solution data, run-time stability diagnostics (ICFL, VCFL, CCFL, Rc), and probe/center-of-mass...
Definition m_data_output.fpp.f90:367
m_data_output::s_open_probe_files
impure subroutine, public s_open_probe_files
Open flow probe data files for writing.
Definition m_data_output.fpp.f90:486
m_data_output::s_report_icfl_violation
impure subroutine s_report_icfl_violation(q_prim_vf)
Locate the grid cell responsible for an ICFL violation on this rank and report its state plus the nea...
Definition m_data_output.fpp.f90:758
m_data_output::s_write_centroid_offsets
impure subroutine s_write_centroid_offsets(step)
Write each global patch's centroid_offset to restart_data/ib_offset_<step>.dat, so a restart continue...
Definition m_data_output.fpp.f90:1933
m_data_output::rc_min
real(wp) rc_min
Rc criterion maximum.
Definition m_data_output.fpp.f90:396
m_data_output::vcfl_max
real(wp) vcfl_max
VCFL criterion maximum.
Definition m_data_output.fpp.f90:393
m_data_output::tcfl_max
real(wp) tcfl_max
TCFL criterion maximum.
Definition m_data_output.fpp.f90:395
m_data_output::s_write_probe_files
impure subroutine, public s_write_probe_files(t_step, q_cons_vf, accel_mag)
Write flow probe data at the current time step.
Definition m_data_output.fpp.f90:1966
m_data_output::s_write_parallel_ib_state
subroutine s_write_parallel_ib_state(t_step)
Writes the IB state information out to file.
Definition m_data_output.fpp.f90:1605
m_data_output::s_open_ib_force_history
impure subroutine s_open_ib_force_history
Open the shared history file. Done once for the run.
Definition m_data_output.fpp.f90:1836
m_data_output::s_initialize_data_output_module
impure subroutine, public s_initialize_data_output_module
Initialize the data output module.
Definition m_data_output.fpp.f90:2569
m_data_output::s_finalize_data_output_module
impure subroutine, public s_finalize_data_output_module
Module deallocation and/or disassociation procedures.
Definition m_data_output.fpp.f90:2601
m_data_output::s_close_ib_force_history
impure subroutine, public s_close_ib_force_history
Close the history file. Nothing is buffered, so there is nothing to flush first.
Definition m_data_output.fpp.f90:1888
m_data_output::s_write_serial_ib_state
subroutine s_write_serial_ib_state(t_step)
Write IB state data to a per-timestep serial (unformatted) file.
Definition m_data_output.fpp.f90:1710
m_data_output::s_write_serial_ib_data
subroutine s_write_serial_ib_data(time_step)
Write immersed boundary marker data to a serial (per-processor) unformatted file.
Definition m_data_output.fpp.f90:1477
m_data_output::s_close_run_time_information_file
impure subroutine, public s_close_run_time_information_file
Write footer with stability criteria extrema and run-time to the information file,...
Definition m_data_output.fpp.f90:2535
m_data_output::s_write_ib_force_history
impure subroutine, public s_write_ib_force_history(t_step)
Record every immersed body's force, torque and kinematics for this step.
Definition m_data_output.fpp.f90:1763
m_data_output::s_write_data_files
impure subroutine, public s_write_data_files(q_cons_vf, q_t_sf, q_prim_vf, t_step, bc_type, beta)
Write data files. Dispatch subroutine that replaces procedure pointer.
Definition m_data_output.fpp.f90:417
m_data_output::ib_rec_len
integer, parameter ib_rec_len
plus the newline
Definition m_data_output.fpp.f90:410
m_data_output::s_write_serial_data_files
impure subroutine, public s_write_serial_data_files(q_cons_vf, q_t_sf, q_prim_vf, t_step, bc_type, beta)
Write grid and conservative variable data files in serial format.
Definition m_data_output.fpp.f90:908
m_data_output::ib_rec_body
integer, parameter ib_rec_body
characters the format emits
Definition m_data_output.fpp.f90:409
m_data_output::f_probe_owned
logical function f_probe_owned(v, lo, hi, dir)
Half-open ownership lo < v <= hi, so a probe on a rank face is sampled once, by the rank holding the ...
Definition m_data_output.fpp.f90:2616
m_data_output::ib_rec_fmt
character(len= *), parameter ib_rec_fmt
One fixed-width text record per immersed body per recorded step, in D/ib_forces.dat.
Definition m_data_output.fpp.f90:408
m_data_output::s_close_probe_files
impure subroutine, public s_close_probe_files
Closes probe files.
Definition m_data_output.fpp.f90:2558
m_data_output::s_write_run_time_information
impure subroutine, public s_write_run_time_information(q_prim_vf, t_step)
Write stability criteria extrema to the run-time information file at the given time step.
Definition m_data_output.fpp.f90:518
m_data_output::s_open_run_time_information_file
impure subroutine, public s_open_run_time_information_file
Open the run-time information file and write the stability criteria table header.
Definition m_data_output.fpp.f90:435
m_data_output::ib_hist_file
integer ib_hist_file
held open for the run; -1 until first write
Definition m_data_output.fpp.f90:411
m_data_output::icfl_max
real(wp) icfl_max
ICFL criterion maximum.
Definition m_data_output.fpp.f90:392
m_data_output::s_write_parallel_data_files
impure subroutine, public s_write_parallel_data_files(q_cons_vf, t_step, bc_type, beta, q_t_sf)
Write grid and conservative variable data files in parallel via MPI I/O.
Definition m_data_output.fpp.f90:1271
m_data_output::s_write_ib_state_file
impure subroutine, public s_write_ib_state_file(time_step)
Writes IB state records to restart_data/ib_state.dat. Called on every rank.
Definition m_data_output.fpp.f90:1905
m_data_output::q_cons_temp_ds
type(scalar_field), dimension(:), allocatable q_cons_temp_ds
Definition m_data_output.fpp.f90:399
m_data_output::s_write_parallel_ib_data
subroutine s_write_parallel_ib_data(time_step)
Write immersed boundary marker data in parallel using MPI I/O.
Definition m_data_output.fpp.f90:1506
m_data_output::s_write_ib_data_file
subroutine, public s_write_ib_data_file(time_step)
Dispatch immersed boundary data output to the serial or parallel writer.
Definition m_data_output.fpp.f90:1579
m_data_output::ccfl_max
real(wp) ccfl_max
CCFL criterion maximum.
Definition m_data_output.fpp.f90:394
m_delay_file_access
Rank-staggered file access delays to prevent I/O contention on parallel file systems.
Definition m_delay_file_access.f90:6
m_delay_file_access::s_delay_file_access
impure subroutine, public s_delay_file_access(process_rank)
Introduce a rank-dependent busy-wait delay to stagger parallel file access and reduce I/O contention.
Definition m_delay_file_access.f90:22
m_derived_types
Shared derived types for field data, patch geometry, bubble dynamics, and MPI I/O structures.
Definition m_derived_types.fpp.f90:356
m_global_parameters
Global parameters for the computational domain, fluid properties, and simulation algorithm configurat...
Definition m_global_parameters.fpp.f90:367
m_global_parameters::mytime
real(wp) mytime
Current simulation time.
Definition m_global_parameters.fpp.f90:1107
m_global_parameters::fluid_inv_re
real(wp), dimension(:), allocatable fluid_inv_re
per-fluid Newtonian inverse-Re
Definition m_global_parameters.fpp.f90:794
m_global_parameters::p_glb
integer p_glb
Definition m_global_parameters.fpp.f90:405
m_global_parameters::idwint
type(int_bounds_info), dimension(1:3) idwint
Definition m_global_parameters.fpp.f90:758
m_global_parameters::z_cb
real(wp), dimension(:), allocatable, target z_cb
Definition m_global_parameters.fpp.f90:427
m_global_parameters::any_non_newtonian
logical any_non_newtonian
.true. if any fluid is non-Newtonian
Definition m_global_parameters.fpp.f90:790
m_global_parameters::qbmm_idx
type(qbmm_idx_info) qbmm_idx
QBMM moment index mappings (allocatable; GPU-managed separately).
Definition m_global_parameters.fpp.f90:755
m_global_parameters::proc_rank
integer proc_rank
Rank of the local processor.
Definition m_global_parameters.fpp.f90:386
m_global_parameters::mpi_io_ib_data
type(mpi_io_ib_var), public mpi_io_ib_data
Definition m_global_parameters.fpp.f90:748
m_global_parameters::weight
real(wp), dimension(:), allocatable weight
Simpson quadrature weights.
Definition m_global_parameters.fpp.f90:951
m_global_parameters::local_ib_patch_ids
integer, dimension(num_local_ibs_max) local_ib_patch_ids
lookup table of IBs in the local compute domain
Definition m_global_parameters.fpp.f90:913
m_global_parameters::n_el_bubs_loc
integer n_el_bubs_loc
Definition m_global_parameters.fpp.f90:1129
m_global_parameters::buff_size
integer buff_size
Number of ghost cells for boundary condition storage.
Definition m_global_parameters.fpp.f90:853
m_global_parameters::num_local_ibs
integer num_local_ibs
Definition m_global_parameters.fpp.f90:712
m_global_parameters::y_cc
real(wp), dimension(:), allocatable, target y_cc
Definition m_global_parameters.fpp.f90:433
m_global_parameters::pb_ts
type(pres_field), dimension(:), allocatable pb_ts
Definition m_global_parameters.fpp.f90:1108
m_global_parameters::n_el_bubs_glb
integer n_el_bubs_glb
Number of Lagrangian bubbles (local and global).
Definition m_global_parameters.fpp.f90:1129
m_global_parameters::mv_ts
type(pres_field), dimension(:), allocatable mv_ts
Definition m_global_parameters.fpp.f90:1109
m_global_parameters::z_cc
real(wp), dimension(:), allocatable, target z_cc
Definition m_global_parameters.fpp.f90:433
m_global_parameters::num_procs
integer num_procs
Number of processors.
Definition m_global_parameters.fpp.f90:384
m_global_parameters::x_cc
real(wp), dimension(:), allocatable, target x_cc
Definition m_global_parameters.fpp.f90:433
m_global_parameters::y_cb
real(wp), dimension(:), allocatable, target y_cb
Definition m_global_parameters.fpp.f90:427
m_global_parameters::ptil
real(wp), dimension(:,:,:), allocatable ptil
Pressure modification.
Definition m_global_parameters.fpp.f90:966
m_global_parameters::m_glb
integer m_glb
Definition m_global_parameters.fpp.f90:405
m_global_parameters::nb
integer nb
Definition m_global_parameters.fpp.f90:934
m_global_parameters::mpi_io_data
type(mpi_io_var), public mpi_io_data
Definition m_global_parameters.fpp.f90:747
m_global_parameters::n_glb
integer n_glb
Definition m_global_parameters.fpp.f90:405
m_global_parameters::num_gbl_ibs
integer num_gbl_ibs
Definition m_global_parameters.fpp.f90:712
m_global_parameters::hb_mu_max
real(wp), dimension(:), allocatable hb_mu_max
Definition m_global_parameters.fpp.f90:793
m_global_parameters::is_non_newtonian
logical, dimension(:), allocatable is_non_newtonian
per-fluid NN flag
Definition m_global_parameters.fpp.f90:791
m_global_parameters::x_cb
real(wp), dimension(:), allocatable, target x_cb
Definition m_global_parameters.fpp.f90:427
m_global_parameters::cfl_dt
logical cfl_dt
Definition m_global_parameters.fpp.f90:467
m_helper_basic
Basic floating-point utilities: approximate equality, default detection, and coordinate bounds.
Definition m_helper_basic.fpp.f90:356
m_helper_basic::f_approx_equal
logical elemental function, public f_approx_equal(a, b, tol_input)
Check if two floating point numbers of wp are within tolerance.
Definition m_helper_basic.fpp.f90:373
m_helper
Utility routines for bubble model setup, coordinate transforms, array sampling, and special functions...
Definition m_helper.fpp.f90:367
m_helper::s_downsample_data
subroutine, public s_downsample_data(q_cons_vf, q_cons_temp, m_ds, n_ds, p_ds, m_glb_ds, n_glb_ds, p_glb_ds)
Downsample conservative variable fields by a factor of 3 in each direction using volume averaging.
Definition m_helper.fpp.f90:1289
m_helper::s_int_to_str
elemental subroutine, public s_int_to_str(i, res)
Convert an integer to its trimmed string representation.
Definition m_helper.fpp.f90:777
m_ibm
Ghost-node immersed boundary method: locates ghost/image points, computes interpolation coefficients,...
Definition m_ibm.fpp.f90:357
m_ibm::ib_markers
type(integer_field), public ib_markers
Definition m_ibm.fpp.f90:392
m_ibm::f_needs_centroid_offset
pure logical function f_needs_centroid_offset(patch)
A moving airfoil or STL (geometries 4, 5, 11, 12) moves its centroid to the centre of mass and keeps ...
Definition m_ibm.fpp.f90:3199
m_ibm::centroid_offsets_active
logical centroid_offsets_active
some patch is a moving airfoil or STL, which carries a centroid_offset
Definition m_ibm.fpp.f90:465
m_mpi_proxy
MPI halo exchange, domain decomposition, and buffer packing/unpacking for the simulation solver.
Definition m_mpi_proxy.fpp.f90:367
m_sim_helpers
Simulation helper routines for cell state, CFL calculation, and stability checks.
Definition m_sim_helpers.fpp.f90:367
m_sim_helpers::s_compute_stability_from_dt
subroutine, public s_compute_stability_from_dt(vel, c, rho, re_l, alpha, alpha_rho, j, k, l, icfl, vcfl, rc, ccfl, tcfl)
Computes stability criterion for a specified dt.
Definition m_sim_helpers.fpp.f90:545
m_sim_helpers::s_compute_cell_state
subroutine, public s_compute_cell_state(q_prim_vf, pres, rho, gamma, pi_inf, re, alpha, alpha_rho, vel, vel_sum, qv, j, k, l)
Computes the mixture coefficients, velocity and pressure of one cell.
Definition m_sim_helpers.fpp.f90:423
m_variables_conversion
Conservative-to-primitive variable conversion, mixture property evaluation, and pressure computation.
Definition m_variables_conversion.fpp.f90:367
m_variables_conversion::f_hypoelastic_energy
real(wp) function, public f_hypoelastic_energy(q_cons_vf, j, k, l, rho, g)
Hypoelastic strain energy at one cell, summed over the stress components.
Definition m_variables_conversion.fpp.f90:2552
m_variables_conversion::s_convert_conservative_to_primitive_variables
subroutine, public s_convert_conservative_to_primitive_variables(qk_cons_vf, q_t_sf, qk_prim_vf, ibounds)
Convert conserved variables (rho*alpha, rho*u, E, alpha) to primitives (rho, u, p,...
Definition m_variables_conversion.fpp.f90:1053
m_variables_conversion::s_compute_pressure
subroutine, public s_compute_pressure(energy, alf, dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t, e_e_in, pres_mag)
Compute the pressure from the appropriate equation of state.
Definition m_variables_conversion.fpp.f90:464
m_variables_conversion::s_convert_to_mixture_variables
subroutine, public s_convert_to_mixture_variables(q_vf, i, j, k, rho, gamma, pi_inf, qv, re_k, g_k, g)
Dispatch to the s_convert_mixture_to_mixture_variables and s_convert_species_to_mixture_variables sub...
Definition m_variables_conversion.fpp.f90:446
m_derived_types::integer_field
Derived type annexing an integer scalar field (SF).
Definition m_derived_types.fpp.f90:380
m_derived_types::scalar_field
Derived type annexing a scalar field (SF).
Definition m_derived_types.fpp.f90:370