hdf5IO.h (3171B)
1 // File : hdf5IO.h 2 // Description: HDF5 parallel MPI dumper 3 // Copyright 2023 Harvard University. All Rights Reserved. 4 #ifndef HDF5IO_H_XUXZ7EB0 5 #define HDF5IO_H_XUXZ7EB0 6 7 #include "xdmf_wrapper.h" 8 #include <hdf5.h> 9 #include <iostream> 10 #include <mpi.h> 11 #include <string> 12 #include <vector> 13 14 void write_hdf5_mpi(const std::string fname, 15 const std::vector<float> &u, 16 const int Nx, 17 const int Ny, 18 const int Nz, 19 const MPI_Comm cart_comm) 20 { 21 int rank, dims[3], periods[3], coords[3]; 22 MPI_Comm_rank(cart_comm, &rank); 23 MPI_Cart_get(cart_comm, 3, dims, periods, coords); 24 25 // open HDF5 file 26 H5open(); 27 hid_t plist_id = H5Pcreate(H5P_FILE_ACCESS); // file access policy 28 H5Pset_fapl_mpio(plist_id, cart_comm, MPI_INFO_NULL); 29 hid_t file_id = H5Fcreate( 30 (fname + ".h5").c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist_id); 31 H5Pclose(plist_id); 32 33 // rank local data and offset and global file dimensions 34 hsize_t file_count[3] = {static_cast<unsigned int>(Nz), 35 static_cast<unsigned int>(Ny), 36 static_cast<unsigned int>(Nx)}; 37 hsize_t file_offset[3] = { 38 static_cast<unsigned int>(coords[2]) * Nz, 39 static_cast<unsigned int>(coords[1]) * Ny, 40 static_cast<unsigned int>(coords[0]) * Nx, 41 }; 42 hsize_t file_dims[3] = { 43 static_cast<unsigned int>(dims[2]) * Nz, 44 static_cast<unsigned int>(dims[1]) * Ny, 45 static_cast<unsigned int>(dims[0]) * Nx, 46 }; 47 if (0 == rank) { 48 std::cout << "HDF5 file size: " 49 << (file_dims[0] * file_dims[1] * file_dims[2] * 50 sizeof(float)) / 51 (1024.0 * 1024.0) 52 << " MB\n"; 53 } 54 55 // dataspace for dataset 56 hid_t fspace_id = H5Screate_simple(3, file_dims, NULL); 57 hid_t mspace_id = H5Screate_simple(3, file_count, NULL); 58 59 // chunked dataset 60 plist_id = H5Pcreate(H5P_DATASET_CREATE); 61 H5Pset_chunk(plist_id, 3, file_count); 62 hid_t dataset_id = H5Dcreate(file_id, 63 "u", 64 H5T_NATIVE_FLOAT, 65 fspace_id, 66 H5P_DEFAULT, 67 plist_id, 68 H5P_DEFAULT); 69 H5Pclose(plist_id); 70 H5Sclose(fspace_id); 71 72 // define the hyperslab 73 fspace_id = H5Dget_space(dataset_id); 74 H5Sselect_hyperslab( 75 fspace_id, H5S_SELECT_SET, file_offset, NULL, file_count, NULL); 76 77 // create property list for collective write 78 plist_id = H5Pcreate(H5P_DATASET_XFER); 79 H5Pset_dxpl_mpio(plist_id, H5FD_MPIO_COLLECTIVE); 80 81 // write the data 82 H5Dwrite( 83 dataset_id, H5T_NATIVE_FLOAT, mspace_id, fspace_id, plist_id, u.data()); 84 85 H5Dclose(dataset_id); 86 H5Sclose(fspace_id); 87 H5Sclose(mspace_id); 88 H5Pclose(plist_id); 89 H5Fclose(file_id); 90 H5close(); 91 92 // root rank writes XDMF wrapper 93 if (0 == rank) { 94 xdmf_wrapper(fname, file_dims); 95 } 96 } 97 98 #endif /* HDF5IO_H_XUXZ7EB0 */