From 82acae136c81132e440db5cfec91c704eb5586c0 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Fri, 23 Aug 2024 07:39:11 -0400 Subject: [PATCH 01/49] build: add .venv to .gitignore --- .gitignore | 1 + 1 file changed, 1 insertion(+) diff --git a/.gitignore b/.gitignore index f8824fa9f4..0c012539a6 100644 --- a/.gitignore +++ b/.gitignore @@ -101,3 +101,4 @@ spack-build-out.txt spack-configure-args.txt .vscode/ **/psc*params.txt +.venv/ \ No newline at end of file From fe283bacb23fb7f2be6ea00a336ba3c296b3ff08 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 29 Aug 2024 11:06:46 -0400 Subject: [PATCH 02/49] build: fix various warnings --- src/include/binary_collision.hxx | 2 +- src/include/ddc_particles.hxx | 4 ++-- src/include/writer_adios2.hxx | 15 +++++++++------ src/libpsc/psc_balance/psc_balance_impl.hxx | 2 +- .../output_particles_ascii_impl.hxx | 7 ++++--- .../output_particles_hdf5_impl.hxx | 7 ++++--- 6 files changed, 21 insertions(+), 16 deletions(-) diff --git a/src/include/binary_collision.hxx b/src/include/binary_collision.hxx index 965b3c728c..d8e995b221 100644 --- a/src/include/binary_collision.hxx +++ b/src/include/binary_collision.hxx @@ -20,7 +20,7 @@ struct RngC { real_t ran; do { - ran = real_t(random()) / RAND_MAX; + ran = real_t(random()) / real_t(RAND_MAX); } while (ran == real_t(0.f)); return ran; diff --git a/src/include/ddc_particles.hxx b/src/include/ddc_particles.hxx index 7f02f6b886..7f3025d4ac 100644 --- a/src/include/ddc_particles.hxx +++ b/src/include/ddc_particles.hxx @@ -246,7 +246,7 @@ inline ddc_particles::ddc_particles(const Grid_t& grid) n_recv_ranks = 0; for (int r = 0; r < size; r++) { if (info[r].n_recv_entries) { - MPI_Irecv(info[r].recv_entry.data(), + MPI_Irecv((int*)info[r].recv_entry.data(), sizeof(drecv_entry) / sizeof(int) * info[r].n_recv_entries, MPI_INT, r, 111, comm, &recv_reqs_[n_recv_ranks++]); } @@ -255,7 +255,7 @@ inline ddc_particles::ddc_particles(const Grid_t& grid) n_send_ranks = 0; for (int r = 0; r < size; r++) { if (info[r].n_send_entries) { - MPI_Isend(info[r].send_entry.data(), + MPI_Isend((int*)info[r].send_entry.data(), sizeof(dsend_entry) / sizeof(int) * info[r].n_send_entries, MPI_INT, r, 111, comm, &send_reqs_[n_send_ranks++]); } diff --git a/src/include/writer_adios2.hxx b/src/include/writer_adios2.hxx index 840e4c10cd..acd4e53bd3 100644 --- a/src/include/writer_adios2.hxx +++ b/src/include/writer_adios2.hxx @@ -74,8 +74,9 @@ public: int step = grid.timestep(); double time = grid.timestep() * grid.dt; - char filename[dir_.size() + pfx_.size() + 20]; - sprintf(filename, "%s/%s.%09d.bp", dir_.c_str(), pfx_.c_str(), step); + int len = dir_.size() + pfx_.size() + 20; + char filename[len]; + snprintf(filename, len, "%s/%s.%09d.bp", dir_.c_str(), pfx_.c_str(), step); file_ = io_.open(filename, kg::io::Mode::Write, comm_, pfx_); file_.beginStep(kg::io::StepMode::Append); file_.put("step", step); @@ -157,9 +158,9 @@ public: Double3 length = grid.domain.length; Double3 corner = grid.domain.corner; - auto write_func = [this, step, time, h_expr = move(h_expr), name, + auto write_func = [this, step, time, h_expr = std::move(h_expr), name, comp_names, ldims, gdims, length, corner, - patch_off = move(patch_off)]() { + patch_off = std::move(patch_off)]() { // std::this_thread::sleep_for(std::chrono::milliseconds(1000)); prof_start(pr_thread); @@ -170,8 +171,10 @@ public: Int3 ib = {-(im[0] - ldims[0]) / 2, -(im[1] - ldims[1]) / 2, -(im[2] - ldims[2]) / 2}; - char filename[dir_.size() + pfx_.size() + 20]; - sprintf(filename, "%s/%s.%09d.bp", dir_.c_str(), pfx_.c_str(), step); + int len = dir_.size() + pfx_.size() + 20; + char filename[len]; + snprintf(filename, len, "%s/%s.%09d.bp", dir_.c_str(), pfx_.c_str(), + step); { auto launch = kg::io::Mode::Blocking; diff --git a/src/libpsc/psc_balance/psc_balance_impl.hxx b/src/libpsc/psc_balance/psc_balance_impl.hxx index 6a29b5f623..97084c4643 100644 --- a/src/libpsc/psc_balance/psc_balance_impl.hxx +++ b/src/libpsc/psc_balance/psc_balance_impl.hxx @@ -204,7 +204,7 @@ inline void write_loads(const input& input, int timestep) { char s[20]; - sprintf(s, "loads2-%06d.asc", timestep); + snprintf(s, 20, "loads2-%06d.asc", timestep); FILE* f = fopen(s, "w"); int gp = 0; diff --git a/src/libpsc/psc_output_particles/output_particles_ascii_impl.hxx b/src/libpsc/psc_output_particles/output_particles_ascii_impl.hxx index ddfe2f2ae8..406825fab3 100644 --- a/src/libpsc/psc_output_particles/output_particles_ascii_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_ascii_impl.hxx @@ -22,9 +22,10 @@ struct OutputParticlesAscii int rank; MPI_Comm_rank(comm_, &rank); - char filename[strlen(data_dir) + strlen(basename) + 21]; - sprintf(filename, "%s/%s.%06d_p%06d.asc", data_dir, basename, - grid.timestep(), rank); + int slen = strlen(data_dir) + strlen(basename) + 21; + char filename[slen]; + snprintf(filename, slen, "%s/%s.%06d_p%06d.asc", data_dir, basename, + grid.timestep(), rank); FILE* file = fopen(filename, "w"); auto accessor = mprts.accessor(); 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 443aa7aafc..9836d837d2 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -653,9 +653,10 @@ public: return; } - char filename[strlen(prm_.data_dir) + strlen(prm_.basename) + 20]; - sprintf(filename, "%s/%s.%06d_p%06d.h5", prm_.data_dir, prm_.basename, - grid.timestep(), 0); + int slen = strlen(prm_.data_dir) + strlen(prm_.basename) + 20; + char filename[slen]; + snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", prm_.data_dir, + prm_.basename, grid.timestep(), 0); detail::OutputParticlesHdf5 impl{grid, lo_, hi_, wdims_, prt_type_}; From 4cfd81881f7da15b9038b614f477732f5d1c2b40 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Fri, 23 Aug 2024 09:39:35 -0400 Subject: [PATCH 03/49] outp: psc_flatfoil_small: output particles every 4 steps; for half of the domain only --- src/psc_flatfoil_yz.cxx | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/psc_flatfoil_yz.cxx b/src/psc_flatfoil_yz.cxx index c288b35931..a853acc059 100644 --- a/src/psc_flatfoil_yz.cxx +++ b/src/psc_flatfoil_yz.cxx @@ -537,9 +537,14 @@ void run() // -- output particles OutputParticlesParams outp_params{}; +#if CASE == CASE_2D_SMALL + outp_params.every_step = 4; +#else outp_params.every_step = -400; +#endif outp_params.data_dir = "."; outp_params.basename = "prt"; + outp_params.lo = {0, 0, grid.domain.gdims[2] / 2}; OutputParticles outp{grid, outp_params}; int oute_interval = -100; From c4ba1b481a5de54d1cb68a0cc0fbb48f9161b6c5 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Mon, 10 Jun 2024 23:04:01 -0400 Subject: [PATCH 04/49] outp: put particles at x=1 and x=40 So that they are in different cells. --- src/libpsc/tests/test_output_particles.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/libpsc/tests/test_output_particles.cxx b/src/libpsc/tests/test_output_particles.cxx index 4fbe3d5fbe..7e8fb3325e 100644 --- a/src/libpsc/tests/test_output_particles.cxx +++ b/src/libpsc/tests/test_output_particles.cxx @@ -93,7 +93,7 @@ TYPED_TEST(OutputParticlesTest, Test1) { auto injector = mprts.injector(); injector[0]({{1., 0., 0.}, {}, 1., 0}); - injector[0]({{2., 0., 0.}, {}, 1., 1}); + injector[0]({{40., 0., 0.}, {}, 1., 1}); } auto params = OutputParticlesParams{}; From 7dbeddc3bf7de138e7aefa911d53ad1458b86649 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 07:07:19 -0400 Subject: [PATCH 05/49] outp/hdf5: use MPI_UNSIGNED_LONG to avoid warning --- .../output_particles_hdf5_impl.hxx | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) 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 9836d837d2..874fc058f0 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -227,8 +227,8 @@ struct OutputParticlesHdf5 assert(sizeof(size_t) == sizeof(unsigned long)); size_t n_total, n_off = 0; - MPI_Allreduce(&n_write, &n_total, 1, MPI_LONG, MPI_SUM, comm_); - MPI_Exscan(&n_write, &n_off, 1, MPI_LONG, MPI_SUM, comm_); + MPI_Allreduce(&n_write, &n_total, 1, MPI_UNSIGNED_LONG, MPI_SUM, comm_); + MPI_Exscan(&n_write, &n_off, 1, MPI_UNSIGNED_LONG, MPI_SUM, comm_); struct hdf5_prt* arr = (struct hdf5_prt*)malloc(n_write * sizeof(*arr)); @@ -470,8 +470,8 @@ struct OutputParticlesHdf5 } recv_buf[r] = (size_t*)malloc(2 * remote_sz[r] * sizeof(*recv_buf[r])); recv_buf_p[r] = recv_buf[r]; - MPI_Irecv(recv_buf[r], 2 * remote_sz[r], MPI_LONG, r, 0x4000, comm_, - &recv_req[r]); + MPI_Irecv(recv_buf[r], 2 * remote_sz[r], MPI_UNSIGNED_LONG, r, 0x4000, + comm_, &recv_req[r]); } MPI_Waitall(size, recv_req, MPI_STATUSES_IGNORE); free(recv_req); @@ -530,7 +530,7 @@ struct OutputParticlesHdf5 l_idx_p += 2 * sz; } - MPI_Send(l_idx, 2 * local_sz, MPI_LONG, 0, 0x4000, comm_); + MPI_Send(l_idx, 2 * local_sz, MPI_UNSIGNED_LONG, 0, 0x4000, comm_); free(l_idx); } From a2bb76b56c54b9700ce85a7c2ff0abbe6b928e55 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 07:08:35 -0400 Subject: [PATCH 06/49] outp/hdf5: compute wdims in impl class --- .../psc_output_particles/output_particles_hdf5_impl.hxx | 9 +++------ 1 file changed, 3 insertions(+), 6 deletions(-) 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 874fc058f0..1ab6c541d2 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -97,10 +97,10 @@ struct OutputParticlesHdf5 using Particles = typename Mparticles::Patch; OutputParticlesHdf5(const Grid_t& grid, const Int3& lo, const Int3& hi, - const Int3& wdims, hid_t prt_type) + hid_t prt_type) : lo_{lo}, hi_{hi}, - wdims_{wdims}, + wdims_{hi - lo}, kinds_{grid.kinds}, comm_{grid.comm()}, prt_type_{prt_type} @@ -641,7 +641,6 @@ public: assert(lo_[d] >= 0); assert(hi_[d] <= grid.domain.gdims[d]); } - wdims_ = hi_ - lo_; } template @@ -658,8 +657,7 @@ public: snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", prm_.data_dir, prm_.basename, grid.timestep(), 0); - detail::OutputParticlesHdf5 impl{grid, lo_, hi_, wdims_, - prt_type_}; + detail::OutputParticlesHdf5 impl{grid, lo_, hi_, prt_type_}; impl(mprts, filename, prm_); } @@ -700,5 +698,4 @@ private: detail::Hdf5ParticleType prt_type_; Int3 lo_; // dimensions of the subdomain we're actually writing Int3 hi_; - Int3 wdims_; }; From 5c9b1116ed5cb6335f5f966869328f5b8265a81b Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 07:08:55 -0400 Subject: [PATCH 07/49] outp/hdf5: use std::vector for off, map, idx --- .../output_particles_hdf5_impl.hxx | 16 ++++++---------- 1 file changed, 6 insertions(+), 10 deletions(-) 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 1ab6c541d2..5df7f76fcd 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -385,12 +385,11 @@ struct OutputParticlesHdf5 struct mrc_patch_info info; int nr_kinds = mprts.grid().kinds.size(); - int** off = (int**)malloc(mprts.n_patches() * sizeof(*off)); - int** map = (int**)malloc(mprts.n_patches() * sizeof(*off)); + std::vector off(mprts.n_patches()); + std::vector map(mprts.n_patches()); + std::vector idx(mprts.n_patches()); - count_sort(mprts, map, off); - - size_t** idx = (size_t**)malloc(mprts.n_patches() * sizeof(*idx)); + count_sort(mprts, map.data(), off.data()); size_t *gidx_begin = NULL, *gidx_end = NULL; if (rank == 0) { @@ -415,8 +414,8 @@ struct OutputParticlesHdf5 // find local particle and idx arrays size_t n_write, n_off, n_total; - hdf5_prt* arr = make_local_particle_array(mprts, map, off, idx, &n_write, - &n_off, &n_total); + hdf5_prt* arr = make_local_particle_array( + mprts, map.data(), off.data(), idx.data(), &n_write, &n_off, &n_total); prof_stop(pr_A); prof_start(pr_B); @@ -611,9 +610,6 @@ struct OutputParticlesHdf5 free(map[p]); free(idx[p]); } - free(off); - free(map); - free(idx); } private: From 19d5e9b6ddea2adb26e2460d84ac8533303175cf Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 09:04:01 -0400 Subject: [PATCH 08/49] outp/hdf5: pass off, map, idx as std::vector --- .../output_particles_hdf5_impl.hxx | 20 +++++++++++-------- 1 file changed, 12 insertions(+), 8 deletions(-) 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 5df7f76fcd..4b1630d8c7 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -143,7 +143,8 @@ struct OutputParticlesHdf5 // ---------------------------------------------------------------------- // count_sort - static void count_sort(Mparticles& mprts, int** off, int** map) + static void count_sort(Mparticles& mprts, std::vector& off, + std::vector& map) { int nr_kinds = mprts.grid().kinds.size(); @@ -199,9 +200,12 @@ struct OutputParticlesHdf5 // ---------------------------------------------------------------------- // make_local_particle_array - hdf5_prt* make_local_particle_array(Mparticles& mprts, int** off, int** map, - size_t** idx, size_t* p_n_write, - size_t* p_n_off, size_t* p_n_total) + hdf5_prt* make_local_particle_array(Mparticles& mprts, + const std::vector& off, + const std::vector& map, + std::vector& idx, + size_t* p_n_write, size_t* p_n_off, + size_t* p_n_total) { const auto& grid = mprts.grid(); int nr_kinds = grid.kinds.size(); @@ -241,7 +245,7 @@ struct OutputParticlesHdf5 const auto& patch = grid.patches[p]; int ilo[3], ihi[3], ld[3]; int sz = find_patch_bounds(grid.ldims, patch.off, ilo, ihi, ld); - idx[p] = (size_t*)malloc(2 * sz * sizeof(*idx)); + idx[p] = (size_t*)malloc(2 * sz * sizeof(size_t)); for (int jz = ilo[2]; jz < ihi[2]; jz++) { for (int jy = ilo[1]; jy < ihi[1]; jy++) { @@ -389,7 +393,7 @@ struct OutputParticlesHdf5 std::vector map(mprts.n_patches()); std::vector idx(mprts.n_patches()); - count_sort(mprts, map.data(), off.data()); + count_sort(mprts, map, off); size_t *gidx_begin = NULL, *gidx_end = NULL; if (rank == 0) { @@ -414,8 +418,8 @@ struct OutputParticlesHdf5 // find local particle and idx arrays size_t n_write, n_off, n_total; - hdf5_prt* arr = make_local_particle_array( - mprts, map.data(), off.data(), idx.data(), &n_write, &n_off, &n_total); + hdf5_prt* arr = make_local_particle_array(mprts, map, off, idx, &n_write, + &n_off, &n_total); prof_stop(pr_A); prof_start(pr_B); From 560443ab2d907da4012470424e4194ee1cc537ca Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 09:10:23 -0400 Subject: [PATCH 09/49] outp/hdf5: std::vector for off2 --- src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) 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 4b1630d8c7..c7db153b57 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -162,7 +162,7 @@ struct OutputParticlesHdf5 } // prefix sum to get offsets int o = 0; - int* off2 = (int*)malloc((nr_indices + 1) * sizeof(*off2)); + std::vector off2(nr_indices + 1); for (int si = 0; si <= nr_indices; si++) { int cnt = off[p][si]; off[p][si] = o; // this will be saved for later @@ -177,7 +177,6 @@ struct OutputParticlesHdf5 int si = get_sort_index(prts, prt); map[p][off2[si]++] = n++; } - free(off2); } } From 2bde15d3dfdfc58a325d3cca6fa7ef8ac83bd42c Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 21:48:13 -0400 Subject: [PATCH 10/49] outp/hdf5: fix swapped off, map args Fortunately they were swapped throughout, so very confusing but not an actual bug... --- .../psc_output_particles/output_particles_hdf5_impl.hxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) 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 c7db153b57..6ae6fb484a 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -392,7 +392,7 @@ struct OutputParticlesHdf5 std::vector map(mprts.n_patches()); std::vector idx(mprts.n_patches()); - count_sort(mprts, map, off); + count_sort(mprts, off, map); size_t *gidx_begin = NULL, *gidx_end = NULL; if (rank == 0) { @@ -417,7 +417,7 @@ struct OutputParticlesHdf5 // find local particle and idx arrays size_t n_write, n_off, n_total; - hdf5_prt* arr = make_local_particle_array(mprts, map, off, idx, &n_write, + hdf5_prt* arr = make_local_particle_array(mprts, off, map, idx, &n_write, &n_off, &n_total); prof_stop(pr_A); From e369068449c9757192a02a3b11306469a384f780 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 21:56:21 -0400 Subject: [PATCH 11/49] outp/hdf5: use std::vector> for off, map --- .../output_particles_hdf5_impl.hxx | 18 ++++++++---------- 1 file changed, 8 insertions(+), 10 deletions(-) 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 6ae6fb484a..01fef0cb0a 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -143,15 +143,15 @@ struct OutputParticlesHdf5 // ---------------------------------------------------------------------- // count_sort - static void count_sort(Mparticles& mprts, std::vector& off, - std::vector& map) + static void count_sort(Mparticles& mprts, std::vector>& off, + std::vector>& map) { int nr_kinds = mprts.grid().kinds.size(); for (int p = 0; p < mprts.n_patches(); p++) { const int* ldims = mprts.grid().ldims; int nr_indices = ldims[0] * ldims[1] * ldims[2] * nr_kinds; - off[p] = (int*)calloc(nr_indices + 1, sizeof(*off[p])); + off[p].resize(nr_indices + 1); auto&& prts = mprts[p]; unsigned int n_prts = prts.size(); @@ -171,7 +171,7 @@ struct OutputParticlesHdf5 } // sort a map only, not the actual particles - map[p] = (int*)malloc(n_prts * sizeof(*map[p])); + map[p].resize(n_prts); int n = 0; for (const auto& prt : prts) { int si = get_sort_index(prts, prt); @@ -200,8 +200,8 @@ struct OutputParticlesHdf5 // make_local_particle_array hdf5_prt* make_local_particle_array(Mparticles& mprts, - const std::vector& off, - const std::vector& map, + const std::vector>& off, + const std::vector>& map, std::vector& idx, size_t* p_n_write, size_t* p_n_off, size_t* p_n_total) @@ -388,8 +388,8 @@ struct OutputParticlesHdf5 struct mrc_patch_info info; int nr_kinds = mprts.grid().kinds.size(); - std::vector off(mprts.n_patches()); - std::vector map(mprts.n_patches()); + std::vector> off(mprts.n_patches()); + std::vector> map(mprts.n_patches()); std::vector idx(mprts.n_patches()); count_sort(mprts, off, map); @@ -609,8 +609,6 @@ struct OutputParticlesHdf5 free(gidx_end); for (int p = 0; p < mprts.n_patches(); p++) { - free(off[p]); - free(map[p]); free(idx[p]); } } From e02627b2b95f73cdaeba365a3fe004045d8d30e2 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 22:00:26 -0400 Subject: [PATCH 12/49] outp/hdf5: use std::vector for gidx_{begin,end} --- .../output_particles_hdf5_impl.hxx | 21 ++++++++----------- 1 file changed, 9 insertions(+), 12 deletions(-) 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 01fef0cb0a..87361d9b61 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -313,7 +313,8 @@ struct OutputParticlesHdf5 CE; } - void write_idx(size_t* gidx_begin, size_t* gidx_end, hid_t group, hid_t dxpl) + void write_idx(const std::vector& gidx_begin, + const std::vector& gidx_end, hid_t group, hid_t dxpl) { herr_t ierr; @@ -341,8 +342,8 @@ struct OutputParticlesHdf5 hid_t dset = H5Dcreate(group, "idx_begin", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(dset); - ierr = - H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, gidx_begin); + ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, + gidx_begin.data()); CE; ierr = H5Dclose(dset); CE; @@ -350,8 +351,8 @@ struct OutputParticlesHdf5 dset = H5Dcreate(group, "idx_end", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(dset); - ierr = - H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, gidx_end); + ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, + gidx_end.data()); CE; ierr = H5Dclose(dset); CE; @@ -394,13 +395,11 @@ struct OutputParticlesHdf5 count_sort(mprts, off, map); - size_t *gidx_begin = NULL, *gidx_end = NULL; + std::vector gidx_begin, gidx_end; if (rank == 0) { // alloc global idx array - gidx_begin = (size_t*)malloc(nr_kinds * wdims_[0] * wdims_[1] * - wdims_[2] * sizeof(*gidx_begin)); - gidx_end = (size_t*)malloc(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2] * - sizeof(*gidx_end)); + gidx_begin.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2]); + gidx_end.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2]); for (int jz = 0; jz < wdims_[2]; jz++) { for (int jy = 0; jy < wdims_[1]; jy++) { for (int jx = 0; jx < wdims_[0]; jx++) { @@ -605,8 +604,6 @@ struct OutputParticlesHdf5 prof_stop(pr_E); free(arr); - free(gidx_begin); - free(gidx_end); for (int p = 0; p < mprts.n_patches(); p++) { free(idx[p]); From 022666300e6df869783f387fddc5a2ff5bf1aace Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 6 Jun 2024 22:03:44 -0400 Subject: [PATCH 13/49] outp/hdf5: much easier to init to -1 now --- .../output_particles_hdf5_impl.hxx | 17 +++-------------- 1 file changed, 3 insertions(+), 14 deletions(-) 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 87361d9b61..6588c88c75 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -398,20 +398,9 @@ struct OutputParticlesHdf5 std::vector gidx_begin, gidx_end; if (rank == 0) { // alloc global idx array - gidx_begin.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2]); - gidx_end.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2]); - for (int jz = 0; jz < wdims_[2]; jz++) { - for (int jy = 0; jy < wdims_[1]; jy++) { - for (int jx = 0; jx < wdims_[0]; jx++) { - for (int kind = 0; kind < nr_kinds; kind++) { - int ii = - ((kind * wdims_[2] + jz) * wdims_[1] + jy) * wdims_[0] + jx; - gidx_begin[ii] = size_t(-1); - gidx_end[ii] = size_t(-1); - } - } - } - } + gidx_begin.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2], + size_t(-1)); + gidx_end.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2], size_t(-1)); } // find local particle and idx arrays From 57eca514bb36b04c03b9f541bbe25578511410ab Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Mon, 10 Jun 2024 23:01:54 -0400 Subject: [PATCH 14/49] outp/hdf5: use std::vector for idx --- .../output_particles_hdf5_impl.hxx | 12 ++++-------- 1 file changed, 4 insertions(+), 8 deletions(-) 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 6588c88c75..8db50ed384 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -202,7 +202,7 @@ struct OutputParticlesHdf5 hdf5_prt* make_local_particle_array(Mparticles& mprts, const std::vector>& off, const std::vector>& map, - std::vector& idx, + std::vector>& idx, size_t* p_n_write, size_t* p_n_off, size_t* p_n_total) { @@ -244,7 +244,7 @@ struct OutputParticlesHdf5 const auto& patch = grid.patches[p]; int ilo[3], ihi[3], ld[3]; int sz = find_patch_bounds(grid.ldims, patch.off, ilo, ihi, ld); - idx[p] = (size_t*)malloc(2 * sz * sizeof(size_t)); + idx[p].resize(2 * sz); for (int jz = ilo[2]; jz < ihi[2]; jz++) { for (int jy = ilo[1]; jy < ihi[1]; jy++) { @@ -391,7 +391,7 @@ struct OutputParticlesHdf5 std::vector> off(mprts.n_patches()); std::vector> map(mprts.n_patches()); - std::vector idx(mprts.n_patches()); + std::vector> idx(mprts.n_patches()); count_sort(mprts, off, map); @@ -516,7 +516,7 @@ struct OutputParticlesHdf5 const auto& patch = grid.patches[p]; int ilo[3], ihi[3], ld[3]; int sz = find_patch_bounds(grid.ldims, patch.off, ilo, ihi, ld); - memcpy(l_idx_p, idx[p], 2 * sz * sizeof(*l_idx_p)); + memcpy(l_idx_p, idx[p].data(), 2 * sz * sizeof(*l_idx_p)); l_idx_p += 2 * sz; } @@ -593,10 +593,6 @@ struct OutputParticlesHdf5 prof_stop(pr_E); free(arr); - - for (int p = 0; p < mprts.n_patches(); p++) { - free(idx[p]); - } } private: From 6b542b4848b44d089cb8942848e539f9d3cf0e1f Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Mon, 10 Jun 2024 23:04:06 -0400 Subject: [PATCH 15/49] outp/hdf5: consolidate sort_index() --- .../output_particles_hdf5_impl.hxx | 21 +++++++++++-------- 1 file changed, 12 insertions(+), 9 deletions(-) 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 8db50ed384..53e4ae4cb5 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -107,13 +107,15 @@ struct OutputParticlesHdf5 {} // ---------------------------------------------------------------------- - // get_cell_index + // get_sort_index // FIXME, lots of stuff here is pretty much duplicated from countsort2 - static inline int cell_index_3_to_1(const int* ldims, int j0, int j1, int j2) + static inline int sort_index(const int* ldims, int nr_kinds, Int3 idx3, + int kind) { - return ((j2)*ldims[1] + j1) * ldims[0] + j0; - } + return ((idx3[2] * ldims[1] + idx3[1]) * ldims[0] + idx3[0]) * nr_kinds + + kind; + }; template static inline int get_sort_index(Particles& prts, const Particle& prt) @@ -137,7 +139,7 @@ struct OutputParticlesHdf5 assert(j2 >= 0 && j2 < ldims[2]); assert(prt.kind < grid.kinds.size()); - return cell_index_3_to_1(ldims, j0, j1, j2) * grid.kinds.size() + prt.kind; + return sort_index(ldims, grid.kinds.size(), {j0, j1, j2}, prt.kind); } // ---------------------------------------------------------------------- @@ -219,8 +221,7 @@ struct OutputParticlesHdf5 for (int jy = ilo[1]; jy < ihi[1]; jy++) { for (int jx = ilo[0]; jx < ihi[0]; jx++) { for (int kind = 0; kind < nr_kinds; kind++) { - int si = - cell_index_3_to_1(grid.ldims, jx, jy, jz) * nr_kinds + kind; + int si = sort_index(grid.ldims, nr_kinds, {jx, jy, jz}, kind); n_write += off[p][si + 1] - off[p][si]; } } @@ -250,8 +251,7 @@ struct OutputParticlesHdf5 for (int jy = ilo[1]; jy < ihi[1]; jy++) { for (int jx = ilo[0]; jx < ihi[0]; jx++) { for (int kind = 0; kind < nr_kinds; kind++) { - int si = - cell_index_3_to_1(grid.ldims, jx, jy, jz) * nr_kinds + kind; + int si = sort_index(grid.ldims, nr_kinds, {jx, jy, jz}, kind); int jj = ((kind * ld[2] + jz - ilo[2]) * ld[1] + jy - ilo[1]) * ld[0] + jx - ilo[0]; @@ -402,6 +402,9 @@ struct OutputParticlesHdf5 size_t(-1)); gidx_end.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2], size_t(-1)); } + auto gidx_shape = gt::shape(wdims_[0], wdims_[1], wdims_[2], nr_kinds); + auto _gidx_begin = gt::adapt(gidx_begin.data(), gidx_shape); + auto _gidx_end = gt::adapt(gidx_end.data(), gidx_shape); // find local particle and idx arrays size_t n_write, n_off, n_total; From ddc12b5d3529f5d80f1c0d43eaeba9fff16068c6 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 00:03:33 -0400 Subject: [PATCH 16/49] outp/hdf5: pass just the data to write_idx() --- .../output_particles_hdf5_impl.hxx | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) 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 53e4ae4cb5..77f6a3ab99 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -313,8 +313,8 @@ struct OutputParticlesHdf5 CE; } - void write_idx(const std::vector& gidx_begin, - const std::vector& gidx_end, hid_t group, hid_t dxpl) + void write_idx(const size_t* gidx_begin, const size_t* gidx_end, hid_t group, + hid_t dxpl) { herr_t ierr; @@ -342,8 +342,8 @@ struct OutputParticlesHdf5 hid_t dset = H5Dcreate(group, "idx_begin", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(dset); - ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, - gidx_begin.data()); + ierr = + H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, gidx_begin); CE; ierr = H5Dclose(dset); CE; @@ -351,8 +351,8 @@ struct OutputParticlesHdf5 dset = H5Dcreate(group, "idx_end", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(dset); - ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, - gidx_end.data()); + ierr = + H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, gidx_end); CE; ierr = H5Dclose(dset); CE; @@ -581,7 +581,7 @@ struct OutputParticlesHdf5 prof_stop(pr_C); prof_start(pr_D); - write_idx(gidx_begin, gidx_end, groupp, dxpl); + write_idx(gidx_begin.data(), gidx_end.data(), groupp, dxpl); prof_stop(pr_D); prof_start(pr_E); From 074b01a2bccb39d02141e1c6a19af09281b532da Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 00:10:58 -0400 Subject: [PATCH 17/49] outp/hdf5: use 4-d xtensor arrays for gidx_{begin,end} --- .../output_particles_hdf5_impl.hxx | 23 +++++++------------ 1 file changed, 8 insertions(+), 15 deletions(-) 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 77f6a3ab99..e7d14b06e7 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -395,16 +395,13 @@ struct OutputParticlesHdf5 count_sort(mprts, off, map); - std::vector gidx_begin, gidx_end; + gt::gtensor gidx_begin, gidx_end; if (rank == 0) { // alloc global idx array - gidx_begin.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2], - size_t(-1)); - gidx_end.resize(nr_kinds * wdims_[0] * wdims_[1] * wdims_[2], size_t(-1)); + gidx_begin = + gt::empty({wdims_[0], wdims_[1], wdims_[2], nr_kinds}); + gidx_end = gt::empty({wdims_[0], wdims_[1], wdims_[2], nr_kinds}); } - auto gidx_shape = gt::shape(wdims_[0], wdims_[1], wdims_[2], nr_kinds); - auto _gidx_begin = gt::adapt(gidx_begin.data(), gidx_shape); - auto _gidx_end = gt::adapt(gidx_end.data(), gidx_shape); // find local particle and idx arrays size_t n_write, n_off, n_total; @@ -443,10 +440,8 @@ struct OutputParticlesHdf5 int jj = ((kind * ld[2] + jz - ilo[2]) * ld[1] + jy - ilo[1]) * ld[0] + jx - ilo[0]; - int ii = - ((kind * wdims_[2] + iz) * wdims_[1] + iy) * wdims_[0] + ix; - gidx_begin[ii] = idx[p][jj]; - gidx_end[ii] = idx[p][jj + sz]; + gidx_begin(ix, iy, iz, kind) = idx[p][jj]; + gidx_end(ix, iy, iz, kind) = idx[p][jj + sz]; } } } @@ -487,10 +482,8 @@ struct OutputParticlesHdf5 int jj = ((kind * ld[2] + jz - ilo[2]) * ld[1] + jy - ilo[1]) * ld[0] + jx - ilo[0]; - int ii = - ((kind * wdims_[2] + iz) * wdims_[1] + iy) * wdims_[0] + ix; - gidx_begin[ii] = recv_buf_p[info.rank][jj]; - gidx_end[ii] = recv_buf_p[info.rank][jj + sz]; + gidx_begin(ix, iy, iz, kind) = recv_buf_p[info.rank][jj]; + gidx_end(ix, iy, iz, kind) = recv_buf_p[info.rank][jj + sz]; } } } From ac9eebe5bbff4539cb8c8850c8aab5c2653142b7 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 22:43:08 -0400 Subject: [PATCH 18/49] outp/hdf5: access idx as 5-d through gtensor --- .../output_particles_hdf5_impl.hxx | 21 ++++++++++--------- 1 file changed, 11 insertions(+), 10 deletions(-) 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 e7d14b06e7..185687ef93 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -246,17 +246,18 @@ struct OutputParticlesHdf5 int ilo[3], ihi[3], ld[3]; int sz = find_patch_bounds(grid.ldims, patch.off, ilo, ihi, ld); idx[p].resize(2 * sz); + auto _idx = gt::adapt<5>(idx[p].data(), + gt::shape(ld[0], ld[1], ld[2], nr_kinds, 2)); for (int jz = ilo[2]; jz < ihi[2]; jz++) { for (int jy = ilo[1]; jy < ihi[1]; jy++) { for (int jx = ilo[0]; jx < ihi[0]; jx++) { for (int kind = 0; kind < nr_kinds; kind++) { int si = sort_index(grid.ldims, nr_kinds, {jx, jy, jz}, kind); - int jj = - ((kind * ld[2] + jz - ilo[2]) * ld[1] + jy - ilo[1]) * ld[0] + - jx - ilo[0]; - idx[p][jj] = nn + n_off; - idx[p][jj + sz] = nn + n_off + off[p][si + 1] - off[p][si]; + _idx(jx - ilo[0], jy - ilo[1], jz - ilo[2], kind, 0) = + nn + n_off; + _idx(jx - ilo[0], jy - ilo[1], jz - ilo[2], kind, 1) = + nn + n_off + off[p][si + 1] - off[p][si]; for (int n = off[p][si]; n < off[p][si + 1]; n++, nn++) { auto prt = prts[map[p][n]]; arr[nn].x = prt.x()[0] + patch.xb[0]; @@ -430,6 +431,9 @@ struct OutputParticlesHdf5 int ilo[3], ihi[3], ld[3]; int sz = find_patch_bounds(info.ldims, info.off, ilo, ihi, ld); + auto _idx = gt::adapt<5>(idx[p].data(), + gt::shape(ld[0], ld[1], ld[2], nr_kinds, 2)); + for (int jz = ilo[2]; jz < ihi[2]; jz++) { for (int jy = ilo[1]; jy < ihi[1]; jy++) { for (int jx = ilo[0]; jx < ihi[0]; jx++) { @@ -437,11 +441,8 @@ struct OutputParticlesHdf5 int ix = jx + info.off[0] - lo_[0]; int iy = jy + info.off[1] - lo_[1]; int iz = jz + info.off[2] - lo_[2]; - int jj = - ((kind * ld[2] + jz - ilo[2]) * ld[1] + jy - ilo[1]) * ld[0] + - jx - ilo[0]; - gidx_begin(ix, iy, iz, kind) = idx[p][jj]; - gidx_end(ix, iy, iz, kind) = idx[p][jj + sz]; + gidx_begin(ix, iy, iz, kind) = _idx(jx, jy, jz, kind, 0); + gidx_end(ix, iy, iz, kind) = _idx(jx, jy, jz, kind, 1); } } } From c2dcdb3ad961570a737e420cc7268e0f65f0a97d Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 22:48:08 -0400 Subject: [PATCH 19/49] outp/hdf5: idx directly gtensor-based --- .../output_particles_hdf5_impl.hxx | 19 +++++++------------ 1 file changed, 7 insertions(+), 12 deletions(-) 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 185687ef93..d6daeeb1bb 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -204,7 +204,7 @@ struct OutputParticlesHdf5 hdf5_prt* make_local_particle_array(Mparticles& mprts, const std::vector>& off, const std::vector>& map, - std::vector>& idx, + std::vector>& idx, size_t* p_n_write, size_t* p_n_off, size_t* p_n_total) { @@ -245,18 +245,16 @@ struct OutputParticlesHdf5 const auto& patch = grid.patches[p]; int ilo[3], ihi[3], ld[3]; int sz = find_patch_bounds(grid.ldims, patch.off, ilo, ihi, ld); - idx[p].resize(2 * sz); - auto _idx = gt::adapt<5>(idx[p].data(), - gt::shape(ld[0], ld[1], ld[2], nr_kinds, 2)); + idx[p] = gt::empty({ld[0], ld[1], ld[2], nr_kinds, 2}); for (int jz = ilo[2]; jz < ihi[2]; jz++) { for (int jy = ilo[1]; jy < ihi[1]; jy++) { for (int jx = ilo[0]; jx < ihi[0]; jx++) { for (int kind = 0; kind < nr_kinds; kind++) { int si = sort_index(grid.ldims, nr_kinds, {jx, jy, jz}, kind); - _idx(jx - ilo[0], jy - ilo[1], jz - ilo[2], kind, 0) = + idx[p](jx - ilo[0], jy - ilo[1], jz - ilo[2], kind, 0) = nn + n_off; - _idx(jx - ilo[0], jy - ilo[1], jz - ilo[2], kind, 1) = + idx[p](jx - ilo[0], jy - ilo[1], jz - ilo[2], kind, 1) = nn + n_off + off[p][si + 1] - off[p][si]; for (int n = off[p][si]; n < off[p][si + 1]; n++, nn++) { auto prt = prts[map[p][n]]; @@ -392,7 +390,7 @@ struct OutputParticlesHdf5 std::vector> off(mprts.n_patches()); std::vector> map(mprts.n_patches()); - std::vector> idx(mprts.n_patches()); + std::vector> idx(mprts.n_patches()); count_sort(mprts, off, map); @@ -431,9 +429,6 @@ struct OutputParticlesHdf5 int ilo[3], ihi[3], ld[3]; int sz = find_patch_bounds(info.ldims, info.off, ilo, ihi, ld); - auto _idx = gt::adapt<5>(idx[p].data(), - gt::shape(ld[0], ld[1], ld[2], nr_kinds, 2)); - for (int jz = ilo[2]; jz < ihi[2]; jz++) { for (int jy = ilo[1]; jy < ihi[1]; jy++) { for (int jx = ilo[0]; jx < ihi[0]; jx++) { @@ -441,8 +436,8 @@ struct OutputParticlesHdf5 int ix = jx + info.off[0] - lo_[0]; int iy = jy + info.off[1] - lo_[1]; int iz = jz + info.off[2] - lo_[2]; - gidx_begin(ix, iy, iz, kind) = _idx(jx, jy, jz, kind, 0); - gidx_end(ix, iy, iz, kind) = _idx(jx, jy, jz, kind, 1); + gidx_begin(ix, iy, iz, kind) = idx[p](jx, jy, jz, kind, 0); + gidx_end(ix, iy, iz, kind) = idx[p](jx, jy, jz, kind, 1); } } } From 8c4c9bdf8c9b6b4385c5e6b6997d5c51d2a3c5a1 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 22:52:21 -0400 Subject: [PATCH 20/49] outp/hdf5: use std::vector some more --- .../output_particles_hdf5_impl.hxx | 12 ++++-------- 1 file changed, 4 insertions(+), 8 deletions(-) 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 d6daeeb1bb..19ee5da7cb 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -444,10 +444,9 @@ struct OutputParticlesHdf5 } } - size_t** recv_buf = (size_t**)calloc(size, sizeof(*recv_buf)); - size_t** recv_buf_p = (size_t**)calloc(size, sizeof(*recv_buf_p)); - MPI_Request* recv_req = (MPI_Request*)calloc(size, sizeof(*recv_req)); - recv_req[rank] = MPI_REQUEST_NULL; + std::vector recv_buf(size); + std::vector recv_buf_p(size); + std::vector recv_req(size, MPI_REQUEST_NULL); for (int r = 0; r < size; r++) { if (r == rank) { // skip this proc continue; @@ -457,8 +456,7 @@ struct OutputParticlesHdf5 MPI_Irecv(recv_buf[r], 2 * remote_sz[r], MPI_UNSIGNED_LONG, r, 0x4000, comm_, &recv_req[r]); } - MPI_Waitall(size, recv_req, MPI_STATUSES_IGNORE); - free(recv_req); + MPI_Waitall(size, recv_req.data(), MPI_STATUSES_IGNORE); // build global idx array, remote part for (int p = 0; p < nr_global_patches; p++) { @@ -490,8 +488,6 @@ struct OutputParticlesHdf5 for (int r = 0; r < size; r++) { free(recv_buf[r]); } - free(recv_buf); - free(recv_buf_p); free(remote_sz); } else { // rank != 0: send to rank 0 // FIXME, alloc idx[] as one array in the first place From d22f04af23456be8f71febe6512ba626b56d1461 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 22:54:03 -0400 Subject: [PATCH 21/49] outp/hdf5: more std::vector --- .../output_particles_hdf5_impl.hxx | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) 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 19ee5da7cb..d775ebf160 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -444,17 +444,17 @@ struct OutputParticlesHdf5 } } - std::vector recv_buf(size); + std::vector> recv_buf(size); std::vector recv_buf_p(size); std::vector recv_req(size, MPI_REQUEST_NULL); for (int r = 0; r < size; r++) { if (r == rank) { // skip this proc continue; } - recv_buf[r] = (size_t*)malloc(2 * remote_sz[r] * sizeof(*recv_buf[r])); - recv_buf_p[r] = recv_buf[r]; - MPI_Irecv(recv_buf[r], 2 * remote_sz[r], MPI_UNSIGNED_LONG, r, 0x4000, - comm_, &recv_req[r]); + recv_buf[r].resize(2 * remote_sz[r]); + recv_buf_p[r] = recv_buf[r].data(); + MPI_Irecv(recv_buf[r].data(), 2 * remote_sz[r], MPI_UNSIGNED_LONG, r, + 0x4000, comm_, &recv_req[r]); } MPI_Waitall(size, recv_req.data(), MPI_STATUSES_IGNORE); @@ -485,9 +485,6 @@ struct OutputParticlesHdf5 recv_buf_p[info.rank] += 2 * sz; } - for (int r = 0; r < size; r++) { - free(recv_buf[r]); - } free(remote_sz); } else { // rank != 0: send to rank 0 // FIXME, alloc idx[] as one array in the first place From f39ac16176304d3075ba511c25eb46e4ebab81b0 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 22:55:27 -0400 Subject: [PATCH 22/49] outp/hdf5: remote_sz as std::vector --- .../psc_output_particles/output_particles_hdf5_impl.hxx | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) 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 d775ebf160..06bc15de69 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -412,7 +412,7 @@ struct OutputParticlesHdf5 if (rank == 0) { int nr_global_patches = grid.nGlobalPatches(); - int* remote_sz = (int*)calloc(size, sizeof(*remote_sz)); + std::vector remote_sz(size); for (int p = 0; p < nr_global_patches; p++) { info = grid.globalPatchInfo(p); if (info.rank == rank) { // skip local patches @@ -453,7 +453,7 @@ struct OutputParticlesHdf5 } recv_buf[r].resize(2 * remote_sz[r]); recv_buf_p[r] = recv_buf[r].data(); - MPI_Irecv(recv_buf[r].data(), 2 * remote_sz[r], MPI_UNSIGNED_LONG, r, + MPI_Irecv(recv_buf[r].data(), recv_buf[r].size(), MPI_UNSIGNED_LONG, r, 0x4000, comm_, &recv_req[r]); } MPI_Waitall(size, recv_req.data(), MPI_STATUSES_IGNORE); @@ -484,8 +484,6 @@ struct OutputParticlesHdf5 } recv_buf_p[info.rank] += 2 * sz; } - - free(remote_sz); } else { // rank != 0: send to rank 0 // FIXME, alloc idx[] as one array in the first place int local_sz = 0; From 8d30beac2801a2f3d3fc7418a409194337359f69 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 22:57:25 -0400 Subject: [PATCH 23/49] outp/hdf5: std::vector for lidx --- .../psc_output_particles/output_particles_hdf5_impl.hxx | 8 +++----- 1 file changed, 3 insertions(+), 5 deletions(-) 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 06bc15de69..d61accbfb4 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -493,8 +493,8 @@ struct OutputParticlesHdf5 int sz = find_patch_bounds(grid.ldims, patch.off, ilo, ihi, ld); local_sz += sz; } - size_t* l_idx = (size_t*)malloc(2 * local_sz * sizeof(*l_idx)); - size_t* l_idx_p = (size_t*)l_idx; + std::vector l_idx(2 * local_sz); + size_t* l_idx_p = l_idx.data(); for (int p = 0; p < mprts.n_patches(); p++) { const auto& patch = grid.patches[p]; int ilo[3], ihi[3], ld[3]; @@ -503,9 +503,7 @@ struct OutputParticlesHdf5 l_idx_p += 2 * sz; } - MPI_Send(l_idx, 2 * local_sz, MPI_UNSIGNED_LONG, 0, 0x4000, comm_); - - free(l_idx); + MPI_Send(l_idx.data(), l_idx.size(), MPI_UNSIGNED_LONG, 0, 0x4000, comm_); } prof_stop(pr_B); From 241d5e04c3a85d6ed5d7f799e82067ff943f2a23 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 23:05:18 -0400 Subject: [PATCH 24/49] outp/hdf5: use loop over dims in get_sort_index() --- .../output_particles_hdf5_impl.hxx | 25 ++++++++----------- 1 file changed, 10 insertions(+), 15 deletions(-) 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 d61accbfb4..51fc5a7f58 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -123,23 +123,18 @@ struct OutputParticlesHdf5 const Grid_t& grid = prts.grid(); const int* ldims = grid.ldims; - int j0 = prts.cellPosition(prt.x[0], 0); - int j1 = prts.cellPosition(prt.x[1], 1); - int j2 = prts.cellPosition(prt.x[2], 2); - // FIXME, this is hoping that reason is that we were just on the right - // bnd... - if (j0 == ldims[0]) - j0--; - if (j1 == ldims[1]) - j1--; - if (j2 == ldims[2]) - j2--; - assert(j0 >= 0 && j0 < ldims[0]); - assert(j1 >= 0 && j1 < ldims[1]); - assert(j2 >= 0 && j2 < ldims[2]); + Int3 pos; + for (int d = 0; d < 3; d++) { + pos[d] = prts.cellPosition(prt.x[d], d); + // FIXME, this is hoping that reason is that we were just on the right + // bnd... + if (pos[d] == ldims[d]) + pos[d]--; + assert(pos[d] >= 0 && pos[d] < ldims[d]); + } assert(prt.kind < grid.kinds.size()); - return sort_index(ldims, grid.kinds.size(), {j0, j1, j2}, prt.kind); + return sort_index(ldims, grid.kinds.size(), pos, prt.kind); } // ---------------------------------------------------------------------- From 6b1e9be8a618af67b7e80dc9388bdd40e4cbcd74 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 23:09:37 -0400 Subject: [PATCH 25/49] outp/hdf5: std::vector for temp array of particles --- .../output_particles_hdf5_impl.hxx | 24 +++++++++---------- 1 file changed, 11 insertions(+), 13 deletions(-) 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 51fc5a7f58..37555e7e26 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -196,12 +196,11 @@ struct OutputParticlesHdf5 // ---------------------------------------------------------------------- // make_local_particle_array - hdf5_prt* make_local_particle_array(Mparticles& mprts, - const std::vector>& off, - const std::vector>& map, - std::vector>& idx, - size_t* p_n_write, size_t* p_n_off, - size_t* p_n_total) + std::vector make_local_particle_array( + Mparticles& mprts, const std::vector>& off, + const std::vector>& map, + std::vector>& idx, size_t* p_n_write, + size_t* p_n_off, size_t* p_n_total) { const auto& grid = mprts.grid(); int nr_kinds = grid.kinds.size(); @@ -229,7 +228,7 @@ struct OutputParticlesHdf5 MPI_Allreduce(&n_write, &n_total, 1, MPI_UNSIGNED_LONG, MPI_SUM, comm_); MPI_Exscan(&n_write, &n_off, 1, MPI_UNSIGNED_LONG, MPI_SUM, comm_); - struct hdf5_prt* arr = (struct hdf5_prt*)malloc(n_write * sizeof(*arr)); + std::vector arr(n_write); // copy particles to be written into temp array int nn = 0; @@ -278,7 +277,8 @@ struct OutputParticlesHdf5 } void write_particles(size_t n_write, size_t n_off, size_t n_total, - hdf5_prt* arr, hid_t group, hid_t dxpl) + const std::vector& arr, hid_t group, + hid_t dxpl) { herr_t ierr; @@ -296,7 +296,7 @@ struct OutputParticlesHdf5 hid_t dset = H5Dcreate(group, "1d", prt_type_, filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(dset); - ierr = H5Dwrite(dset, prt_type_, memspace, filespace, dxpl, arr); + ierr = H5Dwrite(dset, prt_type_, memspace, filespace, dxpl, arr.data()); CE; ierr = H5Dclose(dset); @@ -399,8 +399,8 @@ struct OutputParticlesHdf5 // find local particle and idx arrays size_t n_write, n_off, n_total; - hdf5_prt* arr = make_local_particle_array(mprts, off, map, idx, &n_write, - &n_off, &n_total); + auto arr = make_local_particle_array(mprts, off, map, idx, &n_write, &n_off, + &n_total); prof_stop(pr_A); prof_start(pr_B); @@ -567,8 +567,6 @@ struct OutputParticlesHdf5 H5Gclose(group); H5Fclose(file); prof_stop(pr_E); - - free(arr); } private: From f46adeadb086ca8fdc09f391d7953b01f30114d9 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Tue, 11 Jun 2024 23:31:30 -0400 Subject: [PATCH 26/49] outp: use 2 patches in test, one particle each should work on 2 MPI procs now, too. --- src/libpsc/tests/test_output_particles.cxx | 14 +++++++++++--- 1 file changed, 11 insertions(+), 3 deletions(-) diff --git a/src/libpsc/tests/test_output_particles.cxx b/src/libpsc/tests/test_output_particles.cxx index 7e8fb3325e..8c200dad09 100644 --- a/src/libpsc/tests/test_output_particles.cxx +++ b/src/libpsc/tests/test_output_particles.cxx @@ -38,7 +38,8 @@ struct OutputParticlesTest : ::testing::Test ibn[2] = 0; } - auto grid_domain = Grid_t::Domain{gdims, {L, L, L}}; + auto grid_domain = + Grid_t::Domain{gdims, {L, L, L}, {0., 0., 0.}, {2, 1, 1}}; auto grid_bc = psc::grid::BC{{BND_FLD_PERIODIC, BND_FLD_PERIODIC, BND_FLD_PERIODIC}, {BND_FLD_PERIODIC, BND_FLD_PERIODIC, BND_FLD_PERIODIC}, @@ -92,8 +93,15 @@ TYPED_TEST(OutputParticlesTest, Test1) Mparticles mprts{grid}; { auto injector = mprts.injector(); - injector[0]({{1., 0., 0.}, {}, 1., 0}); - injector[0]({{40., 0., 0.}, {}, 1., 1}); + for (int p = 0; p < grid.n_patches(); p++) { + auto info = grid.mrc_domain().localPatchInfo(p); + if (info.global_patch == 0) { + injector[p]({{1., 0., 0.}, {}, 1., 0}); + } + if (info.global_patch == 1) { + injector[p]({{121., 0., 0.}, {}, 1., 1}); + } + } } auto params = OutputParticlesParams{}; From 59616b5644c8f14bc9e3bda9c7d5be601f70173e Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 22:37:23 -0400 Subject: [PATCH 27/49] outp: pass gidx_{begin,end} as gtensor's to write_idx() --- .../output_particles_hdf5_impl.hxx | 13 +++++++------ 1 file changed, 7 insertions(+), 6 deletions(-) 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 37555e7e26..d07b2f62e9 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -307,7 +307,8 @@ struct OutputParticlesHdf5 CE; } - void write_idx(const size_t* gidx_begin, const size_t* gidx_end, hid_t group, + void write_idx(const gt::gtensor& gidx_begin, + const gt::gtensor& gidx_end, hid_t group, hid_t dxpl) { herr_t ierr; @@ -336,8 +337,8 @@ struct OutputParticlesHdf5 hid_t dset = H5Dcreate(group, "idx_begin", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(dset); - ierr = - H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, gidx_begin); + ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, + gidx_begin.data()); CE; ierr = H5Dclose(dset); CE; @@ -345,8 +346,8 @@ struct OutputParticlesHdf5 dset = H5Dcreate(group, "idx_end", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(dset); - ierr = - H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, gidx_end); + ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, + gidx_end.data()); CE; ierr = H5Dclose(dset); CE; @@ -554,7 +555,7 @@ struct OutputParticlesHdf5 prof_stop(pr_C); prof_start(pr_D); - write_idx(gidx_begin.data(), gidx_end.data(), groupp, dxpl); + write_idx(gidx_begin, gidx_end, groupp, dxpl); prof_stop(pr_D); prof_start(pr_E); From ab60b33bdc98f34b4ce46397d23c605e277660bb Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 22:49:20 -0400 Subject: [PATCH 28/49] outp: start a separate write_hdf5() function The goal being to have the hdf5 specific stuff in one place --- .../output_particles_hdf5_impl.hxx | 34 +++++++++++++------ 1 file changed, 24 insertions(+), 10 deletions(-) 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 d07b2f62e9..b7cf8fcbcd 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -358,6 +358,27 @@ struct OutputParticlesHdf5 CE; } + void write_hdf5(const gt::gtensor& gidx_begin, + const gt::gtensor& gidx_end, size_t n_write, + size_t n_off, size_t n_total, + const std::vector& arr, hid_t groupp, hid_t group, + hid_t dxpl) + { + static int pr_D, pr_E; + if (!pr_D) { + pr_D = prof_register("outp: write idx", 1., 0, 0); + pr_E = prof_register("outp: write prts", 1., 0, 0); + } + + prof_start(pr_D); + write_idx(gidx_begin, gidx_end, group, dxpl); + prof_stop(pr_D); + + prof_start(pr_E); + write_particles(n_write, n_off, n_total, arr, groupp, dxpl); + prof_stop(pr_E); + } + // ---------------------------------------------------------------------- // operator() @@ -367,13 +388,11 @@ struct OutputParticlesHdf5 const auto& grid = mprts.grid(); herr_t ierr; - static int pr_A, pr_B, pr_C, pr_D, pr_E; + static int pr_A, pr_B, pr_C; if (!pr_A) { pr_A = prof_register("outp: local", 1., 0, 0); pr_B = prof_register("outp: comm", 1., 0, 0); pr_C = prof_register("outp: write prep", 1., 0, 0); - pr_D = prof_register("outp: write idx", 1., 0, 0); - pr_E = prof_register("outp: write prts", 1., 0, 0); } prof_start(pr_A); @@ -554,12 +573,8 @@ struct OutputParticlesHdf5 #endif prof_stop(pr_C); - prof_start(pr_D); - write_idx(gidx_begin, gidx_end, groupp, dxpl); - prof_stop(pr_D); - - prof_start(pr_E); - write_particles(n_write, n_off, n_total, arr, groupp, dxpl); + write_hdf5(gidx_begin, gidx_end, n_write, n_off, n_total, arr, groupp, + group, dxpl); ierr = H5Pclose(dxpl); CE; @@ -567,7 +582,6 @@ struct OutputParticlesHdf5 H5Gclose(groupp); H5Gclose(group); H5Fclose(file); - prof_stop(pr_E); } private: From 71fabac7e26a1db1be6afe3765002f3c3612ab10 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:05:43 -0400 Subject: [PATCH 29/49] outp: all of the writing is now in write_hdf5() --- .../output_particles_hdf5_impl.hxx | 134 +++++++++--------- 1 file changed, 70 insertions(+), 64 deletions(-) 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 b7cf8fcbcd..2e64987b5b 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -361,15 +361,70 @@ struct OutputParticlesHdf5 void write_hdf5(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, size_t n_write, size_t n_off, size_t n_total, - const std::vector& arr, hid_t groupp, hid_t group, - hid_t dxpl) + const std::vector& arr, const std::string& filename, + const OutputParticlesParams& params) { - static int pr_D, pr_E; - if (!pr_D) { + int ierr; + + static int pr_C, pr_D, pr_E; + if (!pr_C) { + pr_C = prof_register("outp: write prep", 1., 0, 0); pr_D = prof_register("outp: write idx", 1., 0, 0); pr_E = prof_register("outp: write prts", 1., 0, 0); } + prof_start(pr_C); + + hid_t plist = H5Pcreate(H5P_FILE_ACCESS); + + MPI_Info mpi_info; + MPI_Info_create(&mpi_info); + +#ifdef H5_HAVE_PARALLEL + if (params.romio_cb_write) { + MPI_Info_set(mpi_info, (char*)"romio_cb_write", + (char*)params.romio_cb_write); + } + if (params.romio_ds_write) { + MPI_Info_set(mpi_info, (char*)"romio_ds_write", + (char*)params.romio_ds_write); + } + H5Pset_fapl_mpio(plist, comm_, mpi_info); +#else + mprintf("ERROR: particle output requires parallel hdf5\n"); + abort(); +#endif + + hid_t file = H5Fcreate(filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist); + H5_CHK(file); + H5Pclose(plist); + MPI_Info_free(&mpi_info); + + hid_t group = + H5Gcreate(file, "particles", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(group); + hid_t groupp = + H5Gcreate(group, "p0", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(groupp); + + ierr = H5LTset_attribute_int(group, ".", "lo", lo_, 3); + CE; + ierr = H5LTset_attribute_int(group, ".", "hi", hi_, 3); + CE; + + hid_t dxpl = H5Pcreate(H5P_DATASET_XFER); + H5_CHK(dxpl); +#ifdef H5_HAVE_PARALLEL + if (params.use_independent_io) { + ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_INDEPENDENT); + CE; + } else { + ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_COLLECTIVE); + CE; + } +#endif + prof_stop(pr_C); + prof_start(pr_D); write_idx(gidx_begin, gidx_end, group, dxpl); prof_stop(pr_D); @@ -377,6 +432,15 @@ struct OutputParticlesHdf5 prof_start(pr_E); write_particles(n_write, n_off, n_total, arr, groupp, dxpl); prof_stop(pr_E); + + ierr = H5Pclose(dxpl); + CE; + ierr = H5Gclose(groupp); + CE; + ierr = H5Gclose(group); + CE; + ierr = H5Fclose(file); + CE; } // ---------------------------------------------------------------------- @@ -522,66 +586,8 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - prof_start(pr_C); - - hid_t plist = H5Pcreate(H5P_FILE_ACCESS); - - MPI_Info mpi_info; - MPI_Info_create(&mpi_info); -#ifdef H5_HAVE_PARALLEL - if (params.romio_cb_write) { - MPI_Info_set(mpi_info, (char*)"romio_cb_write", - (char*)params.romio_cb_write); - } - if (params.romio_ds_write) { - MPI_Info_set(mpi_info, (char*)"romio_ds_write", - (char*)params.romio_ds_write); - } - H5Pset_fapl_mpio(plist, comm_, mpi_info); -#else - mprintf("ERROR: particle output requires parallel hdf5\n"); - abort(); -#endif - - hid_t file = H5Fcreate(filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist); - H5_CHK(file); - H5Pclose(plist); - MPI_Info_free(&mpi_info); - - hid_t group = - H5Gcreate(file, "particles", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); - H5_CHK(group); - hid_t groupp = - H5Gcreate(group, "p0", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); - H5_CHK(groupp); - - ierr = H5LTset_attribute_int(group, ".", "lo", lo_, 3); - CE; - ierr = H5LTset_attribute_int(group, ".", "hi", hi_, 3); - CE; - - hid_t dxpl = H5Pcreate(H5P_DATASET_XFER); - H5_CHK(dxpl); -#ifdef H5_HAVE_PARALLEL - if (params.use_independent_io) { - ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_INDEPENDENT); - CE; - } else { - ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_COLLECTIVE); - CE; - } -#endif - prof_stop(pr_C); - - write_hdf5(gidx_begin, gidx_end, n_write, n_off, n_total, arr, groupp, - group, dxpl); - - ierr = H5Pclose(dxpl); - CE; - - H5Gclose(groupp); - H5Gclose(group); - H5Fclose(file); + write_hdf5(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, + params); } private: From 1e9040e05dc662b69ecfe6936524404c4378457d Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:17:15 -0400 Subject: [PATCH 30/49] outp/hdf5: start OutputParticlesWriterHDF5 class That will eventually do all HDF5 specific writing --- .../output_particles_hdf5_impl.hxx | 120 ++++++++++-------- 1 file changed, 68 insertions(+), 52 deletions(-) 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 2e64987b5b..d3ab4f47ca 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -30,6 +30,71 @@ struct hdf5_prt #define MIN(a, b) ((a) < (b) ? (a) : (b)) #define MAX(a, b) ((a) > (b) ? (a) : (b)) +class OutputParticlesWriterHDF5 +{ +public: + explicit OutputParticlesWriterHDF5(const Int3& wdims, + const Grid_t::Kinds& kinds, MPI_Comm comm) + : wdims_{wdims}, kinds_{kinds}, comm_{comm} + {} + + void write_idx(const gt::gtensor& gidx_begin, + const gt::gtensor& gidx_end, hid_t group, + hid_t dxpl) + { + herr_t ierr; + + hsize_t fdims[4]; + fdims[0] = kinds_.size(); + fdims[1] = wdims_[2]; + fdims[2] = wdims_[1]; + fdims[3] = wdims_[0]; + hid_t filespace = H5Screate_simple(4, fdims, NULL); + H5_CHK(filespace); + hid_t memspace; + + assert(sizeof(size_t) == sizeof(hsize_t)); + int rank; + MPI_Comm_rank(comm_, &rank); + if (rank == 0) { + memspace = H5Screate_simple(4, fdims, NULL); + H5_CHK(memspace); + } else { + memspace = H5Screate(H5S_NULL); + H5Sselect_none(memspace); + H5Sselect_none(filespace); + } + + hid_t dset = H5Dcreate(group, "idx_begin", H5T_NATIVE_HSIZE, filespace, + H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(dset); + ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, + gidx_begin.data()); + CE; + ierr = H5Dclose(dset); + CE; + + dset = H5Dcreate(group, "idx_end", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, + H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(dset); + ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, + gidx_end.data()); + CE; + ierr = H5Dclose(dset); + CE; + + ierr = H5Sclose(filespace); + CE; + ierr = H5Sclose(memspace); + CE; + } + +private: + Int3 wdims_; + Grid_t::Kinds kinds_; + MPI_Comm comm_; +}; + // ---------------------------------------------------------------------- template @@ -307,57 +372,6 @@ struct OutputParticlesHdf5 CE; } - void write_idx(const gt::gtensor& gidx_begin, - const gt::gtensor& gidx_end, hid_t group, - hid_t dxpl) - { - herr_t ierr; - - hsize_t fdims[4]; - fdims[0] = kinds_.size(); - fdims[1] = wdims_[2]; - fdims[2] = wdims_[1]; - fdims[3] = wdims_[0]; - hid_t filespace = H5Screate_simple(4, fdims, NULL); - H5_CHK(filespace); - hid_t memspace; - - assert(sizeof(size_t) == sizeof(hsize_t)); - int rank; - MPI_Comm_rank(comm_, &rank); - if (rank == 0) { - memspace = H5Screate_simple(4, fdims, NULL); - H5_CHK(memspace); - } else { - memspace = H5Screate(H5S_NULL); - H5Sselect_none(memspace); - H5Sselect_none(filespace); - } - - hid_t dset = H5Dcreate(group, "idx_begin", H5T_NATIVE_HSIZE, filespace, - H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); - H5_CHK(dset); - ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, - gidx_begin.data()); - CE; - ierr = H5Dclose(dset); - CE; - - dset = H5Dcreate(group, "idx_end", H5T_NATIVE_HSIZE, filespace, H5P_DEFAULT, - H5P_DEFAULT, H5P_DEFAULT); - H5_CHK(dset); - ierr = H5Dwrite(dset, H5T_NATIVE_HSIZE, memspace, filespace, dxpl, - gidx_end.data()); - CE; - ierr = H5Dclose(dset); - CE; - - ierr = H5Sclose(filespace); - CE; - ierr = H5Sclose(memspace); - CE; - } - void write_hdf5(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, size_t n_write, size_t n_off, size_t n_total, @@ -425,8 +439,10 @@ struct OutputParticlesHdf5 #endif prof_stop(pr_C); + OutputParticlesWriterHDF5 writer(wdims_, kinds_, comm_); + prof_start(pr_D); - write_idx(gidx_begin, gidx_end, group, dxpl); + writer.write_idx(gidx_begin, gidx_end, group, dxpl); prof_stop(pr_D); prof_start(pr_E); From 9f7eb4162e5fe4de51674aa2be97772f83a28da5 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:25:44 -0400 Subject: [PATCH 31/49] outp/hdf5: all writing -> OutputParticlesWriterHDF5 class --- .../output_particles_hdf5_impl.hxx | 256 +++++++++--------- 1 file changed, 132 insertions(+), 124 deletions(-) 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 d3ab4f47ca..130da87edc 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -33,11 +33,103 @@ struct hdf5_prt class OutputParticlesWriterHDF5 { public: - explicit OutputParticlesWriterHDF5(const Int3& wdims, - const Grid_t::Kinds& kinds, MPI_Comm comm) - : wdims_{wdims}, kinds_{kinds}, comm_{comm} + explicit OutputParticlesWriterHDF5(const Int3& wdims, const Int3& lo, + const Int3& hi, const Grid_t::Kinds& kinds, + MPI_Comm comm, hid_t prt_type) + : wdims_{wdims}, + lo_{lo}, + hi_{hi}, + kinds_{kinds}, + comm_{comm}, + prt_type_{prt_type} {} + void operator()(const gt::gtensor& gidx_begin, + const gt::gtensor& gidx_end, size_t n_write, + size_t n_off, size_t n_total, + const std::vector& arr, const std::string& filename, + const OutputParticlesParams& params) + { + int ierr; + + static int pr_C, pr_D, pr_E; + if (!pr_C) { + pr_C = prof_register("outp: write prep", 1., 0, 0); + pr_D = prof_register("outp: write idx", 1., 0, 0); + pr_E = prof_register("outp: write prts", 1., 0, 0); + } + + prof_start(pr_C); + + hid_t plist = H5Pcreate(H5P_FILE_ACCESS); + + MPI_Info mpi_info; + MPI_Info_create(&mpi_info); + +#ifdef H5_HAVE_PARALLEL + if (params.romio_cb_write) { + MPI_Info_set(mpi_info, (char*)"romio_cb_write", + (char*)params.romio_cb_write); + } + if (params.romio_ds_write) { + MPI_Info_set(mpi_info, (char*)"romio_ds_write", + (char*)params.romio_ds_write); + } + H5Pset_fapl_mpio(plist, comm_, mpi_info); +#else + mprintf("ERROR: particle output requires parallel hdf5\n"); + abort(); +#endif + + hid_t file = H5Fcreate(filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist); + H5_CHK(file); + H5Pclose(plist); + MPI_Info_free(&mpi_info); + + hid_t group = + H5Gcreate(file, "particles", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(group); + hid_t groupp = + H5Gcreate(group, "p0", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(groupp); + + ierr = H5LTset_attribute_int(group, ".", "lo", lo_, 3); + CE; + ierr = H5LTset_attribute_int(group, ".", "hi", hi_, 3); + CE; + + hid_t dxpl = H5Pcreate(H5P_DATASET_XFER); + H5_CHK(dxpl); +#ifdef H5_HAVE_PARALLEL + if (params.use_independent_io) { + ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_INDEPENDENT); + CE; + } else { + ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_COLLECTIVE); + CE; + } +#endif + prof_stop(pr_C); + + prof_start(pr_D); + write_idx(gidx_begin, gidx_end, group, dxpl); + prof_stop(pr_D); + + prof_start(pr_E); + write_particles(n_write, n_off, n_total, arr, groupp, dxpl); + prof_stop(pr_E); + + ierr = H5Pclose(dxpl); + CE; + ierr = H5Gclose(groupp); + CE; + ierr = H5Gclose(group); + CE; + ierr = H5Fclose(file); + CE; + } + +private: void write_idx(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, hid_t group, hid_t dxpl) @@ -89,10 +181,42 @@ public: CE; } + void write_particles(size_t n_write, size_t n_off, size_t n_total, + const std::vector& arr, hid_t group, + hid_t dxpl) + { + herr_t ierr; + + hsize_t mdims[1] = {n_write}; + hsize_t fdims[1] = {n_total}; + hsize_t foff[1] = {n_off}; + hid_t memspace = H5Screate_simple(1, mdims, NULL); + H5_CHK(memspace); + hid_t filespace = H5Screate_simple(1, fdims, NULL); + H5_CHK(filespace); + ierr = + H5Sselect_hyperslab(filespace, H5S_SELECT_SET, foff, NULL, mdims, NULL); + CE; + + hid_t dset = H5Dcreate(group, "1d", prt_type_, filespace, H5P_DEFAULT, + H5P_DEFAULT, H5P_DEFAULT); + H5_CHK(dset); + ierr = H5Dwrite(dset, prt_type_, memspace, filespace, dxpl, arr.data()); + CE; + + ierr = H5Dclose(dset); + CE; + ierr = H5Sclose(filespace); + CE; + ierr = H5Sclose(memspace); + CE; + } + private: - Int3 wdims_; + Int3 wdims_, lo_, hi_; Grid_t::Kinds kinds_; MPI_Comm comm_; + hid_t prt_type_; }; // ---------------------------------------------------------------------- @@ -341,124 +465,6 @@ struct OutputParticlesHdf5 return arr; } - void write_particles(size_t n_write, size_t n_off, size_t n_total, - const std::vector& arr, hid_t group, - hid_t dxpl) - { - herr_t ierr; - - hsize_t mdims[1] = {n_write}; - hsize_t fdims[1] = {n_total}; - hsize_t foff[1] = {n_off}; - hid_t memspace = H5Screate_simple(1, mdims, NULL); - H5_CHK(memspace); - hid_t filespace = H5Screate_simple(1, fdims, NULL); - H5_CHK(filespace); - ierr = - H5Sselect_hyperslab(filespace, H5S_SELECT_SET, foff, NULL, mdims, NULL); - CE; - - hid_t dset = H5Dcreate(group, "1d", prt_type_, filespace, H5P_DEFAULT, - H5P_DEFAULT, H5P_DEFAULT); - H5_CHK(dset); - ierr = H5Dwrite(dset, prt_type_, memspace, filespace, dxpl, arr.data()); - CE; - - ierr = H5Dclose(dset); - CE; - ierr = H5Sclose(filespace); - CE; - ierr = H5Sclose(memspace); - CE; - } - - void write_hdf5(const gt::gtensor& gidx_begin, - const gt::gtensor& gidx_end, size_t n_write, - size_t n_off, size_t n_total, - const std::vector& arr, const std::string& filename, - const OutputParticlesParams& params) - { - int ierr; - - static int pr_C, pr_D, pr_E; - if (!pr_C) { - pr_C = prof_register("outp: write prep", 1., 0, 0); - pr_D = prof_register("outp: write idx", 1., 0, 0); - pr_E = prof_register("outp: write prts", 1., 0, 0); - } - - prof_start(pr_C); - - hid_t plist = H5Pcreate(H5P_FILE_ACCESS); - - MPI_Info mpi_info; - MPI_Info_create(&mpi_info); - -#ifdef H5_HAVE_PARALLEL - if (params.romio_cb_write) { - MPI_Info_set(mpi_info, (char*)"romio_cb_write", - (char*)params.romio_cb_write); - } - if (params.romio_ds_write) { - MPI_Info_set(mpi_info, (char*)"romio_ds_write", - (char*)params.romio_ds_write); - } - H5Pset_fapl_mpio(plist, comm_, mpi_info); -#else - mprintf("ERROR: particle output requires parallel hdf5\n"); - abort(); -#endif - - hid_t file = H5Fcreate(filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist); - H5_CHK(file); - H5Pclose(plist); - MPI_Info_free(&mpi_info); - - hid_t group = - H5Gcreate(file, "particles", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); - H5_CHK(group); - hid_t groupp = - H5Gcreate(group, "p0", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); - H5_CHK(groupp); - - ierr = H5LTset_attribute_int(group, ".", "lo", lo_, 3); - CE; - ierr = H5LTset_attribute_int(group, ".", "hi", hi_, 3); - CE; - - hid_t dxpl = H5Pcreate(H5P_DATASET_XFER); - H5_CHK(dxpl); -#ifdef H5_HAVE_PARALLEL - if (params.use_independent_io) { - ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_INDEPENDENT); - CE; - } else { - ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_COLLECTIVE); - CE; - } -#endif - prof_stop(pr_C); - - OutputParticlesWriterHDF5 writer(wdims_, kinds_, comm_); - - prof_start(pr_D); - writer.write_idx(gidx_begin, gidx_end, group, dxpl); - prof_stop(pr_D); - - prof_start(pr_E); - write_particles(n_write, n_off, n_total, arr, groupp, dxpl); - prof_stop(pr_E); - - ierr = H5Pclose(dxpl); - CE; - ierr = H5Gclose(groupp); - CE; - ierr = H5Gclose(group); - CE; - ierr = H5Fclose(file); - CE; - } - // ---------------------------------------------------------------------- // operator() @@ -602,8 +608,10 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - write_hdf5(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, - params); + OutputParticlesWriterHDF5 writer(wdims_, lo_, hi_, kinds_, comm_, + prt_type_); + writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, + params); } private: From 7290e33dbe7899a3aa1cc0ca848381b2990c10fd Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:32:38 -0400 Subject: [PATCH 32/49] outp/hdf5: move Hdf5ParticleType class --- .../output_particles_hdf5_impl.hxx | 113 +++++++++--------- 1 file changed, 57 insertions(+), 56 deletions(-) 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 130da87edc..477636c107 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -30,6 +30,62 @@ struct hdf5_prt #define MIN(a, b) ((a) < (b) ? (a) : (b)) #define MAX(a, b) ((a) > (b) ? (a) : (b)) +// ---------------------------------------------------------------------- +// ToHdf5Type helper class + +template +struct ToHdf5Type; + +template <> +struct ToHdf5Type +{ + static hid_t H5Type() { return H5T_NATIVE_INT; } +}; + +template <> +struct ToHdf5Type +{ + static hid_t H5Type() { return H5T_NATIVE_ULONG; } +}; + +template <> +struct ToHdf5Type +{ + static hid_t H5Type() { return H5T_NATIVE_ULLONG; } +}; + +// ====================================================================== +// Hdf5ParticleType + +class Hdf5ParticleType +{ +public: + Hdf5ParticleType() + { + id_ = H5Tcreate(H5T_COMPOUND, sizeof(struct hdf5_prt)); + H5Tinsert(id_, "x", HOFFSET(struct hdf5_prt, x), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "y", HOFFSET(struct hdf5_prt, y), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "z", HOFFSET(struct hdf5_prt, z), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "px", HOFFSET(struct hdf5_prt, px), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "py", HOFFSET(struct hdf5_prt, py), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "pz", HOFFSET(struct hdf5_prt, pz), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "q", HOFFSET(struct hdf5_prt, q), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "m", HOFFSET(struct hdf5_prt, m), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "w", HOFFSET(struct hdf5_prt, w), H5T_NATIVE_FLOAT); + H5Tinsert(id_, "id", HOFFSET(struct hdf5_prt, id), + ToHdf5Type::H5Type()); + H5Tinsert(id_, "tag", HOFFSET(struct hdf5_prt, tag), + ToHdf5Type::H5Type()); + } + + ~Hdf5ParticleType() { H5Tclose(id_); } + + operator hid_t() const { return id_; } + +private: + hid_t id_; +}; + class OutputParticlesWriterHDF5 { public: @@ -219,64 +275,9 @@ private: hid_t prt_type_; }; -// ---------------------------------------------------------------------- - -template -struct ToHdf5Type; - -template <> -struct ToHdf5Type -{ - static hid_t H5Type() { return H5T_NATIVE_INT; } -}; - -template <> -struct ToHdf5Type -{ - static hid_t H5Type() { return H5T_NATIVE_ULONG; } -}; - -template <> -struct ToHdf5Type -{ - static hid_t H5Type() { return H5T_NATIVE_ULLONG; } -}; - namespace detail { -// ====================================================================== -// Hdf5ParticleType - -class Hdf5ParticleType -{ -public: - Hdf5ParticleType() - { - id_ = H5Tcreate(H5T_COMPOUND, sizeof(struct hdf5_prt)); - H5Tinsert(id_, "x", HOFFSET(struct hdf5_prt, x), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "y", HOFFSET(struct hdf5_prt, y), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "z", HOFFSET(struct hdf5_prt, z), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "px", HOFFSET(struct hdf5_prt, px), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "py", HOFFSET(struct hdf5_prt, py), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "pz", HOFFSET(struct hdf5_prt, pz), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "q", HOFFSET(struct hdf5_prt, q), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "m", HOFFSET(struct hdf5_prt, m), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "w", HOFFSET(struct hdf5_prt, w), H5T_NATIVE_FLOAT); - H5Tinsert(id_, "id", HOFFSET(struct hdf5_prt, id), - ToHdf5Type::H5Type()); - H5Tinsert(id_, "tag", HOFFSET(struct hdf5_prt, tag), - ToHdf5Type::H5Type()); - } - - ~Hdf5ParticleType() { H5Tclose(id_); } - - operator hid_t() const { return id_; } - -private: - hid_t id_; -}; - // ====================================================================== // OutputParticlesHdf5 @@ -693,7 +694,7 @@ public: private: const OutputParticlesParams prm_; - detail::Hdf5ParticleType prt_type_; + Hdf5ParticleType prt_type_; Int3 lo_; // dimensions of the subdomain we're actually writing Int3 hi_; }; From 15367f89aa5e79d083dcd9de124cc863a84fd5f7 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:35:37 -0400 Subject: [PATCH 33/49] outp/hdf5: keep writer around through life of class --- .../output_particles_hdf5_impl.hxx | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) 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 477636c107..1fd612c51b 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -293,7 +293,8 @@ struct OutputParticlesHdf5 wdims_{hi - lo}, kinds_{grid.kinds}, comm_{grid.comm()}, - prt_type_{prt_type} + prt_type_{prt_type}, + writer_{wdims_, lo_, hi_, kinds_, comm_, prt_type_} {} // ---------------------------------------------------------------------- @@ -609,10 +610,8 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - OutputParticlesWriterHDF5 writer(wdims_, lo_, hi_, kinds_, comm_, - prt_type_); - writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, - params); + writer_(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, + params); } private: @@ -622,6 +621,7 @@ private: hid_t prt_type_; Grid_t::Kinds kinds_; MPI_Comm comm_; + OutputParticlesWriterHDF5 writer_; }; } // namespace detail From 103ec4ff5927b07cb75be4f2f64739c8273472e4 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:39:52 -0400 Subject: [PATCH 34/49] outp/hdf5: don't keep prt_type_ in intermediate class --- .../output_particles_hdf5_impl.hxx | 14 ++++++-------- 1 file changed, 6 insertions(+), 8 deletions(-) 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 1fd612c51b..c315b83586 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -89,12 +89,12 @@ private: class OutputParticlesWriterHDF5 { public: - explicit OutputParticlesWriterHDF5(const Int3& wdims, const Int3& lo, - const Int3& hi, const Grid_t::Kinds& kinds, - MPI_Comm comm, hid_t prt_type) - : wdims_{wdims}, - lo_{lo}, + explicit OutputParticlesWriterHDF5(const Int3& lo, const Int3& hi, + const Grid_t::Kinds& kinds, MPI_Comm comm, + hid_t prt_type) + : lo_{lo}, hi_{hi}, + wdims_{hi - lo}, kinds_{kinds}, comm_{comm}, prt_type_{prt_type} @@ -293,8 +293,7 @@ struct OutputParticlesHdf5 wdims_{hi - lo}, kinds_{grid.kinds}, comm_{grid.comm()}, - prt_type_{prt_type}, - writer_{wdims_, lo_, hi_, kinds_, comm_, prt_type_} + writer_{lo_, hi_, kinds_, comm_, prt_type} {} // ---------------------------------------------------------------------- @@ -618,7 +617,6 @@ private: Int3 lo_; Int3 hi_; Int3 wdims_; - hid_t prt_type_; Grid_t::Kinds kinds_; MPI_Comm comm_; OutputParticlesWriterHDF5 writer_; From 2f3a810aa45c2fd97f619b995ddb2496ae725dcc Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:53:31 -0400 Subject: [PATCH 35/49] outp/hdf5: keep writer in top-level class so that it's persistent and can, next, take care of creating the hdf5 particle type --- .../output_particles_hdf5_impl.hxx | 41 ++++++++++--------- 1 file changed, 21 insertions(+), 20 deletions(-) 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 c315b83586..b0446e77be 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -288,12 +288,7 @@ struct OutputParticlesHdf5 OutputParticlesHdf5(const Grid_t& grid, const Int3& lo, const Int3& hi, hid_t prt_type) - : lo_{lo}, - hi_{hi}, - wdims_{hi - lo}, - kinds_{grid.kinds}, - comm_{grid.comm()}, - writer_{lo_, hi_, kinds_, comm_, prt_type} + : lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{grid.kinds}, comm_{grid.comm()} {} // ---------------------------------------------------------------------- @@ -470,7 +465,8 @@ struct OutputParticlesHdf5 // operator() void operator()(Mparticles& mprts, const std::string& filename, - const OutputParticlesParams& params) + const OutputParticlesParams& params, + OutputParticlesWriterHDF5& writer) { const auto& grid = mprts.grid(); herr_t ierr; @@ -609,8 +605,8 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - writer_(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, - params); + writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, + params); } private: @@ -619,27 +615,31 @@ private: Int3 wdims_; Grid_t::Kinds kinds_; MPI_Comm comm_; - OutputParticlesWriterHDF5 writer_; }; } // namespace detail class OutputParticlesHdf5 : OutputParticlesBase { -public: - OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& params) - : prm_{params} + static Int3 adjust_hi(const Grid_t& grid, Int3 hi) { - // set hi to gdims by default (if not set differently before) - // and calculate wdims (global dims of region we're writing) - lo_ = prm_.lo; for (int d = 0; d < 3; d++) { - hi_[d] = prm_.hi[d] != 0 ? prm_.hi[d] : grid.domain.gdims[d]; - assert(lo_[d] >= 0); - assert(hi_[d] <= grid.domain.gdims[d]); + if (hi[d] == 0) { + hi[d] = grid.domain.gdims[d]; + } + assert(hi[d] <= grid.domain.gdims[d]); } + return hi; } +public: + OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& params) + : prm_{params}, + lo_{params.lo}, + hi_{adjust_hi(grid, params.hi)}, + writer_{lo_, hi_, grid.kinds, grid.comm(), prt_type_} + {} + template void operator()(Mparticles& mprts) { @@ -655,7 +655,7 @@ public: prm_.basename, grid.timestep(), 0); detail::OutputParticlesHdf5 impl{grid, lo_, hi_, prt_type_}; - impl(mprts, filename, prm_); + impl(mprts, filename, prm_, writer_); } // FIXME, handles MparticlesVpic by conversion for now @@ -695,4 +695,5 @@ private: Hdf5ParticleType prt_type_; Int3 lo_; // dimensions of the subdomain we're actually writing Int3 hi_; + OutputParticlesWriterHDF5 writer_; }; From 5a8c294828b2d1cf0f8c207ccaa919df957ce2e8 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Wed, 21 Aug 2024 23:56:13 -0400 Subject: [PATCH 36/49] outp/hdf5: particle type now within WriterHDF5 --- .../output_particles_hdf5_impl.hxx | 20 ++++++------------- 1 file changed, 6 insertions(+), 14 deletions(-) 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 b0446e77be..d7278a5522 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -90,14 +90,8 @@ class OutputParticlesWriterHDF5 { public: explicit OutputParticlesWriterHDF5(const Int3& lo, const Int3& hi, - const Grid_t::Kinds& kinds, MPI_Comm comm, - hid_t prt_type) - : lo_{lo}, - hi_{hi}, - wdims_{hi - lo}, - kinds_{kinds}, - comm_{comm}, - prt_type_{prt_type} + const Grid_t::Kinds& kinds, MPI_Comm comm) + : lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{kinds}, comm_{comm} {} void operator()(const gt::gtensor& gidx_begin, @@ -272,7 +266,7 @@ private: Int3 wdims_, lo_, hi_; Grid_t::Kinds kinds_; MPI_Comm comm_; - hid_t prt_type_; + Hdf5ParticleType prt_type_; }; namespace detail @@ -286,8 +280,7 @@ struct OutputParticlesHdf5 { using Particles = typename Mparticles::Patch; - OutputParticlesHdf5(const Grid_t& grid, const Int3& lo, const Int3& hi, - hid_t prt_type) + OutputParticlesHdf5(const Grid_t& grid, const Int3& lo, const Int3& hi) : lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{grid.kinds}, comm_{grid.comm()} {} @@ -637,7 +630,7 @@ public: : prm_{params}, lo_{params.lo}, hi_{adjust_hi(grid, params.hi)}, - writer_{lo_, hi_, grid.kinds, grid.comm(), prt_type_} + writer_{lo_, hi_, grid.kinds, grid.comm()} {} template @@ -654,7 +647,7 @@ public: snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", prm_.data_dir, prm_.basename, grid.timestep(), 0); - detail::OutputParticlesHdf5 impl{grid, lo_, hi_, prt_type_}; + detail::OutputParticlesHdf5 impl{grid, lo_, hi_}; impl(mprts, filename, prm_, writer_); } @@ -692,7 +685,6 @@ public: private: const OutputParticlesParams prm_; - Hdf5ParticleType prt_type_; Int3 lo_; // dimensions of the subdomain we're actually writing Int3 hi_; OutputParticlesWriterHDF5 writer_; From 8fa13e56fb4663c6919f90e30e0e4897947611c0 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 08:49:02 -0400 Subject: [PATCH 37/49] outp/hdf5: prettify H5Type helper --- .../output_particles_hdf5_impl.hxx | 23 +++++++++++-------- 1 file changed, 14 insertions(+), 9 deletions(-) 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 d7278a5522..644dfd5765 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -30,30 +30,35 @@ struct hdf5_prt #define MIN(a, b) ((a) < (b) ? (a) : (b)) #define MAX(a, b) ((a) > (b) ? (a) : (b)) +namespace hdf5 +{ + // ---------------------------------------------------------------------- // ToHdf5Type helper class template -struct ToHdf5Type; +struct H5Type; template <> -struct ToHdf5Type +struct H5Type { - static hid_t H5Type() { return H5T_NATIVE_INT; } + static hid_t value() { return H5T_NATIVE_INT; } }; template <> -struct ToHdf5Type +struct H5Type { - static hid_t H5Type() { return H5T_NATIVE_ULONG; } + static hid_t value() { return H5T_NATIVE_ULONG; } }; template <> -struct ToHdf5Type +struct H5Type { - static hid_t H5Type() { return H5T_NATIVE_ULLONG; } + static hid_t value() { return H5T_NATIVE_ULLONG; } }; +} // namespace hdf5 + // ====================================================================== // Hdf5ParticleType @@ -73,9 +78,9 @@ public: H5Tinsert(id_, "m", HOFFSET(struct hdf5_prt, m), H5T_NATIVE_FLOAT); H5Tinsert(id_, "w", HOFFSET(struct hdf5_prt, w), H5T_NATIVE_FLOAT); H5Tinsert(id_, "id", HOFFSET(struct hdf5_prt, id), - ToHdf5Type::H5Type()); + hdf5::H5Type::value()); H5Tinsert(id_, "tag", HOFFSET(struct hdf5_prt, tag), - ToHdf5Type::H5Type()); + hdf5::H5Type::value()); } ~Hdf5ParticleType() { H5Tclose(id_); } From 03224914fa731a804bc2e99a6664d2adb0e158a1 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 09:05:56 -0400 Subject: [PATCH 38/49] outp/hdf5: avoid using hdf5_prt directly --- .../output_particles_hdf5_impl.hxx | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) 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 644dfd5765..79c4c803ba 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -94,6 +94,9 @@ private: class OutputParticlesWriterHDF5 { public: + using particle_type = hdf5_prt; + using particles_type = std::vector; + explicit OutputParticlesWriterHDF5(const Int3& lo, const Int3& hi, const Grid_t::Kinds& kinds, MPI_Comm comm) : lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{kinds}, comm_{comm} @@ -101,8 +104,8 @@ public: void operator()(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, size_t n_write, - size_t n_off, size_t n_total, - const std::vector& arr, const std::string& filename, + size_t n_off, size_t n_total, const particles_type& arr, + const std::string& filename, const OutputParticlesParams& params) { int ierr; @@ -237,8 +240,7 @@ private: } void write_particles(size_t n_write, size_t n_off, size_t n_total, - const std::vector& arr, hid_t group, - hid_t dxpl) + const particles_type& arr, hid_t group, hid_t dxpl) { herr_t ierr; @@ -283,6 +285,8 @@ namespace detail template struct OutputParticlesHdf5 { + using writer_type = OutputParticlesWriterHDF5; + using writer_particles_type = writer_type::particles_type; using Particles = typename Mparticles::Patch; OutputParticlesHdf5(const Grid_t& grid, const Int3& lo, const Int3& hi) @@ -379,7 +383,7 @@ struct OutputParticlesHdf5 // ---------------------------------------------------------------------- // make_local_particle_array - std::vector make_local_particle_array( + writer_particles_type make_local_particle_array( Mparticles& mprts, const std::vector>& off, const std::vector>& map, std::vector>& idx, size_t* p_n_write, @@ -411,7 +415,7 @@ struct OutputParticlesHdf5 MPI_Allreduce(&n_write, &n_total, 1, MPI_UNSIGNED_LONG, MPI_SUM, comm_); MPI_Exscan(&n_write, &n_off, 1, MPI_UNSIGNED_LONG, MPI_SUM, comm_); - std::vector arr(n_write); + writer_particles_type arr(n_write); // copy particles to be written into temp array int nn = 0; From 215beb9ceccc414c64e314f785c9e65433e29d94 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 09:23:04 -0400 Subject: [PATCH 39/49] outp/hdf5: do conversion to hdf5_prt in constructor --- .../output_particles_hdf5_impl.hxx | 33 +++++++++++-------- 1 file changed, 19 insertions(+), 14 deletions(-) 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 79c4c803ba..29468bb695 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -13,10 +13,26 @@ #include "mparticles_cuda.hxx" #endif -// needs to be changed accordingly - struct hdf5_prt { + hdf5_prt() = default; + + template + hdf5_prt(const Particle& prt, const Patch& patch) + { + x = prt.x()[0] + patch.xb[0]; + y = prt.x()[1] + patch.xb[1]; + z = prt.x()[2] + patch.xb[2]; + px = prt.u()[0]; + py = prt.u()[1]; + pz = prt.u()[2]; + q = prt.q(); + m = prt.m(); + w = prt.w(); + id = prt.id(); + tag = prt.tag(); + } + float x, y, z; float px, py, pz; float q, m, w; @@ -438,18 +454,7 @@ struct OutputParticlesHdf5 idx[p](jx - ilo[0], jy - ilo[1], jz - ilo[2], kind, 1) = nn + n_off + off[p][si + 1] - off[p][si]; for (int n = off[p][si]; n < off[p][si + 1]; n++, nn++) { - auto prt = prts[map[p][n]]; - arr[nn].x = prt.x()[0] + patch.xb[0]; - arr[nn].y = prt.x()[1] + patch.xb[1]; - arr[nn].z = prt.x()[2] + patch.xb[2]; - arr[nn].px = prt.u()[0]; - arr[nn].py = prt.u()[1]; - arr[nn].pz = prt.u()[2]; - arr[nn].q = prt.q(); - arr[nn].m = prt.m(); - arr[nn].w = prt.w(); - arr[nn].id = prt.id(); - arr[nn].tag = prt.tag(); + arr[nn] = writer_type::particle_type(prts[map[p][n]], patch); } } } From ac079dba403d74b0b3b4533f275a55684e9e0547 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 09:24:03 -0400 Subject: [PATCH 40/49] outp/hdf5: use std::{min,max} --- .../output_particles_hdf5_impl.hxx | 11 ++++------- 1 file changed, 4 insertions(+), 7 deletions(-) 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 29468bb695..07ac7ea480 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -43,9 +43,6 @@ struct hdf5_prt #define H5_CHK(ierr) assert(ierr >= 0) #define CE assert(ierr == 0) -#define MIN(a, b) ((a) < (b) ? (a) : (b)) -#define MAX(a, b) ((a) > (b) ? (a) : (b)) - namespace hdf5 { @@ -387,10 +384,10 @@ struct OutputParticlesHdf5 int ihi[3], int ld[3]) { for (int d = 0; d < 3; d++) { - ilo[d] = MAX(lo_[d], off[d]) - off[d]; - ihi[d] = MIN(hi_[d], off[d] + ldims[d]) - off[d]; - ilo[d] = MIN(ldims[d], ilo[d]); - ihi[d] = MAX(0, ihi[d]); + ilo[d] = std::max(lo_[d], off[d]) - off[d]; + ihi[d] = std::min(hi_[d], off[d] + ldims[d]) - off[d]; + ilo[d] = std::min(ldims[d], ilo[d]); + ihi[d] = std::max(0, ihi[d]); ld[d] = ihi[d] - ilo[d]; } return kinds_.size() * ld[0] * ld[1] * ld[2]; From 5fed83b92c38114d0a707b107cb84d07c3effc18 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 09:36:18 -0400 Subject: [PATCH 41/49] outp/hdf5: pass/keep Params around --- .../output_particles_hdf5_impl.hxx | 38 +++++++++---------- 1 file changed, 18 insertions(+), 20 deletions(-) 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 07ac7ea480..56b2484108 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -110,16 +110,16 @@ public: using particle_type = hdf5_prt; using particles_type = std::vector; - explicit OutputParticlesWriterHDF5(const Int3& lo, const Int3& hi, - const Grid_t::Kinds& kinds, MPI_Comm comm) - : lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{kinds}, comm_{comm} + explicit OutputParticlesWriterHDF5(const OutputParticlesParams& prm, Int3& lo, + const Int3& hi, const Grid_t::Kinds& kinds, + MPI_Comm comm) + : prm_{prm}, lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{kinds}, comm_{comm} {} void operator()(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, size_t n_write, size_t n_off, size_t n_total, const particles_type& arr, - const std::string& filename, - const OutputParticlesParams& params) + const std::string& filename) { int ierr; @@ -138,13 +138,13 @@ public: MPI_Info_create(&mpi_info); #ifdef H5_HAVE_PARALLEL - if (params.romio_cb_write) { + if (prm_.romio_cb_write) { MPI_Info_set(mpi_info, (char*)"romio_cb_write", - (char*)params.romio_cb_write); + (char*)prm_.romio_cb_write); } - if (params.romio_ds_write) { + if (prm_.romio_ds_write) { MPI_Info_set(mpi_info, (char*)"romio_ds_write", - (char*)params.romio_ds_write); + (char*)prm_.romio_ds_write); } H5Pset_fapl_mpio(plist, comm_, mpi_info); #else @@ -172,7 +172,7 @@ public: hid_t dxpl = H5Pcreate(H5P_DATASET_XFER); H5_CHK(dxpl); #ifdef H5_HAVE_PARALLEL - if (params.use_independent_io) { + if (prm_.use_independent_io) { ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_INDEPENDENT); CE; } else { @@ -283,6 +283,7 @@ private: } private: + OutputParticlesParams prm_; Int3 wdims_, lo_, hi_; Grid_t::Kinds kinds_; MPI_Comm comm_; @@ -302,7 +303,8 @@ struct OutputParticlesHdf5 using writer_particles_type = writer_type::particles_type; using Particles = typename Mparticles::Patch; - OutputParticlesHdf5(const Grid_t& grid, const Int3& lo, const Int3& hi) + OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& prm, + Int3& lo, const Int3& hi) : lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{grid.kinds}, comm_{grid.comm()} {} @@ -469,7 +471,6 @@ struct OutputParticlesHdf5 // operator() void operator()(Mparticles& mprts, const std::string& filename, - const OutputParticlesParams& params, OutputParticlesWriterHDF5& writer) { const auto& grid = mprts.grid(); @@ -609,14 +610,11 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename, - params); + writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename); } private: - Int3 lo_; - Int3 hi_; - Int3 wdims_; + Int3 lo_, hi_, wdims_; Grid_t::Kinds kinds_; MPI_Comm comm_; }; @@ -641,7 +639,7 @@ public: : prm_{params}, lo_{params.lo}, hi_{adjust_hi(grid, params.hi)}, - writer_{lo_, hi_, grid.kinds, grid.comm()} + writer_{params, lo_, hi_, grid.kinds, grid.comm()} {} template @@ -658,8 +656,8 @@ public: snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", prm_.data_dir, prm_.basename, grid.timestep(), 0); - detail::OutputParticlesHdf5 impl{grid, lo_, hi_}; - impl(mprts, filename, prm_, writer_); + detail::OutputParticlesHdf5 impl{grid, prm_, lo_, hi_}; + impl(mprts, filename, writer_); } // FIXME, handles MparticlesVpic by conversion for now From fe6db7ed9c991a13142098006b9596d62a55ffb2 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 09:43:08 -0400 Subject: [PATCH 42/49] outp/hdf5: make filename inside of Writer --- .../output_particles_hdf5_impl.hxx | 21 +++++++++---------- 1 file changed, 10 insertions(+), 11 deletions(-) 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 56b2484108..18816a0651 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -119,7 +119,7 @@ public: void operator()(const gt::gtensor& gidx_begin, const gt::gtensor& gidx_end, size_t n_write, size_t n_off, size_t n_total, const particles_type& arr, - const std::string& filename) + int timestep) { int ierr; @@ -152,7 +152,12 @@ public: abort(); #endif - hid_t file = H5Fcreate(filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist); + int slen = strlen(prm_.data_dir) + strlen(prm_.basename) + 20; + char filename[slen]; + snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", prm_.data_dir, + prm_.basename, timestep, 0); + + hid_t file = H5Fcreate(filename, H5F_ACC_TRUNC, H5P_DEFAULT, plist); H5_CHK(file); H5Pclose(plist); MPI_Info_free(&mpi_info); @@ -470,8 +475,7 @@ struct OutputParticlesHdf5 // ---------------------------------------------------------------------- // operator() - void operator()(Mparticles& mprts, const std::string& filename, - OutputParticlesWriterHDF5& writer) + void operator()(Mparticles& mprts, OutputParticlesWriterHDF5& writer) { const auto& grid = mprts.grid(); herr_t ierr; @@ -610,7 +614,7 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, filename); + writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, grid.timestep()); } private: @@ -651,13 +655,8 @@ public: return; } - int slen = strlen(prm_.data_dir) + strlen(prm_.basename) + 20; - char filename[slen]; - snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", prm_.data_dir, - prm_.basename, grid.timestep(), 0); - detail::OutputParticlesHdf5 impl{grid, prm_, lo_, hi_}; - impl(mprts, filename, writer_); + impl(mprts, writer_); } // FIXME, handles MparticlesVpic by conversion for now From 1cb9467ae8101acacea63083b90653aa35fa0d6c Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 09:53:10 -0400 Subject: [PATCH 43/49] outp/hdf5: adjust hi initially in Params So we don't need to keep passing this otherwise redundant info around. --- .../output_particles_hdf5_impl.hxx | 44 +++++++++---------- 1 file changed, 22 insertions(+), 22 deletions(-) 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 18816a0651..74f4d58315 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -110,10 +110,9 @@ public: using particle_type = hdf5_prt; using particles_type = std::vector; - explicit OutputParticlesWriterHDF5(const OutputParticlesParams& prm, Int3& lo, - const Int3& hi, const Grid_t::Kinds& kinds, - MPI_Comm comm) - : prm_{prm}, lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{kinds}, comm_{comm} + explicit OutputParticlesWriterHDF5(const OutputParticlesParams& prm, + const Grid_t::Kinds& kinds, MPI_Comm comm) + : prm_{prm}, wdims_{prm.hi - prm.lo}, kinds_{kinds}, comm_{comm} {} void operator()(const gt::gtensor& gidx_begin, @@ -169,9 +168,9 @@ public: H5Gcreate(group, "p0", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(groupp); - ierr = H5LTset_attribute_int(group, ".", "lo", lo_, 3); + ierr = H5LTset_attribute_int(group, ".", "lo", prm_.lo, 3); CE; - ierr = H5LTset_attribute_int(group, ".", "hi", hi_, 3); + ierr = H5LTset_attribute_int(group, ".", "hi", prm_.hi, 3); CE; hid_t dxpl = H5Pcreate(H5P_DATASET_XFER); @@ -289,7 +288,7 @@ private: private: OutputParticlesParams prm_; - Int3 wdims_, lo_, hi_; + Int3 wdims_; Grid_t::Kinds kinds_; MPI_Comm comm_; Hdf5ParticleType prt_type_; @@ -308,9 +307,12 @@ struct OutputParticlesHdf5 using writer_particles_type = writer_type::particles_type; using Particles = typename Mparticles::Patch; - OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& prm, - Int3& lo, const Int3& hi) - : lo_{lo}, hi_{hi}, wdims_{hi - lo}, kinds_{grid.kinds}, comm_{grid.comm()} + OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& prm) + : lo_{prm.lo}, + hi_{prm.hi}, + wdims_{prm.hi - prm.lo}, + kinds_{grid.kinds}, + comm_{grid.comm()} {} // ---------------------------------------------------------------------- @@ -627,23 +629,23 @@ private: class OutputParticlesHdf5 : OutputParticlesBase { - static Int3 adjust_hi(const Grid_t& grid, Int3 hi) + static OutputParticlesParams adjust_params( + const OutputParticlesParams& params_in, const Grid_t& grid) { + OutputParticlesParams params = params_in; for (int d = 0; d < 3; d++) { - if (hi[d] == 0) { - hi[d] = grid.domain.gdims[d]; + if (params.hi[d] == 0) { + params.hi[d] = grid.domain.gdims[d]; } - assert(hi[d] <= grid.domain.gdims[d]); + assert(params.lo[d] >= 0); + assert(params.hi[d] <= grid.domain.gdims[d]); } - return hi; + return params; } public: OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& params) - : prm_{params}, - lo_{params.lo}, - hi_{adjust_hi(grid, params.hi)}, - writer_{params, lo_, hi_, grid.kinds, grid.comm()} + : prm_{adjust_params(params, grid)}, writer_{prm_, grid.kinds, grid.comm()} {} template @@ -655,7 +657,7 @@ public: return; } - detail::OutputParticlesHdf5 impl{grid, prm_, lo_, hi_}; + detail::OutputParticlesHdf5 impl{grid, prm_}; impl(mprts, writer_); } @@ -693,7 +695,5 @@ public: private: const OutputParticlesParams prm_; - Int3 lo_; // dimensions of the subdomain we're actually writing - Int3 hi_; OutputParticlesWriterHDF5 writer_; }; From a930803705ecb2bad214dd5a04bfa0d1c8f40cd0 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 09:54:25 -0400 Subject: [PATCH 44/49] outp/hdf5: Consistently use params instead of prm --- .../output_particles_hdf5_impl.hxx | 47 ++++++++++--------- 1 file changed, 24 insertions(+), 23 deletions(-) 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 74f4d58315..0506cd2dda 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -110,9 +110,9 @@ public: using particle_type = hdf5_prt; using particles_type = std::vector; - explicit OutputParticlesWriterHDF5(const OutputParticlesParams& prm, + explicit OutputParticlesWriterHDF5(const OutputParticlesParams& params, const Grid_t::Kinds& kinds, MPI_Comm comm) - : prm_{prm}, wdims_{prm.hi - prm.lo}, kinds_{kinds}, comm_{comm} + : params_{params}, wdims_{params.hi - params.lo}, kinds_{kinds}, comm_{comm} {} void operator()(const gt::gtensor& gidx_begin, @@ -137,13 +137,13 @@ public: MPI_Info_create(&mpi_info); #ifdef H5_HAVE_PARALLEL - if (prm_.romio_cb_write) { + if (params_.romio_cb_write) { MPI_Info_set(mpi_info, (char*)"romio_cb_write", - (char*)prm_.romio_cb_write); + (char*)params_.romio_cb_write); } - if (prm_.romio_ds_write) { + if (params_.romio_ds_write) { MPI_Info_set(mpi_info, (char*)"romio_ds_write", - (char*)prm_.romio_ds_write); + (char*)params_.romio_ds_write); } H5Pset_fapl_mpio(plist, comm_, mpi_info); #else @@ -151,10 +151,10 @@ public: abort(); #endif - int slen = strlen(prm_.data_dir) + strlen(prm_.basename) + 20; + int slen = strlen(params_.data_dir) + strlen(params_.basename) + 20; char filename[slen]; - snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", prm_.data_dir, - prm_.basename, timestep, 0); + snprintf(filename, slen, "%s/%s.%06d_p%06d.h5", params_.data_dir, + params_.basename, timestep, 0); hid_t file = H5Fcreate(filename, H5F_ACC_TRUNC, H5P_DEFAULT, plist); H5_CHK(file); @@ -168,15 +168,15 @@ public: H5Gcreate(group, "p0", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); H5_CHK(groupp); - ierr = H5LTset_attribute_int(group, ".", "lo", prm_.lo, 3); + ierr = H5LTset_attribute_int(group, ".", "lo", params_.lo, 3); CE; - ierr = H5LTset_attribute_int(group, ".", "hi", prm_.hi, 3); + ierr = H5LTset_attribute_int(group, ".", "hi", params_.hi, 3); CE; hid_t dxpl = H5Pcreate(H5P_DATASET_XFER); H5_CHK(dxpl); #ifdef H5_HAVE_PARALLEL - if (prm_.use_independent_io) { + if (params_.use_independent_io) { ierr = H5Pset_dxpl_mpio(dxpl, H5FD_MPIO_INDEPENDENT); CE; } else { @@ -287,7 +287,7 @@ private: } private: - OutputParticlesParams prm_; + OutputParticlesParams params_; Int3 wdims_; Grid_t::Kinds kinds_; MPI_Comm comm_; @@ -307,10 +307,10 @@ struct OutputParticlesHdf5 using writer_particles_type = writer_type::particles_type; using Particles = typename Mparticles::Patch; - OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& prm) - : lo_{prm.lo}, - hi_{prm.hi}, - wdims_{prm.hi - prm.lo}, + OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& params) + : lo_{params.lo}, + hi_{params.hi}, + wdims_{params.hi - params.lo}, kinds_{grid.kinds}, comm_{grid.comm()} {} @@ -645,7 +645,8 @@ class OutputParticlesHdf5 : OutputParticlesBase public: OutputParticlesHdf5(const Grid_t& grid, const OutputParticlesParams& params) - : prm_{adjust_params(params, grid)}, writer_{prm_, grid.kinds, grid.comm()} + : params_{adjust_params(params, grid)}, + writer_{params_, grid.kinds, grid.comm()} {} template @@ -653,11 +654,11 @@ public: { const auto& grid = mprts.grid(); - if (prm_.every_step <= 0 || grid.timestep() % prm_.every_step != 0) { + if (params_.every_step <= 0 || grid.timestep() % params_.every_step != 0) { return; } - detail::OutputParticlesHdf5 impl{grid, prm_}; + detail::OutputParticlesHdf5 impl{grid, params_}; impl(mprts, writer_); } @@ -667,7 +668,7 @@ public: { const auto& grid = mprts.grid(); - if (prm_.every_step <= 0 || grid.timestep() % prm_.every_step != 0) { + if (params_.every_step <= 0 || grid.timestep() % params_.every_step != 0) { return; } @@ -683,7 +684,7 @@ public: { const auto& grid = mprts.grid(); - if (prm_.every_step <= 0 || grid.timestep() % prm_.every_step != 0) { + if (params_.every_step <= 0 || grid.timestep() % params_.every_step != 0) { return; } @@ -694,6 +695,6 @@ public: #endif private: - const OutputParticlesParams prm_; + const OutputParticlesParams params_; OutputParticlesWriterHDF5 writer_; }; From 23d2d902f348f2a02779bf48161111bc363b7683 Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 22 Aug 2024 10:01:29 -0400 Subject: [PATCH 45/49] outp/hdf5: rename {rank,size} -> mpi_{rank,size} --- .../output_particles_hdf5_impl.hxx | 36 +++++++++---------- 1 file changed, 18 insertions(+), 18 deletions(-) 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 0506cd2dda..6a116d7d9a 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -221,9 +221,9 @@ private: hid_t memspace; assert(sizeof(size_t) == sizeof(hsize_t)); - int rank; - MPI_Comm_rank(comm_, &rank); - if (rank == 0) { + int mpi_rank; + MPI_Comm_rank(comm_, &mpi_rank); + if (mpi_rank == 0) { memspace = H5Screate_simple(4, fdims, NULL); H5_CHK(memspace); } else { @@ -490,9 +490,9 @@ struct OutputParticlesHdf5 } prof_start(pr_A); - int rank, size; - MPI_Comm_rank(comm_, &rank); - MPI_Comm_size(comm_, &size); + int mpi_rank, mpi_size; + MPI_Comm_rank(comm_, &mpi_rank); + MPI_Comm_size(comm_, &mpi_size); struct mrc_patch_info info; int nr_kinds = mprts.grid().kinds.size(); @@ -504,7 +504,7 @@ struct OutputParticlesHdf5 count_sort(mprts, off, map); gt::gtensor gidx_begin, gidx_end; - if (rank == 0) { + if (mpi_rank == 0) { // alloc global idx array gidx_begin = gt::empty({wdims_[0], wdims_[1], wdims_[2], nr_kinds}); @@ -518,13 +518,13 @@ struct OutputParticlesHdf5 prof_stop(pr_A); prof_start(pr_B); - if (rank == 0) { + if (mpi_rank == 0) { int nr_global_patches = grid.nGlobalPatches(); - std::vector remote_sz(size); + std::vector remote_sz(mpi_size); for (int p = 0; p < nr_global_patches; p++) { info = grid.globalPatchInfo(p); - if (info.rank == rank) { // skip local patches + if (info.rank == mpi_rank) { // skip local patches continue; } int ilo[3], ihi[3], ld[3]; @@ -553,11 +553,11 @@ struct OutputParticlesHdf5 } } - std::vector> recv_buf(size); - std::vector recv_buf_p(size); - std::vector recv_req(size, MPI_REQUEST_NULL); - for (int r = 0; r < size; r++) { - if (r == rank) { // skip this proc + std::vector> recv_buf(mpi_size); + std::vector recv_buf_p(mpi_size); + std::vector recv_req(mpi_size, MPI_REQUEST_NULL); + for (int r = 0; r < mpi_size; r++) { + if (r == mpi_rank) { // skip this proc continue; } recv_buf[r].resize(2 * remote_sz[r]); @@ -565,12 +565,12 @@ struct OutputParticlesHdf5 MPI_Irecv(recv_buf[r].data(), recv_buf[r].size(), MPI_UNSIGNED_LONG, r, 0x4000, comm_, &recv_req[r]); } - MPI_Waitall(size, recv_req.data(), MPI_STATUSES_IGNORE); + MPI_Waitall(mpi_size, recv_req.data(), MPI_STATUSES_IGNORE); // build global idx array, remote part for (int p = 0; p < nr_global_patches; p++) { info = grid.globalPatchInfo(p); - if (info.rank == rank) { // skip local patches + if (info.rank == mpi_rank) { // skip local patches continue; } int ilo[3], ihi[3], ld[3]; @@ -593,7 +593,7 @@ struct OutputParticlesHdf5 } recv_buf_p[info.rank] += 2 * sz; } - } else { // rank != 0: send to rank 0 + } else { // mpi_rank != 0: send to mpi_rank 0 // FIXME, alloc idx[] as one array in the first place int local_sz = 0; for (int p = 0; p < mprts.n_patches(); p++) { From cf2f0f20edbdd0c850c04368fa2323e1ffd8bf4c Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Fri, 23 Aug 2024 07:34:33 -0400 Subject: [PATCH 46/49] outp/hdf5: don't need to count all particles Since we already know how many particles we have on a proc. --- .../output_particles_hdf5_impl.hxx | 17 +---------------- 1 file changed, 1 insertion(+), 16 deletions(-) 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 6a116d7d9a..920c1bb2e2 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -415,22 +415,7 @@ struct OutputParticlesHdf5 int nr_kinds = grid.kinds.size(); // count all particles to be written locally - size_t n_write = 0; - for (int p = 0; p < mprts.n_patches(); p++) { - auto& patch = grid.patches[p]; - int ilo[3], ihi[3], ld[3]; - int sz = find_patch_bounds(grid.ldims, patch.off, ilo, ihi, ld); - for (int jz = ilo[2]; jz < ihi[2]; jz++) { - for (int jy = ilo[1]; jy < ihi[1]; jy++) { - for (int jx = ilo[0]; jx < ihi[0]; jx++) { - for (int kind = 0; kind < nr_kinds; kind++) { - int si = sort_index(grid.ldims, nr_kinds, {jx, jy, jz}, kind); - n_write += off[p][si + 1] - off[p][si]; - } - } - } - } - } + size_t n_write = mprts.size(); assert(sizeof(size_t) == sizeof(unsigned long)); size_t n_total, n_off = 0; From 6b5170ef653719a3ad6de17bb9dad9fc7e63542a Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Fri, 23 Aug 2024 07:38:46 -0400 Subject: [PATCH 47/49] outp/hdf5: don't pass redundant n_write The number of particles is arr.size() --- .../output_particles_hdf5_impl.hxx | 26 +++++++++---------- 1 file changed, 12 insertions(+), 14 deletions(-) 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 920c1bb2e2..95d591fe4b 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -116,9 +116,8 @@ public: {} void operator()(const gt::gtensor& gidx_begin, - const gt::gtensor& gidx_end, size_t n_write, - size_t n_off, size_t n_total, const particles_type& arr, - int timestep) + const gt::gtensor& gidx_end, size_t n_off, + size_t n_total, const particles_type& arr, int timestep) { int ierr; @@ -191,7 +190,7 @@ public: prof_stop(pr_D); prof_start(pr_E); - write_particles(n_write, n_off, n_total, arr, groupp, dxpl); + write_particles(n_off, n_total, arr, groupp, dxpl); prof_stop(pr_E); ierr = H5Pclose(dxpl); @@ -256,12 +255,12 @@ private: CE; } - void write_particles(size_t n_write, size_t n_off, size_t n_total, - const particles_type& arr, hid_t group, hid_t dxpl) + void write_particles(size_t n_off, size_t n_total, const particles_type& arr, + hid_t group, hid_t dxpl) { herr_t ierr; - hsize_t mdims[1] = {n_write}; + hsize_t mdims[1] = {arr.size()}; hsize_t fdims[1] = {n_total}; hsize_t foff[1] = {n_off}; hid_t memspace = H5Screate_simple(1, mdims, NULL); @@ -408,8 +407,8 @@ struct OutputParticlesHdf5 writer_particles_type make_local_particle_array( Mparticles& mprts, const std::vector>& off, const std::vector>& map, - std::vector>& idx, size_t* p_n_write, - size_t* p_n_off, size_t* p_n_total) + std::vector>& idx, size_t* p_n_off, + size_t* p_n_total) { const auto& grid = mprts.grid(); int nr_kinds = grid.kinds.size(); @@ -453,7 +452,6 @@ struct OutputParticlesHdf5 } } } - *p_n_write = n_write; *p_n_off = n_off; *p_n_total = n_total; return arr; @@ -497,9 +495,9 @@ struct OutputParticlesHdf5 } // find local particle and idx arrays - size_t n_write, n_off, n_total; - auto arr = make_local_particle_array(mprts, off, map, idx, &n_write, &n_off, - &n_total); + size_t n_off, n_total; + auto arr = + make_local_particle_array(mprts, off, map, idx, &n_off, &n_total); prof_stop(pr_A); prof_start(pr_B); @@ -601,7 +599,7 @@ struct OutputParticlesHdf5 } prof_stop(pr_B); - writer(gidx_begin, gidx_end, n_write, n_off, n_total, arr, grid.timestep()); + writer(gidx_begin, gidx_end, n_off, n_total, arr, grid.timestep()); } private: From a1d0010d5f313adf17a7022e5a0451655d3b3a7e Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 29 Aug 2024 11:32:57 -0400 Subject: [PATCH 48/49] outp/hdf5: loop over particle by index --- .../output_particles_hdf5_impl.hxx | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) 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 95d591fe4b..cc2c2ae058 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -325,11 +325,11 @@ struct OutputParticlesHdf5 kind; }; - template - static inline int get_sort_index(Particles& prts, const Particle& prt) + static inline int get_sort_index(Particles& prts, int n) { const Grid_t& grid = prts.grid(); const int* ldims = grid.ldims; + const auto& prt = prts[n]; Int3 pos; for (int d = 0; d < 3; d++) { @@ -361,8 +361,8 @@ struct OutputParticlesHdf5 unsigned int n_prts = prts.size(); // counting sort to get map - for (const auto& prt : prts) { - int si = get_sort_index(prts, prt); + for (int n = 0; n < n_prts; n++) { + int si = get_sort_index(prts, n); off[p][si]++; } // prefix sum to get offsets @@ -377,10 +377,9 @@ struct OutputParticlesHdf5 // sort a map only, not the actual particles map[p].resize(n_prts); - int n = 0; - for (const auto& prt : prts) { - int si = get_sort_index(prts, prt); - map[p][off2[si]++] = n++; + for (int n = 0; n < n_prts; n++) { + int si = get_sort_index(prts, n); + map[p][off2[si]++] = n; } } } From 51224ff34ad062376c3ec11a23dd1fcea99063eb Mon Sep 17 00:00:00 2001 From: Kai Germaschewski Date: Thu, 29 Aug 2024 12:05:04 -0400 Subject: [PATCH 49/49] outp/hdf5: add template argument for particle selector This particle selector decides whether to include a given particle in the output. Other than the trivial ParticleSelectorAll, ParticleSelectorEveryNth is provided ot select every nth particle. --- .../output_particles_hdf5_impl.hxx | 44 ++++++++++++++++--- src/libpsc/tests/test_output_particles.cxx | 12 ++--- src/psc_config.hxx | 4 +- src/psc_flatfoil_yz.cxx | 4 ++ 4 files changed, 51 insertions(+), 13 deletions(-) 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 cc2c2ae058..bf92547279 100644 --- a/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx +++ b/src/libpsc/psc_output_particles/output_particles_hdf5_impl.hxx @@ -293,13 +293,38 @@ private: Hdf5ParticleType prt_type_; }; +template +class ParticleSelectorEveryNth +{ +public: + template + bool operator()(const Particle& prt) + { + // return true for every `every_`th particle + return (cnt_++ % EVERY) == 0; + } + +private: + int cnt_ = 0; +}; + +class ParticleSelectorAll +{ +public: + template + bool operator()(const Particle& prt) + { + return true; + } +}; + namespace detail { // ====================================================================== // OutputParticlesHdf5 -template +template struct OutputParticlesHdf5 { using writer_type = OutputParticlesWriterHDF5; @@ -353,17 +378,24 @@ struct OutputParticlesHdf5 { int nr_kinds = mprts.grid().kinds.size(); + ParticleSelector selector; + for (int p = 0; p < mprts.n_patches(); p++) { const int* ldims = mprts.grid().ldims; int nr_indices = ldims[0] * ldims[1] * ldims[2] * nr_kinds; off[p].resize(nr_indices + 1); auto&& prts = mprts[p]; unsigned int n_prts = prts.size(); + std::vector particle_indices; + particle_indices.reserve(n_prts); // counting sort to get map for (int n = 0; n < n_prts; n++) { - int si = get_sort_index(prts, n); - off[p][si]++; + if (selector(prts[n])) { + particle_indices.push_back(n); + int si = get_sort_index(prts, n); + off[p][si]++; + } } // prefix sum to get offsets int o = 0; @@ -377,7 +409,7 @@ struct OutputParticlesHdf5 // sort a map only, not the actual particles map[p].resize(n_prts); - for (int n = 0; n < n_prts; n++) { + for (auto n : particle_indices) { int si = get_sort_index(prts, n); map[p][off2[si]++] = n; } @@ -609,6 +641,7 @@ private: } // namespace detail +template class OutputParticlesHdf5 : OutputParticlesBase { static OutputParticlesParams adjust_params( @@ -640,7 +673,8 @@ public: return; } - detail::OutputParticlesHdf5 impl{grid, params_}; + detail::OutputParticlesHdf5 impl{grid, + params_}; impl(mprts, writer_); } diff --git a/src/libpsc/tests/test_output_particles.cxx b/src/libpsc/tests/test_output_particles.cxx index 8c200dad09..ff68aa0777 100644 --- a/src/libpsc/tests/test_output_particles.cxx +++ b/src/libpsc/tests/test_output_particles.cxx @@ -64,14 +64,14 @@ struct OutputParticlesTest : ::testing::Test Int3 ibn = {2, 2, 2}; }; -using OutputParticlesTestTypes = - ::testing::Types +using OutputParticlesTestTypes = ::testing::Types< + Config #ifdef H5_HAVE_PARALLEL - , - Config, - Config + , + Config>, + Config> #endif - >; + >; TYPED_TEST_SUITE(OutputParticlesTest, OutputParticlesTestTypes); diff --git a/src/psc_config.hxx b/src/psc_config.hxx index fcf01cd29a..6401a481e0 100644 --- a/src/psc_config.hxx +++ b/src/psc_config.hxx @@ -31,7 +31,7 @@ #include "../libpsc/vpic/fields_item_vpic.hxx" -using OutputParticlesDefault = OutputParticlesHdf5; +using OutputParticlesDefault = OutputParticlesHdf5; struct SimulationNone { @@ -252,7 +252,7 @@ struct PscConfigVpicPsc using BndParticles = BndParticlesVpic; using Checks = ChecksVpic; using Marder = MarderVpic; - using OutputParticles = OutputParticlesHdf5; + using OutputParticles = OutputParticlesHdf5; using OutputHydro = OutputHydroQ; using Dim = dim_xz; diff --git a/src/psc_flatfoil_yz.cxx b/src/psc_flatfoil_yz.cxx index a853acc059..342e520100 100644 --- a/src/psc_flatfoil_yz.cxx +++ b/src/psc_flatfoil_yz.cxx @@ -226,7 +226,11 @@ using Balance = PscConfig::Balance; using Collision = PscConfig::Collision; using Checks = PscConfig::Checks; using Marder = PscConfig::Marder; +#if CASE == CASE_2D_SMALL +using OutputParticles = OutputParticlesHdf5>; +#else using OutputParticles = PscConfig::OutputParticles; +#endif using Moment_n = typename Moment_n_Selector::type; using Heating = typename HeatingSelector::Heating;