diff --git a/.gitignore b/.gitignore index c8b2c016b..8f584b61b 100644 --- a/.gitignore +++ b/.gitignore @@ -92,4 +92,7 @@ benchmarks/*.png *.avi **isolation_rules/ -**.supercode/ \ No newline at end of file +**.supercode/ + +# User-generated sim files +sims/ \ No newline at end of file diff --git a/CMakeLists.txt b/CMakeLists.txt index 751a690b5..246f03ab7 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -460,8 +460,10 @@ function(MFC_SETUP_TARGET) endif() if (ARGS_HDF5) - find_package(HDF5 REQUIRED) + find_package(HDF5 REQUIRED COMPONENTS Fortran) target_link_libraries(${a_target} PRIVATE HDF5::HDF5) + target_include_directories(${a_target} PRIVATE ${HDF5_Fortran_INCLUDE_DIRS}) + target_link_libraries(${a_target} PRIVATE ${HDF5_Fortran_LIBRARIES}) endif() if (ARGS_FFTW) diff --git a/src/post_process/m_data_output.fpp b/src/post_process/m_data_output.fpp index 410ab9ebc..031d687a4 100644 --- a/src/post_process/m_data_output.fpp +++ b/src/post_process/m_data_output.fpp @@ -4,9 +4,10 @@ !> @brief This module enables the restructuring of the raw simulation data !! file(s) into formatted database file(s). The formats that may be -!! chosen from include Silo-HDF5 and Binary. Each of these database -!! structures contains information about the grid as well as each of -!! the flow variable(s) that were chosen by the user to be included. +!! chosen from include Silo-HDF5, Binary, and HDF5+XDMF. Each of +!! these database structures contains information about the grid as +!! well as each of the flow variable(s) that were chosen by the user +!! to be included. module m_data_output use m_derived_types ! Definitions of the derived types @@ -23,6 +24,8 @@ module m_data_output use m_variables_conversion + use m_hdf5_xdmf_output ! Native HDF5+XDMF output support + implicit none private; public :: s_initialize_data_output_module, & @@ -37,6 +40,7 @@ module m_data_output s_write_intf_data_file, & s_write_energy_data_file, & s_close_formatted_database_file, & + s_write_xdmf_metadata_file, & s_close_intf_data_file, & s_close_energy_data_file, & s_finalize_data_output_module @@ -253,7 +257,7 @@ contains ! Generating Binary Directory Tree - else + elseif (format == 2) then ! Creating the directory associated with the local process dbdir = trim(case_dir)//'/binary' @@ -285,6 +289,13 @@ contains end if + ! Generating HDF5+XDMF Directory Tree + + elseif (format == 3) then + + ! Initialize the HDF5+XDMF module (creates directory structure) + call s_initialize_hdf5_xdmf_module() + end if if (bubbles_lagrange) then !Lagrangian solver @@ -539,7 +550,7 @@ contains ! Binary Database Format - else + elseif (format == 2) then ! Generating the relative path to the formatted database slave ! file, that is to be opened for the current time-step, t_step @@ -597,6 +608,13 @@ contains end if + ! HDF5+XDMF Database Format + + elseif (format == 3) then + + ! Open new HDF5 file for this time-step + call s_open_hdf5_file(t_step) + end if end subroutine s_open_formatted_database_file @@ -836,6 +854,13 @@ contains end if + ! HDF5+XDMF Database Format + + elseif (format == 3) then + + ! Write grid coordinates to HDF5 file + call s_write_grid_to_hdf5() + end if end subroutine s_write_grid_to_formatted_database_file @@ -1003,7 +1028,7 @@ contains ! Binary Database Format - else + elseif (format == 2) then ! Writing the name of the flow variable and its data, associated ! with the local processor, to the formatted database slave file @@ -1034,6 +1059,13 @@ contains end if + ! HDF5+XDMF Database Format + + elseif (format == 3) then + + ! Write the flow variable to HDF5 file + call s_write_variable_to_hdf5(varname, q_sf) + end if end subroutine s_write_variable_to_formatted_database_file @@ -1702,10 +1734,14 @@ contains if (proc_rank == 0) ierr = DBCLOSE(dbroot) ! Binary database format - else + elseif (format == 2) then close (dbfile) if (n == 0 .and. proc_rank == 0) close (dbroot) + ! HDF5+XDMF database format + elseif (format == 3) then + call s_close_hdf5_file() + end if end subroutine s_close_formatted_database_file @@ -1722,6 +1758,21 @@ contains end subroutine s_close_energy_data_file + !> Subroutine to write XDMF metadata file for HDF5+XDMF format + !! @param t_step Current time step + !! @param time Physical time value + impure subroutine s_write_xdmf_metadata_file(t_step, time) + + integer, intent(in) :: t_step + real(wp), intent(in) :: time + + ! Only write XDMF for format=3 + if (format == 3) then + call s_write_xdmf_file(t_step, time) + end if + + end subroutine s_write_xdmf_metadata_file + impure subroutine s_finalize_data_output_module() ! Description: Deallocation procedures for the module @@ -1746,6 +1797,11 @@ contains deallocate (dims) end if + ! Finalize HDF5+XDMF module if used + if (format == 3) then + call s_finalize_hdf5_xdmf_module() + end if + end subroutine s_finalize_data_output_module end module m_data_output diff --git a/src/post_process/m_hdf5_xdmf_output.fpp b/src/post_process/m_hdf5_xdmf_output.fpp new file mode 100644 index 000000000..f7a6517cb --- /dev/null +++ b/src/post_process/m_hdf5_xdmf_output.fpp @@ -0,0 +1,354 @@ +!> +!! @file m_hdf5_xdmf_output.fpp +!! @brief Contains module m_hdf5_xdmf_output + +!> @brief This module provides native HDF5 + XDMF output support for MFC. +!! XDMF is an XML-based format that describes HDF5 data for visualization +!! tools like ParaView. This format is more portable than Silo-HDF5 and +!! works reliably with modern ParaView versions (6.0+). +module m_hdf5_xdmf_output + + use m_derived_types + use m_global_parameters + use m_mpi_proxy + use m_compile_specific + use hdf5 + + implicit none + + private; public :: s_initialize_hdf5_xdmf_module, & + s_open_hdf5_file, & + s_write_grid_to_hdf5, & + s_write_variable_to_hdf5, & + s_close_hdf5_file, & + s_write_xdmf_file, & + s_finalize_hdf5_xdmf_module + + ! HDF5 file and group identifiers + integer(HID_T) :: h5file_id !< HDF5 file identifier + integer(HID_T) :: h5dspace_id !< HDF5 dataspace identifier + integer(HID_T) :: h5dset_id !< HDF5 dataset identifier + integer(HID_T) :: h5plist_id !< HDF5 property list identifier + + ! Directory paths for output + character(LEN=path_len + name_len) :: hdf5_dir + character(LEN=path_len + 2*name_len) :: proc_rank_dir_hdf5 + character(LEN=path_len + 2*name_len) :: rootdir_hdf5 + + ! Error flag + integer :: h5err + + ! Track variable names for XDMF generation + character(LEN=name_len), allocatable, dimension(:) :: stored_varnames + integer :: num_stored_vars + +contains + + !> @brief Initialize the HDF5+XDMF output module + subroutine s_initialize_hdf5_xdmf_module() + + character(LEN=len_trim(case_dir) + 2*name_len) :: file_loc + logical :: dir_check + integer :: i + + ! Initialize HDF5 library + call h5open_f(h5err) + if (h5err /= 0) then + call s_mpi_abort('Failed to initialize HDF5 library. Exiting.') + end if + + ! Create HDF5 output directory structure + hdf5_dir = trim(case_dir)//'/hdf5_xdmf' + + write (proc_rank_dir_hdf5, '(A,I0)') '/p', proc_rank + proc_rank_dir_hdf5 = trim(hdf5_dir)//trim(proc_rank_dir_hdf5) + + file_loc = trim(proc_rank_dir_hdf5)//'/.' + call my_inquire(file_loc, dir_check) + if (.not. dir_check) then + call s_create_directory(trim(proc_rank_dir_hdf5)) + end if + + ! Root process creates root directory for XDMF files + if (proc_rank == 0) then + rootdir_hdf5 = trim(hdf5_dir)//'/root' + file_loc = trim(rootdir_hdf5)//'/.' + call my_inquire(file_loc, dir_check) + if (.not. dir_check) then + call s_create_directory(trim(rootdir_hdf5)) + end if + end if + + ! Allocate variable name storage for XDMF generation + allocate (stored_varnames(100)) ! Max 100 variables + num_stored_vars = 0 + + end subroutine s_initialize_hdf5_xdmf_module + + !> @brief Open an HDF5 file for the current time step + !! @param t_step Current time step + subroutine s_open_hdf5_file(t_step) + + integer, intent(in) :: t_step + character(LEN=len_trim(case_dir) + 3*name_len) :: file_loc + + ! Generate file path + write (file_loc, '(A,I0,A)') '/', t_step, '.h5' + file_loc = trim(proc_rank_dir_hdf5)//trim(file_loc) + + ! Create HDF5 file with default properties + call h5fcreate_f(trim(file_loc), H5F_ACC_TRUNC_F, h5file_id, h5err) + if (h5err /= 0) then + call s_mpi_abort('Unable to create HDF5 file '//trim(file_loc)//'. Exiting.') + end if + + ! Reset variable counter for this timestep + num_stored_vars = 0 + + end subroutine s_open_hdf5_file + + !> @brief Write grid coordinates to HDF5 file + subroutine s_write_grid_to_hdf5() + + integer(HSIZE_T), dimension(1) :: dims_1d + integer :: i + + ! Write X coordinates + dims_1d(1) = size(x_cb) + call h5screate_simple_f(1, dims_1d, h5dspace_id, h5err) + call h5dcreate_f(h5file_id, 'x_cb', H5T_NATIVE_DOUBLE, h5dspace_id, h5dset_id, h5err) + call h5dwrite_f(h5dset_id, H5T_NATIVE_DOUBLE, x_cb, dims_1d, h5err) + call h5dclose_f(h5dset_id, h5err) + call h5sclose_f(h5dspace_id, h5err) + + ! Write Y coordinates if 2D or 3D + if (n > 0) then + dims_1d(1) = size(y_cb) + call h5screate_simple_f(1, dims_1d, h5dspace_id, h5err) + call h5dcreate_f(h5file_id, 'y_cb', H5T_NATIVE_DOUBLE, h5dspace_id, h5dset_id, h5err) + call h5dwrite_f(h5dset_id, H5T_NATIVE_DOUBLE, y_cb, dims_1d, h5err) + call h5dclose_f(h5dset_id, h5err) + call h5sclose_f(h5dspace_id, h5err) + end if + + ! Write Z coordinates if 3D + if (p > 0) then + dims_1d(1) = size(z_cb) + call h5screate_simple_f(1, dims_1d, h5dspace_id, h5err) + call h5dcreate_f(h5file_id, 'z_cb', H5T_NATIVE_DOUBLE, h5dspace_id, h5dset_id, h5err) + call h5dwrite_f(h5dset_id, H5T_NATIVE_DOUBLE, z_cb, dims_1d, h5err) + call h5dclose_f(h5dset_id, h5err) + call h5sclose_f(h5dspace_id, h5err) + end if + + end subroutine s_write_grid_to_hdf5 + + !> @brief Write a flow variable to HDF5 file + !! @param varname Name of the variable + !! @param q_sf 3D array of variable data + subroutine s_write_variable_to_hdf5(varname, q_sf) + + character(LEN=*), intent(in) :: varname + real(wp), dimension(-offset_x%beg:, -offset_y%beg:, -offset_z%beg:), intent(in) :: q_sf + + integer(HSIZE_T), dimension(3) :: dims_3d + integer(HSIZE_T), dimension(2) :: dims_2d + integer(HSIZE_T), dimension(1) :: dims_1d + + ! Determine dimensions based on problem dimensionality + if (p > 0) then + ! 3D case + dims_3d(1) = m + offset_x%beg + offset_x%end + 1 + dims_3d(2) = n + offset_y%beg + offset_y%end + 1 + dims_3d(3) = p + offset_z%beg + offset_z%end + 1 + call h5screate_simple_f(3, dims_3d, h5dspace_id, h5err) + call h5dcreate_f(h5file_id, trim(varname), H5T_NATIVE_DOUBLE, h5dspace_id, h5dset_id, h5err) + call h5dwrite_f(h5dset_id, H5T_NATIVE_DOUBLE, q_sf, dims_3d, h5err) + elseif (n > 0) then + ! 2D case + dims_2d(1) = m + offset_x%beg + offset_x%end + 1 + dims_2d(2) = n + offset_y%beg + offset_y%end + 1 + call h5screate_simple_f(2, dims_2d, h5dspace_id, h5err) + call h5dcreate_f(h5file_id, trim(varname), H5T_NATIVE_DOUBLE, h5dspace_id, h5dset_id, h5err) + call h5dwrite_f(h5dset_id, H5T_NATIVE_DOUBLE, q_sf(:, :, 0), dims_2d, h5err) + else + ! 1D case + dims_1d(1) = m + offset_x%beg + offset_x%end + 1 + call h5screate_simple_f(1, dims_1d, h5dspace_id, h5err) + call h5dcreate_f(h5file_id, trim(varname), H5T_NATIVE_DOUBLE, h5dspace_id, h5dset_id, h5err) + call h5dwrite_f(h5dset_id, H5T_NATIVE_DOUBLE, q_sf(:, 0, 0), dims_1d, h5err) + end if + + call h5dclose_f(h5dset_id, h5err) + call h5sclose_f(h5dspace_id, h5err) + + ! Store variable name for XDMF generation + num_stored_vars = num_stored_vars + 1 + stored_varnames(num_stored_vars) = trim(varname) + + end subroutine s_write_variable_to_hdf5 + + !> @brief Close the current HDF5 file + subroutine s_close_hdf5_file() + + call h5fclose_f(h5file_id, h5err) + if (h5err /= 0) then + call s_mpi_abort('Failed to close HDF5 file. Exiting.') + end if + + end subroutine s_close_hdf5_file + + !> @brief Write XDMF metadata file for a time step (root process only) + !! @param t_step Current time step + !! @param time Physical time value + subroutine s_write_xdmf_file(t_step, time) + + integer, intent(in) :: t_step + real(wp), intent(in) :: time + + character(LEN=len_trim(case_dir) + 3*name_len) :: xdmf_file + character(LEN=256) :: dims_str, h5_path + integer :: unit_num, i, ip + integer :: nx, ny, nz + + ! Only root process writes XDMF + if (proc_rank /= 0) return + + ! Calculate grid dimensions (cell-centered values) + nx = m + offset_x%beg + offset_x%end + 1 + if (n > 0) then + ny = n + offset_y%beg + offset_y%end + 1 + else + ny = 1 + end if + if (p > 0) then + nz = p + offset_z%beg + offset_z%end + 1 + else + nz = 1 + end if + + ! Generate XDMF file path + write (xdmf_file, '(A,I0,A)') '/timestep_', t_step, '.xdmf' + xdmf_file = trim(rootdir_hdf5)//trim(xdmf_file) + + unit_num = 100 + open (unit=unit_num, file=trim(xdmf_file), status='replace', action='write') + + ! Write XDMF header + write (unit_num, '(A)') '' + write (unit_num, '(A)') '' + write (unit_num, '(A)') '' + write (unit_num, '(A)') ' ' + + ! Write grid collection for all processors + write (unit_num, '(A)') ' ' + write (unit_num, '(A,E15.8,A)') ' ' + write (unit_num, '(A)') ' ' + write (unit_num, '(A)') '' + + close (unit_num) + + end subroutine s_write_xdmf_file + + !> @brief Finalize HDF5+XDMF module + subroutine s_finalize_hdf5_xdmf_module() + + ! Close HDF5 library + call h5close_f(h5err) + + ! Deallocate variable name storage + if (allocated(stored_varnames)) deallocate (stored_varnames) + + end subroutine s_finalize_hdf5_xdmf_module + +end module m_hdf5_xdmf_output diff --git a/src/post_process/m_start_up.fpp b/src/post_process/m_start_up.fpp index 2b5c3c844..202aa6694 100644 --- a/src/post_process/m_start_up.fpp +++ b/src/post_process/m_start_up.fpp @@ -880,6 +880,9 @@ contains call s_close_energy_data_file() end if + ! Write XDMF metadata file for HDF5+XDMF format (before closing HDF5 file) + call s_write_xdmf_metadata_file(t_step, real(t_step, wp)) + ! Closing the formatted database file call s_close_formatted_database_file() diff --git a/toolchain/mfc/case_validator.py b/toolchain/mfc/case_validator.py index 33d0be5e3..72186e03b 100644 --- a/toolchain/mfc/case_validator.py +++ b/toolchain/mfc/case_validator.py @@ -1322,8 +1322,8 @@ def check_output_format(self): precision = self.get('precision') if format is not None: - self.prohibit(format not in [1, 2], - "format must be 1 or 2") + self.prohibit(format not in [1, 2, 3], + "format must be 1 (Silo-HDF5), 2 (Binary), or 3 (HDF5+XDMF)") if precision is not None: self.prohibit(precision not in [1, 2],