MFC:Pre_process  v1.0
m_data_output.f90
Go to the documentation of this file.
1 !! __ _______________
2 !! / |/ / ____/ ____/
3 !! / /|_/ / /_ / /
4 !! / / / / __/ / /___
5 !! /_/ /_/_/ \____/
6 !!
7 !! This file is part of MFC.
8 !!
9 !! MFC is the legal property of its developers, whose names
10 !! are listed in the copyright file included with this source
11 !! distribution.
12 !!
13 !! MFC is free software: you can redistribute it and/or modify
14 !! it under the terms of the GNU General Public License as published
15 !! by the Free Software Foundation, either version 3 of the license
16 !! or any later version.
17 !!
18 !! MFC is distributed in the hope that it will be useful,
19 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
20 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 !! GNU General Public License for more details.
22 !!
23 !! You should have received a copy of the GNU General Public License
24 !! along with MFC (LICENSE).
25 !! If not, see <http://www.gnu.org/licenses/>.
26 
33 
40 
41  ! Dependencies =============================================================
42  USE m_derived_types !< definitions of the derived types
43 
44  USE m_global_parameters !< global parameters for the code
45 
46  USE m_mpi_proxy !< message passing interface (mpi) module proxy
47 
48  USE mpi !< message passing interface (mpi) module
49 
51 
53  ! ==========================================================================
54 
55 
56  IMPLICIT NONE
57 
58 
59  PRIVATE; PUBLIC :: s_write_serial_data_files, &
64 
65  abstract INTERFACE ! ===================================================
66 
69  SUBROUTINE s_write_abstract_data_files(q_cons_vf)
70 
71  IMPORT :: scalar_field, sys_size
72 
73  ! Conservative variables
74  TYPE(scalar_field), &
75  DIMENSION(sys_size), &
76  INTENT(IN) :: q_cons_vf
77 
78 
79  END SUBROUTINE s_write_abstract_data_files ! -------------------
80  END INTERFACE ! ========================================================
81 
82 
83 
84  CHARACTER(LEN = path_len + 2*name_len), PRIVATE :: t_step_dir
86 
87  CHARACTER(LEN = path_len + 2*name_len), PUBLIC :: restart_dir
89 
90  PROCEDURE(s_write_abstract_data_files), POINTER :: s_write_data_files => null()
91 
92  CONTAINS
93 
97  SUBROUTINE s_write_serial_data_files(q_cons_vf) ! -----------
98  TYPE(scalar_field), &
99  DIMENSION(sys_size), &
100  INTENT(IN) :: q_cons_vf
101 
102  LOGICAL :: file_exist
103 
104  CHARACTER(LEN=15) :: FMT
105 
106  CHARACTER(len = &
107  int(floor(log10(REAL(sys_size, kind(0d0))))) + 1) :: file_num
110 
111  CHARACTER(LEN = LEN_TRIM(t_step_dir) + name_len) :: file_loc
113 
114  INTEGER :: i,j,k,l
115  INTEGER :: t_step
116 
117  REAL(KIND(0d0)), DIMENSION(nb) :: nRtmp
118  REAL(KIND(0d0)) :: nbub
119  REAL(KIND(0d0)) :: gamma, lit_gamma, pi_inf
120  REAL(KIND(0d0)) :: rho
121 
122  t_step = 0
123 
124  ! Outputting the Locations of the Cell-boundaries ==================
125 
126  ! x-coordinate direction
127  file_loc = trim(t_step_dir) // '/x_cb.dat'
128  OPEN(1, file = trim(file_loc), form = 'unformatted', status = 'new')
129  WRITE(1) x_cb
130  CLOSE(1)
131 
132  ! y- and z-coordinate directions
133  IF(n > 0) THEN
134  ! y-coordinate direction
135  file_loc = trim(t_step_dir) // '/y_cb.dat'
136  OPEN(1, file = trim(file_loc), form = 'unformatted', &
137  status = 'new')
138  WRITE(1) y_cb
139  CLOSE(1)
140 
141  ! z-coordinate direction
142  IF(p > 0) THEN
143  file_loc = trim(t_step_dir) // '/z_cb.dat'
144  OPEN(1, file = trim(file_loc), form = 'unformatted', &
145  status = 'new')
146  WRITE(1) z_cb
147  CLOSE(1)
148  END IF
149  END IF
150  ! ==================================================================
151 
152 
153  ! Outputting Conservative Variables ================================
154  DO i = 1, sys_size
155  WRITE(file_num, '(I0)') i
156  file_loc = trim(t_step_dir) // '/q_cons_vf' // trim(file_num) &
157  // '.dat'
158  OPEN(1, file = trim(file_loc), form = 'unformatted', &
159  status = 'new')
160  WRITE(1) q_cons_vf(i)%sf
161  CLOSE(1)
162  END DO
163  ! ==================================================================
164 
165  gamma = fluid_pp(1)%gamma
166  lit_gamma = 1d0/fluid_pp(1)%gamma + 1d0
167  pi_inf = fluid_pp(1)%pi_inf
168 
169  IF (precision==1) THEN
170  fmt="(2F30.7)"
171  ELSE
172  fmt="(2F40.14)"
173  END IF
174 
175  !1D
176  IF (n ==0 .AND. p ==0) THEN
177  ! writting an output directory
178  WRITE(t_step_dir,'(A,I0,A,I0)') trim(case_dir) // '/D'
179  file_loc = trim(t_step_dir) // '/.'
180 
181  INQUIRE( file = trim(file_loc), exist = file_exist )
182 
183  IF(.NOT.file_exist) CALL system('mkdir -p ' // trim(t_step_dir))
184 
185  IF (model_eqns == 2) THEN
186  DO i = 1, sys_size
187  WRITE(file_loc,'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/prim.', i, '.', proc_rank, '.', t_step,'.dat'
188 
189  OPEN(2,file= trim(file_loc) )
190  DO j=0,m
191  CALL s_convert_to_mixture_variables( q_cons_vf, j, 0, 0, rho, gamma, pi_inf)
192 
193  lit_gamma = 1d0/gamma + 1d0
194 
195  IF ( ((i.ge.cont_idx%beg) .AND. (i.le.cont_idx%end)) &
196  .OR. &
197  ((i.ge.adv_idx%beg) .AND. (i.le.adv_idx%end)) &
198  ) THEN
199  WRITE(2,fmt) x_cb(j),q_cons_vf(i)%sf(j,0,0)
200  ELSE IF (i.eq.mom_idx%beg) THEN !u
201  WRITE(2,fmt) x_cb(j),q_cons_vf(mom_idx%beg)%sf(j,0,0)/rho
202  ELSE IF (i.eq.e_idx) THEN !p
203  IF (model_eqns == 4) THEN
204  !Tait pressure from density
205  WRITE(2,fmt) x_cb(j), &
206  (pref + pi_inf) * ( &
207  ( q_cons_vf(1)%sf(j,0,0)/ &
208  (rhoref*(1.d0-q_cons_vf(4)%sf(j,0,0))) &
209  ) ** lit_gamma ) &
210  - pi_inf
211  ELSE IF (model_eqns == 2 .AND. (bubbles .NEQV. .true.)) THEN
212  !Stiffened gas pressure from energy
213  WRITE(2,fmt) x_cb(j), &
214  ( &
215  q_cons_vf(e_idx)%sf(j,0,0) - &
216  0.5d0*(q_cons_vf(mom_idx%beg)%sf(j,0,0)**2.d0)/rho - &
217  pi_inf &
218  ) / gamma
219  ELSE
220  !Stiffened gas pressure from energy with bubbles
221  WRITE(2,fmt) x_cb(j), &
222  ( &
223  (q_cons_vf(e_idx)%sf(j,0,0) - &
224  0.5d0*(q_cons_vf(mom_idx%beg)%sf(j,0,0)**2.d0)/rho) / &
225  (1.d0 - q_cons_vf(alf_idx)%sf(j,0,0)) - &
226  pi_inf &
227  ) / gamma
228  END IF
229  ELSE IF ((i .GT. alf_idx) .AND. bubbles) THEN
230  DO k = 1,nb
231  nrtmp(k) = q_cons_vf(bub_idx%rs(k))%sf(j,0,0)
232  END DO
233  CALL s_comp_n_from_cons( q_cons_vf(alf_idx)%sf(j,0,0), nrtmp, nbub)
234 
235  WRITE(2,fmt) x_cb(j),q_cons_vf(i)%sf(j,0,0)/nbub
236  END IF
237  END DO
238  CLOSE(2)
239  END DO
240  END IF
241 
242  DO i = 1, sys_size
243  WRITE(file_loc,'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/cons.', i, '.', proc_rank, '.', t_step,'.dat'
244 
245  OPEN(2,file= trim(file_loc) )
246  DO j=0,m
247  WRITE(2,fmt) x_cb(j),q_cons_vf(i)%sf(j,0,0)
248  END DO
249  CLOSE(2)
250  END DO
251  END IF
252 
253 
254  !2D
255  IF (n > 0 .AND. p ==0) THEN
256  WRITE(t_step_dir,'(A,I0,A,I0)') trim(case_dir) // '/D'
257  file_loc = trim(t_step_dir) // '/.'
258 
259  INQUIRE( file = trim(file_loc), exist = file_exist )
260 
261  IF(.not.file_exist) CALL system('mkdir -p ' // trim(t_step_dir))
262 
263  DO i = 1, sys_size
264  WRITE(file_loc,'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/prim.', i, '.', proc_rank, '.', t_step,'.dat'
265  OPEN(2,file=trim(file_loc) )
266  DO j=0,m; DO k=0,n
267  CALL s_convert_to_mixture_variables( q_cons_vf, j, k, 0, rho, gamma, pi_inf)
268  lit_gamma = 1d0/gamma + 1d0
269 
270  IF (i==1 .OR. i==5) THEN
271  WRITE(2,fmt) x_cb(j),y_cb(k),q_cons_vf(i)%sf(j,k,0)
272  ELSE IF (i==2 .OR. i==3) THEN !u
273  WRITE(2,fmt) x_cb(j),y_cb(k),q_cons_vf(i)%sf(j,k,0)/rho
274  ELSE IF (i==4) THEN !p
275  IF (model_eqns == 4) THEN
276  WRITE(2,fmt) x_cb(j),y_cb(k), &
277  (pref+fluid_pp(1)%pi_inf) * &
278  (( &
279  q_cons_vf(1)%sf(j,k,0)/ &
280  (rhoref*(1.-q_cons_vf(4)%sf(j,k,0))) &
281  ) ** (1./fluid_pp(1)%gamma + 1.)) - fluid_pp(1)%pi_inf
282  ELSE IF (model_eqns == 2 .AND. (bubbles .NEQV. .true.)) THEN
283  !Stiffened gas pressure from energy
284  WRITE(2,fmt) x_cb(j),y_cb(k), &
285  ( &
286  q_cons_vf(e_idx)%sf(j,k,0) - &
287  0.5d0*(q_cons_vf(2)%sf(j,k,0)**2.d0)/rho - &
288  pi_inf &
289  ) / gamma
290  ELSE
291  !Stiffened gas pressure from energy with bubbles
292  WRITE(2,fmt) x_cb(j),y_cb(k), &
293  ( &
294  (q_cons_vf(e_idx)%sf(j,k,0) - &
295  0.5d0*( q_cons_vf(2)%sf(j,k,0)**2.d0 + q_cons_vf(3)%sf(j,k,0)**2.d0 ) &
296  / rho) / &
297  (1.d0 - q_cons_vf(alf_idx)%sf(j,k,0)) - &
298  pi_inf &
299  ) / gamma
300  END IF
301  ELSE IF ( (i .GT. alf_idx) .AND. bubbles) THEN
302  DO l = 1,nb
303  nrtmp(l) = q_cons_vf(bub_idx%rs(l))%sf(j,k,0)
304  END DO
305  CALL s_comp_n_from_cons( q_cons_vf(alf_idx)%sf(j,k,0), nrtmp, nbub)
306 
307  WRITE(2,fmt) x_cb(j), y_cb(k), q_cons_vf(i)%sf(j,k,0)/nbub
308  END IF
309  END DO
310  WRITE(2,*)
311  END DO
312  CLOSE(2)
313 
314  WRITE(file_loc,'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/cons.', i, '.', proc_rank, '.', t_step,'.dat'
315  OPEN(2,file= trim(file_loc) )
316  DO j=0,m
317  DO k=0,n
318  WRITE(2,fmt) x_cb(j),y_cb(k),q_cons_vf(i)%sf(j,k,0)
319  END DO
320  WRITE(2,*)
321  END DO
322  CLOSE(2)
323  END DO
324  END IF
325 
326  ! 3D
327  IF ( p > 0) THEN
328  DO i = 1,sys_size
329  WRITE(file_loc,'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) // '/cons.', i, '.', proc_rank, '.', t_step,'.dat'
330  OPEN(2,file= trim(file_loc) )
331  DO j=0,m
332  DO k=0,n
333  DO l=0,p
334  WRITE(2,fmt) x_cb(j),y_cb(k),z_cb(l), q_cons_vf(i)%sf(j,k,l)
335  END DO
336  WRITE(2,*)
337  END DO
338  WRITE(2,*)
339  END DO
340  CLOSE(2)
341  END DO
342  END IF
343 
344  END SUBROUTINE s_write_serial_data_files ! ------------------------------------
345 
346 
347 
351  SUBROUTINE s_write_parallel_data_files(q_cons_vf) ! --
353  ! Conservative variables
354  TYPE(scalar_field), &
355  DIMENSION(sys_size), &
356  INTENT(IN) :: q_cons_vf
357 
358  INTEGER :: ifile, ierr, data_size
359  INTEGER, DIMENSION(MPI_STATUS_SIZE) :: status
360  INTEGER(KIND=MPI_OFFSET_KIND) :: disp
361  INTEGER(KIND=MPI_OFFSET_KIND) :: m_MOK, n_MOK, p_MOK
362  INTEGER(KIND=MPI_OFFSET_KIND) :: WP_MOK, var_MOK, str_MOK
363  INTEGER(KIND=MPI_OFFSET_KIND) :: NVARS_MOK
364  INTEGER(KIND=MPI_OFFSET_KIND) :: MOK
365 
366  CHARACTER(LEN=path_len + 2*name_len) :: file_loc
367  LOGICAL :: file_exist
368 
369  ! Generic loop iterator
370  INTEGER :: i
371 
372  ! Initialize MPI data I/O
373  CALL s_initialize_mpi_data(q_cons_vf)
374 
375  ! Open the file to write all flow variables
376  WRITE(file_loc, '(A)') '0.dat'
377  file_loc = trim(restart_dir) // trim(mpiiofs) // trim(file_loc)
378  INQUIRE(file = trim(file_loc),exist = file_exist)
379  IF (file_exist .AND. proc_rank == 0) THEN
380  CALL mpi_file_delete(file_loc,mpi_info_int,ierr)
381  END IF
382  CALL mpi_file_open(mpi_comm_world,file_loc,ior(mpi_mode_wronly,mpi_mode_create),&
383  mpi_info_int,ifile,ierr)
384 
385  ! Size of local arrays
386  data_size = (m+1)*(n+1)*(p+1)
387 
388  ! Resize some integers so MPI can write even the biggest files
389  m_mok = int(m_glb+1, mpi_offset_kind)
390  n_mok = int(n_glb+1, mpi_offset_kind)
391  p_mok = int(p_glb+1, mpi_offset_kind)
392  wp_mok = int(8d0, mpi_offset_kind)
393  mok = int(1d0, mpi_offset_kind)
394  str_mok = int(name_len, mpi_offset_kind)
395  nvars_mok = int(sys_size, mpi_offset_kind)
396 
397  ! Write the data for each variable
398  IF (bubbles) THEN
399  DO i = 1, sys_size! adv_idx%end
400  var_mok = int(i, mpi_offset_kind)
401 
402  ! Initial displacement to skip at beginning of file
403  disp = m_mok*max(mok,n_mok)*max(mok,p_mok)*wp_mok*(var_mok-1)
404 
405  CALL mpi_file_set_view(ifile,disp,mpi_double_precision,mpi_io_data%view(i), &
406  'native',mpi_info_int,ierr)
407  CALL mpi_file_write_all(ifile,mpi_io_data%var(i)%sf,data_size, &
408  mpi_double_precision,status,ierr)
409  END DO
410  ELSE
411  DO i = 1, adv_idx%end
412  var_mok = int(i, mpi_offset_kind)
413 
414  ! Initial displacement to skip at beginning of file
415  disp = m_mok*max(mok,n_mok)*max(mok,p_mok)*wp_mok*(var_mok-1)
416 
417  CALL mpi_file_set_view(ifile,disp,mpi_double_precision,mpi_io_data%view(i), &
418  'native',mpi_info_int,ierr)
419  CALL mpi_file_write_all(ifile,mpi_io_data%var(i)%sf,data_size, &
420  mpi_double_precision,status,ierr)
421  END DO
422  END IF
423 
424  CALL mpi_file_close(ifile,ierr)
425 
426  END SUBROUTINE s_write_parallel_data_files ! ---------------------------
427 
428 
431  SUBROUTINE s_initialize_data_output_module() ! ----------------------------
432  ! Generic string used to store the address of a particular file
433  CHARACTER(LEN = LEN_TRIM(case_dir) + 2*name_len) :: file_loc
434 
435  ! Generic logical used to check the existence of directories
436  LOGICAL :: dir_check
437 
438 
439  IF (parallel_io .NEQV. .true.) THEN
440  ! Setting the address of the time-step directory
441  WRITE(t_step_dir, '(A,I0,A)') '/p_all/p', proc_rank, '/0'
442  t_step_dir = trim(case_dir) // trim(t_step_dir)
443 
444 
445  ! Checking the existence of the time-step directory, removing it, if
446  ! it exists, and creating a new copy. Note that if preexisting grid
447  ! and/or initial condition data are to be read in from the very same
448  ! location, then the above described steps are not executed here but
449  ! rather in the module m_start_up.f90.
450  IF(old_grid .NEQV. .true.) THEN
451 
452  file_loc = trim(t_step_dir) // '/'
453 
454  CALL my_inquire(file_loc,dir_check)
455 
456  IF(dir_check) CALL system('rm -rf ' // trim(t_step_dir))
457  CALL system('mkdir -p ' // trim(t_step_dir))
458 
459  END IF
460 
462  ELSE
463  WRITE(restart_dir, '(A)') '/restart_data'
464  restart_dir = trim(case_dir) // trim(restart_dir)
465 
466  IF ((old_grid .NEQV. .true.) .AND. (proc_rank == 0)) THEN
467 
468  file_loc = trim(restart_dir) // '/'
469  CALL my_inquire(file_loc,dir_check)
470 
471  IF (dir_check) CALL system('rm -rf ' // trim(restart_dir))
472  CALL system('mkdir -p ' // trim(restart_dir))
473  END IF
474 
475  CALL s_mpi_barrier()
476 
478 
479  END IF
480 
481  END SUBROUTINE s_initialize_data_output_module ! --------------------------
482 
483 
484 
486  SUBROUTINE s_finalize_data_output_module() ! ---------------------------
488  s_write_data_files => null()
489 
490  END SUBROUTINE s_finalize_data_output_module ! -------------------------
491 
492 
493 
494 
495 
496 END MODULE m_data_output
character(len=path_len+2 *name_len), private t_step_dir
Time-step folder into which grid and initial condition data will be placed.
This module takes care of writing the grid and initial condition data files into the "0" time-step di...
subroutine s_mpi_barrier()
Halts all processes until all have reached barrier.
type(physical_parameters), dimension(num_fluids_max) fluid_pp
Database of the physical parameters of each of the fluids that is present in the flow. These include the stiffened gas equation of state parameters, the Reynolds numbers and the Weber numbers.
integer mpi_info_int
MPI info for parallel IO with Lustre file systems.
This module contains all of the parameters characterizing the computational domain, simulation algorithm, initial condition and the stiffened equation of state.
Derived type annexing a scalar field (SF)
logical old_grid
Use existing grid data.
procedure(s_convert_xxxxx_to_mixture_variables), pointer s_convert_to_mixture_variables
real(kind(0d0)), dimension(:), allocatable y_cb
real(kind(0d0)), dimension(:), allocatable x_cb
subroutine, public s_initialize_data_output_module()
Computation of parameters, allocation procedures, and/or any other tasks needed to properly setup the...
type(int_bounds_info) cont_idx
Indexes of first & last continuity eqns.
This module contains subroutines that are compiler specific.
integer p_glb
Global number of cells in each direction.
character(len=path_len+2 *name_len), public restart_dir
Restart data folder.
This file contains the definitions of all of the custom-defined types used in the pre-process code...
type(int_bounds_info) adv_idx
Indexes of first & last advection eqns.
real(kind(0d0)), dimension(:), allocatable z_cb
Locations of cell-boundaries (cb) in x-, y- and z-directions, respectively.
integer p
Number of cells in the x-, y- and z-coordinate directions.
integer alf_idx
Index of void fraction.
subroutine, public s_finalize_data_output_module()
Resets s_write_data_files pointer.
integer sys_size
Number of unknowns in the system of equations.
integer model_eqns
Multicomponent flow model.
type(bub_bounds_info) bub_idx
Indexes of first & last bubble variable eqns.
integer precision
Precision of output files.
type(mpi_io_var), public mpi_io_data
character(len=path_len) case_dir
Case folder location.
character(len=name_len) mpiiofs
Interface for the conservative data.
integer e_idx
Index of total energy equation.
real(kind(0d0)) pref
Reference parameters for Tait EOS.
type(int_bounds_info) mom_idx
Indexes of first & last momentum eqns.
subroutine, public s_write_parallel_data_files(q_cons_vf)
Writes grid and initial condition data files in parallel to the "0" time-step directory in the local ...
subroutine, public s_write_serial_data_files(q_cons_vf)
Writes grid and initial condition data files to the "0" time-step directory in the local processor ra...
subroutine my_inquire(fileloc, dircheck)
Inquires on the existence of a directory.
subroutine s_initialize_mpi_data(q_cons_vf)
logical parallel_io
Format of the data files.
subroutine s_comp_n_from_cons(vftmp, nRtmp, ntmp)
Computes the bubble number density n from the conservative variables.
procedure(s_write_abstract_data_files), pointer, public s_write_data_files
integer, parameter name_len
Maximum name length.
This module consists of subroutines used in the conversion of the conservative variables into the pri...
integer proc_rank
Rank of the local processor.
This module serves as a proxy to the parameters and subroutines available in the MPI implementation&#39;s...
Definition: m_mpi_proxy.f90:39