From 460ff4bd049dd0849a393aa560e5b0b8bfb265e2 Mon Sep 17 00:00:00 2001 From: pdziekan Date: Thu, 9 Feb 2023 15:21:33 +0100 Subject: [PATCH 1/4] recrod halo hdf5 --- libmpdata++/output/hdf5.hpp | 59 ++++++++++++++++++++++++++++++------- 1 file changed, 49 insertions(+), 10 deletions(-) diff --git a/libmpdata++/output/hdf5.hpp b/libmpdata++/output/hdf5.hpp index 0391e729..0bdb3e9c 100644 --- a/libmpdata++/output/hdf5.hpp +++ b/libmpdata++/output/hdf5.hpp @@ -55,10 +55,10 @@ namespace libmpdataxx H5::PredType::NATIVE_FLOAT, flttype_output = H5::PredType::NATIVE_FLOAT; // using floats not to waste disk space - blitz::TinyVector cshape, shape, chunk, srfcshape, srfcchunk, offst; + blitz::TinyVector cshape, shape, chunk, srfcshape, srfcchunk, offst, shape_h, chunk_h, offst_h; H5::DSetCreatPropList params; - H5::DataSpace sspace, cspace, srfcspace; + H5::DataSpace sspace, cspace, srfcspace, sspace_h; #if defined(USE_MPI) hid_t fapl_id; #endif @@ -95,21 +95,28 @@ namespace libmpdataxx { // creating the dimensions // x,y,z - offst = 0; + offst = 0; + offst_h = 0; for (int d = 0; d < parent_t::n_dims; ++d) - shape[d] = this->mem->distmem.grid_size[d]; + { + shape[d] = this->mem->distmem.grid_size[d]; + shape_h[d] = this->mem->distmem.grid_size[d] + 2 * this->halo; + } - chunk = shape; + chunk = shape; + chunk_h = shape_h; // there is one more coordinate than cell index in each dimension cshape = shape + 1; srfcshape = shape; *(srfcshape.end()-1) = 1; - sspace = H5::DataSpace(parent_t::n_dims, shape.data()); + + sspace = H5::DataSpace(parent_t::n_dims, shape.data()); + sspace_h = H5::DataSpace(parent_t::n_dims, shape_h.data()); srfcspace = H5::DataSpace(parent_t::n_dims, srfcshape.data()); - cspace = H5::DataSpace(parent_t::n_dims, cshape.data()); + cspace = H5::DataSpace(parent_t::n_dims, cshape.data()); #if defined(USE_MPI) if (this->mem->distmem.size() > 1) @@ -120,10 +127,11 @@ namespace libmpdataxx if (this->mem->distmem.rank() == this->mem->distmem.size() - 1) cshape[0] += 1; - offst[0] = this->mem->grid_size[0].first(); + offst[0] = this->mem->grid_size[0].first(); // TODO: same for offst_h // chunk size has to be common to all processes ! // TODO: something better ? + // TODO: same for chunk_h chunk[0] = ( (typename solver_t::real_t) (this->mem->distmem.grid_size[0])) / this->mem->distmem.size() + 0.5 ; } #endif @@ -135,6 +143,7 @@ namespace libmpdataxx *(srfcchunk.end()-1) = 1; params.setChunk(parent_t::n_dims, chunk.data()); + #if !defined(USE_MPI) params.setDeflate(5); // TODO: move such constant to the header #endif @@ -203,10 +212,10 @@ namespace libmpdataxx } } - std::string base_name() + std::string base_name(const std::string &name = "timestep") { std::stringstream ss; - ss << "timestep" << std::setw(10) << std::setfill('0') << this->timestep; + ss << name << std::setw(10) << std::setfill('0') << this->timestep; return ss.str(); } @@ -216,6 +225,12 @@ namespace libmpdataxx return base_name() + ".h5"; } + std::string hdf_name(const std::string &base_name) + { + // TODO: add option of .nc extension for Paraview sake ? + return base_name + ".h5"; + } + void record_all() { // in concurrent setup only the first solver does output @@ -362,6 +377,30 @@ namespace libmpdataxx record_dsc_helper(aux, arr); } + // for array + halo + void record_aux_halo_hlpr(const std::string &name, const typename solver_t::arr_t &arr, H5::H5File hdf) + { + assert(this->rank == 0); + + params.setChunk(parent_t::n_dims, chunk_h.data()); + + auto aux = hdf.createDataSet( + name, + flttype_output, + sspace_h, + params + ); + + // revert to default chunk + params.setChunk(parent_t::n_dims, chunk.data()); + + record_dsc_helper(aux, arr); + + auto space = aux.getSpace(); + space.selectHyperslab(H5S_SELECT_SET, shape_h.data(), offst_h.data()); + aux.write(arr.data(), flttype_solver, H5::DataSpace(parent_t::n_dims, shape_h.data()), space, dxpl_id); + } + void record_scalar_hlpr(const std::string &name, const std::string &group_name, typename solver_t::real_t data, H5::H5File hdf) { assert(this->rank == 0); From 3cdcc484827c819c3c1aa66df7ae587ff0935e3f Mon Sep 17 00:00:00 2001 From: Piotr Dziekan Date: Fri, 10 Feb 2023 11:24:43 +0100 Subject: [PATCH 2/4] hdf5 recorc halo: fix for MPI --- libmpdata++/output/hdf5.hpp | 13 ++++++++++--- 1 file changed, 10 insertions(+), 3 deletions(-) diff --git a/libmpdata++/output/hdf5.hpp b/libmpdata++/output/hdf5.hpp index 0bdb3e9c..74a6908e 100644 --- a/libmpdata++/output/hdf5.hpp +++ b/libmpdata++/output/hdf5.hpp @@ -124,15 +124,22 @@ namespace libmpdataxx shape[0] = this->mem->grid_size[0].length(); cshape[0] = this->mem->grid_size[0].length(); + shape_h[0] = + this->mem->distmem.rank() == 0 || this->mem->distmem.rank() == this->mem->distmem.size()-1 ? + this->mem->grid_size[0].length() + this->halo : + this->mem->grid_size[0].length(); + + if (this->mem->distmem.rank() == this->mem->distmem.size() - 1) cshape[0] += 1; - offst[0] = this->mem->grid_size[0].first(); // TODO: same for offst_h + offst[0] = this->mem->grid_size[0].first(); + offst_h[0] = this->mem->distmem.rank() == 0 ? 0 : this->mem->grid_size[0].first() + this->halo; // chunk size has to be common to all processes ! // TODO: something better ? - // TODO: same for chunk_h - chunk[0] = ( (typename solver_t::real_t) (this->mem->distmem.grid_size[0])) / this->mem->distmem.size() + 0.5 ; + chunk[0] = ( (typename solver_t::real_t) (this->mem->distmem.grid_size[0])) / this->mem->distmem.size() + 0.5 ; + chunk_h[0] = 1;//chunk[0]; } #endif From 8b21e920d87a81041afa18367b7ea54946338db9 Mon Sep 17 00:00:00 2001 From: Piotr Dziekan Date: Fri, 10 Feb 2023 11:29:24 +0100 Subject: [PATCH 3/4] hdf5: some comments --- libmpdata++/output/hdf5.hpp | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/libmpdata++/output/hdf5.hpp b/libmpdata++/output/hdf5.hpp index 74a6908e..02e3a742 100644 --- a/libmpdata++/output/hdf5.hpp +++ b/libmpdata++/output/hdf5.hpp @@ -137,7 +137,7 @@ namespace libmpdataxx offst_h[0] = this->mem->distmem.rank() == 0 ? 0 : this->mem->grid_size[0].first() + this->halo; // chunk size has to be common to all processes ! - // TODO: something better ? + // TODO: set to 1? Test performance... chunk[0] = ( (typename solver_t::real_t) (this->mem->distmem.grid_size[0])) / this->mem->distmem.size() + 0.5 ; chunk_h[0] = 1;//chunk[0]; } @@ -152,7 +152,8 @@ namespace libmpdataxx params.setChunk(parent_t::n_dims, chunk.data()); #if !defined(USE_MPI) - params.setDeflate(5); // TODO: move such constant to the header + params.setDeflate(5); // TODO: move such constant to the header + // TODO2: why not deflate without MPI? #endif // creating variables From 484eaf3caa18fbacd899372eca5babf99c1ed8ba Mon Sep 17 00:00:00 2001 From: Piotr Dziekan Date: Fri, 10 Feb 2023 14:50:51 +0100 Subject: [PATCH 4/4] hdf5 halo fix for mpi --- libmpdata++/output/hdf5.hpp | 31 ++++++++++++++++++------------- 1 file changed, 18 insertions(+), 13 deletions(-) diff --git a/libmpdata++/output/hdf5.hpp b/libmpdata++/output/hdf5.hpp index 02e3a742..2d6bf8eb 100644 --- a/libmpdata++/output/hdf5.hpp +++ b/libmpdata++/output/hdf5.hpp @@ -55,10 +55,10 @@ namespace libmpdataxx H5::PredType::NATIVE_FLOAT, flttype_output = H5::PredType::NATIVE_FLOAT; // using floats not to waste disk space - blitz::TinyVector cshape, shape, chunk, srfcshape, srfcchunk, offst, shape_h, chunk_h, offst_h; + blitz::TinyVector cshape, shape, chunk, srfcshape, srfcchunk, offst, shape_h, chunk_h, offst_h, shape_mem_h, offst_mem_h; H5::DSetCreatPropList params; - H5::DataSpace sspace, cspace, srfcspace, sspace_h; + H5::DataSpace sspace, cspace, srfcspace, sspace_h, sspace_mem_h; #if defined(USE_MPI) hid_t fapl_id; #endif @@ -97,11 +97,13 @@ namespace libmpdataxx // x,y,z offst = 0; offst_h = 0; + offst_mem_h = 0; for (int d = 0; d < parent_t::n_dims; ++d) { - shape[d] = this->mem->distmem.grid_size[d]; - shape_h[d] = this->mem->distmem.grid_size[d] + 2 * this->halo; + shape[d] = this->mem->distmem.grid_size[d]; // shape of arrays stored in file + shape_h[d] = this->mem->distmem.grid_size[d] + 2 * this->halo; // shape of arrays with halos stored in files + shape_mem_h[d] = this->mem->grid_size[d].length() + 2 * this->halo; // shape of the array with halo stored in memory of given MPI rank } chunk = shape; @@ -113,10 +115,11 @@ namespace libmpdataxx srfcshape = shape; *(srfcshape.end()-1) = 1; - sspace = H5::DataSpace(parent_t::n_dims, shape.data()); - sspace_h = H5::DataSpace(parent_t::n_dims, shape_h.data()); - srfcspace = H5::DataSpace(parent_t::n_dims, srfcshape.data()); - cspace = H5::DataSpace(parent_t::n_dims, cshape.data()); + sspace = H5::DataSpace(parent_t::n_dims, shape.data()); + sspace_h = H5::DataSpace(parent_t::n_dims, shape_h.data()); + sspace_mem_h = H5::DataSpace(parent_t::n_dims, shape_mem_h.data()); + srfcspace = H5::DataSpace(parent_t::n_dims, srfcshape.data()); + cspace = H5::DataSpace(parent_t::n_dims, cshape.data()); #if defined(USE_MPI) if (this->mem->distmem.size() > 1) @@ -133,8 +136,10 @@ namespace libmpdataxx if (this->mem->distmem.rank() == this->mem->distmem.size() - 1) cshape[0] += 1; - offst[0] = this->mem->grid_size[0].first(); - offst_h[0] = this->mem->distmem.rank() == 0 ? 0 : this->mem->grid_size[0].first() + this->halo; + offst[0] = this->mem->grid_size[0].first(); + offst_h[0] = this->mem->distmem.rank() == 0 ? 0 : this->mem->grid_size[0].first() + this->halo; + if (this->mem->distmem.rank() > 0) + offst_mem_h[0] = this->halo; // chunk size has to be common to all processes ! // TODO: set to 1? Test performance... @@ -402,11 +407,11 @@ namespace libmpdataxx // revert to default chunk params.setChunk(parent_t::n_dims, chunk.data()); - record_dsc_helper(aux, arr); - auto space = aux.getSpace(); space.selectHyperslab(H5S_SELECT_SET, shape_h.data(), offst_h.data()); - aux.write(arr.data(), flttype_solver, H5::DataSpace(parent_t::n_dims, shape_h.data()), space, dxpl_id); + sspace_mem_h.selectHyperslab(H5S_SELECT_SET, shape_h.data(), offst_mem_h.data()); + + aux.write(arr.data(), flttype_solver, sspace_mem_h, space, dxpl_id); } void record_scalar_hlpr(const std::string &name, const std::string &group_name, typename solver_t::real_t data, H5::H5File hdf)