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/post_process/m_data_output.fpp"
2!>
3!! @file
4!! @brief Contains module m_data_output
5
6!> @brief Writes post-processed grid and flow-variable data to Silo-HDF5 or binary database files
7module m_data_output
8
9 use m_derived_types
10 use m_global_parameters
11 use m_derived_variables
12 use m_mpi_proxy
13 use m_compile_specific
14 use m_helper
15 use m_variables_conversion
16 use m_eos
17 use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, format_silo, format_binary, precision_single
18
19 implicit none
20
21 private; public :: s_initialize_data_output_module, s_define_output_region, s_open_formatted_database_file, &
22 & s_open_intf_data_file, s_open_energy_data_file, s_write_grid_to_formatted_database_file, &
23 & s_write_variable_to_formatted_database_file, s_write_lag_bubbles_results_to_text, &
24 & s_write_lag_bubbles_to_formatted_database_file, s_write_ib_state_files, s_write_intf_data_file, &
25 & s_write_energy_data_file, s_write_ib_bodies_to_formatted_database_file, s_close_formatted_database_file, &
26 & s_close_intf_data_file, s_close_energy_data_file, s_finalize_data_output_module, out
27
28 ! Include Silo-HDF5 interface library
29 include 'silo_f9x.inc'
30
31 !> Output workspace: flow variable buffers, VisIt extents/offsets, directory paths, file handles, and variable count.
32 type(output_context) :: out
33
34 ! Generic error flag for Silo-HDF5 and Binary I/O operations
35 integer, private :: err
36
37contains
38
39 !> Allocate storage arrays, configure output directories, and count flow variables for formatted database output.
40 impure subroutine s_initialize_data_output_module()
41
42 character(LEN=len_trim(case_dir) + 2*name_len) :: file_loc
43 logical :: dir_check
44 integer :: i
45
46 allocate (out%q_sf(-offset_x%beg:m + offset_x%end,-offset_y%beg:n + offset_y%end,-offset_z%beg:p + offset_z%end))
47 if (grid_geometry == 3) then
48 allocate (out%cyl_q_sf(-offset_y%beg:n + offset_y%end,-offset_z%beg:p + offset_z%end,-offset_x%beg:m + offset_x%end))
49 end if
50
51 if (precision == precision_single) then
52 allocate (out%q_sf_s(-offset_x%beg:m + offset_x%end,-offset_y%beg:n + offset_y%end,-offset_z%beg:p + offset_z%end))
53 if (grid_geometry == 3) then
54 allocate (out%cyl_q_sf_s(-offset_y%beg:n + offset_y%end,-offset_z%beg:p + offset_z%end, &
55 & -offset_x%beg:m + offset_x%end))
56 end if
57 end if
58
59 if (n == 0) then
60 allocate (out%q_root_sf(0:m_root,0:0,0:0))
61 if (precision == precision_single) then
62 allocate (out%q_root_sf_s(0:m_root,0:0,0:0))
63 end if
64 end if
65
66 ! Allocating the spatial and data extents and also the variables for the offsets and the one bookkeeping the number of
67 ! cell-boundaries in each active coordinate direction. Note that all these variables are only needed by the Silo-HDF5 format
68 ! for multidimensional data.
69 if (format == format_silo) then
70 allocate (out%data_extents(1:2,0:num_procs - 1))
71
72 if (p > 0) then
73 allocate (out%spatial_extents(1:6,0:num_procs - 1))
74 allocate (out%lo_offset(1:3))
75 allocate (out%hi_offset(1:3))
76 allocate (out%dims(1:3))
77 else if (n > 0) then
78 allocate (out%spatial_extents(1:4,0:num_procs - 1))
79 allocate (out%lo_offset(1:2))
80 allocate (out%hi_offset(1:2))
81 allocate (out%dims(1:2))
82 else
83 allocate (out%spatial_extents(1:2,0:num_procs - 1))
84 allocate (out%lo_offset(1:1))
85 allocate (out%hi_offset(1:1))
86 allocate (out%dims(1:1))
87 end if
88 end if
89
90 ! The size of the ghost zone layer in each of the active coordinate directions was set in the module m_mpi_proxy.f90. The
91 ! results are now transferred to the local variables of this module when they are required by the Silo-HDF5 format, for
92 ! multidimensional data sets. With the same, latter, requirements, the variables bookkeeping the number of cell-boundaries
93 ! in each active coordinate direction are also set here.
94 if (format == format_silo) then
95 if (p > 0) then
96 if (grid_geometry == 3) then
97 out%lo_offset(:) = (/offset_y%beg, offset_z%beg, offset_x%beg/)
98 out%hi_offset(:) = (/offset_y%end, offset_z%end, offset_x%end/)
99 else
100 out%lo_offset(:) = (/offset_x%beg, offset_y%beg, offset_z%beg/)
101 out%hi_offset(:) = (/offset_x%end, offset_y%end, offset_z%end/)
102 end if
103
104 if (grid_geometry == 3) then
105 out%dims(:) = (/n + offset_y%beg + offset_y%end + 2, p + offset_z%beg + offset_z%end + 2, &
106 & m + offset_x%beg + offset_x%end + 2/)
107 else
108 out%dims(:) = (/m + offset_x%beg + offset_x%end + 2, n + offset_y%beg + offset_y%end + 2, &
109 & p + offset_z%beg + offset_z%end + 2/)
110 end if
111 else if (n > 0) then
112 out%lo_offset(:) = (/offset_x%beg, offset_y%beg/)
113 out%hi_offset(:) = (/offset_x%end, offset_y%end/)
114
115 out%dims(:) = (/m + offset_x%beg + offset_x%end + 2, n + offset_y%beg + offset_y%end + 2/)
116 else
117 out%lo_offset(:) = (/offset_x%beg/)
118 out%hi_offset(:) = (/offset_x%end/)
119 out%dims(:) = (/m + offset_x%beg + offset_x%end + 2/)
120 end if
121 end if
122
123 if (format == format_silo) then
124 out%dbdir = trim(case_dir) // '/silo_hdf5'
125
126 write (out%proc_rank_dir, '(A,I0)') '/p', proc_rank
127
128 out%proc_rank_dir = trim(out%dbdir) // trim(out%proc_rank_dir)
129
130 file_loc = trim(out%proc_rank_dir) // '/.'
131
132 call my_inquire(file_loc, dir_check)
133 if (dir_check .neqv. .true.) then
134 call s_create_directory(trim(out%proc_rank_dir))
135 end if
136
137 if (proc_rank == 0) then
138 out%rootdir = trim(out%dbdir) // '/root'
139
140 file_loc = trim(out%rootdir) // '/.'
141
142 call my_inquire(file_loc, dir_check)
143 if (dir_check .neqv. .true.) then
144 call s_create_directory(trim(out%rootdir))
145 end if
146 end if
147 else
148 out%dbdir = trim(case_dir) // '/binary'
149
150 write (out%proc_rank_dir, '(A,I0)') '/p', proc_rank
151
152 out%proc_rank_dir = trim(out%dbdir) // trim(out%proc_rank_dir)
153
154 file_loc = trim(out%proc_rank_dir) // '/.'
155
156 call my_inquire(file_loc, dir_check)
157
158 if (dir_check .neqv. .true.) then
159 call s_create_directory(trim(out%proc_rank_dir))
160 end if
161
162 if (n == 0 .and. proc_rank == 0) then
163 out%rootdir = trim(out%dbdir) // '/root'
164
165 file_loc = trim(out%rootdir) // '/.'
166
167 call my_inquire(file_loc, dir_check)
168
169 if (dir_check .neqv. .true.) then
170 call s_create_directory(trim(out%rootdir))
171 end if
172 end if
173 end if
174
175 if (bubbles_lagrange) then ! Lagrangian solver
176 if (lag_txt_wrt) then
177 out%dbdir = trim(case_dir) // '/lag_bubbles_post_process'
178 file_loc = trim(out%dbdir) // '/.'
179 call my_inquire(file_loc, dir_check)
180
181 if (dir_check .neqv. .true.) then
182 call s_create_directory(trim(out%dbdir))
183 end if
184 end if
185 end if
186
187 ! Contrary to the Silo-HDF5 database format, handles of the Binary database master/root and slave/local process files are
188 ! perfectly static throughout post-process. Hence, they are set here so that they do not have to be repetitively computed in
189 ! later procedures.
190 if (format == format_binary) then
191 if (n == 0 .and. proc_rank == 0) out%dbroot = 2
192 out%dbfile = 1
193 end if
194
195 if (format == format_binary) then
196 out%dbvars = 0
197
198 if ((model_eqns == model_eqns_5eq) .or. (model_eqns == model_eqns_6eq)) then
199 do i = 1, num_fluids
200 if (alpha_rho_wrt(i) .or. (cons_vars_wrt .or. prim_vars_wrt)) then
201 out%dbvars = out%dbvars + 1
202 end if
203 end do
204 end if
205
206 if ((rho_wrt .or. (model_eqns == model_eqns_gamma_law .and. (cons_vars_wrt .or. prim_vars_wrt))) &
207 & .and. (.not. relativity)) then
208 out%dbvars = out%dbvars + 1
209 end if
210
211 if (relativity .and. (rho_wrt .or. prim_vars_wrt)) out%dbvars = out%dbvars + 1
212 if (relativity .and. (rho_wrt .or. cons_vars_wrt)) out%dbvars = out%dbvars + 1
213
214 do i = 1, eqn_idx%E - eqn_idx%mom%beg
215 if (mom_wrt(i) .or. cons_vars_wrt) out%dbvars = out%dbvars + 1
216 end do
217
218 do i = 1, eqn_idx%E - eqn_idx%mom%beg
219 if (vel_wrt(i) .or. prim_vars_wrt) out%dbvars = out%dbvars + 1
220 end do
221
222 do i = 1, eqn_idx%E - eqn_idx%mom%beg
223 if (flux_wrt(i)) out%dbvars = out%dbvars + 1
224 end do
225
226 if (e_wrt .or. cons_vars_wrt) out%dbvars = out%dbvars + 1
227 if (pres_wrt .or. prim_vars_wrt) out%dbvars = out%dbvars + 1
228 if (hypoelasticity) out%dbvars = out%dbvars + (num_dims*(num_dims + 1))/2
229 if (cont_damage) out%dbvars = out%dbvars + 1
230 if (hyper_cleaning) out%dbvars = out%dbvars + 1
231
232 if (mhd) then
233 if (n == 0) then
234 out%dbvars = out%dbvars + 2
235 else
236 out%dbvars = out%dbvars + 3
237 end if
238 end if
239
240 if ((model_eqns == model_eqns_5eq) .or. (model_eqns == model_eqns_6eq)) then
241 do i = 1, num_fluids - 1
242 if (alpha_wrt(i) .or. (cons_vars_wrt .or. prim_vars_wrt)) then
243 out%dbvars = out%dbvars + 1
244 end if
245 end do
246
247 if (alpha_wrt(num_fluids) .or. (cons_vars_wrt .or. prim_vars_wrt)) then
248 out%dbvars = out%dbvars + 1
249 end if
250 end if
251
252 if (gamma_wrt .or. (model_eqns == model_eqns_gamma_law .and. (cons_vars_wrt .or. prim_vars_wrt))) then
253 out%dbvars = out%dbvars + 1
254 end if
255
256 if (heat_ratio_wrt) out%dbvars = out%dbvars + 1
257
258 if (pi_inf_wrt .or. (model_eqns == model_eqns_gamma_law .and. (cons_vars_wrt .or. prim_vars_wrt))) then
259 out%dbvars = out%dbvars + 1
260 end if
261
262 if (pres_inf_wrt) out%dbvars = out%dbvars + 1
263 if (c_wrt) out%dbvars = out%dbvars + 1
264
265 if (p > 0) then
266 do i = 1, num_vels
267 if (omega_wrt(i)) out%dbvars = out%dbvars + 1
268 end do
269 else if (n > 0) then
270 do i = 1, num_vels
271 if (omega_wrt(i)) out%dbvars = out%dbvars + 1
272 end do
273 end if
274
275 if (schlieren_wrt) out%dbvars = out%dbvars + 1
276 end if
277
278 end subroutine s_initialize_data_output_module
279
280 !> Compute the cell-index bounds for the user-specified partial output domain in each coordinate direction.
281 impure subroutine s_define_output_region
282
283 integer :: i
284 integer :: lower_bound, upper_bound
285
286# 286 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
287 if (m == 0) return ! Early return for y or z if simulation is 1D or 2D
288# 287 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
289
290 lower_bound = -offset_x%beg
291# 289 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
292 upper_bound = m + offset_x%end
293# 290 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
294
295 do i = lower_bound, upper_bound
296 if (x_cc(i) > x_output%beg) then
297# 293 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
298 x_output_idx%beg = i + offset_x%beg
299# 294 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
300 exit
301 end if
302 end do
303
304 do i = upper_bound, lower_bound, -1
305 if (x_cc(i) < x_output%end) then
306# 300 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
307 x_output_idx%end = i + offset_x%beg
308# 301 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
309 exit
310 end if
311 end do
312
313 ! If no grid points are within the output region
314 if ((x_cc(lower_bound) > x_output%end) .or. (x_cc(upper_bound) < x_output%beg)) then
315# 307 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
316 x_output_idx%beg = 0
317# 308 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
318 x_output_idx%end = 0
319# 309 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
320 end if
321# 286 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
322 if (n == 0) return ! Early return for y or z if simulation is 1D or 2D
323# 287 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
324
325 lower_bound = -offset_y%beg
326# 289 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
327 upper_bound = n + offset_y%end
328# 290 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
329
330 do i = lower_bound, upper_bound
331 if (y_cc(i) > y_output%beg) then
332# 293 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
333 y_output_idx%beg = i + offset_y%beg
334# 294 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
335 exit
336 end if
337 end do
338
339 do i = upper_bound, lower_bound, -1
340 if (y_cc(i) < y_output%end) then
341# 300 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
342 y_output_idx%end = i + offset_y%beg
343# 301 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
344 exit
345 end if
346 end do
347
348 ! If no grid points are within the output region
349 if ((y_cc(lower_bound) > y_output%end) .or. (y_cc(upper_bound) < y_output%beg)) then
350# 307 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
351 y_output_idx%beg = 0
352# 308 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
353 y_output_idx%end = 0
354# 309 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
355 end if
356# 286 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
357 if (p == 0) return ! Early return for y or z if simulation is 1D or 2D
358# 287 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
359
360 lower_bound = -offset_z%beg
361# 289 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
362 upper_bound = p + offset_z%end
363# 290 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
364
365 do i = lower_bound, upper_bound
366 if (z_cc(i) > z_output%beg) then
367# 293 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
368 z_output_idx%beg = i + offset_z%beg
369# 294 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
370 exit
371 end if
372 end do
373
374 do i = upper_bound, lower_bound, -1
375 if (z_cc(i) < z_output%end) then
376# 300 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
377 z_output_idx%end = i + offset_z%beg
378# 301 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
379 exit
380 end if
381 end do
382
383 ! If no grid points are within the output region
384 if ((z_cc(lower_bound) > z_output%end) .or. (z_cc(upper_bound) < z_output%beg)) then
385# 307 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
386 z_output_idx%beg = 0
387# 308 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
388 z_output_idx%end = 0
389# 309 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
390 end if
391# 311 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
392
393 end subroutine s_define_output_region
394
395 !> Open (or create) the Silo-HDF5 or Binary formatted database slave and master files for a given time step.
396 impure subroutine s_open_formatted_database_file(t_step)
397
398 integer, intent(in) :: t_step
399 character(LEN=len_trim(case_dir) + 3*name_len) :: file_loc
400 integer :: ierr
401
402 if (format == format_silo) then
403 write (file_loc, '(A,I0,A)') '/', t_step, '.silo'
404 file_loc = trim(out%proc_rank_dir) // trim(file_loc)
405
406 ierr = dbcreate(trim(file_loc), len_trim(file_loc), db_clobber, db_local, 'MFC v3.0', 8, db_hdf5, out%dbfile)
407
408 if (out%dbfile == -1) then
409 call s_mpi_abort('Unable to create Silo-HDF5 database ' // 'slave file ' // trim(file_loc) // '. ' // 'Exiting.')
410 end if
411
412 if (proc_rank == 0) then
413 write (file_loc, '(A,I0,A)') '/collection_', t_step, '.silo'
414 file_loc = trim(out%rootdir) // trim(file_loc)
415
416 ierr = dbcreate(trim(file_loc), len_trim(file_loc), db_clobber, db_local, 'MFC v3.0', 8, db_hdf5, out%dbroot)
417
418 if (out%dbroot == -1) then
419 call s_mpi_abort('Unable to create Silo-HDF5 database ' // 'master file ' // trim(file_loc) // '. ' &
420 & // 'Exiting.')
421 end if
422 end if
423 else
424 write (file_loc, '(A,I0,A)') '/', t_step, '.dat'
425 file_loc = trim(out%proc_rank_dir) // trim(file_loc)
426
427 open (out%dbfile, iostat=err, file=trim(file_loc), form='unformatted', status='replace')
428
429 if (err /= 0) then
430 call s_mpi_abort('Unable to create Binary database slave ' // 'file ' // trim(file_loc) // '. Exiting.')
431 end if
432
433 if (output_partial_domain) then
434 write (out%dbfile) x_output_idx%end - x_output_idx%beg, y_output_idx%end - y_output_idx%beg, &
435 & z_output_idx%end - z_output_idx%beg, out%dbvars
436 else
437 write (out%dbfile) m, n, p, out%dbvars
438 end if
439
440 if (n == 0 .and. proc_rank == 0) then
441 write (file_loc, '(A,I0,A)') '/', t_step, '.dat'
442 file_loc = trim(out%rootdir) // trim(file_loc)
443
444 open (out%dbroot, iostat=err, file=trim(file_loc), form='unformatted', status='replace')
445
446 if (err /= 0) then
447 call s_mpi_abort('Unable to create Binary database ' // 'master file ' // trim(file_loc) // '. Exiting.')
448 end if
449
450 if (output_partial_domain) then
451 write (out%dbroot) x_output_idx%end - x_output_idx%beg, 0, 0, out%dbvars
452 else
453 write (out%dbroot) m_root, 0, 0, out%dbvars
454 end if
455 end if
456 end if
457
458 end subroutine s_open_formatted_database_file
459
460 !> Open the interface data file for appending extracted interface coordinates.
461 impure subroutine s_open_intf_data_file()
462
463 character(LEN=path_len + 3*name_len) :: file_path
464
465 write (file_path, '(A)') '/intf_data.dat'
466 file_path = trim(case_dir) // trim(file_path)
467
468 open (211, file=trim(file_path), form='formatted', position='append', status='unknown')
469
470 end subroutine s_open_intf_data_file
471
472 !> Open the energy data file for appending volume-integrated energy budget quantities.
473 impure subroutine s_open_energy_data_file()
474
475 character(LEN=path_len + 3*name_len) :: file_path
476
477 write (file_path, '(A)') '/eng_data.dat'
478 file_path = trim(case_dir) // trim(file_path)
479
480 open (251, file=trim(file_path), form='formatted', position='append', status='unknown')
481
482 end subroutine s_open_energy_data_file
483
484 !> Write the computational grid (cell-boundary coordinates) to the formatted database slave and master files.
485 impure subroutine s_write_grid_to_formatted_database_file(t_step)
486
487 integer, intent(in) :: t_step
488
489 ! NAG compiler requires these to be statically sized
490 character(LEN=4*name_len), dimension(num_procs) :: meshnames
491 integer, dimension(num_procs) :: meshtypes
492 integer :: i
493 integer :: ierr
494 integer :: extents_size
495
496 if (format == format_silo) then
497 ! For multidimensional data sets, the spatial extents of all of the grid(s) handled by the local processor(s) are
498 ! recorded so that they may be written, by root processor, to the formatted database master file.
499 if (num_procs > 1) then
500 call s_mpi_gather_spatial_extents(out%spatial_extents)
501 else if (p > 0) then
502 if (grid_geometry == 3) then
503 out%spatial_extents(:,0) = (/minval(y_cb), minval(z_cb), minval(x_cb), maxval(y_cb), maxval(z_cb), &
504 & maxval(x_cb)/)
505 else
506 out%spatial_extents(:,0) = (/minval(x_cb), minval(y_cb), minval(z_cb), maxval(x_cb), maxval(y_cb), &
507 & maxval(z_cb)/)
508 end if
509 else if (n > 0) then
510 out%spatial_extents(:,0) = (/minval(x_cb), minval(y_cb), maxval(x_cb), maxval(y_cb)/)
511 else
512 out%spatial_extents(:,0) = (/minval(x_cb), maxval(x_cb)/)
513 end if
514
515 ! Next, the root processor proceeds to record all of the spatial extents in the formatted database master file. In
516 ! addition, it also records a sub-domain connectivity map so that the entire grid may be reassembled by looking at the
517 ! master file.
518 if (proc_rank == 0) then
519 do i = 1, num_procs
520 write (meshnames(i), '(A,I0,A,I0,A)') '../p', i - 1, '/', t_step, '.silo:rectilinear_grid'
521 end do
522
523 meshtypes = db_quad_rect
524
525 err = dbset2dstrlen(len(meshnames(1)))
526 err = dbmkoptlist(2, out%optlist)
527 extents_size = size(out%spatial_extents, 1)
528 err = dbaddiopt(out%optlist, dbopt_extents_size, extents_size)
529 err = dbadddopt(out%optlist, dbopt_extents, out%spatial_extents)
530 err = dbputmmesh(out%dbroot, 'rectilinear_grid', 16, num_procs, meshnames, len_trim(meshnames), meshtypes, &
531 & out%optlist, ierr)
532 err = dbfreeoptlist(out%optlist)
533 end if
534
535 ! Finally, the local quadrilateral mesh, either 2D or 3D, along with its offsets that indicate the presence and size of
536 ! ghost zone layer(s), are put in the formatted database slave file.
537
538 ! Silo carries the mesh coordinates in their own datatype, separate from the flow variables, so the cell boundaries
539 ! are copied down to single precision here when that is what was asked for. Without this the mesh is always written as
540 ! DB_DOUBLE, which keeps downstream readers on a double-precision path regardless of `precision`.
541 if (precision == precision_single) then
542 x_cb_s = real(x_cb, sp)
543 if (n > 0) then
544 y_cb_s = real(y_cb, sp)
545 if (p > 0) z_cb_s = real(z_cb, sp)
546 end if
547 end if
548
549 if (p > 0) then
550 err = dbmkoptlist(2, out%optlist)
551 err = dbaddiaopt(out%optlist, dbopt_lo_offset, size(out%lo_offset), out%lo_offset)
552 err = dbaddiaopt(out%optlist, dbopt_hi_offset, size(out%hi_offset), out%hi_offset)
553# 473 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
554 if (precision == 1) then
555# 474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
556 if (grid_geometry == 3) then
557 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, y_cb_s, z_cb_s, &
558# 476 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
559 & x_cb_s, out%dims, 3, db_float, db_collinear, out%optlist, ierr)
560# 477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
561 else
562 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, x_cb_s, y_cb_s, &
563# 479 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
564 & z_cb_s, out%dims, 3, db_float, db_collinear, out%optlist, ierr)
565# 480 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
566 end if
567 end if
568# 473 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
569 if (precision == 2) then
570# 474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
571 if (grid_geometry == 3) then
572 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, y_cb, z_cb, &
573# 476 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
574 & x_cb, out%dims, 3, db_double, db_collinear, out%optlist, ierr)
575# 477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
576 else
577 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, x_cb, y_cb, &
578# 479 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
579 & z_cb, out%dims, 3, db_double, db_collinear, out%optlist, ierr)
580# 480 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
581 end if
582 end if
583# 483 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
584 err = dbfreeoptlist(out%optlist)
585 else if (n > 0) then
586 err = dbmkoptlist(2, out%optlist)
587 err = dbaddiaopt(out%optlist, dbopt_lo_offset, size(out%lo_offset), out%lo_offset)
588 err = dbaddiaopt(out%optlist, dbopt_hi_offset, size(out%hi_offset), out%hi_offset)
589# 489 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
590 if (precision == 1) then
591# 490 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
592 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, x_cb_s, y_cb_s, &
593# 491 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
594 & db_f77null, out%dims, 2, db_float, db_collinear, out%optlist, ierr)
595# 492 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
596 end if
597# 489 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
598 if (precision == 2) then
599# 490 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
600 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, x_cb, y_cb, &
601# 491 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
602 & db_f77null, out%dims, 2, db_double, db_collinear, out%optlist, ierr)
603# 492 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
604 end if
605# 494 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
606 err = dbfreeoptlist(out%optlist)
607 else
608 err = dbmkoptlist(2, out%optlist)
609 err = dbaddiaopt(out%optlist, dbopt_lo_offset, size(out%lo_offset), out%lo_offset)
610 err = dbaddiaopt(out%optlist, dbopt_hi_offset, size(out%hi_offset), out%hi_offset)
611# 500 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
612 if (precision == 1) then
613# 501 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
614 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, x_cb_s, db_f77null, &
615# 502 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
616 & db_f77null, out%dims, 1, db_float, db_collinear, out%optlist, ierr)
617# 503 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
618 end if
619# 500 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
620 if (precision == 2) then
621# 501 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
622 err = dbputqm(out%dbfile, 'rectilinear_grid', 16, 'x', 1, 'y', 1, 'z', 1, x_cb, db_f77null, &
623# 502 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
624 & db_f77null, out%dims, 1, db_double, db_collinear, out%optlist, ierr)
625# 503 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
626 end if
627# 505 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
628 err = dbfreeoptlist(out%optlist)
629 end if
630 else if (format == format_binary) then
631 ! Multidimensional local grid data is written to the formatted database slave file. Recall that no master file to
632 ! maintained in multidimensions.
633 if (p > 0) then
634 if (precision == precision_single) then
635 write (out%dbfile) real(x_cb, sp), real(y_cb, sp), real(z_cb, sp)
636 else
637 if (output_partial_domain) then
638 write (out%dbfile) x_cb(x_output_idx%beg - 1:x_output_idx%end), &
639 & y_cb(y_output_idx%beg - 1:y_output_idx%end), z_cb(z_output_idx%beg - 1:z_output_idx%end)
640 else
641 write (out%dbfile) x_cb, y_cb, z_cb
642 end if
643 end if
644 else if (n > 0) then
645 if (precision == precision_single) then
646 write (out%dbfile) real(x_cb, sp), real(y_cb, sp)
647 else
648 if (output_partial_domain) then
649 write (out%dbfile) x_cb(x_output_idx%beg - 1:x_output_idx%end), y_cb(y_output_idx%beg - 1:y_output_idx%end)
650 else
651 write (out%dbfile) x_cb, y_cb
652 end if
653 end if
654
655 ! One-dimensional local grid data is written to the formatted database slave file. In addition, the local grid data
656 ! is put together by the root process and written to the master file.
657 else
658 if (precision == precision_single) then
659 write (out%dbfile) real(x_cb, sp)
660 else
661 if (output_partial_domain) then
662 write (out%dbfile) x_cb(x_output_idx%beg - 1:x_output_idx%end)
663 else
664 write (out%dbfile) x_cb
665 end if
666 end if
667
668 if (num_procs > 1) then
669 call s_mpi_defragment_1d_grid_variable()
670 else
671 x_root_cb(:) = x_cb(:)
672 end if
673
674 if (proc_rank == 0) then
675 if (precision == precision_single) then
676 write (out%dbroot) real(x_root_cb, wp)
677 else
678 if (output_partial_domain) then
679 write (out%dbroot) x_root_cb(x_output_idx%beg - 1:x_output_idx%end)
680 else
681 write (out%dbroot) x_root_cb
682 end if
683 end if
684 end if
685 end if
686 end if
687
688 end subroutine s_write_grid_to_formatted_database_file
689
690 !> Write a single flow variable field to the formatted database slave and master files for a given time step.
691 impure subroutine s_write_variable_to_formatted_database_file(varname, t_step)
692
693 character(LEN=*), intent(in) :: varname
694 integer, intent(in) :: t_step
695
696 ! NAG compiler requires these to be statically sized
697 character(LEN=4*name_len), dimension(num_procs) :: varnames
698 integer, dimension(num_procs) :: vartypes
699 integer :: i, j, k
700 integer :: ierr
701 integer :: extents_size
702
703 if (format == format_silo) then
704 ! Determining the extents of the flow variable on each local process and gathering all this information on root process
705 if (num_procs > 1) then
706 call s_mpi_gather_data_extents(out%q_sf, out%data_extents)
707 else
708 out%data_extents(:,0) = (/minval(out%q_sf), maxval(out%q_sf)/)
709 end if
710
711 if (proc_rank == 0) then
712 do i = 1, num_procs
713 write (varnames(i), '(A,I0,A,I0,A)') '../p', i - 1, '/', t_step, '.silo:' // trim(varname)
714 end do
715
716 vartypes = db_quadvar
717
718 err = dbset2dstrlen(len(varnames(1)))
719 err = dbmkoptlist(2, out%optlist)
720 extents_size = size(out%data_extents, 1)
721 err = dbaddiopt(out%optlist, dbopt_extents_size, extents_size)
722 err = dbadddopt(out%optlist, dbopt_extents, out%data_extents)
723 err = dbputmvar(out%dbroot, trim(varname), len_trim(varname), num_procs, varnames, len_trim(varnames), vartypes, &
724 & out%optlist, ierr)
725 err = dbfreeoptlist(out%optlist)
726 end if
727
728 if (wp == dp) then
729 if (precision == precision_single) then
730 do i = -offset_x%beg, m + offset_x%end
731 do j = -offset_y%beg, n + offset_y%end
732 do k = -offset_z%beg, p + offset_z%end
733 out%q_sf_s(i, j, k) = real(out%q_sf(i, j, k), sp)
734 end do
735 end do
736 end do
737 if (grid_geometry == 3) then
738 do i = -offset_x%beg, m + offset_x%end
739 do j = -offset_y%beg, n + offset_y%end
740 do k = -offset_z%beg, p + offset_z%end
741 out%cyl_q_sf_s(j, k, i) = out%q_sf_s(i, j, k)
742 end do
743 end do
744 end do
745 end if
746 else
747 if (grid_geometry == 3) then
748 do i = -offset_x%beg, m + offset_x%end
749 do j = -offset_y%beg, n + offset_y%end
750 do k = -offset_z%beg, p + offset_z%end
751 out%cyl_q_sf(j, k, i) = out%q_sf(i, j, k)
752 end do
753 end do
754 end do
755 end if
756 end if
757 else if (wp == sp) then
758 do i = -offset_x%beg, m + offset_x%end
759 do j = -offset_y%beg, n + offset_y%end
760 do k = -offset_z%beg, p + offset_z%end
761 out%q_sf_s(i, j, k) = out%q_sf(i, j, k)
762 end do
763 end do
764 end do
765 if (grid_geometry == 3) then
766 do i = -offset_x%beg, m + offset_x%end
767 do j = -offset_y%beg, n + offset_y%end
768 do k = -offset_z%beg, p + offset_z%end
769 out%cyl_q_sf_s(j, k, i) = out%q_sf_s(i, j, k)
770 end do
771 end do
772 end do
773 end if
774 end if
775
776# 654 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
777 if (precision == 1) then
778# 655 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
779 if (p > 0) then
780 if (grid_geometry == 3) then
781 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, &
782 & out%cyl_q_sf_s, out%dims - 1, 3, db_f77null, 0, db_float, db_zonecent, &
783# 659 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
784 & db_f77null, ierr)
785 else
786 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, out%q_sf_s, &
787# 662 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
788 & out%dims - 1, 3, db_f77null, 0, db_float, db_zonecent, db_f77null, ierr)
789# 663 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
790 end if
791 else if (n > 0) then
792 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, out%q_sf_s, &
793# 666 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
794 & out%dims - 1, 2, db_f77null, 0, db_float, db_zonecent, db_f77null, ierr)
795# 667 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
796 else
797 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, out%q_sf_s, &
798# 669 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
799 & out%dims - 1, 1, db_f77null, 0, db_float, db_zonecent, db_f77null, ierr)
800# 670 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
801 end if
802 end if
803# 654 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
804 if (precision == 2) then
805# 655 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
806 if (p > 0) then
807 if (grid_geometry == 3) then
808 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, &
809 & out%cyl_q_sf, out%dims - 1, 3, db_f77null, 0, db_double, db_zonecent, &
810# 659 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
811 & db_f77null, ierr)
812 else
813 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, out%q_sf, &
814# 662 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
815 & out%dims - 1, 3, db_f77null, 0, db_double, db_zonecent, db_f77null, ierr)
816# 663 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
817 end if
818 else if (n > 0) then
819 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, out%q_sf, &
820# 666 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
821 & out%dims - 1, 2, db_f77null, 0, db_double, db_zonecent, db_f77null, ierr)
822# 667 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
823 else
824 err = dbputqv1(out%dbfile, trim(varname), len_trim(varname), 'rectilinear_grid', 16, out%q_sf, &
825# 669 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
826 & out%dims - 1, 1, db_f77null, 0, db_double, db_zonecent, db_f77null, ierr)
827# 670 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
828 end if
829 end if
830# 673 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
831 else
832 ! Writing the name of the flow variable and its data, associated with the local processor, to the formatted database
833 ! slave file
834 if (precision == precision_single) then
835 write (out%dbfile) varname, real(out%q_sf, wp)
836 else
837 write (out%dbfile) varname, out%q_sf
838 end if
839
840 ! In 1D, the root process also takes care of gathering the flow variable data from all of the local processor(s) and
841 ! writes it to the formatted database master file.
842 if (n == 0) then
843 if (num_procs > 1) then
844 call s_mpi_defragment_1d_flow_variable(out%q_sf, out%q_root_sf)
845 else
846 out%q_root_sf(:,:,:) = out%q_sf(:,:,:)
847 end if
848
849 if (proc_rank == 0) then
850 if (precision == precision_single) then
851 write (out%dbroot) varname, real(out%q_root_sf, wp)
852 else
853 write (out%dbroot) varname, out%q_root_sf
854 end if
855 end if
856 end if
857 end if
858
859 end subroutine s_write_variable_to_formatted_database_file
860
861 !> Write the post-processed results in the folder 'lag_bubbles_data'
862 impure subroutine s_write_lag_bubbles_results_to_text(t_step)
863
864 integer, intent(in) :: t_step
865 character(len=len_trim(case_dir) + 3*name_len) :: file_loc
866 integer :: id
867
868#ifdef MFC_MPI
869 real(wp), dimension(20) :: inputvals
870 real(wp) :: time_real
871 integer, dimension(MPI_STATUS_SIZE) :: status
872 integer(KIND=MPI_OFFSET_KIND) :: disp
873 integer :: view
874 logical :: file_exist
875 integer, dimension(2) :: gsizes, lsizes, start_idx_part
876 integer :: ifile
877 integer :: ierr
878 real(wp) :: file_time, file_dt
879 integer :: file_num_procs, file_tot_part
880 integer :: i
881 integer, dimension(:), allocatable :: proc_bubble_counts
882 real(wp), dimension(1:1,1:lag_io_vars) :: lag_io_null
883
884 lag_io_null = 0._wp
885
886 ! Construct file path
887 write (file_loc, '(A,I0,A)') 'lag_bubbles_', t_step, '.dat'
888 file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // trim(file_loc)
889
890 ! Check if file exists
891 inquire (file=trim(file_loc), exist=file_exist)
892 if (.not. file_exist) then
893 call s_mpi_abort('Restart file ' // trim(file_loc) // ' does not exist!')
894 end if
895
896 if (.not. parallel_io) return
897
898 if (proc_rank == 0) then
899 call mpi_file_open(mpi_comm_self, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
900 call s_check_mpi_file_open(ierr, file_loc)
901
902 call mpi_file_read(ifile, file_tot_part, 1, mpi_integer, status, ierr)
903 call mpi_file_read(ifile, file_time, 1, mpi_p, status, ierr)
904 call mpi_file_read(ifile, file_dt, 1, mpi_p, status, ierr)
905 call mpi_file_read(ifile, file_num_procs, 1, mpi_integer, status, ierr)
906
907 call mpi_file_close(ifile, ierr)
908 end if
909
910 call mpi_bcast(file_tot_part, 1, mpi_integer, 0, mpi_comm_world, ierr)
911 call mpi_bcast(file_time, 1, mpi_p, 0, mpi_comm_world, ierr)
912 call mpi_bcast(file_dt, 1, mpi_p, 0, mpi_comm_world, ierr)
913 call mpi_bcast(file_num_procs, 1, mpi_integer, 0, mpi_comm_world, ierr)
914 time_real = file_time
915
916 allocate (proc_bubble_counts(file_num_procs))
917
918 if (proc_rank == 0) then
919 call mpi_file_open(mpi_comm_self, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
920 call s_check_mpi_file_open(ierr, file_loc)
921
922 ! Skip to processor counts position
923 disp = int(sizeof(file_tot_part) + 2*sizeof(file_time) + sizeof(file_num_procs), mpi_offset_kind)
924 call mpi_file_seek(ifile, disp, mpi_seek_set, ierr)
925 call mpi_file_read(ifile, proc_bubble_counts, file_num_procs, mpi_integer, status, ierr)
926
927 call mpi_file_close(ifile, ierr)
928 end if
929
930 call mpi_bcast(proc_bubble_counts, file_num_procs, mpi_integer, 0, mpi_comm_world, ierr)
931
932 if (file_tot_part > 0) then
933 gsizes(1) = file_tot_part
934 gsizes(2) = lag_io_vars
935 lsizes(1) = file_tot_part
936 lsizes(2) = lag_io_vars
937 start_idx_part(1) = 0
938 start_idx_part(2) = 0
939
940 call mpi_type_create_subarray(2, gsizes, lsizes, start_idx_part, mpi_order_fortran, mpi_p, view, ierr)
941 call mpi_type_commit(view, ierr)
942
943 call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
944 call s_check_mpi_file_open(ierr, file_loc)
945
946 disp = int(sizeof(file_tot_part) + 2*sizeof(file_time) + sizeof(file_num_procs) &
947 & + file_num_procs*sizeof(proc_bubble_counts(1)), mpi_offset_kind)
948 call mpi_file_set_view(ifile, disp, mpi_p, view, 'native', mpi_info_null, ierr)
949
950 allocate (mpi_io_data_lg_bubbles(file_tot_part,1:lag_io_vars))
951
952 call mpi_file_read_all(ifile, mpi_io_data_lg_bubbles, lag_io_vars*file_tot_part, mpi_p, status, ierr)
953
954 write (file_loc, '(A,I0,A)') 'lag_bubbles_post_process_', t_step, '.dat'
955 file_loc = trim(case_dir) // '/lag_bubbles_post_process/' // trim(file_loc)
956
957 if (proc_rank == 0) then
958 open (unit=29, file=file_loc, form='formatted', position='rewind')
959
960 if (lag_header) then
961 write (29, '(A)', advance='no')
962 if (lag_id_wrt) write (29, '(A8)', advance='no') 'id, '
963 if (lag_pos_wrt) write (29, '(3(A17))', advance='no') 'px, ', 'py, ', 'pz, '
964 if (lag_pos_prev_wrt) write (29, '(3(A17))', advance='no') 'pvx, ', 'pvy, ', 'pvz, '
965 if (lag_vel_wrt) write (29, '(3(A17))', advance='no') 'vx, ', 'vy, ', 'vz, '
966 if (lag_rad_wrt) write (29, '(A17)', advance='no') 'radius, '
967 if (lag_rvel_wrt) write (29, '(A17)', advance='no') 'rvel, '
968 if (lag_r0_wrt) write (29, '(A17)', advance='no') 'r0, '
969 if (lag_rmax_wrt) write (29, '(A17)', advance='no') 'rmax, '
970 if (lag_rmin_wrt) write (29, '(A17)', advance='no') 'rmin, '
971 if (lag_dphidt_wrt) write (29, '(A17)', advance='no') 'dphidt, '
972 if (lag_pres_wrt) write (29, '(A17)', advance='no') 'pressure, '
973 if (lag_mv_wrt) write (29, '(A17)', advance='no') 'mv, '
974 if (lag_mg_wrt) write (29, '(A17)', advance='no') 'mg, '
975 if (lag_betat_wrt) write (29, '(A17)', advance='no') 'betaT, '
976 if (lag_betac_wrt) write (29, '(A17)', advance='no') 'betaC, '
977 write (29, '(A15)') 'time'
978 end if
979
980 do i = 1, file_tot_part
981 id = int(mpi_io_data_lg_bubbles(i, 1))
982 inputvals(1:20) = mpi_io_data_lg_bubbles(i,2:21)
983 if (id > 0) then
984 write (29, '(100(A))', advance='no') ''
985 if (lag_id_wrt) write (29, '(I6, A)', advance='no') id, ', '
986 if (lag_pos_wrt) write (29, '(3(E15.7, A))', advance='no') inputvals(1), ', ', inputvals(2), ', ', &
987 & inputvals(3), ', '
988 if (lag_pos_prev_wrt) write (29, '(3(E15.7, A))', advance='no') inputvals(4), ', ', inputvals(5), ', ', &
989 & inputvals(6), ', '
990 if (lag_vel_wrt) write (29, '(3(E15.7, A))', advance='no') inputvals(7), ', ', inputvals(8), ', ', &
991 & inputvals(9), ', '
992 if (lag_rad_wrt) write (29, '(E15.7, A)', advance='no') inputvals(10), ', '
993 if (lag_rvel_wrt) write (29, '(E15.7, A)', advance='no') inputvals(11), ', '
994 if (lag_r0_wrt) write (29, '(E15.7, A)', advance='no') inputvals(12), ', '
995 if (lag_rmax_wrt) write (29, '(E15.7, A)', advance='no') inputvals(13), ', '
996 if (lag_rmin_wrt) write (29, '(E15.7, A)', advance='no') inputvals(14), ', '
997 if (lag_dphidt_wrt) write (29, '(E15.7, A)', advance='no') inputvals(15), ', '
998 if (lag_pres_wrt) write (29, '(E15.7, A)', advance='no') inputvals(16), ', '
999 if (lag_mv_wrt) write (29, '(E15.7, A)', advance='no') inputvals(17), ', '
1000 if (lag_mg_wrt) write (29, '(E15.7, A)', advance='no') inputvals(18), ', '
1001 if (lag_betat_wrt) write (29, '(E15.7, A)', advance='no') inputvals(19), ', '
1002 if (lag_betac_wrt) write (29, '(E15.7, A)', advance='no') inputvals(20), ', '
1003 write (29, '(E15.7)') time_real
1004 end if
1005 end do
1006 close (29)
1007 end if
1008
1009 deallocate (mpi_io_data_lg_bubbles)
1010
1011 call s_mpi_barrier()
1012
1013 call mpi_file_close(ifile, ierr)
1014 else
1015 call mpi_type_contiguous(0, mpi_p, view, ierr)
1016 call mpi_type_commit(view, ierr)
1017
1018 call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
1019 call s_check_mpi_file_open(ierr, file_loc)
1020
1021 disp = int(sizeof(file_tot_part) + 2*sizeof(file_time) + sizeof(file_num_procs) &
1022 & + file_num_procs*sizeof(proc_bubble_counts(1)), mpi_offset_kind)
1023 call mpi_file_set_view(ifile, disp, mpi_p, view, 'native', mpi_info_null, ierr)
1024
1025 call mpi_file_read_all(ifile, lag_io_null, 0, mpi_p, status, ierr)
1026
1027 call mpi_file_close(ifile, ierr)
1028 call mpi_type_free(view, ierr)
1029 end if
1030#endif
1031
1032 end subroutine s_write_lag_bubbles_results_to_text
1033
1034 !> Read Lagrangian bubble restart data and write bubble positions and scalar fields to the Silo database.
1035 impure subroutine s_write_lag_bubbles_to_formatted_database_file(t_step)
1036
1037 integer, intent(in) :: t_step
1038 character(len=len_trim(case_dir) + 3*name_len) :: file_loc
1039 integer :: id
1040
1041#ifdef MFC_MPI
1042 real(wp) :: time_real
1043 integer, dimension(MPI_STATUS_SIZE) :: status
1044 integer(KIND=MPI_OFFSET_KIND) :: disp
1045 integer :: view
1046 logical :: file_exist
1047 integer, dimension(2) :: gsizes, lsizes, start_idx_part
1048 integer :: ifile, ierr, nbub
1049 real(wp) :: file_time, file_dt
1050 integer :: file_num_procs, file_tot_part
1051 integer, dimension(:), allocatable :: proc_bubble_counts
1052 real(wp), dimension(1:1,1:lag_io_vars) :: dummy
1053 character(LEN=4*name_len), dimension(num_procs) :: meshnames
1054 integer, dimension(num_procs) :: meshtypes
1055 real(wp) :: dummy_data
1056 integer :: i
1057 real(wp), dimension(:), allocatable :: bub_id
1058 real(wp), dimension(:), allocatable :: px, py, pz, ppx, ppy, ppz, vx, vy, vz
1059 real(wp), dimension(:), allocatable :: radius, rvel, rnot, rmax, rmin, dphidt
1060 real(wp), dimension(:), allocatable :: pressure, mv, mg, betat, betac
1061
1062 dummy = 0._wp
1063 dummy_data = 0._wp
1064
1065 ! Construct file path
1066 write (file_loc, '(A,I0,A)') 'lag_bubbles_', t_step, '.dat'
1067 file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // trim(file_loc)
1068
1069 ! Check if file exists
1070 inquire (file=trim(file_loc), exist=file_exist)
1071 if (.not. file_exist) then
1072 call s_mpi_abort('Restart file ' // trim(file_loc) // ' does not exist!')
1073 end if
1074
1075 if (.not. parallel_io) return
1076
1077 if (proc_rank == 0) then
1078 call mpi_file_open(mpi_comm_self, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
1079 call s_check_mpi_file_open(ierr, file_loc)
1080
1081 call mpi_file_read(ifile, file_tot_part, 1, mpi_integer, status, ierr)
1082 call mpi_file_read(ifile, file_time, 1, mpi_p, status, ierr)
1083 call mpi_file_read(ifile, file_dt, 1, mpi_p, status, ierr)
1084 call mpi_file_read(ifile, file_num_procs, 1, mpi_integer, status, ierr)
1085
1086 call mpi_file_close(ifile, ierr)
1087 end if
1088
1089 call mpi_bcast(file_tot_part, 1, mpi_integer, 0, mpi_comm_world, ierr)
1090 call mpi_bcast(file_time, 1, mpi_p, 0, mpi_comm_world, ierr)
1091 call mpi_bcast(file_dt, 1, mpi_p, 0, mpi_comm_world, ierr)
1092 call mpi_bcast(file_num_procs, 1, mpi_integer, 0, mpi_comm_world, ierr)
1093 time_real = file_time
1094
1095 allocate (proc_bubble_counts(file_num_procs))
1096
1097 if (proc_rank == 0) then
1098 call mpi_file_open(mpi_comm_self, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
1099 call s_check_mpi_file_open(ierr, file_loc)
1100
1101 ! Skip to processor counts position
1102 disp = int(sizeof(file_tot_part) + 2*sizeof(file_time) + sizeof(file_num_procs), mpi_offset_kind)
1103 call mpi_file_seek(ifile, disp, mpi_seek_set, ierr)
1104 call mpi_file_read(ifile, proc_bubble_counts, file_num_procs, mpi_integer, status, ierr)
1105
1106 call mpi_file_close(ifile, ierr)
1107 end if
1108
1109 call mpi_bcast(proc_bubble_counts, file_num_procs, mpi_integer, 0, mpi_comm_world, ierr)
1110
1111 ! Set time variables from file
1112
1113 nbub = proc_bubble_counts(proc_rank + 1)
1114
1115 start_idx_part(1) = 0
1116 do i = 1, proc_rank
1117 start_idx_part(1) = start_idx_part(1) + proc_bubble_counts(i)
1118 end do
1119
1120 start_idx_part(2) = 0
1121 lsizes(1) = nbub
1122 lsizes(2) = lag_io_vars
1123
1124 gsizes(1) = file_tot_part
1125 gsizes(2) = lag_io_vars
1126
1127 if (nbub > 0) then
1128# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1129 allocate (bub_id(nbub))
1130# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1131# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1132 allocate (px(nbub))
1133# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1134# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1135 allocate (py(nbub))
1136# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1137# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1138 allocate (pz(nbub))
1139# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1140# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1141 allocate (ppx(nbub))
1142# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1143# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1144 allocate (ppy(nbub))
1145# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1146# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1147 allocate (ppz(nbub))
1148# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1149# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1150 allocate (vx(nbub))
1151# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1152# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1153 allocate (vy(nbub))
1154# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1155# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1156 allocate (vz(nbub))
1157# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1158# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1159 allocate (radius(nbub))
1160# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1161# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1162 allocate (rvel(nbub))
1163# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1164# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1165 allocate (rnot(nbub))
1166# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1167# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1168 allocate (rmax(nbub))
1169# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1170# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1171 allocate (rmin(nbub))
1172# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1173# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1174 allocate (dphidt(nbub))
1175# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1176# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1177 allocate (pressure(nbub))
1178# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1179# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1180 allocate (mv(nbub))
1181# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1182# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1183 allocate (mg(nbub))
1184# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1185# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1186 allocate (betat(nbub))
1187# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1188# 973 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1189 allocate (betac(nbub))
1190# 974 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1191# 975 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1192 allocate (mpi_io_data_lg_bubbles(nbub,1:lag_io_vars))
1193
1194 call mpi_type_create_subarray(2, gsizes, lsizes, start_idx_part, mpi_order_fortran, mpi_p, view, ierr)
1195 call mpi_type_commit(view, ierr)
1196
1197 call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
1198 call s_check_mpi_file_open(ierr, file_loc)
1199
1200 ! Skip extended header
1201 disp = int(sizeof(file_tot_part) + 2*sizeof(file_time) + sizeof(file_num_procs) &
1202 & + file_num_procs*sizeof(proc_bubble_counts(1)), mpi_offset_kind)
1203 call mpi_file_set_view(ifile, disp, mpi_p, view, 'native', mpi_info_int, ierr)
1204
1205 call mpi_file_read_all(ifile, mpi_io_data_lg_bubbles, lag_io_vars*nbub, mpi_p, status, ierr)
1206
1207 call mpi_file_close(ifile, ierr)
1208 call mpi_type_free(view, ierr)
1209
1210 ! Extract data from MPI_IO_DATA_lg_bubbles array Adjust these indices based on your actual data layout
1211# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1212 bub_id(:) = mpi_io_data_lg_bubbles(:,1)
1213# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1214# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1215 px(:) = mpi_io_data_lg_bubbles(:,2)
1216# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1217# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1218 py(:) = mpi_io_data_lg_bubbles(:,3)
1219# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1220# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1221 pz(:) = mpi_io_data_lg_bubbles(:,4)
1222# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1223# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1224 ppx(:) = mpi_io_data_lg_bubbles(:,5)
1225# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1226# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1227 ppy(:) = mpi_io_data_lg_bubbles(:,6)
1228# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1229# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1230 ppz(:) = mpi_io_data_lg_bubbles(:,7)
1231# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1232# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1233 vx(:) = mpi_io_data_lg_bubbles(:,8)
1234# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1235# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1236 vy(:) = mpi_io_data_lg_bubbles(:,9)
1237# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1238# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1239 vz(:) = mpi_io_data_lg_bubbles(:,10)
1240# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1241# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1242 radius(:) = mpi_io_data_lg_bubbles(:,11)
1243# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1244# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1245 rvel(:) = mpi_io_data_lg_bubbles(:,12)
1246# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1247# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1248 rnot(:) = mpi_io_data_lg_bubbles(:,13)
1249# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1250# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1251 rmax(:) = mpi_io_data_lg_bubbles(:,14)
1252# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1253# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1254 rmin(:) = mpi_io_data_lg_bubbles(:,15)
1255# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1256# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1257 dphidt(:) = mpi_io_data_lg_bubbles(:,16)
1258# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1259# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1260 pressure(:) = mpi_io_data_lg_bubbles(:,17)
1261# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1262# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1263 mv(:) = mpi_io_data_lg_bubbles(:,18)
1264# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1265# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1266 mg(:) = mpi_io_data_lg_bubbles(:,19)
1267# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1268# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1269 betat(:) = mpi_io_data_lg_bubbles(:,20)
1270# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1271# 998 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1272 betac(:) = mpi_io_data_lg_bubbles(:,21)
1273# 999 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1274# 1000 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1275
1276 ! Next, the root processor proceeds to record all of the spatial extents in the formatted database master file. In
1277 ! addition, it also records a sub-domain connectivity map so that the entire grid may be reassembled by looking at the
1278 ! master file.
1279 if (proc_rank == 0) then
1280 do i = 1, num_procs
1281 write (meshnames(i), '(A,I0,A,I0,A)') '../p', i - 1, '/', t_step, '.silo:lag_bubbles'
1282 meshtypes(i) = db_pointmesh
1283 end do
1284 err = dbset2dstrlen(len(meshnames(1)))
1285 err = dbputmmesh(out%dbroot, 'lag_bubbles', 16, num_procs, meshnames, len_trim(meshnames), meshtypes, db_f77null, &
1286 & ierr)
1287 end if
1288
1289 err = dbputpm(out%dbfile, 'lag_bubbles', 11, 3, px, py, pz, nbub, db_double, db_f77null, ierr)
1290
1291 if (lag_id_wrt) call s_write_lag_variable_to_formatted_database_file('part_id', t_step, bub_id, nbub)
1292 if (lag_vel_wrt) then
1293 call s_write_lag_variable_to_formatted_database_file('part_vel1', t_step, vx, nbub)
1294 call s_write_lag_variable_to_formatted_database_file('part_vel2', t_step, vy, nbub)
1295 if (p > 0) call s_write_lag_variable_to_formatted_database_file('part_vel3', t_step, vz, nbub)
1296 end if
1297 if (lag_rad_wrt) call s_write_lag_variable_to_formatted_database_file('part_radius', t_step, radius, nbub)
1298 if (lag_rvel_wrt) call s_write_lag_variable_to_formatted_database_file('part_rdot', t_step, rvel, nbub)
1299 if (lag_r0_wrt) call s_write_lag_variable_to_formatted_database_file('part_r0', t_step, rnot, nbub)
1300 if (lag_rmax_wrt) call s_write_lag_variable_to_formatted_database_file('part_rmax', t_step, rmax, nbub)
1301 if (lag_rmin_wrt) call s_write_lag_variable_to_formatted_database_file('part_rmin', t_step, rmin, nbub)
1302 if (lag_dphidt_wrt) call s_write_lag_variable_to_formatted_database_file('part_dphidt', t_step, dphidt, nbub)
1303 if (lag_pres_wrt) call s_write_lag_variable_to_formatted_database_file('part_pressure', t_step, pressure, nbub)
1304 if (lag_mv_wrt) call s_write_lag_variable_to_formatted_database_file('part_mv', t_step, mv, nbub)
1305 if (lag_mg_wrt) call s_write_lag_variable_to_formatted_database_file('part_mg', t_step, mg, nbub)
1306 if (lag_betat_wrt) call s_write_lag_variable_to_formatted_database_file('part_betaT', t_step, betat, nbub)
1307 if (lag_betac_wrt) call s_write_lag_variable_to_formatted_database_file('part_betaC', t_step, betac, nbub)
1308
1309 deallocate (bub_id, px, py, pz, ppx, ppy, ppz, vx, vy, vz, radius, rvel, rnot, rmax, rmin, dphidt, pressure, mv, mg, &
1310 & betat, betac)
1311 deallocate (mpi_io_data_lg_bubbles)
1312 else
1313 call mpi_type_contiguous(0, mpi_p, view, ierr)
1314 call mpi_type_commit(view, ierr)
1315
1316 call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
1317 call s_check_mpi_file_open(ierr, file_loc)
1318
1319 ! Skip extended header
1320 disp = int(sizeof(file_tot_part) + 2*sizeof(file_time) + sizeof(file_num_procs) &
1321 & + file_num_procs*sizeof(proc_bubble_counts(1)), mpi_offset_kind)
1322 call mpi_file_set_view(ifile, disp, mpi_p, view, 'native', mpi_info_int, ierr)
1323
1324 call mpi_file_read_all(ifile, dummy, 0, mpi_p, status, ierr)
1325
1326 call mpi_file_close(ifile, ierr)
1327 call mpi_type_free(view, ierr)
1328
1329 if (proc_rank == 0) then
1330 do i = 1, num_procs
1331 write (meshnames(i), '(A,I0,A,I0,A)') '../p', i - 1, '/', t_step, '.silo:lag_bubbles'
1332 meshtypes(i) = db_pointmesh
1333 end do
1334 err = dbset2dstrlen(len(meshnames(1)))
1335 err = dbputmmesh(out%dbroot, 'lag_bubbles', 16, num_procs, meshnames, len_trim(meshnames), meshtypes, db_f77null, &
1336 & ierr)
1337 end if
1338
1339 err = dbsetemptyok(1)
1340 err = dbputpm(out%dbfile, 'lag_bubbles', 11, 3, dummy_data, dummy_data, dummy_data, 0, db_double, db_f77null, ierr)
1341
1342 if (lag_id_wrt) call s_write_lag_variable_to_formatted_database_file('part_id', t_step)
1343 if (lag_vel_wrt) then
1344 call s_write_lag_variable_to_formatted_database_file('part_vel1', t_step)
1345 call s_write_lag_variable_to_formatted_database_file('part_vel2', t_step)
1346 if (p > 0) call s_write_lag_variable_to_formatted_database_file('part_vel3', t_step)
1347 end if
1348 if (lag_rad_wrt) call s_write_lag_variable_to_formatted_database_file('part_radius', t_step)
1349 if (lag_rvel_wrt) call s_write_lag_variable_to_formatted_database_file('part_rdot', t_step)
1350 if (lag_r0_wrt) call s_write_lag_variable_to_formatted_database_file('part_r0', t_step)
1351 if (lag_rmax_wrt) call s_write_lag_variable_to_formatted_database_file('part_rmax', t_step)
1352 if (lag_rmin_wrt) call s_write_lag_variable_to_formatted_database_file('part_rmin', t_step)
1353 if (lag_dphidt_wrt) call s_write_lag_variable_to_formatted_database_file('part_dphidt', t_step)
1354 if (lag_pres_wrt) call s_write_lag_variable_to_formatted_database_file('part_pressure', t_step)
1355 if (lag_mv_wrt) call s_write_lag_variable_to_formatted_database_file('part_mv', t_step)
1356 if (lag_mg_wrt) call s_write_lag_variable_to_formatted_database_file('part_mg', t_step)
1357 if (lag_betat_wrt) call s_write_lag_variable_to_formatted_database_file('part_betaT', t_step)
1358 if (lag_betac_wrt) call s_write_lag_variable_to_formatted_database_file('part_betaC', t_step)
1359 end if
1360#endif
1361
1362 end subroutine s_write_lag_bubbles_to_formatted_database_file
1363
1364 !> Write a single Lagrangian bubble point-variable to the Silo database slave and master files.
1365 subroutine s_write_lag_variable_to_formatted_database_file(varname, t_step, data, nBubs)
1366
1367 character(len=*), intent(in) :: varname
1368 integer, intent(in) :: t_step
1369 real(wp), dimension(1:), intent(in), optional :: data
1370 integer, intent(in), optional :: nBubs
1371 character(len=64), dimension(num_procs) :: var_names
1372 integer, dimension(num_procs) :: var_types
1373 real(wp) :: dummy_data
1374 integer :: ierr
1375 integer :: i
1376
1377 dummy_data = 0._wp
1378
1379 if (present(nbubs) .and. present(data)) then
1380 if (proc_rank == 0) then
1381 do i = 1, num_procs
1382 write (var_names(i), '(A,I0,A,I0,A)') '../p', i - 1, '/', t_step, '.silo:' // trim(varname)
1383 var_types(i) = db_pointvar
1384 end do
1385 err = dbset2dstrlen(len(var_names(1)))
1386 err = dbputmvar(out%dbroot, trim(varname), len_trim(varname), num_procs, var_names, len_trim(var_names), &
1387 & var_types, db_f77null, ierr)
1388 end if
1389
1390 err = dbputpv1(out%dbfile, trim(varname), len_trim(varname), 'lag_bubbles', 11, data, nbubs, db_double, db_f77null, &
1391 & ierr)
1392 else
1393 if (proc_rank == 0) then
1394 do i = 1, num_procs
1395 write (var_names(i), '(A,I0,A,I0,A)') '../p', i - 1, '/', t_step, '.silo:' // trim(varname)
1396 var_types(i) = db_pointvar
1397 end do
1398 err = dbset2dstrlen(len(var_names(1)))
1399 err = dbsetemptyok(1)
1400 err = dbputmvar(out%dbroot, trim(varname), len_trim(varname), num_procs, var_names, len_trim(var_names), &
1401 & var_types, db_f77null, ierr)
1402 end if
1403
1404 err = dbsetemptyok(1)
1405 err = dbputpv1(out%dbfile, trim(varname), len_trim(varname), 'lag_bubbles', 11, dummy_data, 0, db_double, db_f77null, &
1406 & ierr)
1407 end if
1408
1409 end subroutine s_write_lag_variable_to_formatted_database_file
1410
1411 !> Convert the binary immersed-boundary state file to per-body formatted text files
1412 impure subroutine s_write_ib_state_files()
1413
1414 character(len=len_trim(case_dir) + 4*name_len) :: in_file, out_file, file_loc
1415 integer :: iu_in, ios, i, rec_id
1416 integer, allocatable, dimension(:) :: iu_out
1417 real(wp) :: rec_time
1418 real(wp), dimension(3) :: rec_force, rec_torque
1419 real(wp), dimension(3) :: rec_vel, rec_angular_vel
1420 real(wp), dimension(3) :: rec_angles, rec_centroid
1421
1422 file_loc = trim(case_dir) // '/D'
1423
1424 in_file = trim(file_loc) // '/ib_state.dat'
1425 open (newunit=iu_in, file=trim(in_file), form='unformatted', access='stream', status='old', action='read', iostat=ios)
1426 if (ios /= 0) then
1427 call s_mpi_abort('Cannot open IB state input file: ' // trim(in_file))
1428 end if
1429
1430 allocate (iu_out(num_ibs))
1431 do i = 1, num_ibs
1432 write (out_file, '(A,I0,A)') trim(file_loc) // '/ib_', i, '.txt'
1433 open (newunit=iu_out(i), file=trim(out_file), form='formatted', status='replace', action='write', iostat=ios)
1434 if (ios /= 0) then
1435 call s_mpi_abort('Cannot open IB state output file: ' // trim(out_file))
1436 end if
1437 write (iu_out(i), &
1438 & '(A)') 'mytime fx fy fz Tau_x Tau_y Tau_z vx vy vz omega_x omega_y omega_z angle_x angle_y angle_z x_c y_c z_c'
1439 end do
1440
1441 do
1442 read (iu_in, iostat=ios) rec_time, rec_id, rec_force, rec_torque, rec_vel, rec_angular_vel, rec_angles, &
1443 & rec_centroid(1), rec_centroid(2), rec_centroid(3)
1444 if (ios /= 0) exit
1445
1446 if (rec_id >= 1 .and. rec_id <= num_ibs) then
1447 write (iu_out(rec_id), '(19(ES24.16E3,1X))') rec_time, rec_force(1), rec_force(2), rec_force(3), rec_torque(1), &
1448 & rec_torque(2), rec_torque(3), rec_vel(1), rec_vel(2), rec_vel(3), rec_angular_vel(1), &
1449 & rec_angular_vel(2), rec_angular_vel(3), rec_angles(1), rec_angles(2), rec_angles(3), rec_centroid(1), &
1450 & rec_centroid(2), rec_centroid(3)
1451 end if
1452 end do
1453
1454 close (iu_in)
1455 do i = 1, num_ibs
1456 close (iu_out(i))
1457 end do
1458 deallocate (iu_out)
1459
1460 end subroutine s_write_ib_state_files
1461
1462 !> Extract the volume-fraction interface contour from primitive fields and write the coordinates to the interface data file.
1463 impure subroutine s_write_intf_data_file(q_prim_vf)
1464
1465 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf
1466 integer :: i, j, k, l, cent
1467 integer :: counter, root !< number of data points extracted to fit shape to SH perturbations
1468 real(wp), allocatable :: x_td(:), y_td(:), x_d1(:), y_d1(:), y_d(:), x_d(:)
1469 real(wp) :: axp, axm, ayp, aym, tgp, euc_d, thres, maxalph_loc, maxalph_glb
1470
1471 allocate (x_d1(m*n))
1472 allocate (y_d1(m*n))
1473 counter = 0
1474 maxalph_loc = 0._wp
1475 do k = 0, p
1476 do j = 0, n
1477 do i = 0, m
1478 if (q_prim_vf(eqn_idx%E + 2)%sf(i, j, k) > maxalph_loc) then
1479 maxalph_loc = q_prim_vf(eqn_idx%E + 2)%sf(i, j, k)
1480 end if
1481 end do
1482 end do
1483 end do
1484
1485 call s_mpi_allreduce_max(maxalph_loc, maxalph_glb)
1486 if (p > 0) then
1487 do l = 0, p
1488 if (z_cc(l) < dz(l) .and. z_cc(l) > 0) then
1489 cent = l
1490 end if
1491 end do
1492 else
1493 cent = 0
1494 end if
1495
1496 thres = 0.9_wp*maxalph_glb
1497 do k = 0, n
1498 do j = 0, m
1499 axp = q_prim_vf(eqn_idx%E + 2)%sf(j + 1, k, cent)
1500 axm = q_prim_vf(eqn_idx%E + 2)%sf(j, k, cent)
1501 ayp = q_prim_vf(eqn_idx%E + 2)%sf(j, k + 1, cent)
1502 aym = q_prim_vf(eqn_idx%E + 2)%sf(j, k, cent)
1503 if ((axp > thres .and. axm < thres) .or. (axp < thres .and. axm > thres) .or. (ayp > thres .and. aym < thres) &
1504 & .or. (ayp < thres .and. aym > thres)) then
1505 if (counter == 0) then
1506 counter = counter + 1
1507 x_d1(counter) = x_cc(j)
1508 y_d1(counter) = y_cc(k)
1509 else
1510 tgp = sqrt(dx(j)**2 + dy(k)**2)
1511 do i = 1, counter
1512 euc_d = sqrt((x_cc(j) - x_d1(i))**2 + (y_cc(k) - y_d1(i))**2)
1513 if (euc_d < tgp) then
1514 exit
1515 else if (i == counter) then
1516 counter = counter + 1
1517 x_d1(counter) = x_cc(j)
1518 y_d1(counter) = y_cc(k)
1519 end if
1520 end do
1521 end if
1522 end if
1523 end do
1524 end do
1525
1526 allocate (x_d(counter), y_d(counter))
1527
1528 do i = 1, counter
1529 y_d(i) = y_d1(i)
1530 x_d(i) = x_d1(i)
1531 end do
1532 root = 0
1533
1534 call s_mpi_gather_data(x_d, counter, x_td, root)
1535 call s_mpi_gather_data(y_d, counter, y_td, root)
1536 if (proc_rank == 0) then
1537 do i = 1, size(x_td)
1538 if (i == size(x_td)) then
1539 write (211, '(F12.9,1X,F12.9,1X,I4)') x_td(i), y_td(i), size(x_td)
1540 else
1541 write (211, '(F12.9,1X,F12.9,1X,F3.1)') x_td(i), y_td(i), 0._wp
1542 end if
1543 end do
1544 end if
1545
1546 end subroutine s_write_intf_data_file
1547
1548 !> Compute volume-integrated kinetic, potential, and internal energies and write the energy budget to the energy data file.
1549 impure subroutine s_write_energy_data_file(q_prim_vf, q_cons_vf)
1550
1551 type(scalar_field), dimension(sys_size), intent(in) :: q_prim_vf, q_cons_vf
1552 real(wp) :: elk, egk, elp, egint, eg_phase, vb, vl, pres_av, et
1553 real(wp) :: rho, pres, dv, tmp, gamma, pi_inf, qv, maxma, maxma_glb, maxvel, c, ma
1554 real(wp), dimension(num_vels) :: vel
1555 real(wp), dimension(num_fluids) :: adv, alpha_rho
1556 integer :: i, j, k, l, s !< looping indices
1557
1558 egk = 0._wp
1559 elp = 0._wp
1560 egint = 0._wp
1561 vb = 0._wp
1562 maxvel = 0._wp
1563 maxma = 0._wp
1564 vl = 0._wp
1565 elk = 0._wp
1566 et = 0._wp
1567 vb = 0._wp
1568 dv = 0._wp
1569 pres_av = 0._wp
1570 pres = 0._wp
1571 c = 0._wp
1572
1573 do k = 0, p
1574 do j = 0, n
1575 do i = 0, m
1576 dv = dx(i)*dy(j)*dz(k)
1577 pres = q_prim_vf(eqn_idx%E)%sf(i, j, k)
1578 do s = 1, num_vels
1579 vel(s) = q_prim_vf(num_fluids + s)%sf(i, j, k)
1580 egk = egk + 0.5_wp*q_prim_vf(eqn_idx%E + 2)%sf(i, j, k)*q_prim_vf(2)%sf(i, j, k)*vel(s)*vel(s)*dv
1581 elk = elk + 0.5_wp*q_prim_vf(eqn_idx%E + 1)%sf(i, j, k)*q_prim_vf(1)%sf(i, j, k)*vel(s)*vel(s)*dv
1582 if (abs(vel(s)) > maxvel) then
1583 maxvel = abs(vel(s))
1584 end if
1585 end do
1586 do l = 1, num_fluids
1587 adv(l) = q_prim_vf(eqn_idx%E + l)%sf(i, j, k)
1588 alpha_rho(l) = q_prim_vf(l)%sf(i, j, k)
1589 end do
1590
1591 call s_phase_internal_energy(pres, adv(2), alpha_rho(2), 2, eg_phase)
1592 egint = egint + eg_phase*dv
1593
1594 call s_compute_mixture_coefficients(alpha_rho, adv, rho, gamma, pi_inf, qv)
1595
1596 call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, adv, c, alpha_rho)
1597
1598 ma = maxvel/c
1599 if (ma > maxma .and. (adv(1) > (1.0_wp - 1.0e-10_wp))) then
1600 maxma = ma
1601 end if
1602 vl = vl + adv(1)*dv
1603 vb = vb + adv(2)*dv
1604 pres_av = pres_av + adv(1)*pres*dv
1605 et = et + q_cons_vf(eqn_idx%E)%sf(i, j, k)*dv
1606 end do
1607 end do
1608 end do
1609
1610 tmp = pres_av
1611 call s_mpi_allreduce_sum(tmp, pres_av)
1612 tmp = vl
1613 call s_mpi_allreduce_sum(tmp, vl)
1614
1615 call s_mpi_allreduce_max(maxma, maxma_glb)
1616 tmp = elk
1617 call s_mpi_allreduce_sum(tmp, elk)
1618 tmp = egint
1619 call s_mpi_allreduce_sum(tmp, egint)
1620 tmp = egk
1621 call s_mpi_allreduce_sum(tmp, egk)
1622 tmp = vb
1623 call s_mpi_allreduce_sum(tmp, vb)
1624 tmp = et
1625 call s_mpi_allreduce_sum(tmp, et)
1626
1627 elp = pres_av/vl*vb
1628 if (proc_rank == 0) then
1629 write (251, '(10X, 8F24.8)') elp, egint, elk, egk, et, vb, vl, maxma_glb
1630 end if
1631
1632 end subroutine s_write_energy_data_file
1633
1634 !> Read IB state and write a Silo point mesh with per-body scalar fields.
1635 impure subroutine s_write_ib_bodies_to_formatted_database_file(t_step)
1636
1637 integer, intent(in) :: t_step
1638 character(len=len_trim(case_dir) + 3*name_len) :: file_loc
1639
1640#ifdef MFC_MPI
1641 integer, parameter :: nfields_per_ib = 20
1642 real(wp) :: ib_buf(nfields_per_ib)
1643 real(wp), dimension(:,:), allocatable :: ib_data
1644 logical :: file_exist
1645 character(LEN=4*name_len), dimension(num_procs) :: meshnames
1646 integer, dimension(num_procs) :: meshtypes
1647 integer :: i, ios, file_unit
1648 integer :: ierr, nbodies
1649 integer :: r, nlocal, gbl_id
1650 character(len=10) :: t_step_string
1651 real(wp), dimension(:), allocatable :: px, py, pz
1652 real(wp), dimension(:), allocatable :: force_x, force_y, force_z
1653 real(wp), dimension(:), allocatable :: torque_x, torque_y, torque_z
1654 real(wp), dimension(:), allocatable :: vel_x, vel_y, vel_z
1655 real(wp), dimension(:), allocatable :: omega_x, omega_y, omega_z
1656 real(wp), dimension(:), allocatable :: angle_x, angle_y, angle_z
1657 real(wp), dimension(:), allocatable :: ib_diameter
1658 real(sp), dimension(:), allocatable :: px_s, py_s, pz_s
1659 logical, dimension(:), allocatable :: keep
1660 integer :: nkept
1661 real(wp) :: r_ib
1662
1663 if (proc_rank == 0) then
1664 nbodies = num_ibs
1665
1666 if (nbodies > 0) then
1667 allocate (ib_data(nbodies, nfields_per_ib))
1668 allocate (px(nbodies), py(nbodies), pz(nbodies))
1669 allocate (force_x(nbodies), force_y(nbodies), force_z(nbodies))
1670 allocate (torque_x(nbodies), torque_y(nbodies), torque_z(nbodies))
1671 allocate (vel_x(nbodies), vel_y(nbodies), vel_z(nbodies))
1672 allocate (omega_x(nbodies), omega_y(nbodies), omega_z(nbodies))
1673 allocate (angle_x(nbodies), angle_y(nbodies), angle_z(nbodies))
1674 allocate (ib_diameter(nbodies))
1675
1676 if (file_per_process) then
1677 call s_int_to_str(t_step, t_step_string)
1678 ib_data = 0._wp
1679 do r = 0, num_procs - 1
1680 write (file_loc, '(A,I0,A,i7.7,A)') 'ib_state_', t_step, '_', r, '.dat'
1681 file_loc = trim(case_dir) // '/restart_data/lustre_' // trim(t_step_string) // '/' // trim(file_loc)
1682
1683 inquire (file=trim(file_loc), exist=file_exist)
1684 if (.not. file_exist) call s_mpi_abort('Restart file ' // trim(file_loc) // ' does not exist!')
1685
1686 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='old', iostat=ios)
1687 if (ios /= 0) call s_mpi_abort('Cannot open IB state file: ' // trim(file_loc))
1688
1689 read (file_unit, iostat=ios) nlocal
1690 if (ios /= 0) call s_mpi_abort('Error reading IB state file header: ' // trim(file_loc))
1691
1692 do i = 1, nlocal
1693 read (file_unit, iostat=ios) gbl_id
1694 if (ios /= 0) call s_mpi_abort('Error reading IB patch ID: ' // trim(file_loc))
1695 read (file_unit, iostat=ios) ib_buf
1696 if (ios /= 0) call s_mpi_abort('Error reading IB state data: ' // trim(file_loc))
1697 ib_data(gbl_id,:) = ib_buf(:)
1698 end do
1699
1700 close (file_unit)
1701 end do
1702 else
1703 ! Build path to per-timestep IB state file
1704 write (file_loc, '(A,I0,A)') '/restart_data/ib_state_', t_step, '.dat'
1705 file_loc = trim(case_dir) // trim(file_loc)
1706
1707 inquire (file=trim(file_loc), exist=file_exist)
1708 if (.not. file_exist) call s_mpi_abort('Restart file ' // trim(file_loc) // ' does not exist!')
1709
1710 open (newunit=file_unit, file=trim(file_loc), form='unformatted', access='stream', status='old', iostat=ios)
1711 if (ios /= 0) call s_mpi_abort('Cannot open IB state file: ' // trim(file_loc))
1712
1713 do i = 1, nbodies
1714 read (file_unit, iostat=ios) ib_buf
1715 if (ios /= 0) call s_mpi_abort('Error reading IB state file')
1716 ib_data(i,:) = ib_buf(:)
1717 end do
1718
1719 close (file_unit)
1720 end if
1721
1722 do i = 1, nbodies
1723 force_x(i) = ib_data(i, 2); force_y(i) = ib_data(i, 3); force_z(i) = ib_data(i, 4)
1724 torque_x(i) = ib_data(i, 5); torque_y(i) = ib_data(i, 6); torque_z(i) = ib_data(i, 7)
1725 vel_x(i) = ib_data(i, 8); vel_y(i) = ib_data(i, 9); vel_z(i) = ib_data(i, 10)
1726 omega_x(i) = ib_data(i, 11); omega_y(i) = ib_data(i, 12); omega_z(i) = ib_data(i, 13)
1727 angle_x(i) = ib_data(i, 14); angle_y(i) = ib_data(i, 15); angle_z(i) = ib_data(i, 16)
1728 px(i) = ib_data(i, 17); py(i) = ib_data(i, 18); pz(i) = ib_data(i, 19)
1729 ib_diameter(i) = ib_data(i, 20)*2.0_wp
1730 end do
1731
1732 ! When only part of the domain is written, the bodies outside that window are dropped so the point mesh matches
1733 ! the cropped grid. A body is kept when its bounding sphere overlaps the window rather than when its centroid is
1734 ! inside it, so one straddling the boundary still appears instead of vanishing at the edge.
1735 if (output_partial_domain) then
1736 allocate (keep(nbodies))
1737
1738 do i = 1, nbodies
1739 r_ib = 0.5_wp*ib_diameter(i)
1740 keep(i) = (px(i) + r_ib >= x_output%beg) .and. (px(i) - r_ib <= x_output%end)
1741 if (n > 0) keep(i) = keep(i) .and. (py(i) + r_ib >= y_output%beg) .and. (py(i) - r_ib <= y_output%end)
1742 if (p > 0) keep(i) = keep(i) .and. (pz(i) + r_ib >= z_output%beg) .and. (pz(i) - r_ib <= z_output%end)
1743 end do
1744
1745 nkept = count(keep)
1746
1747 if (nkept < nbodies) then
1748# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1749 px(1:nkept) = pack(px(1:nbodies), keep)
1750# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1751# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1752 py(1:nkept) = pack(py(1:nbodies), keep)
1753# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1754# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1755 pz(1:nkept) = pack(pz(1:nbodies), keep)
1756# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1757# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1758 ib_diameter(1:nkept) = pack(ib_diameter(1:nbodies), keep)
1759# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1760# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1761 force_x(1:nkept) = pack(force_x(1:nbodies), keep)
1762# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1763# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1764 force_y(1:nkept) = pack(force_y(1:nbodies), keep)
1765# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1766# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1767 force_z(1:nkept) = pack(force_z(1:nbodies), keep)
1768# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1769# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1770 torque_x(1:nkept) = pack(torque_x(1:nbodies), keep)
1771# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1772# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1773 torque_y(1:nkept) = pack(torque_y(1:nbodies), keep)
1774# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1775# 1474 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1776 torque_z(1:nkept) = pack(torque_z(1:nbodies), keep)
1777# 1475 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1778# 1476 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1779# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1780 vel_x(1:nkept) = pack(vel_x(1:nbodies), keep)
1781# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1782# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1783 vel_y(1:nkept) = pack(vel_y(1:nbodies), keep)
1784# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1785# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1786 vel_z(1:nkept) = pack(vel_z(1:nbodies), keep)
1787# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1788# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1789 omega_x(1:nkept) = pack(omega_x(1:nbodies), keep)
1790# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1791# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1792 omega_y(1:nkept) = pack(omega_y(1:nbodies), keep)
1793# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1794# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1795 omega_z(1:nkept) = pack(omega_z(1:nbodies), keep)
1796# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1797# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1798 angle_x(1:nkept) = pack(angle_x(1:nbodies), keep)
1799# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1800# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1801 angle_y(1:nkept) = pack(angle_y(1:nbodies), keep)
1802# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1803# 1477 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1804 angle_z(1:nkept) = pack(angle_z(1:nbodies), keep)
1805# 1478 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1806# 1479 "/home/runner/work/MFC/MFC/src/post_process/m_data_output.fpp"
1807 end if
1808
1809 nbodies = nkept
1810 deallocate (keep)
1811 end if
1812
1813 write (meshnames(1), '(A,I0,A)') '../p0/', t_step, '.silo:ib_bodies'
1814 meshtypes(1) = db_pointmesh
1815 err = dbset2dstrlen(len(meshnames(1)))
1816 err = dbputmmesh(out%dbroot, 'ib_bodies', 16, 1, meshnames, len_trim(meshnames), meshtypes, db_f77null, ierr)
1817
1818 ! Silo carries the point-mesh coordinates in their own datatype, so they need the same single-precision
1819 ! treatment as the rectilinear mesh.
1820 if (precision == precision_single) then
1821 allocate (px_s(nbodies), py_s(nbodies), pz_s(nbodies))
1822 px_s(1:nbodies) = real(px(1:nbodies), sp)
1823 py_s(1:nbodies) = real(py(1:nbodies), sp)
1824 pz_s(1:nbodies) = real(pz(1:nbodies), sp)
1825 err = dbputpm(out%dbfile, 'ib_bodies', 9, 3, px_s, py_s, pz_s, nbodies, db_float, db_f77null, ierr)
1826 deallocate (px_s, py_s, pz_s)
1827 else
1828 err = dbputpm(out%dbfile, 'ib_bodies', 9, 3, px, py, pz, nbodies, db_double, db_f77null, ierr)
1829 end if
1830
1831 call s_write_ib_variable('ib_force_x', t_step, force_x, nbodies)
1832 call s_write_ib_variable('ib_force_y', t_step, force_y, nbodies)
1833 call s_write_ib_variable('ib_force_z', t_step, force_z, nbodies)
1834 call s_write_ib_variable('ib_torque_x', t_step, torque_x, nbodies)
1835 call s_write_ib_variable('ib_torque_y', t_step, torque_y, nbodies)
1836 call s_write_ib_variable('ib_torque_z', t_step, torque_z, nbodies)
1837 call s_write_ib_variable('ib_vel_x', t_step, vel_x, nbodies)
1838 call s_write_ib_variable('ib_vel_y', t_step, vel_y, nbodies)
1839 call s_write_ib_variable('ib_vel_z', t_step, vel_z, nbodies)
1840 call s_write_ib_variable('ib_omega_x', t_step, omega_x, nbodies)
1841 call s_write_ib_variable('ib_omega_y', t_step, omega_y, nbodies)
1842 call s_write_ib_variable('ib_omega_z', t_step, omega_z, nbodies)
1843 call s_write_ib_variable('ib_angle_x', t_step, angle_x, nbodies)
1844 call s_write_ib_variable('ib_angle_y', t_step, angle_y, nbodies)
1845 call s_write_ib_variable('ib_angle_z', t_step, angle_z, nbodies)
1846 call s_write_ib_variable('ib_diameter', t_step, ib_diameter, nbodies)
1847
1848 deallocate (ib_data, px, py, pz, force_x, force_y, force_z)
1849 deallocate (torque_x, torque_y, torque_z, vel_x, vel_y, vel_z)
1850 deallocate (omega_x, omega_y, omega_z, angle_x, angle_y, angle_z)
1851 deallocate (ib_diameter)
1852 end if
1853 end if
1854#endif
1855
1856 end subroutine s_write_ib_bodies_to_formatted_database_file
1857
1858 !> Write a single IB point-variable to the Silo database slave and master files.
1859 subroutine s_write_ib_variable(varname, t_step, data, nBodies)
1860
1861 character(len=*), intent(in) :: varname
1862 integer, intent(in) :: t_step
1863 real(wp), dimension(:), intent(in) :: data
1864 integer, intent(in) :: nBodies
1865 character(len=4*name_len) :: var_name_entry
1866 integer :: var_type_entry, ierr
1867 real(sp), dimension(:), allocatable :: data_s
1868
1869 write (var_name_entry, '(A,I0,A)') '../p0/', t_step, '.silo:' // trim(varname)
1870 var_type_entry = db_pointvar
1871 err = dbset2dstrlen(len(var_name_entry))
1872 err = dbputmvar(out%dbroot, trim(varname), len_trim(varname), 1, var_name_entry, len_trim(var_name_entry), &
1873 & var_type_entry, db_f77null, ierr)
1874
1875 if (precision == precision_single) then
1876 allocate (data_s(nbodies))
1877 data_s(1:nbodies) = real(data(1:nbodies), sp)
1878 err = dbputpv1(out%dbfile, trim(varname), len_trim(varname), 'ib_bodies', 9, data_s, nbodies, db_float, db_f77null, &
1879 & ierr)
1880 deallocate (data_s)
1881 else
1882 err = dbputpv1(out%dbfile, trim(varname), len_trim(varname), 'ib_bodies', 9, data, nbodies, db_double, db_f77null, ierr)
1883 end if
1884
1885 end subroutine s_write_ib_variable
1886
1887 !> Close the formatted database slave file and, for the root process, the master file.
1888 impure subroutine s_close_formatted_database_file()
1889
1890 integer :: ierr
1891
1892 if (format == format_silo) then
1893 ierr = dbclose(out%dbfile)
1894 if (proc_rank == 0) ierr = dbclose(out%dbroot)
1895 else
1896 close (out%dbfile)
1897 if (n == 0 .and. proc_rank == 0) close (out%dbroot)
1898 end if
1899
1900 end subroutine s_close_formatted_database_file
1901
1902 !> Close the interface data file.
1903 impure subroutine s_close_intf_data_file()
1904
1905 close (211)
1906
1907 end subroutine s_close_intf_data_file
1908
1909 !> Close the energy data file.
1910 impure subroutine s_close_energy_data_file()
1911
1912 close (251)
1913
1914 end subroutine s_close_energy_data_file
1915
1916 !> Deallocate module arrays and release all data-output resources.
1917 impure subroutine s_finalize_data_output_module()
1918
1919 deallocate (out%q_sf)
1920 if (n == 0) deallocate (out%q_root_sf)
1921 if (grid_geometry == 3) then
1922 deallocate (out%cyl_q_sf)
1923 end if
1924
1925 ! Deallocating spatial and data extents and also the variables for the offsets and the one bookkeeping the number of
1926 ! cell-boundaries in each active coordinate direction. Note that all these variables were only needed by Silo-HDF5 format
1927 ! for multidimensional data.
1928 if (format == format_silo) then
1929 deallocate (out%spatial_extents)
1930 deallocate (out%data_extents)
1931 deallocate (out%lo_offset)
1932 deallocate (out%hi_offset)
1933 deallocate (out%dims)
1934 end if
1935
1936 end subroutine s_finalize_data_output_module
1937
1938end 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_compile_specific
Platform-specific file and directory operations: create, delete, inquire, getcwd, and basename.
Definition m_compile_specific.f90:6
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::lag_io_vars
integer, parameter lag_io_vars
Number of variables per particle for MPI_IO.
Definition m_constants.fpp.f90:77
m_constants::model_eqns_5eq
integer, parameter model_eqns_5eq
Definition m_constants.fpp.f90:155
m_constants::format_silo
integer, parameter format_silo
Definition m_constants.fpp.f90:149
m_constants::format_binary
integer, parameter format_binary
Definition m_constants.fpp.f90:150
m_constants::precision_single
integer, parameter precision_single
Definition m_constants.fpp.f90:165
m_constants::model_eqns_6eq
integer, parameter model_eqns_6eq
Definition m_constants.fpp.f90:156
m_constants::model_eqns_gamma_law
integer, parameter model_eqns_gamma_law
Definition m_constants.fpp.f90:154
m_data_output
Writes post-processed grid and flow-variable data to Silo-HDF5 or binary database files.
Definition m_data_output.fpp.f90:7
m_data_output::s_write_grid_to_formatted_database_file
impure subroutine, public s_write_grid_to_formatted_database_file(t_step)
Write the computational grid (cell-boundary coordinates) to the formatted database slave and master f...
Definition m_data_output.fpp.f90:486
m_data_output::s_write_variable_to_formatted_database_file
impure subroutine, public s_write_variable_to_formatted_database_file(varname, t_step)
Write a single flow variable field to the formatted database slave and master files for a given time ...
Definition m_data_output.fpp.f90:692
m_data_output::s_open_energy_data_file
impure subroutine, public s_open_energy_data_file()
Open the energy data file for appending volume-integrated energy budget quantities.
Definition m_data_output.fpp.f90:474
m_data_output::s_open_intf_data_file
impure subroutine, public s_open_intf_data_file()
Open the interface data file for appending extracted interface coordinates.
Definition m_data_output.fpp.f90:462
m_data_output::s_write_ib_state_files
impure subroutine, public s_write_ib_state_files()
Convert the binary immersed-boundary state file to per-body formatted text files.
Definition m_data_output.fpp.f90:1413
m_data_output::s_write_energy_data_file
impure subroutine, public s_write_energy_data_file(q_prim_vf, q_cons_vf)
Compute volume-integrated kinetic, potential, and internal energies and write the energy budget to th...
Definition m_data_output.fpp.f90:1550
m_data_output::s_write_lag_bubbles_to_formatted_database_file
impure subroutine, public s_write_lag_bubbles_to_formatted_database_file(t_step)
Read Lagrangian bubble restart data and write bubble positions and scalar fields to the Silo database...
Definition m_data_output.fpp.f90:1036
m_data_output::s_write_intf_data_file
impure subroutine, public s_write_intf_data_file(q_prim_vf)
Extract the volume-fraction interface contour from primitive fields and write the coordinates to the ...
Definition m_data_output.fpp.f90:1464
m_data_output::s_write_ib_variable
subroutine s_write_ib_variable(varname, t_step, data, nbodies)
Write a single IB point-variable to the Silo database slave and master files.
Definition m_data_output.fpp.f90:1860
m_data_output::s_initialize_data_output_module
impure subroutine, public s_initialize_data_output_module()
Allocate storage arrays, configure output directories, and count flow variables for formatted databas...
Definition m_data_output.fpp.f90:41
m_data_output::s_close_energy_data_file
impure subroutine, public s_close_energy_data_file()
Close the energy data file.
Definition m_data_output.fpp.f90:1911
m_data_output::s_close_formatted_database_file
impure subroutine, public s_close_formatted_database_file()
Close the formatted database slave file and, for the root process, the master file.
Definition m_data_output.fpp.f90:1889
m_data_output::s_open_formatted_database_file
impure subroutine, public s_open_formatted_database_file(t_step)
Open (or create) the Silo-HDF5 or Binary formatted database slave and master files for a given time s...
Definition m_data_output.fpp.f90:397
m_data_output::s_close_intf_data_file
impure subroutine, public s_close_intf_data_file()
Close the interface data file.
Definition m_data_output.fpp.f90:1904
m_data_output::s_write_lag_bubbles_results_to_text
impure subroutine, public s_write_lag_bubbles_results_to_text(t_step)
Write the post-processed results in the folder 'lag_bubbles_data'.
Definition m_data_output.fpp.f90:863
m_data_output::s_finalize_data_output_module
impure subroutine, public s_finalize_data_output_module()
Deallocate module arrays and release all data-output resources.
Definition m_data_output.fpp.f90:1918
m_data_output::err
integer, private err
Definition m_data_output.fpp.f90:35
m_data_output::s_define_output_region
impure subroutine, public s_define_output_region
Compute the cell-index bounds for the user-specified partial output domain in each coordinate directi...
Definition m_data_output.fpp.f90:282
m_data_output::s_write_ib_bodies_to_formatted_database_file
impure subroutine, public s_write_ib_bodies_to_formatted_database_file(t_step)
Read IB state and write a Silo point mesh with per-body scalar fields.
Definition m_data_output.fpp.f90:1636
m_data_output::out
type(output_context), public out
Output workspace: flow variable buffers, VisIt extents/offsets, directory paths, file handles,...
Definition m_data_output.fpp.f90:32
m_data_output::s_write_lag_variable_to_formatted_database_file
subroutine s_write_lag_variable_to_formatted_database_file(varname, t_step, data, nbubs)
Write a single Lagrangian bubble point-variable to the Silo database slave and master files.
Definition m_data_output.fpp.f90:1366
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_derived_variables
Computes derived flow quantities (sound speed, vorticity, Schlieren, etc.) from conservative and prim...
Definition m_derived_variables.fpp.f90:8
m_eos
Equations of state in Gamma/Pi form, rho e = Gamma(rho) p + Pi(rho).
Definition m_eos.fpp.f90:402
m_eos::s_phase_internal_energy
subroutine, public s_phase_internal_energy(pres, alpha, alpha_rho, i, e_phase)
Internal energy per unit volume of phase i at pressure pres: alpha (Gamma p + Pi) + alpha_rho qv,...
Definition m_eos.fpp.f90:1733
m_eos::s_compute_mixture_coefficients
subroutine, public s_compute_mixture_coefficients(alpha_rho_k, alpha_k, rho_k, gamma_k, pi_inf_k, qv_k)
Mixture coefficients of one state. Under bubbles_euler with num_fluids == 1 the sole advection slot a...
Definition m_eos.fpp.f90:1980
m_eos::s_compute_speed_of_sound
subroutine, public s_compute_speed_of_sound(pres, rho, gamma, pi_inf, adv, c, alpha_rho)
Speed of sound of a thermodynamic state. Enthalpy is not an argument: for a real state H,...
Definition m_eos.fpp.f90:2153
m_global_parameters
Global parameters for the post-process: domain geometry, equation of state, and output database setti...
Definition m_global_parameters.fpp.f90:19
m_global_parameters::grid_geometry
integer grid_geometry
Definition m_global_parameters.fpp.f90:53
m_global_parameters::z_cb_s
real(sp), dimension(:), allocatable z_cb_s
Definition m_global_parameters.fpp.f90:69
m_global_parameters::x_cb_s
real(sp), dimension(:), allocatable x_cb_s
Definition m_global_parameters.fpp.f90:69
m_global_parameters::offset_y
type(int_bounds_info) offset_y
Definition m_global_parameters.fpp.f90:142
m_global_parameters::mpi_io_data_lg_bubbles
real(wp), dimension(:,:), allocatable, public mpi_io_data_lg_bubbles
Definition m_global_parameters.fpp.f90:129
m_global_parameters::y_cc
real(wp), dimension(:), allocatable y_cc
Definition m_global_parameters.fpp.f90:74
m_global_parameters::proc_rank
integer proc_rank
Rank of the local processor.
Definition m_global_parameters.fpp.f90:41
m_global_parameters::adv
real(wp), dimension(:), allocatable adv
Advection variables.
Definition m_global_parameters.fpp.f90:133
m_global_parameters::z_output_idx
type(int_bounds_info) z_output_idx
Indices of domain to output for post-processing.
Definition m_global_parameters.fpp.f90:137
m_global_parameters::y_cb
real(wp), dimension(:), allocatable y_cb
Definition m_global_parameters.fpp.f90:64
m_global_parameters::x_output
type(bounds_info) x_output
Definition m_global_parameters.fpp.f90:136
m_global_parameters::dz
real(wp), dimension(:), allocatable dz
Definition m_global_parameters.fpp.f90:80
m_global_parameters::x_root_cb
real(wp), dimension(:), allocatable x_root_cb
Definition m_global_parameters.fpp.f90:64
m_global_parameters::z_cb
real(wp), dimension(:), allocatable z_cb
Definition m_global_parameters.fpp.f90:64
m_global_parameters::z_output
type(bounds_info) z_output
Portion of domain to output for post-processing.
Definition m_global_parameters.fpp.f90:136
m_global_parameters::x_output_idx
type(int_bounds_info) x_output_idx
Definition m_global_parameters.fpp.f90:137
m_global_parameters::x_cc
real(wp), dimension(:), allocatable x_cc
Definition m_global_parameters.fpp.f90:74
m_global_parameters::x_cb
real(wp), dimension(:), allocatable x_cb
Definition m_global_parameters.fpp.f90:64
m_global_parameters::dy
real(wp), dimension(:), allocatable dy
Definition m_global_parameters.fpp.f90:80
m_global_parameters::offset_x
type(int_bounds_info) offset_x
Definition m_global_parameters.fpp.f90:142
m_global_parameters::z_cc
real(wp), dimension(:), allocatable z_cc
Definition m_global_parameters.fpp.f90:74
m_global_parameters::y_cb_s
real(sp), dimension(:), allocatable y_cb_s
Definition m_global_parameters.fpp.f90:69
m_global_parameters::num_procs
integer num_procs
Number of processors.
Definition m_global_parameters.fpp.f90:36
m_global_parameters::y_output_idx
type(int_bounds_info) y_output_idx
Definition m_global_parameters.fpp.f90:137
m_global_parameters::offset_z
type(int_bounds_info) offset_z
Definition m_global_parameters.fpp.f90:142
m_global_parameters::dx
real(wp), dimension(:), allocatable dx
Cell-width distributions in the x-, y- and z-coordinate directions.
Definition m_global_parameters.fpp.f90:80
m_global_parameters::m_root
integer m_root
Definition m_global_parameters.fpp.f90:44
m_global_parameters::y_output
type(bounds_info) y_output
Definition m_global_parameters.fpp.f90:136
m_helper
Utility routines for bubble model setup, coordinate transforms, array sampling, and special functions...
Definition m_helper.fpp.f90:367
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_mpi_proxy
MPI gather and scatter operations for distributing post-process grid and flow-variable data.
Definition m_mpi_proxy.fpp.f90:7
m_mpi_proxy::s_mpi_defragment_1d_grid_variable
impure subroutine s_mpi_defragment_1d_grid_variable
Collect the sub-domain cell-boundary or cell-center location data from all processors and put back to...
Definition m_mpi_proxy.fpp.f90:533
m_mpi_proxy::s_mpi_defragment_1d_flow_variable
impure subroutine s_mpi_defragment_1d_flow_variable(q_sf, q_root_sf)
Gather the sub-domain flow variable data from all processors and reassemble it for the entire computa...
Definition m_mpi_proxy.fpp.f90:589
m_mpi_proxy::s_mpi_gather_spatial_extents
impure subroutine s_mpi_gather_spatial_extents(spatial_extents)
Gather spatial extents from all ranks for Silo database metadata.
Definition m_mpi_proxy.fpp.f90:441
m_mpi_proxy::s_mpi_gather_data_extents
impure subroutine s_mpi_gather_data_extents(q_sf, data_extents)
Gather the Silo database metadata for the flow variable's extents to boost performance of the multidi...
Definition m_mpi_proxy.fpp.f90:554
m_variables_conversion
Conservative-to-primitive variable conversion, mixture property evaluation, and pressure computation.
Definition m_variables_conversion.fpp.f90:367
m_derived_types::output_context
Output workspace for post_process: flow variable buffers, VisIt extents/offsets, directory paths,...
Definition m_derived_types.fpp.f90:534
m_derived_types::scalar_field
Derived type annexing a scalar field (SF).
Definition m_derived_types.fpp.f90:370