From bce9c4c262e2ab7e64997425ba6639d618d5f234 Mon Sep 17 00:00:00 2001 From: Haocheng Yu <98496194+haochey@users.noreply.github.com> Date: Mon, 1 Jul 2024 13:10:59 -0400 Subject: [PATCH] Add ib_markers to Silo Output (#496) --- src/common/m_mpi_common.fpp | 2 - src/post_process/m_data_input.f90 | 110 ++++++++++++++++++++++- src/post_process/m_global_parameters.fpp | 6 ++ src/post_process/m_mpi_proxy.fpp | 2 +- src/post_process/m_start_up.f90 | 8 +- src/simulation/m_data_output.fpp | 6 +- toolchain/mfc/run/case_dicts.py | 3 +- 7 files changed, 129 insertions(+), 8 deletions(-) diff --git a/src/common/m_mpi_common.fpp b/src/common/m_mpi_common.fpp index ab86ce52cf..78fe8df374 100644 --- a/src/common/m_mpi_common.fpp +++ b/src/common/m_mpi_common.fpp @@ -133,7 +133,6 @@ contains end if #endif -#ifndef MFC_POST_PROCESS if (present(ib_markers)) then #ifdef MFC_PRE_PROCESS @@ -146,7 +145,6 @@ contains call MPI_TYPE_COMMIT(MPI_IO_IB_DATA%view, ierr) end if -#endif #ifndef MFC_POST_PROCESS if (present(ib_markers)) then diff --git a/src/post_process/m_data_input.f90 b/src/post_process/m_data_input.f90 index be3eadd8ec..567816f73b 100644 --- a/src/post_process/m_data_input.f90 +++ b/src/post_process/m_data_input.f90 @@ -51,6 +51,9 @@ end subroutine s_read_abstract_data_files type(scalar_field), allocatable, dimension(:), public :: q_prim_vf !< !! Primitive variables + ! type(scalar_field), public :: ib_markers !< + type(integer_field), public :: ib_markers + procedure(s_read_abstract_data_files), pointer :: s_read_data_files => null() contains @@ -75,6 +78,11 @@ subroutine s_read_serial_data_files(t_step) !! Used to store the variable position, in character form, of the !! currently manipulated conservative variable file + character(LEN=len_trim(case_dir) + 2*name_len) :: t_step_ib_dir !< + !! Location of the time-step directory associated with t_step + + character(LEN=len_trim(case_dir) + 3*name_len) :: file_loc_ib !< + logical :: dir_check !< !! Generic logical used to test the existence of a particular folder @@ -87,9 +95,14 @@ subroutine s_read_serial_data_files(t_step) write (t_step_dir, '(A,I0,A,I0)') '/p_all/p', proc_rank, '/', t_step t_step_dir = trim(case_dir)//trim(t_step_dir) + write (t_step_ib_dir, '(A,I0,A,I0)') '/p_all/p', proc_rank, '/', 0 + t_step_ib_dir = trim(case_dir)//trim(t_step_ib_dir) + ! Inquiring as to the existence of the time-step directory file_loc = trim(t_step_dir)//'/.' + file_loc_ib = trim(t_step_ib_dir)//'/.' + call my_inquire(file_loc, dir_check) ! If the time-step directory is missing, the post-process exits. @@ -98,6 +111,14 @@ subroutine s_read_serial_data_files(t_step) ' is missing. Exiting ...') end if + call my_inquire(file_loc_ib, dir_check) + + ! If the time-step directory is missing, the post-process exits. + if (dir_check .neqv. .true.) then + call s_mpi_abort('Time-step folder '//trim(t_step_ib_dir)// & + ' is missing. Exiting ...') + end if + ! Reading the Grid Data File for the x-direction =================== ! Checking whether x_cb.dat exists @@ -205,6 +226,19 @@ subroutine s_read_serial_data_files(t_step) end do + if (ib) then + write (file_loc_ib, '(A,I0,A)') & + trim(t_step_ib_dir)//'/ib.dat' + inquire (FILE=trim(file_loc_ib), EXIST=file_check) + if (file_check) then + open (2, FILE=trim(file_loc_ib), & + FORM='unformatted', & + ACTION='read', & + STATUS='old') + call s_mpi_abort(trim(file_loc)//' is missing. Exiting ...') + end if + end if + ! ================================================================== end subroutine s_read_serial_data_files @@ -317,7 +351,11 @@ subroutine s_read_parallel_data_files(t_step) call MPI_FILE_OPEN(MPI_COMM_SELF, file_loc, MPI_MODE_RDONLY, mpi_info_int, ifile, ierr) ! Initialize MPI data I/O - call s_initialize_mpi_data(q_cons_vf) + if (ib) then + call s_initialize_mpi_data(q_cons_vf, ib_markers) + else + call s_initialize_mpi_data(q_cons_vf) + end if ! Size of local arrays data_size = (m + 1)*(n + 1)*(p + 1) @@ -351,6 +389,29 @@ subroutine s_read_parallel_data_files(t_step) call s_mpi_barrier() call MPI_FILE_CLOSE(ifile, ierr) + + if (ib) then + + write (file_loc, '(A)') 'ib.dat' + file_loc = trim(case_dir)//'/restart_data'//trim(mpiiofs)//trim(file_loc) + inquire (FILE=trim(file_loc), EXIST=file_exist) + + if (file_exist) then + + call MPI_FILE_OPEN(MPI_COMM_WORLD, file_loc, MPI_MODE_RDONLY, mpi_info_int, ifile, ierr) + + disp = 0 + + call MPI_FILE_SET_VIEW(ifile, disp, MPI_INTEGER, MPI_IO_IB_DATA%view, & + 'native', mpi_info_int, ierr) + call MPI_FILE_READ(ifile, MPI_IO_IB_DATA%var%sf, data_size, & + MPI_INTEGER, status, ierr) + + else + call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting...') + end if + + end if else call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting...') end if @@ -364,7 +425,11 @@ subroutine s_read_parallel_data_files(t_step) call MPI_FILE_OPEN(MPI_COMM_WORLD, file_loc, MPI_MODE_RDONLY, mpi_info_int, ifile, ierr) ! Initialize MPI data I/O - call s_initialize_mpi_data(q_cons_vf) + if (ib) then + call s_initialize_mpi_data(q_cons_vf, ib_markers) + else + call s_initialize_mpi_data(q_cons_vf) + end if ! Size of local arrays data_size = (m + 1)*(n + 1)*(p + 1) @@ -408,6 +473,28 @@ subroutine s_read_parallel_data_files(t_step) call s_mpi_barrier() call MPI_FILE_CLOSE(ifile, ierr) + + if (ib) then + + write (file_loc, '(A)') 'ib.dat' + file_loc = trim(case_dir)//'/restart_data'//trim(mpiiofs)//trim(file_loc) + inquire (FILE=trim(file_loc), EXIST=file_exist) + + if (file_exist) then + + call MPI_FILE_OPEN(MPI_COMM_WORLD, file_loc, MPI_MODE_RDONLY, mpi_info_int, ifile, ierr) + + disp = 0 + + call MPI_FILE_SET_VIEW(ifile, disp, MPI_INTEGER, MPI_IO_IB_DATA%view, & + 'native', mpi_info_int, ierr) + call MPI_FILE_READ(ifile, MPI_IO_IB_DATA%var%sf, data_size, & + MPI_INTEGER, status, ierr) + + else + call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting...') + end if + end if else call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting...') end if @@ -1079,6 +1166,12 @@ subroutine s_initialize_data_input_module -buff_size:p + buff_size)) end do + if (ib) then + allocate (ib_markers%sf(-buff_size:m + buff_size, & + -buff_size:n + buff_size, & + -buff_size:p + buff_size)) + end if + ! Simulation is 2D else @@ -1091,6 +1184,11 @@ subroutine s_initialize_data_input_module 0:0)) end do + if (ib) then + allocate (ib_markers%sf(-buff_size:m + buff_size, & + -buff_size:n + buff_size, & + 0:0)) + end if end if ! Simulation is 1D @@ -1105,6 +1203,10 @@ subroutine s_initialize_data_input_module 0:0)) end do + if (ib) then + allocate (ib_markers%sf(-buff_size:m + buff_size, 0:0, 0:0)) + end if + end if if (parallel_io .neqv. .true.) then @@ -1129,6 +1231,10 @@ subroutine s_finalize_data_input_module deallocate (q_cons_vf) deallocate (q_prim_vf) + if (ib) then + deallocate (ib_markers%sf) + end if + s_read_data_files => null() end subroutine s_finalize_data_input_module diff --git a/src/post_process/m_global_parameters.fpp b/src/post_process/m_global_parameters.fpp index 5c5ed3ef1f..5575672652 100644 --- a/src/post_process/m_global_parameters.fpp +++ b/src/post_process/m_global_parameters.fpp @@ -128,6 +128,7 @@ module m_global_parameters #ifdef MFC_MPI type(mpi_io_var), public :: MPI_IO_DATA + type(mpi_io_ib_var), public :: MPI_IO_IB_DATA #endif @@ -191,6 +192,7 @@ module m_global_parameters logical :: qm_wrt logical :: schlieren_wrt logical :: cf_wrt + logical :: ib !> @} real(kind(0d0)), dimension(num_fluids_max) :: schlieren_alpha !< @@ -334,6 +336,7 @@ contains qm_wrt = .false. schlieren_wrt = .false. cf_wrt = .false. + ib = .false. schlieren_alpha = dflt_real @@ -598,6 +601,8 @@ contains allocate (MPI_IO_DATA%var(i)%sf(0:m, 0:n, 0:p)) MPI_IO_DATA%var(i)%sf => null() end do + + if (ib) allocate (MPI_IO_IB_DATA%var%sf(0:m, 0:n, 0:p)) #endif ! Size of the ghost zone layer is non-zero only when post-processing @@ -752,6 +757,7 @@ contains deallocate (MPI_IO_DATA%view) end if + if (ib) MPI_IO_IB_DATA%var%sf => null() #endif end subroutine s_finalize_global_parameters_module diff --git a/src/post_process/m_mpi_proxy.fpp b/src/post_process/m_mpi_proxy.fpp index adf9a8085e..1da63aceef 100644 --- a/src/post_process/m_mpi_proxy.fpp +++ b/src/post_process/m_mpi_proxy.fpp @@ -168,7 +168,7 @@ contains & 'heat_ratio_wrt', 'pi_inf_wrt', 'pres_inf_wrt', 'cons_vars_wrt', & & 'prim_vars_wrt', 'c_wrt', 'qm_wrt','schlieren_wrt', 'bubbles', 'qbmm', & & 'polytropic', 'polydisperse', 'file_per_process', 'relax', 'cf_wrt', & - & 'adv_n' ] + & 'adv_n', 'ib' ] call MPI_BCAST(${VAR}$, 1, MPI_LOGICAL, 0, MPI_COMM_WORLD, ierr) #:endfor diff --git a/src/post_process/m_start_up.f90 b/src/post_process/m_start_up.f90 index 55659e2547..b20bc4df8a 100644 --- a/src/post_process/m_start_up.f90 +++ b/src/post_process/m_start_up.f90 @@ -74,7 +74,7 @@ subroutine s_read_input_file parallel_io, rhoref, pref, bubbles, qbmm, sigR, & R0ref, nb, polytropic, thermal, Ca, Web, Re_inv, & polydisperse, poly_sigma, file_per_process, relax, & - relax_model, cf_wrt, sigma, adv_n + relax_model, cf_wrt, sigma, adv_n, ib ! Inquiring the status of the post_process.inp file file_loc = 'post_process.inp' @@ -496,6 +496,12 @@ subroutine s_save_data(t_step, varname, pres, c, H) end if ! ---------------------------------------------------------------------- + if (ib) then + q_sf = real(ib_markers%sf(-offset_x%beg:m + offset_x%end, -offset_y%beg:n + offset_y%end, -offset_z%beg:p + offset_z%end)) + varname = 'ib_markers' + call s_write_variable_to_formatted_database_file(varname, t_step) + end if + ! Adding Q_M to the formatted database file ------------------ if (p > 0 .and. qm_wrt) then call s_derive_qm(q_prim_vf, q_sf) diff --git a/src/simulation/m_data_output.fpp b/src/simulation/m_data_output.fpp index f1b9d8ba76..f630073fda 100644 --- a/src/simulation/m_data_output.fpp +++ b/src/simulation/m_data_output.fpp @@ -928,7 +928,11 @@ contains else ! Initialize MPI data I/O - call s_initialize_mpi_data(q_cons_vf) + if (ib) then + call s_initialize_mpi_data(q_cons_vf, ib_markers) + else + call s_initialize_mpi_data(q_cons_vf) + end if ! Open the file to write all flow variables write (file_loc, '(I0,A)') t_step, '.dat' diff --git a/toolchain/mfc/run/case_dicts.py b/toolchain/mfc/run/case_dicts.py index 542eb7767b..b7f4360e84 100644 --- a/toolchain/mfc/run/case_dicts.py +++ b/toolchain/mfc/run/case_dicts.py @@ -308,7 +308,8 @@ def analytic(self): 'omega_wrt': ParamType.LOG, 'qbmm': ParamType.LOG, 'qm_wrt': ParamType.LOG, - 'cf_wrt': ParamType.LOG + 'cf_wrt': ParamType.LOG, + 'ib': ParamType.LOG }) for cmp_id in range(1,3+1):