diff --git a/src/include/DiagEnergies.inl b/src/include/DiagEnergies.inl index bd139746e..7a81a0ecd 100644 --- a/src/include/DiagEnergies.inl +++ b/src/include/DiagEnergies.inl @@ -41,7 +41,7 @@ inline void DiagEnergies::operator()(Mparticles& mprts, MfieldsState& mflds) if (rank_ == 0) { assert(file_); - fprintf(file_.get(), "%g", grid.timestep() * grid.dt); + fprintf(file_.get(), "%g", grid.time()); } write_one(ef_, mprts, mflds); diff --git a/src/include/grid.hxx b/src/include/grid.hxx index d3e3c60e1..ca0881c11 100644 --- a/src/include/grid.hxx +++ b/src/include/grid.hxx @@ -119,6 +119,8 @@ struct Grid_ int timestep() const { return timestep_; } + real_t time() const { return timestep_ * dt; } + template void Foreach_3d(int l, int r, FUNC F) const { diff --git a/src/include/writer_adios2.hxx b/src/include/writer_adios2.hxx index acd4e53bd..9eaf69aaa 100644 --- a/src/include/writer_adios2.hxx +++ b/src/include/writer_adios2.hxx @@ -72,7 +72,7 @@ public: void begin_step(const Grid_t& grid) { int step = grid.timestep(); - double time = grid.timestep() * grid.dt; + double time = grid.time(); int len = dir_.size() + pfx_.size() + 20; char filename[len]; diff --git a/src/include/writer_mrc.hxx b/src/include/writer_mrc.hxx index b6c5c6dca..f29109514 100644 --- a/src/include/writer_mrc.hxx +++ b/src/include/writer_mrc.hxx @@ -29,13 +29,13 @@ public: void begin_step(const Grid_t& grid) { - mrc_io_open(io_.get(), "w", grid.timestep(), grid.timestep() * grid.dt); + mrc_io_open(io_.get(), "w", grid.timestep(), grid.time()); // save some basic info about the run in the output file struct mrc_obj* obj = mrc_obj_create(mrc_io_comm(io_.get())); mrc_obj_set_name(obj, "psc"); mrc_obj_dict_add_int(obj, "timestep", grid.timestep()); - mrc_obj_dict_add_float(obj, "time", grid.timestep() * grid.dt); + mrc_obj_dict_add_float(obj, "time", grid.time()); mrc_obj_dict_add_float(obj, "cc", grid.norm.cc); mrc_obj_dict_add_float(obj, "dt", grid.dt); mrc_obj_write(obj, io_.get()); diff --git a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx index 710e2210a..af1798b7a 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -4,7 +4,9 @@ #include #include +// FIXME it would be much nicer to use +#include "grid.hxx" #include "output_particles.hxx" #include "psc_particles_single.h" @@ -117,7 +119,7 @@ public: void operator()(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, size_t n_off, - size_t n_total, const particles_type& arr, int timestep) + size_t n_total, const particles_type& arr, const Grid_t& grid) { int ierr; @@ -158,13 +160,13 @@ public: // This is probably a problem with other outputs, too int slen = strlen(params_.data_dir) + strlen(params_.basename) + 15; char filename[slen]; - if (timestep > 999999999) { + if (grid.timestep() > 999999999) { LOG_WARN("%s step index exceeds digit limit in file name. Data still " "written, but file name is truncated.\n", params_.basename); } snprintf(filename, slen, "%s/%s.%09d.h5", params_.data_dir, - params_.basename, timestep); + params_.basename, grid.timestep()); hid_t file = H5Fcreate(filename, H5F_ACC_TRUNC, H5P_DEFAULT, plist); H5_CHK(file); @@ -196,6 +198,9 @@ public: #endif prof_stop(pr_C); + write_time(grid.time(), file, dxpl); + write_domain(grid.domain, file, dxpl); + prof_start(pr_D); write_idx(gidx_begin, gidx_end, group, dxpl); prof_stop(pr_D); @@ -215,6 +220,94 @@ public: } private: + void write_time(double time, hid_t group, hid_t dxpl) + { + herr_t ierr; + + hid_t filespace = H5Screate(H5S_SCALAR); + H5_CHK(filespace); + hid_t memspace; + + int mpi_rank; + MPI_Comm_rank(comm_, &mpi_rank); + if (mpi_rank == 0) { + memspace = H5Screate(H5S_SCALAR); + H5_CHK(memspace); + } else { + memspace = H5Screate(H5S_NULL); + H5Sselect_none(memspace); + H5Sselect_none(filespace); + } + + hid_t dset = H5Dcreate(group, "time", H5T_NATIVE_DOUBLE, filespace, + H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(dset); + ierr = H5Dwrite(dset, H5T_NATIVE_DOUBLE, memspace, filespace, dxpl, &time); + CE; + ierr = H5Dclose(dset); + CE; + + ierr = H5Sclose(filespace); + CE; + ierr = H5Sclose(memspace); + CE; + } + +private: + void write_domain(const Grid_t::Domain& domain, hid_t group, hid_t dxpl) + { + herr_t ierr; + + hsize_t data_len = 3; // 3 spatial dimensions + hid_t filespace = H5Screate_simple(1, &data_len, NULL); + H5_CHK(filespace); + hid_t memspace; + + assert(sizeof(size_t) == sizeof(hsize_t)); + int mpi_rank; + MPI_Comm_rank(comm_, &mpi_rank); + if (mpi_rank == 0) { + memspace = H5Screate_simple(1, &data_len, NULL); + H5_CHK(memspace); + } else { + memspace = H5Screate(H5S_NULL); + H5Sselect_none(memspace); + H5Sselect_none(filespace); + } + + hid_t dset = H5Dcreate(group, "gdims", H5T_NATIVE_INT32, filespace, + H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(dset); + ierr = + H5Dwrite(dset, H5T_NATIVE_INT32, memspace, filespace, dxpl, domain.gdims); + CE; + ierr = H5Dclose(dset); + CE; + + dset = H5Dcreate(group, "corner", H5T_NATIVE_DOUBLE, filespace, H5P_DEFAULT, + H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(dset); + ierr = H5Dwrite(dset, H5T_NATIVE_DOUBLE, memspace, filespace, dxpl, + domain.corner); + CE; + ierr = H5Dclose(dset); + CE; + + dset = H5Dcreate(group, "length", H5T_NATIVE_DOUBLE, filespace, H5P_DEFAULT, + H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(dset); + ierr = H5Dwrite(dset, H5T_NATIVE_DOUBLE, memspace, filespace, dxpl, + domain.length); + CE; + ierr = H5Dclose(dset); + CE; + + ierr = H5Sclose(filespace); + CE; + ierr = H5Sclose(memspace); + CE; + } + void write_idx(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, hid_t group, hid_t dxpl) @@ -519,7 +612,7 @@ struct OutputParticlesHdf5 void operator()(Mparticles& mprts, OutputParticlesWriterHDF5& writer) { - const auto& grid = mprts.grid(); + const Grid_t& grid = mprts.grid(); herr_t ierr; static int pr_A, pr_B, pr_C; @@ -656,7 +749,7 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - writer(gidx_begin, gidx_end, n_off, n_total, arr, grid.timestep()); + writer(gidx_begin, gidx_end, n_off, n_total, arr, grid); } private: diff --git a/src/psc_bgk.cxx b/src/psc_bgk.cxx index da9fdea3d..b58190eb6 100644 --- a/src/psc_bgk.cxx +++ b/src/psc_bgk.cxx @@ -144,7 +144,7 @@ void writeGT(const GT& gt, const Grid_t& grid, const std::string& name, { WriterMRC writer; writer.open(name); - writer.begin_step(grid.timestep(), grid.timestep() * grid.dt); + writer.begin_step(grid.timestep(), grid.time()); writer.write(gt, grid, name, compNames); writer.end_step(); writer.close(); diff --git a/src/psc_radiation.cxx b/src/psc_radiation.cxx index 369b76a94..1f94a27f2 100644 --- a/src/psc_radiation.cxx +++ b/src/psc_radiation.cxx @@ -235,7 +235,7 @@ void run() Marder marder(grid, marder_diffusion, marder_loop, marder_dump); auto lf_ext_current = [&](const Grid_t& grid, MfieldsState& mflds) { - double time = grid.timestep() * grid.dt; + double time = grid.time(); auto& gdims = grid.domain.gdims; for (int p = 0; p < mflds.n_patches(); ++p) { auto& patch = grid.patches[p];