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, ierr)
1462 end if
1463
1464 call mpi_file_close(ifile, ierr)
1465
1466 if (ib) then
1467 call s_write_parallel_ib_data(t_step)
1468 end if
1469 end if
1470#endif
1471
1472 end subroutine s_write_parallel_data_files
1473
1474 !> Write immersed boundary marker data to a serial (per-processor) unformatted file
1475 subroutine s_write_serial_ib_data(time_step)
1476
1477 integer, intent(in) :: time_step
1478 character(LEN=path_len + 2*name_len) :: file_path
1479 character(LEN=path_len + 2*name_len) :: t_step_dir
1480
1481 write (t_step_dir, '(A,I0,A,I0)') trim(case_dir) // '/p_all'
1482 write (t_step_dir, '(a,i0,a,i0)') trim(case_dir) // '/p_all/p', proc_rank, '/', time_step
1483 write (file_path, '(A,I0,A)') trim(t_step_dir) // '/ib_data.dat'
1484
1485 open (2, file=trim(file_path), form='unformatted', status='new')
1486
1487
1488# 1050 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1489#if defined(MFC_OpenACC)
1490# 1050 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1491!$acc update host(ib_markers%sf)
1492# 1050 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1493#elif defined(MFC_OpenMP)
1494# 1050 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1495!$omp target update from(ib_markers%sf)
1496# 1050 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1497#endif
1498# 1051 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1499 write (2) ib_markers%sf(0:m,0:n,0:p); close (2)
1500
1501 end subroutine s_write_serial_ib_data
1502
1503 !> Write immersed boundary marker data in parallel using MPI I/O
1504 subroutine s_write_parallel_ib_data(time_step)
1505
1506 integer, intent(in) :: time_step
1507
1508#ifdef MFC_MPI
1509 character(LEN=path_len + 2*name_len) :: file_loc
1510 integer(kind=MPI_OFFSET_kind) :: disp
1511 integer(kind=MPI_OFFSET_kind) :: m_MOK, n_MOK, p_MOK
1512 integer(kind=MPI_OFFSET_kind) :: WP_MOK, var_MOK, MOK
1513 integer :: ifile, ierr, data_size
1514 integer, dimension(MPI_STATUS_SIZE) :: status
1515 character(len=10) :: t_step_string
1516
1517
1518# 1069 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1519#if defined(MFC_OpenACC)
1520# 1069 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1521!$acc update host(ib_markers%sf)
1522# 1069 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1523#elif defined(MFC_OpenMP)
1524# 1069 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1525!$omp target update from(ib_markers%sf)
1526# 1069 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1527#endif
1528# 1070 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1529
1530 data_size = (m + 1)*(n + 1)*(p + 1)
1531 m_mok = int(m_glb + 1, mpi_offset_kind)
1532 n_mok = int(n_glb + 1, mpi_offset_kind)
1533 p_mok = int(p_glb + 1, mpi_offset_kind)
1534 wp_mok = int(storage_size(0._stp)/8, mpi_offset_kind)
1535 mok = int(1._wp, mpi_offset_kind)
1536
1537 if (file_per_process) then
1538 call s_int_to_str(time_step, t_step_string)
1539
1540 if (proc_rank == 0) then
1541 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string)
1542 call s_create_directory(trim(file_loc))
1543 end if
1544 call s_mpi_barrier()
1545 call s_delay_file_access(proc_rank)
1546
1547 write (file_loc, '(A,I0,A,i7.7,A)') 'ib_markers_', time_step, '_', proc_rank, '.dat'
1548 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string) // '/' // trim(file_loc)
1549
1550 call mpi_file_open(mpi_comm_self, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1551 call s_check_mpi_file_open(ierr, file_loc)
1552 call mpi_file_write_all(ifile, mpi_io_ib_data%var%sf, data_size, mpi_integer, status, ierr)
1553 call mpi_file_close(ifile, ierr)
1554 else
1555 write (file_loc, '(A)') 'ib.dat'
1556 file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // trim(file_loc)
1557
1558 call s_mpi_barrier()
1559 call s_delay_file_access(proc_rank)
1560
1561 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1562 call s_check_mpi_file_open(ierr, file_loc)
1563
1564 var_mok = int(sys_size + 1, mpi_offset_kind)
1565 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1 + int(time_step/t_step_save))
1566 if (time_step == 0) disp = 0
1567
1568 call mpi_file_set_view(ifile, disp, mpi_integer, mpi_io_ib_data%view, 'native', mpi_info_int, ierr)
1569 call mpi_file_write_all(ifile, mpi_io_ib_data%var%sf, data_size, mpi_integer, status, ierr)
1570 call mpi_file_close(ifile, ierr)
1571 end if
1572#endif
1573
1574 end subroutine s_write_parallel_ib_data
1575
1576 !> Dispatch immersed boundary data output to the serial or parallel writer
1577 subroutine s_write_ib_data_file(time_step)
1578
1579 integer, intent(in) :: time_step
1580
1581
1582# 1122 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1583#if defined(MFC_OpenACC)
1584# 1122 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1585!$acc update host(patch_ib(1:num_ibs))
1586# 1122 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1587#elif defined(MFC_OpenMP)
1588# 1122 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1589!$omp target update from(patch_ib(1:num_ibs))
1590# 1122 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1591#endif
1592# 1123 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1593
1594 if (parallel_io) then
1595 call s_write_parallel_ib_data(time_step)
1596 else
1597 call s_write_serial_ib_data(time_step)
1598 end if
1599
1600 end subroutine s_write_ib_data_file
1601
1602 !> Writes the IB state information out to file
1603 subroutine s_write_parallel_ib_state(t_step)
1604
1605 integer, intent(in) :: t_step
1606
1607#ifdef MFC_MPI
1608 character(LEN=path_len + 2*name_len) :: file_loc
1609 integer(kind=MPI_OFFSET_KIND) :: disp
1610 integer(kind=MPI_OFFSET_KIND) :: WP_MOK
1611 integer :: ifile, ierr
1612 integer, dimension(MPI_STATUS_SIZE) :: status
1613 logical :: file_exist, dir_check
1614 integer :: i, ib_idx
1615 integer, parameter :: NFIELDS_PER_IB = 20
1616 real(wp) :: ib_buf(NFIELDS_PER_IB)
1617 integer :: file_unit
1618 character(len=10) :: t_step_string
1619
1620 ! Partition IBs across ranks round-robin style
1621 integer :: ib_start, ib_end, nibs_per_rank, remainder
1622
1623 wp_mok = int(storage_size(0._wp)/8, mpi_offset_kind)
1624
1625 if (file_per_process) then
1626 call s_int_to_str(t_step, t_step_string)
1627
1628 if (proc_rank == 0) then
1629 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string)
1630 call s_create_directory(trim(file_loc))
1631 end if
1632 call s_mpi_barrier()
1633 call s_delay_file_access(proc_rank)
1634
1635 write (file_loc, '(A,I0,A,i7.7,A)') 'ib_state_', t_step, '_', proc_rank, '.dat'
1636 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string) // '/' // trim(file_loc)
1637
1638 inquire (file=trim(file_loc), exist=file_exist)
1639 if (file_exist) then
1640 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='replace')
1641 else
1642 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='new')
1643 end if
1644
1645 write (file_unit) num_local_ibs
1646 do i = 1, num_local_ibs
1647 ib_idx = local_ib_patch_ids(i)
1648 ib_buf(1) = mytime
1649 ib_buf(2:4) = patch_ib(ib_idx)%force(1:3)
1650 ib_buf(5:7) = patch_ib(ib_idx)%torque(1:3)
1651 ib_buf(8:10) = patch_ib(ib_idx)%vel(1:3)
1652 ib_buf(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1653 ib_buf(14:16) = patch_ib(ib_idx)%angles(1:3)
1654 ib_buf(17) = patch_ib(ib_idx)%x_centroid
1655 ib_buf(18) = patch_ib(ib_idx)%y_centroid
1656 ib_buf(19) = patch_ib(ib_idx)%z_centroid
1657 ib_buf(20) = patch_ib(ib_idx)%radius
1658
1659 write (file_unit) patch_ib(ib_idx)%gbl_patch_id
1660 write (file_unit) ib_buf
1661 end do
1662
1663 close (file_unit)
1664 else
1665 if (proc_rank == 0) then
1666 call s_create_directory(trim(case_dir) // '/restart_data')
1667 end if
1668 call s_mpi_barrier()
1669
1670 write (file_loc, '(A,I0,A)') '/restart_data/ib_state_', t_step, '.dat'
1671 file_loc = trim(case_dir) // trim(file_loc)
1672
1673 inquire (file=trim(file_loc), exist=file_exist)
1674 if (file_exist .and. proc_rank == 0) then
1675 call mpi_file_delete(file_loc, mpi_info_int, ierr)
1676 end if
1677 call s_mpi_barrier()
1678
1679 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1680 call s_check_mpi_file_open(ierr, file_loc)
1681
1682 do i = 1, num_local_ibs
1683 ib_idx = local_ib_patch_ids(i)
1684 ib_buf(1) = mytime
1685 ib_buf(2:4) = patch_ib(ib_idx)%force(1:3)
1686 ib_buf(5:7) = patch_ib(ib_idx)%torque(1:3)
1687 ib_buf(8:10) = patch_ib(ib_idx)%vel(1:3)
1688 ib_buf(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1689 ib_buf(14:16) = patch_ib(ib_idx)%angles(1:3)
1690 ib_buf(17) = patch_ib(ib_idx)%x_centroid
1691 ib_buf(18) = patch_ib(ib_idx)%y_centroid
1692 ib_buf(19) = patch_ib(ib_idx)%z_centroid
1693 ib_buf(20) = patch_ib(ib_idx)%radius
1694
1695 ! Global IB index determines position in file
1696 disp = int(patch_ib(ib_idx)%gbl_patch_id - 1, mpi_offset_kind)*int(nfields_per_ib, mpi_offset_kind)*wp_mok
1697
1698 call mpi_file_write_at(ifile, disp, ib_buf, nfields_per_ib, mpi_p, status, ierr)
1699 end do
1700
1701 call mpi_file_close(ifile, ierr)
1702 end if
1703#endif
1704
1705 end subroutine s_write_parallel_ib_state
1706
1707 !> Write IB state data to a per-timestep serial (unformatted) file
1708 subroutine s_write_serial_ib_state(t_step)
1709
1710 integer, intent(in) :: t_step
1711 character(LEN=path_len + 2*name_len) :: file_loc
1712 integer :: i, ios, file_unit
1713 integer, parameter :: NFIELDS_PER_IB = 20
1714 real(wp) :: ib_buf(NFIELDS_PER_IB)
1715
1716 call s_create_directory(trim(case_dir) // '/restart_data')
1717
1718 write (file_loc, '(A,I0,A)') '/restart_data/ib_state_', t_step, '.dat'
1719 file_loc = trim(case_dir) // trim(file_loc)
1720
1721 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='replace', iostat=ios)
1722 if (ios /= 0) call s_mpi_abort('Cannot open IB state output file: ' // trim(file_loc))
1723
1724 do i = 1, num_ibs
1725 ib_buf(1) = mytime
1726 ib_buf(2:4) = patch_ib(i)%force(1:3)
1727 ib_buf(5:7) = patch_ib(i)%torque(1:3)
1728 ib_buf(8:10) = patch_ib(i)%vel(1:3)
1729 ib_buf(11:13) = patch_ib(i)%angular_vel(1:3)
1730 ib_buf(14:16) = patch_ib(i)%angles(1:3)
1731 ib_buf(17) = patch_ib(i)%x_centroid
1732 ib_buf(18) = patch_ib(i)%y_centroid
1733 ib_buf(19) = patch_ib(i)%z_centroid
1734 ib_buf(20) = patch_ib(i)%radius
1735
1736 write (file_unit) ib_buf
1737 end do
1738
1739 close (file_unit)
1740
1741 end subroutine s_write_serial_ib_state
1742
1743 !> Record every immersed body's force, torque and kinematics for this step.
1744 !!
1745 !! One shared text file, D/ib_forces.dat, opened once for the run. Each rank writes only the
1746 !! bodies it owns, at a byte offset computed from the step and the global body id, so the file
1747 !! is byte-identical however the domain is decomposed and needs no merge step. Writing a file
1748 !! per body instead costs an inquire, open and close per body per rank per step, which is
1749 !! 3e5 filesystem metadata operations per step at 1000 ranks holding 100 bodies each.
1750 !!
1751 !! Layout: row r holds all num_gbl_ibs bodies in global id order, so body g occupies bytes
1752 !! ((r*num_gbl_ibs) + g - 1)*IB_REC_LEN. The file carries no header line, which would shift
1753 !! every offset after it; the columns are listed in docs/documentation/case.md.
1754 !!
1755 !! Rows are numbered from the first step this run records, not from t_step, so that row 0 is
1756 !! always written. Nothing pre-fills the file -- both open paths create it empty and every
1757 !! write lands at a computed offset -- so a row no rank ever writes stays a hole, and a hole
1758 !! reads back as NUL bytes rather than blanks. Counting from t_step would leave exactly such a
1759 !! hole wherever the run starts: the skip below means step t_step_start is never written, and
1760 !! on a restart every row beneath it would be missing as well.
1761 impure subroutine s_write_ib_force_history(t_step)
1762
1763 integer, intent(in) :: t_step
1764 character(LEN=IB_REC_LEN) :: rec
1765 character(LEN=path_len + 2*name_len) :: file_loc
1766 real(wp) :: fields(19)
1767 integer :: i, ib_idx, n_write, row
1768
1769#ifdef MFC_MPI
1770 integer(kind=MPI_OFFSET_KIND) :: disp
1771 integer :: ierr, status(mpi_status_size)
1772#endif
1773
1774 if (.not. ib_force_wrt) return
1775 if (mod(t_step, max(ib_force_stride, 1)) /= 0) return
1776 ! This runs at RK stage 1, before the step's force has been computed, so the row for step N carries the
1777 ! 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
1778 ! the row would record identically zero force. That is not a measurement, and on a run chained across a
1779 ! queue's walltime limit it lands once per restart -- in a six-wingbeat case, zeros at steps 20649,
1780 ! 34649 and 48649 sitting among neighbours of -0.134, +0.474 and -0.475, corrupting every per-beat
1781 ! trough taken over the joined trace.
1782 if (t_step == t_step_start) return
1783
1784 n_write = num_local_ibs
1785 if (num_procs == 1) n_write = num_ibs
1786 ! Relative to the first recorded step, so the first one written is row 0 and the file is
1787 ! dense. See the hole discussion above.
1788 row = t_step/max(ib_force_stride, 1) - t_step_start/max(ib_force_stride, 1) - 1
1789
1790
1791# 1320 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1792#if defined(MFC_OpenACC)
1793# 1320 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1794!$acc update host(patch_ib(1:num_ibs))
1795# 1320 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1796#elif defined(MFC_OpenMP)
1797# 1320 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1798!$omp target update from(patch_ib(1:num_ibs))
1799# 1320 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1800#endif
1801# 1321 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1802
1803 call s_open_ib_force_history()
1804
1805 do i = 1, n_write
1806 ib_idx = i
1807 if (num_procs > 1) ib_idx = local_ib_patch_ids(i)
1808
1809 fields(1) = mytime
1810 fields(2:4) = patch_ib(ib_idx)%force(1:3)
1811 fields(5:7) = patch_ib(ib_idx)%torque(1:3)
1812 fields(8:10) = patch_ib(ib_idx)%vel(1:3)
1813 fields(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1814 fields(14:16) = patch_ib(ib_idx)%angles(1:3)
1815 fields(17) = patch_ib(ib_idx)%x_centroid
1816 fields(18) = patch_ib(ib_idx)%y_centroid
1817 fields(19) = patch_ib(ib_idx)%z_centroid
1818
1819 write (rec, ib_rec_fmt) patch_ib(ib_idx)%gbl_patch_id, fields
1820 rec(ib_rec_len:ib_rec_len) = new_line('a')
1821
1822#ifdef MFC_MPI
1823 disp = (int(row, mpi_offset_kind)*int(num_gbl_ibs, mpi_offset_kind) + int(patch_ib(ib_idx)%gbl_patch_id - 1, &
1824 & mpi_offset_kind))*int(ib_rec_len, mpi_offset_kind)
1825 call mpi_file_write_at(ib_hist_file, disp, rec, ib_rec_len, mpi_character, status, ierr)
1826#else
1827 write (ib_hist_file, rec=row*num_gbl_ibs + patch_ib(ib_idx)%gbl_patch_id) rec
1828#endif
1829 end do
1830
1831 end subroutine s_write_ib_force_history
1832
1833 !> Open the shared history file. Done once for the run.
1834 impure subroutine s_open_ib_force_history
1835
1836 character(LEN=path_len + 2*name_len) :: file_loc
1837 character(LEN=IB_REC_LEN) :: probe
1838 integer :: i
1839
1840#ifdef MFC_MPI
1841 integer :: ierr
1842 logical :: file_exist
1843#endif
1844
1845 if (ib_hist_file /= -1) return
1846
1847 ! Every offset below assumes the format emits exactly IB_REC_BODY characters. Measure it once
1848 ! rather than trusting that the format and the constant were edited together: a format one
1849 ! character wider would shear every record past the first without any other symptom.
1850 write (probe, ib_rec_fmt) 0, [(0._wp, i=1, 19)]
1851 if (len_trim(probe) /= ib_rec_body) then
1852# 1370 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1853 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")
1854# 1370 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1855 end if
1856# 1372 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1857
1858 file_loc = trim(case_dir) // '/D/ib_forces.dat'
1859#ifdef MFC_MPI
1860 ! MPI_MODE_CREATE does not truncate, so a shorter run following a longer one in the same
1861 ! directory would keep the old tail past its last record. Delete first, as the ib_state
1862 ! writer does, then barrier so no rank opens before the delete lands.
1863 inquire (file=trim(file_loc), exist=file_exist)
1864 if (file_exist .and. proc_rank == 0) call mpi_file_delete(file_loc, mpi_info_null, ierr)
1865
1866 ! MPI_INFO_NULL, not mpi_info_int: the latter is only created when parallel_io is on
1867 ! (m_global_parameters_common.fpp returns before MPI_INFO_CREATE otherwise), and the
1868 ! force history is written whatever parallel_io is set to. Passing the uninitialised
1869 ! handle aborted every IBM case that runs the solver with parallel_io = F, in
1870 ! MPI_Info_dup, at the first recorded step. The hint it carries only disables ROMIO
1871 ! write data sieving, which this writer does not depend on.
1872 ! Collective: every rank opens, including one holding no body this step.
1873 call s_mpi_barrier()
1874 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_null, ib_hist_file, ierr)
1875 call s_check_mpi_file_open(ierr, file_loc)
1876#else
1877 ! Unformatted: the record is already a formatted string, so this writes its bytes verbatim and
1878 ! produces the same file the MPI branch does. A formatted direct-access write would need a
1879 ! format and would pad rather than emit the string as-is.
1880 open (newunit=ib_hist_file, file=trim(file_loc), form='unformatted', access='direct', recl=ib_rec_len, status='replace')
1881#endif
1882
1883 end subroutine s_open_ib_force_history
1884
1885 !> Close the history file. Nothing is buffered, so there is nothing to flush first.
1886 impure subroutine s_close_ib_force_history
1887
1888#ifdef MFC_MPI
1889 integer :: ierr
1890#endif
1891
1892 if (ib_hist_file == -1) return
1893#ifdef MFC_MPI
1894 call mpi_file_close(ib_hist_file, ierr)
1895#else
1896 close (ib_hist_file)
1897#endif
1898 ib_hist_file = -1
1899
1900 end subroutine s_close_ib_force_history
1901
1902 !> @brief Writes IB state records to restart_data/ib_state.dat. Called on every rank.
1903 impure subroutine s_write_ib_state_file(time_step)
1904
1905 integer, intent(in) :: time_step
1906
1907
1908# 1422 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1909#if defined(MFC_OpenACC)
1910# 1422 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1911!$acc update host(patch_ib(1:num_ibs))
1912# 1422 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1913#elif defined(MFC_OpenMP)
1914# 1422 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1915!$omp target update from(patch_ib(1:num_ibs))
1916# 1422 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1917#endif
1918# 1423 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1919
1920 if (parallel_io) then
1921 call s_write_parallel_ib_state(time_step)
1922 else
1923 call s_write_serial_ib_state(time_step)
1924 end if
1925 call s_write_centroid_offsets(time_step)
1926
1927 end subroutine s_write_ib_state_file
1928
1929 !> Write each global patch's centroid_offset to restart_data/ib_offset_<step>.dat, so a restart continues about the same centre
1930 !! of mass instead of re-measuring it from the body voxelised at the restart attitude. Collective.
1931 impure subroutine s_write_centroid_offsets(step)
1932
1933 integer, intent(in) :: step
1934 character(len=path_len + 2*name_len) :: file_loc
1935 real(wp), dimension(:,:), allocatable :: off_loc, off_glb
1936 integer :: gid, i, ib_idx, n_own, file_unit
1937
1938 if (.not. centroid_offsets_active) return
1939 ! column 4 flags a patch that has an offset; only the owning rank contributes, so the sum is the value
1940 allocate (off_loc(num_gbl_ibs, 4), off_glb(num_gbl_ibs, 4))
1941 off_loc = 0._wp
1942 n_own = num_local_ibs
1943 if (num_procs == 1) n_own = num_ibs
1944 do i = 1, n_own
1945 ib_idx = i
1946 if (num_procs > 1) ib_idx = local_ib_patch_ids(i)
1947 if (.not. f_needs_centroid_offset(patch_ib(ib_idx))) cycle
1948 off_loc(patch_ib(ib_idx)%gbl_patch_id,:) = [patch_ib(ib_idx)%centroid_offset, 1._wp]
1949 end do
1950 call s_mpi_allreduce_vectors_sum(off_loc, off_glb, num_gbl_ibs, 4)
1951 if (proc_rank == 0) then
1952 write (file_loc, '(A,I0,A)') trim(case_dir) // '/restart_data/ib_offset_', step, '.dat'
1953 open (newunit=file_unit, file=trim(file_loc), status='replace', action='write')
1954 do gid = 1, num_gbl_ibs
1955 if (off_glb(gid, 4) > 0.5_wp) write (file_unit, '(I0,3(1X,ES24.16))') gid, off_glb(gid,1:3)
1956 end do
1957 close (file_unit)
1958 end if
1959 deallocate (off_loc, off_glb)
1960
1961 end subroutine s_write_centroid_offsets
1962
1963 !> Write flow probe data at the current time step
1964 impure subroutine s_write_probe_files(t_step, q_cons_vf, accel_mag)
1965
1966 integer, intent(in) :: t_step
1967 type(scalar_field), dimension(sys_size), intent(in) :: q_cons_vf
1968 real(wp), dimension(0:m,0:n,0:p), intent(in) :: accel_mag
1969 real(wp), dimension(-1:m) :: distx
1970 real(wp), dimension(-1:n) :: disty
1971 real(wp), dimension(-1:p) :: distz
1972
1973 ! The cell-averaged partial densities, density, velocity, pressure, volume fractions, specific heat ratio function, liquid
1974 ! stiffness function, and sound speed.
1975 real(wp) :: lit_gamma, nbub
1976 real(wp) :: rho
1977 real(wp), dimension(num_vels) :: vel
1978 real(wp) :: pres
1979 real(wp) :: ptilde
1980 real(wp) :: ptot
1981 real(wp) :: alf
1982 real(wp) :: alfgr
1983 real(wp), dimension(num_fluids) :: alpha, alpha_rho
1984 real(wp) :: gamma
1985 real(wp) :: pi_inf
1986 real(wp) :: qv
1987 real(wp) :: c
1988 real(wp) :: m00, m10, m01, m20, m02
1989 real(wp) :: varr, varv
1990 real(wp), dimension(Nb) :: nr, r, nrdot, rdot
1991 real(wp) :: nr3
1992 real(wp) :: accel
1993 real(wp) :: int_pres
1994 real(wp) :: max_pres
1995 real(wp), dimension(2) :: re
1996 real(wp), dimension(6) :: tau_e
1997 real(wp) :: g_undamaged, g_damaged
1998 real(wp) :: dyn_p, t
1999 real(wp) :: damage_state
2000 real(wp) :: solid_partial_density !< damageable-solid partial density at the probe cell
2001 integer :: i, j, k, l, s, d !< Generic loop iterator
2002 real(wp) :: nondim_time !< Non-dimensional time
2003 real(wp) :: tmp !< Temporary variable to store quantity for mpi_allreduce
2004 real(wp) :: rhoyks(1:num_species)
2005
2006 t = dflt_t_guess
2007
2008 if (time_stepper == 23) then
2009 nondim_time = mytime
2010 else
2011 if (t_step_old /= dflt_int) then
2012 nondim_time = real(t_step + t_step_old, wp)*dt
2013 else
2014 nondim_time = real(t_step, wp)*dt
2015 end if
2016 end if
2017
2018 do i = 1, num_probes
2019 rho = 0._wp
2020 do s = 1, num_vels
2021 vel(s) = 0._wp
2022 end do
2023 pres = 0._wp
2024 gamma = 0._wp
2025 pi_inf = 0._wp
2026 qv = 0._wp
2027 c = 0._wp
2028 accel = 0._wp
2029 nr = 0._wp; r = 0._wp
2030 nrdot = 0._wp; rdot = 0._wp
2031 nbub = 0._wp
2032 m00 = 0._wp
2033 m10 = 0._wp
2034 m01 = 0._wp
2035 m20 = 0._wp
2036 m02 = 0._wp
2037 varr = 0._wp; varv = 0._wp
2038 alf = 0._wp
2039 do s = 1, (num_dims*(num_dims + 1))/2
2040 tau_e(s) = 0._wp
2041 end do
2042 damage_state = 0._wp
2043 g_damaged = 0._wp
2044
2045 if (n == 0) then
2046 if (f_probe_owned(probe(i)%x, x_cb(-1), x_cb(m), 1)) then
2047 do s = -1, m
2048 distx(s) = x_cb(s) - probe(i)%x
2049 if (distx(s) < 0._wp) distx(s) = 1000._wp
2050 end do
2051 j = minloc(distx, 1)
2052 if (j == 1) j = 2 ! Pick first point if probe is at edge
2053 k = 0
2054 l = 0
2055
2056 if (chemistry) then
2057 do d = 1, num_species
2058 rhoyks(d) = q_cons_vf(eqn_idx%species%beg + d - 1)%sf(j - 2, k, l)
2059 end do
2060 end if
2061
2062 ! Computing/Sharing necessary state variables
2063 if (hypoelasticity) then
2064 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k, l, rho, gamma, pi_inf, qv, re, g_undamaged, &
2065 & fluid_pp(:)%G)
2066 else
2067 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k, l, rho, gamma, pi_inf, qv)
2068 end if
2069 do s = 1, num_vels
2070 vel(s) = q_cons_vf(eqn_idx%cont%end + s)%sf(j - 2, k, l)/rho
2071 end do
2072 do s = 1, num_fluids
2073 alpha(s) = q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(j - 2, k, l)
2074 alpha_rho(s) = q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, k, l)
2075 end do
2076
2077 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
2078
2079 if (hypoelasticity) then
2080 if (cont_damage) then
2081 ! Recover D = U_D/m_s, clamped to [0, 1]
2082 solid_partial_density = 0._wp
2083 do s = 1, num_fluids
2084 if (fluid_pp(s)%G > verysmall) then
2085 solid_partial_density = solid_partial_density + q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, &
2086 & k, l)
2087 end if
2088 end do
2089 damage_state = min(max(q_cons_vf(eqn_idx%damage)%sf(j - 2, k, l)/max(solid_partial_density, &
2090 & verysmall), 0._wp), 1._wp)
2091 end if
2092 g_damaged = g_undamaged*max(1._wp - damage_state, 0._wp)
2093
2094 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), &
2095 & dyn_p, pi_inf, gamma, rho, qv, rhoyks(:), pres, t, &
2096 & f_hypoelastic_energy(q_cons_vf, j - 2, k, l, rho, g_undamaged))
2097 else
2098 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), &
2099 & dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t)
2100 end if
2101
2102 if (hypoelasticity) then
2103 tau_e(1) = q_cons_vf(eqn_idx%stress%end)%sf(j - 2, k, l)/rho
2104 end if
2105
2106 if (bubbles_euler) then
2107 alf = q_cons_vf(eqn_idx%alf)%sf(j - 2, k, l)
2108 if (num_fluids == 3) then
2109 alfgr = q_cons_vf(eqn_idx%alf - 1)%sf(j - 2, k, l)
2110 end if
2111 do s = 1, nb
2112 nr(s) = q_cons_vf(qbmm_idx%rs(s))%sf(j - 2, k, l)
2113 nrdot(s) = q_cons_vf(qbmm_idx%vs(s))%sf(j - 2, k, l)
2114 end do
2115
2116 if (adv_n) then
2117 nbub = q_cons_vf(eqn_idx%n)%sf(j - 2, k, l)
2118 else
2119 nr3 = 0._wp
2120 do s = 1, nb
2121 nr3 = nr3 + weight(s)*(nr(s)**3._wp)
2122 end do
2123
2124 nbub = sqrt((4._wp*pi/3._wp)*nr3/alf)
2125 end if
2126#ifdef MFC_DEBUG
2127 print *, 'In probe, nbub: ', nbub
2128#endif
2129 if (qbmm) then
2130 m00 = q_cons_vf(qbmm_idx%moms(1, 1))%sf(j - 2, k, l)/nbub
2131 m10 = q_cons_vf(qbmm_idx%moms(1, 2))%sf(j - 2, k, l)/nbub
2132 m01 = q_cons_vf(qbmm_idx%moms(1, 3))%sf(j - 2, k, l)/nbub
2133 m20 = q_cons_vf(qbmm_idx%moms(1, 4))%sf(j - 2, k, l)/nbub
2134 m02 = q_cons_vf(qbmm_idx%moms(1, 6))%sf(j - 2, k, l)/nbub
2135
2136 m10 = m10/m00
2137 m01 = m01/m00
2138 m20 = m20/m00
2139 m02 = m02/m00
2140
2141 varr = m20 - m10**2._wp
2142 varv = m02 - m01**2._wp
2143 end if
2144 r(:) = nr(:)/nbub
2145 rdot(:) = nrdot(:)/nbub
2146
2147 ptilde = ptil(j - 2, k, l)
2148 ptot = pres - ptilde
2149 end if
2150
2151 ! Compute mixture sound Speed
2152 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
2153 if (hypoelasticity) c = sqrt(c*c + (4._wp/3._wp)*g_damaged/rho)
2154
2155 accel = accel_mag(j - 2, k, l)
2156 end if
2157 else if (p == 0) then
2158 if (chemistry) then
2159 do d = 1, num_species
2160 rhoyks(d) = q_cons_vf(eqn_idx%species%beg + d - 1)%sf(j - 2, k - 2, l)
2161 end do
2162 end if
2163
2164 if (f_probe_owned(probe(i)%x, x_cb(-1), x_cb(m), 1)) then
2165 if (f_probe_owned(probe(i)%y, y_cb(-1), y_cb(n), 2)) then
2166 do s = -1, m
2167 distx(s) = x_cb(s) - probe(i)%x
2168 if (distx(s) < 0._wp) distx(s) = 1000._wp
2169 end do
2170 do s = -1, n
2171 disty(s) = y_cb(s) - probe(i)%y
2172 if (disty(s) < 0._wp) disty(s) = 1000._wp
2173 end do
2174 j = minloc(distx, 1)
2175 k = minloc(disty, 1)
2176 if (j == 1) j = 2 ! Pick first point if probe is at edge
2177 if (k == 1) k = 2 ! Pick first point if probe is at edge
2178 l = 0
2179
2180 ! Computing/Sharing necessary state variables
2181 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k - 2, l, rho, gamma, pi_inf, qv, re, g_undamaged, &
2182 & fluid_pp(:)%G)
2183 do s = 1, num_vels
2184 vel(s) = q_cons_vf(eqn_idx%cont%end + s)%sf(j - 2, k - 2, l)/rho
2185 end do
2186 do s = 1, num_fluids
2187 alpha(s) = q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(j - 2, k - 2, l)
2188 alpha_rho(s) = q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, k - 2, l)
2189 end do
2190
2191 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
2192
2193 if (hypoelasticity) then
2194 if (cont_damage) then
2195 ! Recover D = U_D/m_s, clamped to [0, 1]
2196 solid_partial_density = 0._wp
2197 do s = 1, num_fluids
2198 if (fluid_pp(s)%G > verysmall) then
2199 solid_partial_density = solid_partial_density + q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j &
2200 & - 2, k - 2, l)
2201 end if
2202 end do
2203 damage_state = min(max(q_cons_vf(eqn_idx%damage)%sf(j - 2, k - 2, l)/max(solid_partial_density, &
2204 & verysmall), 0._wp), 1._wp)
2205 end if
2206 g_damaged = g_undamaged*max(1._wp - damage_state, 0._wp)
2207
2208 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, &
2209 & k - 2, l), dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t, &
2210 & f_hypoelastic_energy(q_cons_vf, j - 2, k - 2, l, rho, g_undamaged))
2211 else
2212 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, &
2213 & k - 2, l), dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t)
2214 end if
2215
2216 if (hypoelasticity) then
2217 do s = 1, 3
2218 tau_e(s) = q_cons_vf(eqn_idx%stress%beg + s - 1)%sf(j - 2, k - 2, l)/rho
2219 end do
2220 end if
2221
2222 if (bubbles_euler) then
2223 alf = q_cons_vf(eqn_idx%alf)%sf(j - 2, k - 2, l)
2224 do s = 1, nb
2225 nr(s) = q_cons_vf(qbmm_idx%rs(s))%sf(j - 2, k - 2, l)
2226 nrdot(s) = q_cons_vf(qbmm_idx%vs(s))%sf(j - 2, k - 2, l)
2227 end do
2228
2229 if (adv_n) then
2230 nbub = q_cons_vf(eqn_idx%n)%sf(j - 2, k - 2, l)
2231 else
2232 nr3 = 0._wp
2233 do s = 1, nb
2234 nr3 = nr3 + weight(s)*(nr(s)**3._wp)
2235 end do
2236
2237 nbub = sqrt((4._wp*pi/3._wp)*nr3/alf)
2238 end if
2239
2240 r(:) = nr(:)/nbub
2241 rdot(:) = nrdot(:)/nbub
2242 end if
2243 ! Compute mixture sound speed
2244 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
2245 if (hypoelasticity) c = sqrt(c*c + (4._wp/3._wp)*g_damaged/rho)
2246 end if
2247 end if
2248 else
2249 if (f_probe_owned(probe(i)%x, x_cb(-1), x_cb(m), 1)) then
2250 if (f_probe_owned(probe(i)%y, y_cb(-1), y_cb(n), 2)) then
2251 if (f_probe_owned(probe(i)%z, z_cb(-1), z_cb(p), 3)) then
2252 do s = -1, m
2253 distx(s) = x_cb(s) - probe(i)%x
2254 if (distx(s) < 0._wp) distx(s) = 1000._wp
2255 end do
2256 do s = -1, n
2257 disty(s) = y_cb(s) - probe(i)%y
2258 if (disty(s) < 0._wp) disty(s) = 1000._wp
2259 end do
2260 do s = -1, p
2261 distz(s) = z_cb(s) - probe(i)%z
2262 if (distz(s) < 0._wp) distz(s) = 1000._wp
2263 end do
2264 j = minloc(distx, 1)
2265 k = minloc(disty, 1)
2266 l = minloc(distz, 1)
2267 if (j == 1) j = 2 ! Pick first point if probe is at edge
2268 if (k == 1) k = 2 ! Pick first point if probe is at edge
2269 if (l == 1) l = 2 ! Pick first point if probe is at edge
2270
2271 ! Computing/Sharing necessary state variables
2272 call s_convert_to_mixture_variables(q_cons_vf, j - 2, k - 2, l - 2, rho, gamma, pi_inf, qv, re, &
2273 & g_undamaged, fluid_pp(:)%G)
2274 do s = 1, num_vels
2275 vel(s) = q_cons_vf(eqn_idx%cont%end + s)%sf(j - 2, k - 2, l - 2)/rho
2276 end do
2277 do s = 1, num_fluids
2278 alpha(s) = q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(j - 2, k - 2, l - 2)
2279 alpha_rho(s) = q_cons_vf(eqn_idx%cont%beg + s - 1)%sf(j - 2, k - 2, l - 2)
2280 end do
2281
2282 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
2283
2284 if (chemistry) then
2285 do d = 1, num_species
2286 rhoyks(d) = q_cons_vf(eqn_idx%species%beg + d - 1)%sf(j - 2, k - 2, l - 2)
2287 end do
2288 end if
2289
2290 if (hypoelasticity) then
2291 if (cont_damage) then
2292 ! Recover D = U_D/m_s, clamped to [0, 1]
2293 solid_partial_density = 0._wp
2294 do s = 1, num_fluids
2295 if (fluid_pp(s)%G > verysmall) then
2296 solid_partial_density = solid_partial_density + q_cons_vf(eqn_idx%cont%beg + s &
2297 & - 1)%sf(j - 2, k - 2, l - 2)
2298 end if
2299 end do
2300 damage_state = min(max(q_cons_vf(eqn_idx%damage)%sf(j - 2, k - 2, &
2301 & l - 2)/max(solid_partial_density, verysmall), 0._wp), 1._wp)
2302 end if
2303 g_damaged = g_undamaged*max(1._wp - damage_state, 0._wp)
2304
2305 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k - 2, l - 2), &
2306 & q_cons_vf(eqn_idx%alf)%sf(j - 2, k - 2, l - 2), dyn_p, pi_inf, gamma, &
2307 & rho, qv, rhoyks, pres, t, f_hypoelastic_energy(q_cons_vf, j - 2, k - 2, &
2308 & l - 2, rho, g_undamaged))
2309 else
2310 call s_compute_pressure(q_cons_vf(eqn_idx%E)%sf(j - 2, k - 2, l - 2), &
2311 & q_cons_vf(eqn_idx%alf)%sf(j - 2, k - 2, l - 2), dyn_p, pi_inf, gamma, &
2312 & rho, qv, rhoyks, pres, t)
2313 end if
2314
2315 if (hypoelasticity) then
2316 do s = 1, 6
2317 tau_e(s) = q_cons_vf(eqn_idx%stress%beg + s - 1)%sf(j - 2, k - 2, l - 2)/rho
2318 end do
2319 end if
2320
2321 ! Compute mixture sound speed
2322 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, alpha, c, alpha_rho)
2323 if (hypoelasticity) c = sqrt(c*c + (4._wp/3._wp)*g_damaged/rho)
2324
2325 accel = accel_mag(j - 2, k - 2, l - 2)
2326 end if
2327 end if
2328 end if
2329 end if
2330 if (num_procs > 1) then
2331# 1836 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2332 tmp = rho
2333# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2334 call s_mpi_allreduce_sum(tmp, rho)
2335# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2336# 1836 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2337 tmp = pres
2338# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2339 call s_mpi_allreduce_sum(tmp, pres)
2340# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2341# 1836 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2342 tmp = gamma
2343# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2344 call s_mpi_allreduce_sum(tmp, gamma)
2345# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2346# 1836 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2347 tmp = pi_inf
2348# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2349 call s_mpi_allreduce_sum(tmp, pi_inf)
2350# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2351# 1836 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2352 tmp = qv
2353# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2354 call s_mpi_allreduce_sum(tmp, qv)
2355# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2356# 1836 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2357 tmp = c
2358# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2359 call s_mpi_allreduce_sum(tmp, c)
2360# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2361# 1836 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2362 tmp = accel
2363# 1837 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2364 call s_mpi_allreduce_sum(tmp, accel)
2365# 1838 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2366# 1839 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2367
2368 do s = 1, num_vels
2369 tmp = vel(s)
2370 call s_mpi_allreduce_sum(tmp, vel(s))
2371 end do
2372
2373 if (bubbles_euler) then
2374# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2375 tmp = alf
2376# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2377 call s_mpi_allreduce_sum(tmp, alf)
2378# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2379# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2380 tmp = alfgr
2381# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2382 call s_mpi_allreduce_sum(tmp, alfgr)
2383# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2384# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2385 tmp = nbub
2386# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2387 call s_mpi_allreduce_sum(tmp, nbub)
2388# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2389# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2390 tmp = nr(1)
2391# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2392 call s_mpi_allreduce_sum(tmp, nr(1))
2393# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2394# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2395 tmp = nrdot(1)
2396# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2397 call s_mpi_allreduce_sum(tmp, nrdot(1))
2398# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2399# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2400 tmp = m00
2401# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2402 call s_mpi_allreduce_sum(tmp, m00)
2403# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2404# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2405 tmp = r(1)
2406# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2407 call s_mpi_allreduce_sum(tmp, r(1))
2408# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2409# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2410 tmp = rdot(1)
2411# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2412 call s_mpi_allreduce_sum(tmp, rdot(1))
2413# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2414# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2415 tmp = ptilde
2416# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2417 call s_mpi_allreduce_sum(tmp, ptilde)
2418# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2419# 1847 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2420 tmp = ptot
2421# 1848 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2422 call s_mpi_allreduce_sum(tmp, ptot)
2423# 1849 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2424# 1850 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2425
2426 if (qbmm) then
2427# 1853 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2428 tmp = varr
2429# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2430 call s_mpi_allreduce_sum(tmp, varr)
2431# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2432# 1853 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2433 tmp = varv
2434# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2435 call s_mpi_allreduce_sum(tmp, varv)
2436# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2437# 1853 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2438 tmp = m10
2439# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2440 call s_mpi_allreduce_sum(tmp, m10)
2441# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2442# 1853 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2443 tmp = m01
2444# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2445 call s_mpi_allreduce_sum(tmp, m01)
2446# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2447# 1853 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2448 tmp = m20
2449# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2450 call s_mpi_allreduce_sum(tmp, m20)
2451# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2452# 1853 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2453 tmp = m02
2454# 1854 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2455 call s_mpi_allreduce_sum(tmp, m02)
2456# 1855 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2457# 1856 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2458 end if
2459 end if
2460
2461 if (hypoelasticity) then
2462 do s = 1, (num_dims*(num_dims + 1))/2
2463 tmp = tau_e(s)
2464 call s_mpi_allreduce_sum(tmp, tau_e(s))
2465 end do
2466 end if
2467
2468 if (cont_damage) then
2469 tmp = damage_state
2470 call s_mpi_allreduce_sum(tmp, damage_state)
2471 end if
2472 end if
2473 if (proc_rank == 0) then
2474 if (n == 0) then
2475 if (bubbles_euler .and. (num_fluids <= 2)) then
2476 if (qbmm) then
2477 write (i + 30, '(6x,f12.6,14f28.16)') nondim_time, rho, vel(1), pres, alf, r(1), rdot(1), nr(1), &
2478 & nrdot(1), varr, varv, m10, m01, m20, m02
2479 else
2480 write (i + 30, '(6x,f12.6,8f24.8)') nondim_time, rho, vel(1), pres, alf, r(1), rdot(1), nr(1), nrdot(1)
2481 ! ptilde, & ptot
2482 end if
2483 else if (bubbles_euler .and. (num_fluids == 3)) then
2484 write (i + 30, &
2485 & '(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)') &
2486 & nondim_time, rho, vel(1), pres, alf, alfgr, nr(1), nrdot(1), r(1), rdot(1), ptilde, ptot
2487 else if (bubbles_euler .and. num_fluids == 4) then
2488 write (i + 30, &
2489 & '(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)') &
2490 & 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, &
2491 & 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), &
2492 & 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), &
2493 & q_cons_vf(10)%sf(j - 2, 0, 0), nbub, r(1), rdot(1)
2494 else if (hypoelasticity) then
2495 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), pres, tau_e(1)
2496 else
2497 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), pres
2498 end if
2499 else if (p == 0) then
2500 if (bubbles_euler) then
2501# 1900 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2502 write (i + 30, '(6X,10F24.8)') nondim_time, rho, vel(1), vel(2), pres, alf, nr(1), nrdot(1), r(1), &
2503 & rdot(1)
2504# 1903 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2505 else if (hypoelasticity) then
2506# 1905 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2507 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8,F24.8,' // 'F24.8,F24.8,F24.8)') nondim_time, rho, &
2508 & vel(1), vel(2), pres, tau_e(1), tau_e(2), tau_e(3)
2509# 1908 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2510 else
2511# 1910 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2512 write (i + 30, '(6X,F12.6,F24.8,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), vel(2), pres
2513# 1912 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2514 end if
2515 else
2516# 1915 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2517 if (hypoelasticity) then
2518 write (i + 30, '(6X,F12.6,16F24.8)') nondim_time, rho, vel(1), vel(2), vel(3), pres, gamma, pi_inf, &
2519 & qv, c, accel, tau_e(1), tau_e(2), tau_e(3), tau_e(4), tau_e(5), tau_e(6)
2520 else
2521 write (i + 30, &
2522 & '(6X,F12.6,F24.8,F24.8,F24.8,F24.8,' // 'F24.8,F24.8,F24.8,F24.8,F24.8,' // 'F24.8)') &
2523 & nondim_time, rho, vel(1), vel(2), vel(3), pres, gamma, pi_inf, qv, c, accel
2524 end if
2525# 1924 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2526 end if
2527 end if
2528 end do
2529
2530 end subroutine s_write_probe_files
2531
2532 !> Write footer with stability criteria extrema and run-time to the information file, then close it
2533 impure subroutine s_close_run_time_information_file
2534
2535 real(wp) :: run_time !< Run-time of the simulation
2536
2537 write (3, '(A)') ' '
2538 write (3, '(A)') ''
2539
2540 write (3, '(A,F9.6)') 'ICFL Max: ', icfl_max
2541 if (surface_tension) write (3, '(A,F9.6)') 'CCFL Max: ', ccfl_max
2542 if (heat_conduction) write (3, '(A,F9.6)') 'TCFL Max: ', tcfl_max
2543 if (viscous) write (3, '(A,F9.6)') 'VCFL Max: ', vcfl_max
2544 if (viscous) write (3, '(A,ES16.6)') 'Rc Min: ', rc_min
2545
2546 call cpu_time(run_time)
2547
2548 write (3, '(A)') ''
2549 write (3, '(A,I0,A)') 'Run-time: ', int(anint(run_time)), 's'
2550 write (3, '(A)') ' '
2551 close (3)
2552
2553 end subroutine s_close_run_time_information_file
2554
2555 !> Closes probe files
2556 impure subroutine s_close_probe_files
2557
2558 integer :: i !< Generic loop iterator
2559
2560 do i = 1, num_probes
2561 close (i + 30)
2562 end do
2563
2564 end subroutine s_close_probe_files
2565
2566 !> Initialize the data output module
2567 impure subroutine s_initialize_data_output_module
2568
2569 integer :: i, m_ds, n_ds, p_ds
2570
2571 if (run_time_info) then
2572 icfl_max = 0._wp
2573 if (surface_tension) then
2574 ccfl_max = 0._wp
2575 end if
2576 if (heat_conduction) then
2577 tcfl_max = 0._wp
2578 end if
2579 if (viscous) then
2580 vcfl_max = 0._wp
2581 rc_min = 1.e12_wp
2582 end if
2583 end if
2584
2585 if (down_sample) then
2586 m_ds = int((m + 1)/3) - 1
2587 n_ds = int((n + 1)/3) - 1
2588 p_ds = int((p + 1)/3) - 1
2589
2590 allocate (q_cons_temp_ds(1:sys_size))
2591 do i = 1, sys_size
2592 allocate (q_cons_temp_ds(i)%sf(-1:m_ds + 1,-1:n_ds + 1,-1:p_ds + 1))
2593 end do
2594 end if
2595
2596 end subroutine s_initialize_data_output_module
2597
2598 !> Module deallocation and/or disassociation procedures
2599 impure subroutine s_finalize_data_output_module
2600
2601 integer :: i
2602
2603 if (down_sample) then
2604 do i = 1, sys_size
2605 deallocate (q_cons_temp_ds(i)%sf)
2606 end do
2607 deallocate (q_cons_temp_ds)
2608 end if
2609
2610 end subroutine s_finalize_data_output_module
2611
2612 !> 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
2613 !! samples (left of the face); the first rank in a direction also owns its lower face.
2614 logical function f_probe_owned(v, lo, hi, dir)
2615
2616 real(wp), intent(in) :: v, lo, hi
2617 integer, intent(in) :: dir
2618 logical :: first
2619
2620 first = num_procs == 1
2621 if (.not. first) first = proc_coords(dir) == 0
2622 f_probe_owned = v <= hi .and. (v > lo .or. (first .and. v >= lo))
2623
2624 end function f_probe_owned
2625
2626end module m_data_output
q_cons_vf
type(scalar_field), dimension(sys_size), intent(inout) q_cons_vf
Definition m_phase_change.fpp.f90:1041
k
integer, intent(in) k
Definition m_phase_change.fpp.f90:1043
j
integer, intent(in) j
Definition m_phase_change.fpp.f90:1043
l
integer, intent(in) l
Definition m_phase_change.fpp.f90:1043
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:1932
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:1965
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:1604
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:1835
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:2568
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:2600
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:1887
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:1709
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:1476
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:2534
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:1762
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:2615
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:2557
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:1904
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:1505
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:1578
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:3141
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:2534
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