From 5ed7c6b7abef0effb95937b6edc2e7db51810ce9 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 13:41:16 -0400 Subject: [PATCH 01/29] +test_bnd_prt_inflow: copy from test_mparticles for now --- src/libpsc/tests/CMakeLists.txt | 1 + src/libpsc/tests/test_bnd_prt_inflow.cxx | 243 +++++++++++++++++++++++ 2 files changed, 244 insertions(+) create mode 100644 src/libpsc/tests/test_bnd_prt_inflow.cxx diff --git a/src/libpsc/tests/CMakeLists.txt b/src/libpsc/tests/CMakeLists.txt index 913d29da9e..1c29b8f494 100644 --- a/src/libpsc/tests/CMakeLists.txt +++ b/src/libpsc/tests/CMakeLists.txt @@ -24,6 +24,7 @@ add_psc_test(test_output_particles) add_psc_test(test_mfields) add_psc_test(test_mfields_cuda) add_psc_test(test_bnd) +add_psc_test(test_bnd_prt_inflow) add_psc_test(test_deposit) add_psc_test(test_current_deposition) add_psc_test(test_push_particles) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx new file mode 100644 index 0000000000..57e8260a1e --- /dev/null +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -0,0 +1,243 @@ + +#include + +#include "test_common.hxx" + +#include "../libpsc/vpic/PscGridBase.h" +#include "../libpsc/vpic/PscParticleBc.h" +#include "../libpsc/vpic/PscParticlesBase.h" +#include "../libpsc/vpic/mparticles_vpic.hxx" +#include "../libpsc/vpic/vpic_config.h" +#include "psc_particles_double.h" +#include "psc_particles_single.h" +#include "particle_with_id.h" +#include "setup_particles.hxx" +#ifdef USE_CUDA +#include "../libpsc/cuda/mparticles_cuda.hxx" +#include "../libpsc/cuda/mparticles_cuda.inl" +#endif +#include "particles_simple.inl" +#include + +#ifdef DO_VPIC +using VpicConfig = VpicConfigWrap; +#else +using VpicConfig = VpicConfigPsc; +#endif + +using Grid = VpicConfig::Grid; +using MparticlesVpic = VpicConfig::Mparticles; + +template +struct Config +{ + using Mparticles = _Mparticles; + using MakeGrid = _MakeGrid; +}; + +using MparticlesInflowTestTypes = ::testing::Types< + Config, Config, + Config, Config +#ifdef USE_CUDA + , + Config, MakeTestGridYZ1>, + Config, MakeTestGridYZ> +#endif + >; + +TYPED_TEST_SUITE(MparticlesInflowTest, MparticlesInflowTestTypes); + +// ====================================================================== +// MparticlesInflowTest + +template +struct MparticlesInflowTest : ::testing::Test +{ + using Mparticles = typename T::Mparticles; + using Particle = typename Mparticles::Particle; + using MakeGrid = typename T::MakeGrid; + + MparticlesInflowTest() : grid_{MakeGrid{}()} + { + grid_.kinds.emplace_back(Grid_t::Kind(1., 1., "test_species")); + } + + template + Mparticles mk_mprts(tag dummy) + { + Mparticles mprts(grid_); + mprts.define_species("test_species", 1., 1., 100, 10, 10, 0); + return mprts; + } + + Mparticles mk_mprts(MparticlesVpic* dummy) + { + // FIXME, vgrid_ is bad, and this kinda belongs to where grid_ is set up + + // Setup basic grid parameters + auto& domain = grid_.domain; + double dx[3], xl[3], xh[3]; + for (int d = 0; d < 3; d++) { + dx[d] = domain.length[d] / domain.gdims[d]; + xl[d] = domain.corner[d]; + xh[d] = xl[d] + domain.length[d]; + } + vgrid_.setup(dx, grid_.dt, 1., 1.); + + // Define the grid + vgrid_.partition_periodic_box(xl, xh, domain.gdims, domain.np); + + Mparticles mprts(grid_, &vgrid_); + mprts.define_species("test_species", 1., 1., 100, 10, 10, 0); + return mprts; + } + + Mparticles mk_mprts() { return mk_mprts(static_cast(nullptr)); } + + template + void inject_test_particles(_Mparticles& mprts, int n_prts) + { + auto inj = mprts.injector(); + for (int p = 0; p < mprts.n_patches(); ++p) { + auto injector = inj[p]; + auto& patch = mprts.grid().patches[p]; + for (int n = 0; n < n_prts; n++) { + double nn = double(n) / n_prts; + auto L = patch.xe - patch.xb; + psc::particle::Inject prt = {{patch.xb[0] + nn * L[0], + patch.xb[1] + nn * L[1], + patch.xb[2] + nn * L[2]}, + {}, + 1., + 0}; + injector(prt); + } + } + } + + const Grid_t& grid() { return grid_; } + +private: + Grid_t grid_; + Grid vgrid_; // FIXME +}; + +// ---------------------------------------------------------------------- +// Inject + +TYPED_TEST(MparticlesInflowTest, Inject) +{ + const int n_prts = 4; + + auto mprts = this->mk_mprts(); + + this->inject_test_particles(mprts, n_prts); +} + +// ----------------------------------------------------------------------- +// Inject2 + +TYPED_TEST(MparticlesInflowTest, Inject2) +{ + using Mparticles = typename TypeParam::Mparticles; + const int n_prts = 1; + + auto mprts = this->mk_mprts(); + + // FIXME, kinda the cheap way out + if (mprts.n_patches() != 4) { + return; + } + ASSERT_EQ(mprts.n_patches(), 4); + + int nn = 0; + { + auto inj = mprts.injector(); + for (int p = 0; p < mprts.n_patches(); ++p) { + auto injector = inj[p]; + auto& patch = mprts.grid().patches[p]; + for (int n = 0; n < n_prts; n++) { + psc::particle::Inject prt{ + .5 * (patch.xb + patch.xe), {}, double(nn), 0}; + // use weight to store particle number for testing + injector(prt); + nn++; + } + } + } + + EXPECT_EQ(mprts.size(), 4); + + auto n_prts_by_patch = mprts.sizeByPatch(); + EXPECT_EQ(n_prts_by_patch, + std::vector({n_prts, n_prts, n_prts, n_prts})); + + nn = 0; + auto accessor = mprts.accessor(); + for (int p = 0; p < mprts.n_patches(); ++p) { + auto& patch = mprts.grid().patches[p]; + for (auto prt : accessor[p]) { + // xm is patch-relative position + auto xm = .5 * (patch.xe - patch.xb); + EXPECT_EQ(prt.x()[0], xm[0]); + EXPECT_EQ(prt.x()[1], xm[1]); + EXPECT_EQ(prt.x()[2], xm[2]); + EXPECT_EQ(prt.qni_wni(), nn); + EXPECT_EQ(prt.w(), nn); + EXPECT_EQ(prt.kind(), 0); + + auto x = .5 * (patch.xb + patch.xe); + EXPECT_EQ(prt.position()[0], x[0]); + EXPECT_EQ(prt.position()[1], x[1]); + EXPECT_EQ(prt.position()[2], x[2]); + nn++; + } + } +} + +// ---------------------------------------------------------------------- +// Injector + +TYPED_TEST(MparticlesInflowTest, Injector) +{ + const int n_prts = 4; + + auto mprts = this->mk_mprts(); + + this->inject_test_particles(mprts, n_prts); + + auto accessor = mprts.accessor(); + for (int p = 0; p < mprts.n_patches(); ++p) { + auto prts = accessor[p]; + auto& patch = mprts.grid().patches[p]; + EXPECT_EQ(prts.size(), n_prts); + int n = 0; + for (auto prt : prts) { + double nn = double(n) / n_prts; + auto L = patch.xe - patch.xb; + auto x = prt.position(); + EXPECT_EQ(x[0], patch.xb[0] + nn * L[0]); + EXPECT_EQ(x[1], patch.xb[1] + nn * L[1]); + EXPECT_EQ(x[2], patch.xb[2] + nn * L[2]); + EXPECT_EQ(prt.w(), 1.); + EXPECT_EQ(prt.kind(), 0); + ++n; + } + } +} + +int main(int argc, char** argv) +{ + MPI_Init(&argc, &argv); + // #ifdef USE_VPIC FIXME + MPI_Comm_dup(MPI_COMM_WORLD, &psc_comm_world); + MPI_Comm_rank(psc_comm_world, &psc_world_rank); + MPI_Comm_size(psc_comm_world, &psc_world_size); + // #endif + + ::testing::InitGoogleTest(&argc, argv); + int rc = RUN_ALL_TESTS(); + + MPI_Finalize(); + return rc; +} From 48368e14d5d3fceaf60814046a14efc81c566f92 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 13:44:35 -0400 Subject: [PATCH 02/29] test_bnd_prt_inflow: rm vpic stuff --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 56 +++--------------------- 1 file changed, 7 insertions(+), 49 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 57e8260a1e..6f03aa51be 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -3,11 +3,6 @@ #include "test_common.hxx" -#include "../libpsc/vpic/PscGridBase.h" -#include "../libpsc/vpic/PscParticleBc.h" -#include "../libpsc/vpic/PscParticlesBase.h" -#include "../libpsc/vpic/mparticles_vpic.hxx" -#include "../libpsc/vpic/vpic_config.h" #include "psc_particles_double.h" #include "psc_particles_single.h" #include "particle_with_id.h" @@ -19,15 +14,6 @@ #include "particles_simple.inl" #include -#ifdef DO_VPIC -using VpicConfig = VpicConfigWrap; -#else -using VpicConfig = VpicConfigPsc; -#endif - -using Grid = VpicConfig::Grid; -using MparticlesVpic = VpicConfig::Mparticles; - template struct Config { @@ -35,15 +21,15 @@ struct Config using MakeGrid = _MakeGrid; }; -using MparticlesInflowTestTypes = ::testing::Types< - Config, Config, - Config, Config +using MparticlesInflowTestTypes = + ::testing::Types, + Config #ifdef USE_CUDA - , - Config, MakeTestGridYZ1>, - Config, MakeTestGridYZ> + , + Config, MakeTestGridYZ1>, + Config, MakeTestGridYZ> #endif - >; + >; TYPED_TEST_SUITE(MparticlesInflowTest, MparticlesInflowTestTypes); @@ -70,28 +56,6 @@ struct MparticlesInflowTest : ::testing::Test return mprts; } - Mparticles mk_mprts(MparticlesVpic* dummy) - { - // FIXME, vgrid_ is bad, and this kinda belongs to where grid_ is set up - - // Setup basic grid parameters - auto& domain = grid_.domain; - double dx[3], xl[3], xh[3]; - for (int d = 0; d < 3; d++) { - dx[d] = domain.length[d] / domain.gdims[d]; - xl[d] = domain.corner[d]; - xh[d] = xl[d] + domain.length[d]; - } - vgrid_.setup(dx, grid_.dt, 1., 1.); - - // Define the grid - vgrid_.partition_periodic_box(xl, xh, domain.gdims, domain.np); - - Mparticles mprts(grid_, &vgrid_); - mprts.define_species("test_species", 1., 1., 100, 10, 10, 0); - return mprts; - } - Mparticles mk_mprts() { return mk_mprts(static_cast(nullptr)); } template @@ -119,7 +83,6 @@ struct MparticlesInflowTest : ::testing::Test private: Grid_t grid_; - Grid vgrid_; // FIXME }; // ---------------------------------------------------------------------- @@ -229,11 +192,6 @@ TYPED_TEST(MparticlesInflowTest, Injector) int main(int argc, char** argv) { MPI_Init(&argc, &argv); - // #ifdef USE_VPIC FIXME - MPI_Comm_dup(MPI_COMM_WORLD, &psc_comm_world); - MPI_Comm_rank(psc_comm_world, &psc_world_rank); - MPI_Comm_size(psc_comm_world, &psc_world_size); - // #endif ::testing::InitGoogleTest(&argc, argv); int rc = RUN_ALL_TESTS(); From 84fb632b034c890c1a6b6d9c8b108df768f7d4b9 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 13:46:35 -0400 Subject: [PATCH 03/29] test_bnd_prt_inflow: bring back TestSetupParticles --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 6f03aa51be..7b9579d464 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -189,6 +189,26 @@ TYPED_TEST(MparticlesInflowTest, Injector) } } +TEST(TestSetupParticlesInflow, Simple) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 1, 1}}; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, .1}; + Mparticles mprts{grid}; + + SetupParticles setup_particles(grid); + setup_particles( + mprts, [&](int kind, Double3 crd, psc_particle_npt& npt) { npt.n = 1; }); + + auto n_cells = + grid.domain.gdims[0] * grid.domain.gdims[1] * grid.domain.gdims[2]; + EXPECT_EQ(mprts.size(), n_cells * kinds.size() * prm.nicell); +} + int main(int argc, char** argv) { MPI_Init(&argc, &argv); From 5557cca987c7a0d47c4ebd7e377626c96fcef1d2 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 14:00:27 -0400 Subject: [PATCH 04/29] particle: impl operator<< for Inject --- src/include/particle.h | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/src/include/particle.h b/src/include/particle.h index 225addd9ac..d07b1e918e 100644 --- a/src/include/particle.h +++ b/src/include/particle.h @@ -37,3 +37,11 @@ struct Inject } // namespace particle } // namespace psc + +inline std::ostream& operator<<(std::ostream& os, + const psc::particle::Inject& inj) +{ + os << "Inject{x=" << inj.x << ", u=" << inj.u << ", w=" << inj.w + << ", kind=" << inj.kind << ", tag=" << inj.tag << "}"; + return os; +} \ No newline at end of file From 5bab7aac411c319b1196f29e009d5d3afb76671b Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 14:07:37 -0400 Subject: [PATCH 05/29] test_bnd_prt_inflow: test prt from setupParticle --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 17 ++++++++++++----- 1 file changed, 12 insertions(+), 5 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 7b9579d464..94435f8405 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -201,12 +201,19 @@ TEST(TestSetupParticlesInflow, Simple) Mparticles mprts{grid}; SetupParticles setup_particles(grid); - setup_particles( - mprts, [&](int kind, Double3 crd, psc_particle_npt& npt) { npt.n = 1; }); - auto n_cells = - grid.domain.gdims[0] * grid.domain.gdims[1] * grid.domain.gdims[2]; - EXPECT_EQ(mprts.size(), n_cells * kinds.size() * prm.nicell); + Double3 pos = {5., 5., 5.}; + Double3 u = Double3{0.0, 0.0, 0.0}; + int tag = 0; + psc_particle_np np = {0, 1.0, [u]() { return u; }, tag}; + double wni = 1.0; + auto prt = setup_particles.setupParticle(np, pos, wni); + + EXPECT_EQ(prt.x, pos); + EXPECT_EQ(prt.u, u); + EXPECT_EQ(prt.w, wni); + EXPECT_EQ(prt.kind, 0); + EXPECT_EQ(prt.tag, tag); } int main(int argc, char** argv) From da2d632ee1fa5fc404afc48adcb790755b0a0def Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 14:16:21 -0400 Subject: [PATCH 06/29] setup_particles: use Double3 in npt --- src/include/setup_particles.hxx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/include/setup_particles.hxx b/src/include/setup_particles.hxx index bad72d12c0..632e4a63d0 100644 --- a/src/include/setup_particles.hxx +++ b/src/include/setup_particles.hxx @@ -9,10 +9,10 @@ struct psc_particle_npt { - int kind; ///< particle kind - double n; ///< density - double p[3]; ///< momentum - double T[3]; ///< temperature + int kind; ///< particle kind + double n; ///< density + Double3 p; ///< momentum + Double3 T; ///< temperature psc::particle::Tag tag; }; From cd8bb5271342eac09b156e2d8ee6610c94957b8c Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 14:18:15 -0400 Subject: [PATCH 07/29] test_bnd_prt_inflow: test createMaxwellian --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 31 ++++++++++++++++++++++++ 1 file changed, 31 insertions(+) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 94435f8405..7115c69d13 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -216,6 +216,37 @@ TEST(TestSetupParticlesInflow, Simple) EXPECT_EQ(prt.tag, tag); } +TEST(TestSetupParticlesInflow, Maxwellian) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 1, 1}}; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, .1}; + Mparticles mprts{grid}; + + SetupParticles setup_particles(grid); + + Double3 pos = {5., 5., 5.}; + Double3 u = {0.0, 0.5, 0.0}; + Double3 T = {0.0, 0.0, 0.0}; + int tag = 0; + + psc_particle_npt npt = {0, 1.0, u, T, tag}; + auto p = setup_particles.createMaxwellian(npt); + psc_particle_np np = {0, 1.0, p, tag}; + double wni = 1.0; + auto prt = setup_particles.setupParticle(np, pos, wni); + + EXPECT_EQ(prt.x, pos); + EXPECT_EQ(prt.u, u); + EXPECT_EQ(prt.w, wni); + EXPECT_EQ(prt.kind, 0); + EXPECT_EQ(prt.tag, tag); +} + int main(int argc, char** argv) { MPI_Init(&argc, &argv); From 1357f4eda6030eb0fa8adc04709c794f07f5ecca Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 14:33:31 -0400 Subject: [PATCH 08/29] test_bnd_prt_inflow: test AdvanceParticle.push_x --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 37 ++++++++++++++++++++++++ 1 file changed, 37 insertions(+) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 7115c69d13..8ac62785a9 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -7,6 +7,8 @@ #include "psc_particles_single.h" #include "particle_with_id.h" #include "setup_particles.hxx" +#include "pushp.hxx" +#include "dim.hxx" #ifdef USE_CUDA #include "../libpsc/cuda/mparticles_cuda.hxx" #include "../libpsc/cuda/mparticles_cuda.inl" @@ -247,6 +249,41 @@ TEST(TestSetupParticlesInflow, Maxwellian) EXPECT_EQ(prt.tag, tag); } +TEST(TestSetupParticlesInflow, Advance) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 1, 1}}; + double dt = 10.; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, dt}; + Mparticles mprts{grid}; + + SetupParticles setup_particles(grid); + + Double3 pos = {5., -2., 5.}; + Double3 u = {0.0, 0.5, 0.0}; + Double3 T = {0.0, 0.0, 0.0}; + int tag = 0; + + psc_particle_npt npt = {0, 1.0, u, T, tag}; + auto p = setup_particles.createMaxwellian(npt); + psc_particle_np np = {0, 1.0, p, tag}; + double wni = 1.0; + auto prt = setup_particles.setupParticle(np, pos, wni); + + auto advance = AdvanceParticle(dt); + + auto v = advance.calc_v(prt.u); + advance.push_x(prt.x, v, 1.0); + + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], 2.47214, 1e-5); + EXPECT_NEAR(prt.x[2], 5., 1e-5); +} + int main(int argc, char** argv) { MPI_Init(&argc, &argv); From 7e4b07dc3e97abed10f9d1a3bd2a0bebe321dd60 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 14:54:11 -0400 Subject: [PATCH 09/29] test_bnd_prt_inflow: factor out Inflow class --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 46 +++++++++++++++++------- 1 file changed, 34 insertions(+), 12 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 8ac62785a9..4d091260ad 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -249,6 +249,37 @@ TEST(TestSetupParticlesInflow, Maxwellian) EXPECT_EQ(prt.tag, tag); } +template +class Inflow +{ +public: + using Mparticles = MP; + using real_t = typename Mparticles::real_t; + + Inflow(const Grid_t& grid, psc_particle_npt npt) + : grid_(grid), setup_particles_(grid), npt_(npt) + {} + + auto get_advanced_prt(Double3 pos) + { + auto p = setup_particles_.createMaxwellian(npt_); + psc_particle_np np = {npt_.kind, npt_.n, p}; + double wni = 1.0; + auto prt = setup_particles_.setupParticle(np, pos, wni); + + auto advance = AdvanceParticle(grid_.dt); + + auto v = advance.calc_v(prt.u); + advance.push_x(prt.x, v, 1.0); + + return prt; + } + + const Grid_t& grid_; + SetupParticles setup_particles_; + psc_particle_npt npt_; +}; + TEST(TestSetupParticlesInflow, Advance) { using Mparticles = MparticlesDouble; @@ -261,23 +292,14 @@ TEST(TestSetupParticlesInflow, Advance) Grid_t grid{domain, {}, kinds, {prm}, dt}; Mparticles mprts{grid}; - SetupParticles setup_particles(grid); - Double3 pos = {5., -2., 5.}; Double3 u = {0.0, 0.5, 0.0}; Double3 T = {0.0, 0.0, 0.0}; - int tag = 0; - - psc_particle_npt npt = {0, 1.0, u, T, tag}; - auto p = setup_particles.createMaxwellian(npt); - psc_particle_np np = {0, 1.0, p, tag}; - double wni = 1.0; - auto prt = setup_particles.setupParticle(np, pos, wni); + psc_particle_npt npt = {0, 1.0, u, T}; - auto advance = AdvanceParticle(dt); + Inflow inflow(grid, npt); - auto v = advance.calc_v(prt.u); - advance.push_x(prt.x, v, 1.0); + auto prt = inflow.get_advanced_prt(pos); EXPECT_NEAR(prt.x[0], 5., 1e-5); EXPECT_NEAR(prt.x[1], 2.47214, 1e-5); From b4e0eff58ff19eada807c539ebf1d40c2260e8c8 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 14:56:47 -0400 Subject: [PATCH 10/29] test_bnd_prt_inflow: optimizations --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 4d091260ad..fd9e882131 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -257,27 +257,27 @@ class Inflow using real_t = typename Mparticles::real_t; Inflow(const Grid_t& grid, psc_particle_npt npt) - : grid_(grid), setup_particles_(grid), npt_(npt) + : grid_(grid), + setup_particles_(grid), + advance_(grid.dt), + np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)} {} auto get_advanced_prt(Double3 pos) { - auto p = setup_particles_.createMaxwellian(npt_); - psc_particle_np np = {npt_.kind, npt_.n, p}; double wni = 1.0; - auto prt = setup_particles_.setupParticle(np, pos, wni); + auto prt = setup_particles_.setupParticle(np_, pos, wni); - auto advance = AdvanceParticle(grid_.dt); - - auto v = advance.calc_v(prt.u); - advance.push_x(prt.x, v, 1.0); + auto v = advance_.calc_v(prt.u); + advance_.push_x(prt.x, v, 1.0); return prt; } const Grid_t& grid_; + AdvanceParticle advance_; SetupParticles setup_particles_; - psc_particle_npt npt_; + psc_particle_np np_; }; TEST(TestSetupParticlesInflow, Advance) From ac56f54a08dab3d0aaabb6bd7f8e3b4fbe0b007d Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:05:53 -0400 Subject: [PATCH 11/29] test_bnd_prt_inflow: +inject_into and test --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 35 ++++++++++++++++++++++++ 1 file changed, 35 insertions(+) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index fd9e882131..06ed74dd15 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -274,6 +274,14 @@ class Inflow return prt; } + // TODO: make this signature inject_into(injector, cell, np) + auto inject_into(Double3 pos) + { + auto prt = get_advanced_prt(pos); + + return prt; + } + const Grid_t& grid_; AdvanceParticle advance_; SetupParticles setup_particles_; @@ -306,6 +314,33 @@ TEST(TestSetupParticlesInflow, Advance) EXPECT_NEAR(prt.x[2], 5., 1e-5); } +TEST(TestSetupParticlesInflow, InjectInto) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 1, 1}}; + double dt = 10.; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, dt}; + Mparticles mprts{grid}; + + Double3 pos = {5., -2., 5.}; + Double3 u = {0.0, 0.5, 0.0}; + Double3 T = {0.0, 0.0, 0.0}; + psc_particle_npt npt = {0, 1.0, u, T}; + + Inflow inflow(grid, npt); + + auto prt = inflow.inject_into(pos); + + // TODO write actual tests + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], 2.47214, 1e-5); + EXPECT_NEAR(prt.x[2], 5., 1e-5); +} + int main(int argc, char** argv) { MPI_Init(&argc, &argv); From c384c494d21ccbe6bbba0e2ba3b71cb85b68ffba Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:11:26 -0400 Subject: [PATCH 12/29] test_bnd_prt_inflow: use TestInjector --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 23 +++++++++++++++++------ 1 file changed, 17 insertions(+), 6 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 06ed74dd15..c231289b04 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -249,11 +249,12 @@ TEST(TestSetupParticlesInflow, Maxwellian) EXPECT_EQ(prt.tag, tag); } -template +template class Inflow { public: using Mparticles = MP; + using Injector = INJ; using real_t = typename Mparticles::real_t; Inflow(const Grid_t& grid, psc_particle_npt npt) @@ -275,11 +276,11 @@ class Inflow } // TODO: make this signature inject_into(injector, cell, np) - auto inject_into(Double3 pos) + void inject_into(Injector& injector, Double3 pos) { auto prt = get_advanced_prt(pos); - return prt; + injector(prt); } const Grid_t& grid_; @@ -288,6 +289,14 @@ class Inflow psc_particle_np np_; }; +class TestInjector +{ +public: + void operator()(psc::particle::Inject prt) { prts.push_back(prt); } + + std::vector prts; +}; + TEST(TestSetupParticlesInflow, Advance) { using Mparticles = MparticlesDouble; @@ -305,7 +314,7 @@ TEST(TestSetupParticlesInflow, Advance) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt); auto prt = inflow.get_advanced_prt(pos); @@ -331,9 +340,11 @@ TEST(TestSetupParticlesInflow, InjectInto) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt); - auto prt = inflow.inject_into(pos); + TestInjector injector; + inflow.inject_into(injector, pos); + auto prt = injector.prts[0]; // TODO write actual tests EXPECT_NEAR(prt.x[0], 5., 1e-5); From 2c00b528a0cbba07046770fc20a8e597aacb8f86 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:22:14 -0400 Subject: [PATCH 13/29] test_bnd_prt_inflow: take cell_idx --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index c231289b04..91488291b6 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -276,8 +276,11 @@ class Inflow } // TODO: make this signature inject_into(injector, cell, np) - void inject_into(Injector& injector, Double3 pos) + void inject_into(Injector& injector, Int3 cell_idx) { + cell_idx[1] = -1; + auto pos = (Double3(cell_idx) + Double3(.5, .5, .5)) * grid_.domain.dx + + grid_.domain.corner; auto prt = get_advanced_prt(pos); injector(prt); @@ -343,12 +346,12 @@ TEST(TestSetupParticlesInflow, InjectInto) Inflow inflow(grid, npt); TestInjector injector; - inflow.inject_into(injector, pos); + inflow.inject_into(injector, {0, 0, 0}); auto prt = injector.prts[0]; // TODO write actual tests EXPECT_NEAR(prt.x[0], 5., 1e-5); - EXPECT_NEAR(prt.x[1], 2.47214, 1e-5); + EXPECT_NEAR(prt.x[1], -0.5278640, 1e-5); EXPECT_NEAR(prt.x[2], 5., 1e-5); } From 8a34d88a33d7e043a1f0b2654b44f7cd8a016045 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:29:10 -0400 Subject: [PATCH 14/29] test_bnd_prt_inflow: take offset_in_cell_dist --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 91488291b6..7fcf0568c1 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -276,11 +276,14 @@ class Inflow } // TODO: make this signature inject_into(injector, cell, np) - void inject_into(Injector& injector, Int3 cell_idx) + void inject_into(Injector& injector, Int3 cell_idx, + std::function offset_in_cell_dist) { cell_idx[1] = -1; - auto pos = (Double3(cell_idx) + Double3(.5, .5, .5)) * grid_.domain.dx + - grid_.domain.corner; + Double3 offset = {offset_in_cell_dist(), offset_in_cell_dist(), + offset_in_cell_dist()}; + auto pos = + (Double3(cell_idx) + offset) * grid_.domain.dx + grid_.domain.corner; auto prt = get_advanced_prt(pos); injector(prt); @@ -346,7 +349,7 @@ TEST(TestSetupParticlesInflow, InjectInto) Inflow inflow(grid, npt); TestInjector injector; - inflow.inject_into(injector, {0, 0, 0}); + inflow.inject_into(injector, {0, 0, 0}, []() { return .5; }); auto prt = injector.prts[0]; // TODO write actual tests From 6a5974c6b98a8cbae83de8e0fd8da671dea799a7 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:35:39 -0400 Subject: [PATCH 15/29] setup_particles: +getWeight --- src/include/setup_particles.hxx | 19 +++++++++++++------ 1 file changed, 13 insertions(+), 6 deletions(-) diff --git a/src/include/setup_particles.hxx b/src/include/setup_particles.hxx index 632e4a63d0..c522f34549 100644 --- a/src/include/setup_particles.hxx +++ b/src/include/setup_particles.hxx @@ -241,6 +241,18 @@ struct SetupParticles }); } + // ---------------------------------------------------------------------- + // getWeight + + real_t getWeight(const psc_particle_np& np, int n_in_cell) + { + if (fractional_n_particles_per_cell) { + return 1.; + } else { + return np.n / (n_in_cell * norm_.cori); + } + } + // ---------------------------------------------------------------------- // setupParticles @@ -269,12 +281,7 @@ struct SetupParticles op_cellwise(grid, p, init_np, [&](int n_in_cell, psc_particle_np& np, Double3& pos) { for (int cnt = 0; cnt < n_in_cell; cnt++) { - real_t wni; - if (fractional_n_particles_per_cell) { - wni = 1.; - } else { - wni = np.n / (n_in_cell * norm_.cori); - } + real_t wni = getWeight(np, n_in_cell); auto prt = setupParticle(np, pos, wni); injector(prt); } From eca97459b24720ae245a61eac7e5c3765c83d78f Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:40:23 -0400 Subject: [PATCH 16/29] test_bnd_prt_inflow: loop over n_in_cell --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 39 +++++++++++++++--------- 1 file changed, 25 insertions(+), 14 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 7fcf0568c1..8edbbfb25a 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -264,9 +264,8 @@ class Inflow np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)} {} - auto get_advanced_prt(Double3 pos) + auto get_advanced_prt(Double3 pos, real_t wni) { - double wni = 1.0; auto prt = setup_particles_.setupParticle(np_, pos, wni); auto v = advance_.calc_v(prt.u); @@ -279,14 +278,22 @@ class Inflow void inject_into(Injector& injector, Int3 cell_idx, std::function offset_in_cell_dist) { + // TODO don't hardcode this cell_idx[1] = -1; - Double3 offset = {offset_in_cell_dist(), offset_in_cell_dist(), - offset_in_cell_dist()}; - auto pos = - (Double3(cell_idx) + offset) * grid_.domain.dx + grid_.domain.corner; - auto prt = get_advanced_prt(pos); - injector(prt); + int n_in_cell = setup_particles_.get_n_in_cell(np_); + + for (int cnt = 0; cnt < n_in_cell; cnt++) { + double wni = setup_particles_.getWeight(np_, n_in_cell); + + Double3 offset = {offset_in_cell_dist(), offset_in_cell_dist(), + offset_in_cell_dist()}; + auto pos = + (Double3(cell_idx) + offset) * grid_.domain.dx + grid_.domain.corner; + auto prt = get_advanced_prt(pos, wni); + + injector(prt); + } } const Grid_t& grid_; @@ -322,7 +329,7 @@ TEST(TestSetupParticlesInflow, Advance) Inflow inflow(grid, npt); - auto prt = inflow.get_advanced_prt(pos); + auto prt = inflow.get_advanced_prt(pos, 1.0); EXPECT_NEAR(prt.x[0], 5., 1e-5); EXPECT_NEAR(prt.x[1], 2.47214, 1e-5); @@ -350,12 +357,16 @@ TEST(TestSetupParticlesInflow, InjectInto) TestInjector injector; inflow.inject_into(injector, {0, 0, 0}, []() { return .5; }); - auto prt = injector.prts[0]; - // TODO write actual tests - EXPECT_NEAR(prt.x[0], 5., 1e-5); - EXPECT_NEAR(prt.x[1], -0.5278640, 1e-5); - EXPECT_NEAR(prt.x[2], 5., 1e-5); + EXPECT_EQ(injector.prts.size(), prm.nicell); + + for (int i = 0; i < prm.nicell; i++) { + auto prt = injector.prts[i]; + + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], -0.5278640, 1e-5); + EXPECT_NEAR(prt.x[2], 5., 1e-5); + } } int main(int argc, char** argv) From 40ae798e73b0d8f883f4e74e507805d00e24b7bc Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:42:11 -0400 Subject: [PATCH 17/29] test_bnd_prt_inflow: adjust numbers --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 8edbbfb25a..6349852951 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -348,8 +348,7 @@ TEST(TestSetupParticlesInflow, InjectInto) Grid_t grid{domain, {}, kinds, {prm}, dt}; Mparticles mprts{grid}; - Double3 pos = {5., -2., 5.}; - Double3 u = {0.0, 0.5, 0.0}; + Double3 u = {0.0, 0.8, 0.0}; Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; @@ -364,7 +363,7 @@ TEST(TestSetupParticlesInflow, InjectInto) auto prt = injector.prts[i]; EXPECT_NEAR(prt.x[0], 5., 1e-5); - EXPECT_NEAR(prt.x[1], -0.5278640, 1e-5); + EXPECT_NEAR(prt.x[1], 1.246950, 1e-5); EXPECT_NEAR(prt.x[2], 5., 1e-5); } } From 1b25481f2b648a433f188d1e4696a8ec9ff6d506 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:46:23 -0400 Subject: [PATCH 18/29] test_bnd_prt_inflow: filter out prts not in domain --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 30 ++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 6349852951..b4406d3816 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -292,6 +292,11 @@ class Inflow (Double3(cell_idx) + offset) * grid_.domain.dx + grid_.domain.corner; auto prt = get_advanced_prt(pos, wni); + // TODO make this work for other boundaries and dimensions + if (prt.x[1] < grid_.domain.corner[1]) { + continue; + } + injector(prt); } } @@ -368,6 +373,31 @@ TEST(TestSetupParticlesInflow, InjectInto) } } +TEST(TestSetupParticlesInflow, InjectIntoFilter) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 1, 1}}; + double dt = 10.; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, dt}; + Mparticles mprts{grid}; + + // too slow to enter domain + Double3 u = {0.0, 0.5, 0.0}; + Double3 T = {0.0, 0.0, 0.0}; + psc_particle_npt npt = {0, 1.0, u, T}; + + Inflow inflow(grid, npt); + + TestInjector injector; + inflow.inject_into(injector, {0, 0, 0}, []() { return .5; }); + + EXPECT_EQ(injector.prts.size(), 0); +} + int main(int argc, char** argv) { MPI_Init(&argc, &argv); From e6b513e0c0f8a111eb85b17d8b72e06d64b43f25 Mon Sep 17 00:00:00 2001 From: James McClung Date: Wed, 26 Mar 2025 15:55:50 -0400 Subject: [PATCH 19/29] test_bnd_prt_inflow: +inject_into_patch --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 77 +++++++++++++++++++----- 1 file changed, 63 insertions(+), 14 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index b4406d3816..b2dae18cca 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -249,12 +249,11 @@ TEST(TestSetupParticlesInflow, Maxwellian) EXPECT_EQ(prt.tag, tag); } -template +template class Inflow { public: using Mparticles = MP; - using Injector = INJ; using real_t = typename Mparticles::real_t; Inflow(const Grid_t& grid, psc_particle_npt npt) @@ -275,8 +274,9 @@ class Inflow } // TODO: make this signature inject_into(injector, cell, np) - void inject_into(Injector& injector, Int3 cell_idx, - std::function offset_in_cell_dist) + template + void inject_into_cell(Injector& injector, Int3 cell_idx, + std::function offset_in_cell_dist) { // TODO don't hardcode this cell_idx[1] = -1; @@ -301,6 +301,19 @@ class Inflow } } + template + void inject_into_patch(Injector& injector, + std::function offset_in_cell_dist) + { + Int3 cell_idx = {0, 0, 0}; + for (cell_idx[0] = 0; cell_idx[0] < grid_.domain.ldims[0]; cell_idx[0]++) { + for (cell_idx[2] = 0; cell_idx[2] < grid_.domain.ldims[2]; + cell_idx[2]++) { + inject_into_cell(injector, cell_idx, []() { return .5; }); + } + } + } + const Grid_t& grid_; AdvanceParticle advance_; SetupParticles setup_particles_; @@ -332,7 +345,7 @@ TEST(TestSetupParticlesInflow, Advance) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt); auto prt = inflow.get_advanced_prt(pos, 1.0); @@ -341,7 +354,7 @@ TEST(TestSetupParticlesInflow, Advance) EXPECT_NEAR(prt.x[2], 5., 1e-5); } -TEST(TestSetupParticlesInflow, InjectInto) +TEST(TestSetupParticlesInflow, InjectIntoCell) { using Mparticles = MparticlesDouble; @@ -357,23 +370,21 @@ TEST(TestSetupParticlesInflow, InjectInto) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt); TestInjector injector; - inflow.inject_into(injector, {0, 0, 0}, []() { return .5; }); + inflow.inject_into_cell(injector, {0, 0, 0}, []() { return .5; }); EXPECT_EQ(injector.prts.size(), prm.nicell); - for (int i = 0; i < prm.nicell; i++) { - auto prt = injector.prts[i]; - + for (auto prt : injector.prts) { EXPECT_NEAR(prt.x[0], 5., 1e-5); EXPECT_NEAR(prt.x[1], 1.246950, 1e-5); EXPECT_NEAR(prt.x[2], 5., 1e-5); } } -TEST(TestSetupParticlesInflow, InjectIntoFilter) +TEST(TestSetupParticlesInflow, InjectIntoCellFilter) { using Mparticles = MparticlesDouble; @@ -390,14 +401,52 @@ TEST(TestSetupParticlesInflow, InjectIntoFilter) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt); TestInjector injector; - inflow.inject_into(injector, {0, 0, 0}, []() { return .5; }); + inflow.inject_into_cell(injector, {0, 0, 0}, []() { return .5; }); EXPECT_EQ(injector.prts.size(), 0); } +TEST(TestSetupParticlesInflow, InjectIntoPatch) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 1, 1}}; + double dt = 10.; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, dt}; + Mparticles mprts{grid}; + + Double3 u = {0.0, 0.8, 0.0}; + Double3 T = {0.0, 0.0, 0.0}; + psc_particle_npt npt = {0, 1.0, u, T}; + + Inflow inflow(grid, npt); + + TestInjector injector; + inflow.inject_into_patch(injector, []() { return .5; }); + + EXPECT_EQ(injector.prts.size(), + domain.ldims[0] * domain.ldims[2] * prm.nicell); + + for (int i = 0; i < injector.prts.size(); i++) { + auto prt = injector.prts[i]; + + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], 1.246950, 1e-5); + + if (i < prm.nicell) { + EXPECT_NEAR(prt.x[2], 5., 1e-5); + } else { + EXPECT_NEAR(prt.x[2], 15., 1e-5); + } + } +} + int main(int argc, char** argv) { MPI_Init(&argc, &argv); From 646605727d4b72941f2b2c879eea3c0340a82dd7 Mon Sep 17 00:00:00 2001 From: James McClung Date: Thu, 27 Mar 2025 13:20:00 -0400 Subject: [PATCH 20/29] test_bnd_prt_inflow: +inject_into_boundary --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 71 ++++++++++++++++++++++-- 1 file changed, 65 insertions(+), 6 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index b2dae18cca..8ab044033f 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -302,18 +302,35 @@ class Inflow } template - void inject_into_patch(Injector& injector, + void inject_into_patch(Injector& injector, const Grid_t::Patch& patch, std::function offset_in_cell_dist) { - Int3 cell_idx = {0, 0, 0}; - for (cell_idx[0] = 0; cell_idx[0] < grid_.domain.ldims[0]; cell_idx[0]++) { - for (cell_idx[2] = 0; cell_idx[2] < grid_.domain.ldims[2]; - cell_idx[2]++) { + Int3 ilo = patch.off; + Int3 ihi = ilo + grid_.ldims; + for (Int3 cell_idx = ilo; cell_idx[0] < ihi[0]; cell_idx[0]++) { + for (cell_idx[2] = ilo[2]; cell_idx[2] < ihi[2]; cell_idx[2]++) { inject_into_cell(injector, cell_idx, []() { return .5; }); } } } + template + void inject_into_boundary(Injectors& injectors_by_patch, + std::function offset_in_cell_dist) + { + for (int patch_idx = 0; patch_idx < grid_.n_patches(); patch_idx++) { + const auto& patch = grid_.patches[patch_idx]; + + if (patch.off[1] != 0) { + // TODO make work for other dims + continue; + } + + auto& injector = injectors_by_patch[patch_idx]; + inject_into_patch(injector, patch, offset_in_cell_dist); + } + } + const Grid_t& grid_; AdvanceParticle advance_; SetupParticles setup_particles_; @@ -428,7 +445,7 @@ TEST(TestSetupParticlesInflow, InjectIntoPatch) Inflow inflow(grid, npt); TestInjector injector; - inflow.inject_into_patch(injector, []() { return .5; }); + inflow.inject_into_patch(injector, grid.patches[0], []() { return .5; }); EXPECT_EQ(injector.prts.size(), domain.ldims[0] * domain.ldims[2] * prm.nicell); @@ -446,6 +463,48 @@ TEST(TestSetupParticlesInflow, InjectIntoPatch) } } } +TEST(TestSetupParticlesInflow, InjectIntoBoundary) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 2, 2}}; + double dt = 10.; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, dt}; + Mparticles mprts{grid}; + + Double3 u = {0.0, 0.8, 0.0}; + Double3 T = {0.0, 0.0, 0.0}; + psc_particle_npt npt = {0, 1.0, u, T}; + + Inflow inflow(grid, npt); + + std::vector injectors = {TestInjector(), TestInjector(), + TestInjector(), TestInjector()}; + inflow.inject_into_boundary(injectors, []() { return .5; }); + + EXPECT_EQ(injectors[0].prts.size(), + domain.ldims[0] * domain.ldims[2] * prm.nicell); + for (auto prt : injectors[0].prts) { + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], 1.246950, 1e-5); + EXPECT_NEAR(prt.x[2], 5., 1e-5); + } + + EXPECT_EQ(injectors[1].prts.size(), 0); + + EXPECT_EQ(injectors[2].prts.size(), + domain.ldims[0] * domain.ldims[2] * prm.nicell); + for (auto prt : injectors[2].prts) { + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], 1.246950, 1e-5); + EXPECT_NEAR(prt.x[2], 15., 1e-5); + } + + EXPECT_EQ(injectors[3].prts.size(), 0); +} int main(int argc, char** argv) { From 2a1894482267ccab0292fd7e87d7cf587ba71971 Mon Sep 17 00:00:00 2001 From: James McClung Date: Thu, 27 Mar 2025 13:44:24 -0400 Subject: [PATCH 21/29] test_bnd_prt_inflow: put offset_in_cell_dist in ctor and make it a function pointer. Lambdas would require templating the type (and thus obfuscate the code), and avoid std::function's virtual function calls, since this could be called many times per step --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 48 +++++++++++++----------- 1 file changed, 26 insertions(+), 22 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 8ab044033f..dc30081926 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -256,11 +256,16 @@ class Inflow using Mparticles = MP; using real_t = typename Mparticles::real_t; - Inflow(const Grid_t& grid, psc_particle_npt npt) + // offset_in_cell: () -> double + Inflow(const Grid_t& grid, psc_particle_npt npt, + double (*offset_in_cell_dist)()) : grid_(grid), setup_particles_(grid), advance_(grid.dt), - np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)} + // FIXME np.p is a std::function and is called many times; better to use a + // lambda + np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)}, + offset_in_cell_dist_(offset_in_cell_dist) {} auto get_advanced_prt(Double3 pos, real_t wni) @@ -273,10 +278,8 @@ class Inflow return prt; } - // TODO: make this signature inject_into(injector, cell, np) template - void inject_into_cell(Injector& injector, Int3 cell_idx, - std::function offset_in_cell_dist) + void inject_into_cell(Injector& injector, Int3 cell_idx) { // TODO don't hardcode this cell_idx[1] = -1; @@ -286,8 +289,8 @@ class Inflow for (int cnt = 0; cnt < n_in_cell; cnt++) { double wni = setup_particles_.getWeight(np_, n_in_cell); - Double3 offset = {offset_in_cell_dist(), offset_in_cell_dist(), - offset_in_cell_dist()}; + Double3 offset = {offset_in_cell_dist_(), offset_in_cell_dist_(), + offset_in_cell_dist_()}; auto pos = (Double3(cell_idx) + offset) * grid_.domain.dx + grid_.domain.corner; auto prt = get_advanced_prt(pos, wni); @@ -302,21 +305,19 @@ class Inflow } template - void inject_into_patch(Injector& injector, const Grid_t::Patch& patch, - std::function offset_in_cell_dist) + void inject_into_patch(Injector& injector, const Grid_t::Patch& patch) { Int3 ilo = patch.off; Int3 ihi = ilo + grid_.ldims; for (Int3 cell_idx = ilo; cell_idx[0] < ihi[0]; cell_idx[0]++) { for (cell_idx[2] = ilo[2]; cell_idx[2] < ihi[2]; cell_idx[2]++) { - inject_into_cell(injector, cell_idx, []() { return .5; }); + inject_into_cell(injector, cell_idx); } } } template - void inject_into_boundary(Injectors& injectors_by_patch, - std::function offset_in_cell_dist) + void inject_into_boundary(Injectors& injectors_by_patch) { for (int patch_idx = 0; patch_idx < grid_.n_patches(); patch_idx++) { const auto& patch = grid_.patches[patch_idx]; @@ -327,7 +328,7 @@ class Inflow } auto& injector = injectors_by_patch[patch_idx]; - inject_into_patch(injector, patch, offset_in_cell_dist); + inject_into_patch(injector, patch); } } @@ -335,6 +336,7 @@ class Inflow AdvanceParticle advance_; SetupParticles setup_particles_; psc_particle_np np_; + double (*offset_in_cell_dist_)(); }; class TestInjector @@ -345,6 +347,8 @@ class TestInjector std::vector prts; }; +double half() { return 0.5; } + TEST(TestSetupParticlesInflow, Advance) { using Mparticles = MparticlesDouble; @@ -362,7 +366,7 @@ TEST(TestSetupParticlesInflow, Advance) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt, *half); auto prt = inflow.get_advanced_prt(pos, 1.0); @@ -387,10 +391,10 @@ TEST(TestSetupParticlesInflow, InjectIntoCell) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt, *half); TestInjector injector; - inflow.inject_into_cell(injector, {0, 0, 0}, []() { return .5; }); + inflow.inject_into_cell(injector, {0, 0, 0}); EXPECT_EQ(injector.prts.size(), prm.nicell); @@ -418,10 +422,10 @@ TEST(TestSetupParticlesInflow, InjectIntoCellFilter) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt, *half); TestInjector injector; - inflow.inject_into_cell(injector, {0, 0, 0}, []() { return .5; }); + inflow.inject_into_cell(injector, {0, 0, 0}); EXPECT_EQ(injector.prts.size(), 0); } @@ -442,10 +446,10 @@ TEST(TestSetupParticlesInflow, InjectIntoPatch) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt, *half); TestInjector injector; - inflow.inject_into_patch(injector, grid.patches[0], []() { return .5; }); + inflow.inject_into_patch(injector, grid.patches[0]); EXPECT_EQ(injector.prts.size(), domain.ldims[0] * domain.ldims[2] * prm.nicell); @@ -479,11 +483,11 @@ TEST(TestSetupParticlesInflow, InjectIntoBoundary) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt); + Inflow inflow(grid, npt, *half); std::vector injectors = {TestInjector(), TestInjector(), TestInjector(), TestInjector()}; - inflow.inject_into_boundary(injectors, []() { return .5; }); + inflow.inject_into_boundary(injectors); EXPECT_EQ(injectors[0].prts.size(), domain.ldims[0] * domain.ldims[2] * prm.nicell); From cba843ac63cbb5c0778d499982059d39091d3430 Mon Sep 17 00:00:00 2001 From: James McClung Date: Thu, 27 Mar 2025 14:08:44 -0400 Subject: [PATCH 22/29] test_bnd_prt_inflow: mv getWeight call out of loop --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index dc30081926..cf882f84d2 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -285,10 +285,9 @@ class Inflow cell_idx[1] = -1; int n_in_cell = setup_particles_.get_n_in_cell(np_); + double wni = setup_particles_.getWeight(np_, n_in_cell); for (int cnt = 0; cnt < n_in_cell; cnt++) { - double wni = setup_particles_.getWeight(np_, n_in_cell); - Double3 offset = {offset_in_cell_dist_(), offset_in_cell_dist_(), offset_in_cell_dist_()}; auto pos = From 22cd653c70d7fe77277377c17795bc328200cf36 Mon Sep 17 00:00:00 2001 From: James McClung Date: Thu, 27 Mar 2025 14:27:45 -0400 Subject: [PATCH 23/29] test_bnd_prt_inflow: support other inj dims --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 58 ++++++++++++++---------- 1 file changed, 34 insertions(+), 24 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index cf882f84d2..dab933e9da 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -258,14 +258,15 @@ class Inflow // offset_in_cell: () -> double Inflow(const Grid_t& grid, psc_particle_npt npt, - double (*offset_in_cell_dist)()) + double (*offset_in_cell_dist)(), int inj_dim_idx) : grid_(grid), setup_particles_(grid), advance_(grid.dt), // FIXME np.p is a std::function and is called many times; better to use a // lambda np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)}, - offset_in_cell_dist_(offset_in_cell_dist) + offset_in_cell_dist_(offset_in_cell_dist), + inj_dim_idx_(inj_dim_idx) {} auto get_advanced_prt(Double3 pos, real_t wni) @@ -279,10 +280,11 @@ class Inflow } template - void inject_into_cell(Injector& injector, Int3 cell_idx) + void inject_into_boundary_cell(Injector& injector, + Int3 boundary_cell_global_idx) { - // TODO don't hardcode this - cell_idx[1] = -1; + assert(boundary_cell_global_idx[inj_dim_idx_] == 0); + boundary_cell_global_idx[inj_dim_idx_] = -1; int n_in_cell = setup_particles_.get_n_in_cell(np_); double wni = setup_particles_.getWeight(np_, n_in_cell); @@ -291,11 +293,11 @@ class Inflow Double3 offset = {offset_in_cell_dist_(), offset_in_cell_dist_(), offset_in_cell_dist_()}; auto pos = - (Double3(cell_idx) + offset) * grid_.domain.dx + grid_.domain.corner; + (Double3(boundary_cell_global_idx) + offset) * grid_.domain.dx + + grid_.domain.corner; auto prt = get_advanced_prt(pos, wni); - // TODO make this work for other boundaries and dimensions - if (prt.x[1] < grid_.domain.corner[1]) { + if (prt.x[inj_dim_idx_] < grid_.domain.corner[inj_dim_idx_]) { continue; } @@ -304,13 +306,21 @@ class Inflow } template - void inject_into_patch(Injector& injector, const Grid_t::Patch& patch) + void inject_into_boundary_patch(Injector& injector, + const Grid_t::Patch& boundary_patch) { - Int3 ilo = patch.off; + Int3 ilo = boundary_patch.off; Int3 ihi = ilo + grid_.ldims; - for (Int3 cell_idx = ilo; cell_idx[0] < ihi[0]; cell_idx[0]++) { - for (cell_idx[2] = ilo[2]; cell_idx[2] < ihi[2]; cell_idx[2]++) { - inject_into_cell(injector, cell_idx); + + assert(ilo[inj_dim_idx_] == 0); + + int dim1 = std::min((inj_dim_idx_ + 1) % 3, (inj_dim_idx_ + 2) % 3); + int dim2 = std::max((inj_dim_idx_ + 1) % 3, (inj_dim_idx_ + 2) % 3); + + for (Int3 cell_idx = ilo; cell_idx[dim1] < ihi[dim1]; cell_idx[dim1]++) { + for (cell_idx[dim2] = ilo[dim2]; cell_idx[dim2] < ihi[dim2]; + cell_idx[dim2]++) { + inject_into_boundary_cell(injector, cell_idx); } } } @@ -321,13 +331,12 @@ class Inflow for (int patch_idx = 0; patch_idx < grid_.n_patches(); patch_idx++) { const auto& patch = grid_.patches[patch_idx]; - if (patch.off[1] != 0) { - // TODO make work for other dims + if (patch.off[inj_dim_idx_] != 0) { continue; } auto& injector = injectors_by_patch[patch_idx]; - inject_into_patch(injector, patch); + inject_into_boundary_patch(injector, patch); } } @@ -336,6 +345,7 @@ class Inflow SetupParticles setup_particles_; psc_particle_np np_; double (*offset_in_cell_dist_)(); + int inj_dim_idx_; }; class TestInjector @@ -365,7 +375,7 @@ TEST(TestSetupParticlesInflow, Advance) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half); + Inflow inflow(grid, npt, *half, 1); auto prt = inflow.get_advanced_prt(pos, 1.0); @@ -390,10 +400,10 @@ TEST(TestSetupParticlesInflow, InjectIntoCell) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half); + Inflow inflow(grid, npt, *half, 1); TestInjector injector; - inflow.inject_into_cell(injector, {0, 0, 0}); + inflow.inject_into_boundary_cell(injector, {0, 0, 0}); EXPECT_EQ(injector.prts.size(), prm.nicell); @@ -421,10 +431,10 @@ TEST(TestSetupParticlesInflow, InjectIntoCellFilter) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half); + Inflow inflow(grid, npt, *half, 1); TestInjector injector; - inflow.inject_into_cell(injector, {0, 0, 0}); + inflow.inject_into_boundary_cell(injector, {0, 0, 0}); EXPECT_EQ(injector.prts.size(), 0); } @@ -445,10 +455,10 @@ TEST(TestSetupParticlesInflow, InjectIntoPatch) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half); + Inflow inflow(grid, npt, *half, 1); TestInjector injector; - inflow.inject_into_patch(injector, grid.patches[0]); + inflow.inject_into_boundary_patch(injector, grid.patches[0]); EXPECT_EQ(injector.prts.size(), domain.ldims[0] * domain.ldims[2] * prm.nicell); @@ -482,7 +492,7 @@ TEST(TestSetupParticlesInflow, InjectIntoBoundary) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half); + Inflow inflow(grid, npt, *half, 1); std::vector injectors = {TestInjector(), TestInjector(), TestInjector(), TestInjector()}; From d75f2af47b613eb7dc1b0f5d8ba1a61f03ed908b Mon Sep 17 00:00:00 2001 From: James McClung Date: Thu, 27 Mar 2025 14:53:39 -0400 Subject: [PATCH 24/29] test_bnd_prt_inflow: +InjectIntoBoundaryZ it fails --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 46 +++++++++++++++++++++++- 1 file changed, 45 insertions(+), 1 deletion(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index dab933e9da..c19350ad46 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -476,7 +476,8 @@ TEST(TestSetupParticlesInflow, InjectIntoPatch) } } } -TEST(TestSetupParticlesInflow, InjectIntoBoundary) + +TEST(TestSetupParticlesInflow, InjectIntoBoundaryY) { using Mparticles = MparticlesDouble; @@ -519,6 +520,49 @@ TEST(TestSetupParticlesInflow, InjectIntoBoundary) EXPECT_EQ(injectors[3].prts.size(), 0); } +TEST(TestSetupParticlesInflow, InjectIntoBoundaryZ) +{ + using Mparticles = MparticlesDouble; + + auto domain = Grid_t::Domain{{1, 2, 2}, {10., 20., 20.}, {}, {1, 2, 2}}; + double dt = 10.; + auto kinds = Grid_t::Kinds{{1., 100., "i"}}; + auto prm = Grid_t::NormalizationParams::dimensionless(); + prm.nicell = 2; + Grid_t grid{domain, {}, kinds, {prm}, dt}; + Mparticles mprts{grid}; + + Double3 u = {0.0, 0.0, 0.8}; + Double3 T = {0.0, 0.0, 0.0}; + psc_particle_npt npt = {0, 1.0, u, T}; + + Inflow inflow(grid, npt, *half, 2); + + std::vector injectors = {TestInjector(), TestInjector(), + TestInjector(), TestInjector()}; + inflow.inject_into_boundary(injectors); + + EXPECT_EQ(injectors[0].prts.size(), + domain.ldims[0] * domain.ldims[1] * prm.nicell); + for (auto prt : injectors[0].prts) { + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], 5., 1e-5); + EXPECT_NEAR(prt.x[2], 1.246950, 1e-5); + } + + EXPECT_EQ(injectors[1].prts.size(), + domain.ldims[0] * domain.ldims[1] * prm.nicell); + for (auto prt : injectors[1].prts) { + EXPECT_NEAR(prt.x[0], 5., 1e-5); + EXPECT_NEAR(prt.x[1], 15., 1e-5); + EXPECT_NEAR(prt.x[2], 1.246950, 1e-5); + } + + EXPECT_EQ(injectors[2].prts.size(), 0); + + EXPECT_EQ(injectors[3].prts.size(), 0); +} + int main(int argc, char** argv) { MPI_Init(&argc, &argv); From 2d5b5e1e3381a1bf2c6cd4b9e683838712b8b9ce Mon Sep 17 00:00:00 2001 From: James McClung Date: Thu, 27 Mar 2025 15:04:32 -0400 Subject: [PATCH 25/29] test_bnd_prt_inflow: infer inj_dim from template --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 28 ++++++++++++++---------- 1 file changed, 17 insertions(+), 11 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index c19350ad46..05caa38922 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -249,16 +249,17 @@ TEST(TestSetupParticlesInflow, Maxwellian) EXPECT_EQ(prt.tag, tag); } -template +template class Inflow { public: using Mparticles = MP; using real_t = typename Mparticles::real_t; + using Dim = DIM; // offset_in_cell: () -> double Inflow(const Grid_t& grid, psc_particle_npt npt, - double (*offset_in_cell_dist)(), int inj_dim_idx) + double (*offset_in_cell_dist)()) : grid_(grid), setup_particles_(grid), advance_(grid.dt), @@ -266,8 +267,13 @@ class Inflow // lambda np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)}, offset_in_cell_dist_(offset_in_cell_dist), - inj_dim_idx_(inj_dim_idx) - {} + inj_dim_idx_(!Dim::InvarX::value ? 0 + : !Dim::InvarY::value ? 1 + : !Dim::InvarZ::value ? 2 + : -1) + { + assert(inj_dim_idx_ >= 0); + } auto get_advanced_prt(Double3 pos, real_t wni) { @@ -341,7 +347,7 @@ class Inflow } const Grid_t& grid_; - AdvanceParticle advance_; + AdvanceParticle advance_; SetupParticles setup_particles_; psc_particle_np np_; double (*offset_in_cell_dist_)(); @@ -375,7 +381,7 @@ TEST(TestSetupParticlesInflow, Advance) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half, 1); + Inflow inflow(grid, npt, *half); auto prt = inflow.get_advanced_prt(pos, 1.0); @@ -400,7 +406,7 @@ TEST(TestSetupParticlesInflow, InjectIntoCell) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half, 1); + Inflow inflow(grid, npt, *half); TestInjector injector; inflow.inject_into_boundary_cell(injector, {0, 0, 0}); @@ -431,7 +437,7 @@ TEST(TestSetupParticlesInflow, InjectIntoCellFilter) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half, 1); + Inflow inflow(grid, npt, *half); TestInjector injector; inflow.inject_into_boundary_cell(injector, {0, 0, 0}); @@ -455,7 +461,7 @@ TEST(TestSetupParticlesInflow, InjectIntoPatch) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half, 1); + Inflow inflow(grid, npt, *half); TestInjector injector; inflow.inject_into_boundary_patch(injector, grid.patches[0]); @@ -493,7 +499,7 @@ TEST(TestSetupParticlesInflow, InjectIntoBoundaryY) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half, 1); + Inflow inflow(grid, npt, *half); std::vector injectors = {TestInjector(), TestInjector(), TestInjector(), TestInjector()}; @@ -536,7 +542,7 @@ TEST(TestSetupParticlesInflow, InjectIntoBoundaryZ) Double3 T = {0.0, 0.0, 0.0}; psc_particle_npt npt = {0, 1.0, u, T}; - Inflow inflow(grid, npt, *half, 2); + Inflow inflow(grid, npt, *half); std::vector injectors = {TestInjector(), TestInjector(), TestInjector(), TestInjector()}; From 14b76dfdac5f60cd7a27b016442682797907ddb5 Mon Sep 17 00:00:00 2001 From: James McClung Date: Thu, 27 Mar 2025 15:14:53 -0400 Subject: [PATCH 26/29] test_bnd_prt_inflow: rm/rename some tests the removed ones were copied and unchanged --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 160 +---------------------- 1 file changed, 2 insertions(+), 158 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 05caa38922..0ddc808e8d 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -35,163 +35,7 @@ using MparticlesInflowTestTypes = TYPED_TEST_SUITE(MparticlesInflowTest, MparticlesInflowTestTypes); -// ====================================================================== -// MparticlesInflowTest - -template -struct MparticlesInflowTest : ::testing::Test -{ - using Mparticles = typename T::Mparticles; - using Particle = typename Mparticles::Particle; - using MakeGrid = typename T::MakeGrid; - - MparticlesInflowTest() : grid_{MakeGrid{}()} - { - grid_.kinds.emplace_back(Grid_t::Kind(1., 1., "test_species")); - } - - template - Mparticles mk_mprts(tag dummy) - { - Mparticles mprts(grid_); - mprts.define_species("test_species", 1., 1., 100, 10, 10, 0); - return mprts; - } - - Mparticles mk_mprts() { return mk_mprts(static_cast(nullptr)); } - - template - void inject_test_particles(_Mparticles& mprts, int n_prts) - { - auto inj = mprts.injector(); - for (int p = 0; p < mprts.n_patches(); ++p) { - auto injector = inj[p]; - auto& patch = mprts.grid().patches[p]; - for (int n = 0; n < n_prts; n++) { - double nn = double(n) / n_prts; - auto L = patch.xe - patch.xb; - psc::particle::Inject prt = {{patch.xb[0] + nn * L[0], - patch.xb[1] + nn * L[1], - patch.xb[2] + nn * L[2]}, - {}, - 1., - 0}; - injector(prt); - } - } - } - - const Grid_t& grid() { return grid_; } - -private: - Grid_t grid_; -}; - -// ---------------------------------------------------------------------- -// Inject - -TYPED_TEST(MparticlesInflowTest, Inject) -{ - const int n_prts = 4; - - auto mprts = this->mk_mprts(); - - this->inject_test_particles(mprts, n_prts); -} - -// ----------------------------------------------------------------------- -// Inject2 - -TYPED_TEST(MparticlesInflowTest, Inject2) -{ - using Mparticles = typename TypeParam::Mparticles; - const int n_prts = 1; - - auto mprts = this->mk_mprts(); - - // FIXME, kinda the cheap way out - if (mprts.n_patches() != 4) { - return; - } - ASSERT_EQ(mprts.n_patches(), 4); - - int nn = 0; - { - auto inj = mprts.injector(); - for (int p = 0; p < mprts.n_patches(); ++p) { - auto injector = inj[p]; - auto& patch = mprts.grid().patches[p]; - for (int n = 0; n < n_prts; n++) { - psc::particle::Inject prt{ - .5 * (patch.xb + patch.xe), {}, double(nn), 0}; - // use weight to store particle number for testing - injector(prt); - nn++; - } - } - } - - EXPECT_EQ(mprts.size(), 4); - - auto n_prts_by_patch = mprts.sizeByPatch(); - EXPECT_EQ(n_prts_by_patch, - std::vector({n_prts, n_prts, n_prts, n_prts})); - - nn = 0; - auto accessor = mprts.accessor(); - for (int p = 0; p < mprts.n_patches(); ++p) { - auto& patch = mprts.grid().patches[p]; - for (auto prt : accessor[p]) { - // xm is patch-relative position - auto xm = .5 * (patch.xe - patch.xb); - EXPECT_EQ(prt.x()[0], xm[0]); - EXPECT_EQ(prt.x()[1], xm[1]); - EXPECT_EQ(prt.x()[2], xm[2]); - EXPECT_EQ(prt.qni_wni(), nn); - EXPECT_EQ(prt.w(), nn); - EXPECT_EQ(prt.kind(), 0); - - auto x = .5 * (patch.xb + patch.xe); - EXPECT_EQ(prt.position()[0], x[0]); - EXPECT_EQ(prt.position()[1], x[1]); - EXPECT_EQ(prt.position()[2], x[2]); - nn++; - } - } -} - -// ---------------------------------------------------------------------- -// Injector - -TYPED_TEST(MparticlesInflowTest, Injector) -{ - const int n_prts = 4; - - auto mprts = this->mk_mprts(); - - this->inject_test_particles(mprts, n_prts); - - auto accessor = mprts.accessor(); - for (int p = 0; p < mprts.n_patches(); ++p) { - auto prts = accessor[p]; - auto& patch = mprts.grid().patches[p]; - EXPECT_EQ(prts.size(), n_prts); - int n = 0; - for (auto prt : prts) { - double nn = double(n) / n_prts; - auto L = patch.xe - patch.xb; - auto x = prt.position(); - EXPECT_EQ(x[0], patch.xb[0] + nn * L[0]); - EXPECT_EQ(x[1], patch.xb[1] + nn * L[1]); - EXPECT_EQ(x[2], patch.xb[2] + nn * L[2]); - EXPECT_EQ(prt.w(), 1.); - EXPECT_EQ(prt.kind(), 0); - ++n; - } - } -} - -TEST(TestSetupParticlesInflow, Simple) +TEST(TestSetupParticlesInflow, SetupSimple) { using Mparticles = MparticlesDouble; @@ -218,7 +62,7 @@ TEST(TestSetupParticlesInflow, Simple) EXPECT_EQ(prt.tag, tag); } -TEST(TestSetupParticlesInflow, Maxwellian) +TEST(TestSetupParticlesInflow, SetupMaxwellian) { using Mparticles = MparticlesDouble; From 207af5e23afa9034ace1852aac2a123c4dabfcc4 Mon Sep 17 00:00:00 2001 From: James McClung Date: Fri, 28 Mar 2025 11:24:58 -0400 Subject: [PATCH 27/29] test_bnd_prt_inflow: make inj_dim_idx_ static --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 11 +++++------ 1 file changed, 5 insertions(+), 6 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 0ddc808e8d..b82e1af29a 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -110,11 +110,7 @@ class Inflow // FIXME np.p is a std::function and is called many times; better to use a // lambda np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)}, - offset_in_cell_dist_(offset_in_cell_dist), - inj_dim_idx_(!Dim::InvarX::value ? 0 - : !Dim::InvarY::value ? 1 - : !Dim::InvarZ::value ? 2 - : -1) + offset_in_cell_dist_(offset_in_cell_dist) { assert(inj_dim_idx_ >= 0); } @@ -195,7 +191,10 @@ class Inflow SetupParticles setup_particles_; psc_particle_np np_; double (*offset_in_cell_dist_)(); - int inj_dim_idx_; + const static int inj_dim_idx_ = (!Dim::InvarX::value ? 0 + : !Dim::InvarY::value ? 1 + : !Dim::InvarZ::value ? 2 + : -1); }; class TestInjector From 06761a15d9e85eed6bb06a922ca190e462e0b395 Mon Sep 17 00:00:00 2001 From: James McClung Date: Fri, 28 Mar 2025 11:26:44 -0400 Subject: [PATCH 28/29] test_bnd_prt_inflow: rename to INJECT_DIM_IDX_ --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 24 ++++++++++++------------ 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index b82e1af29a..82a2d02c8c 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -112,7 +112,7 @@ class Inflow np_{npt.kind, npt.n, setup_particles_.createMaxwellian(npt)}, offset_in_cell_dist_(offset_in_cell_dist) { - assert(inj_dim_idx_ >= 0); + assert(INJECT_DIM_IDX_ >= 0); } auto get_advanced_prt(Double3 pos, real_t wni) @@ -129,8 +129,8 @@ class Inflow void inject_into_boundary_cell(Injector& injector, Int3 boundary_cell_global_idx) { - assert(boundary_cell_global_idx[inj_dim_idx_] == 0); - boundary_cell_global_idx[inj_dim_idx_] = -1; + assert(boundary_cell_global_idx[INJECT_DIM_IDX_] == 0); + boundary_cell_global_idx[INJECT_DIM_IDX_] = -1; int n_in_cell = setup_particles_.get_n_in_cell(np_); double wni = setup_particles_.getWeight(np_, n_in_cell); @@ -143,7 +143,7 @@ class Inflow grid_.domain.corner; auto prt = get_advanced_prt(pos, wni); - if (prt.x[inj_dim_idx_] < grid_.domain.corner[inj_dim_idx_]) { + if (prt.x[INJECT_DIM_IDX_] < grid_.domain.corner[INJECT_DIM_IDX_]) { continue; } @@ -158,10 +158,10 @@ class Inflow Int3 ilo = boundary_patch.off; Int3 ihi = ilo + grid_.ldims; - assert(ilo[inj_dim_idx_] == 0); + assert(ilo[INJECT_DIM_IDX_] == 0); - int dim1 = std::min((inj_dim_idx_ + 1) % 3, (inj_dim_idx_ + 2) % 3); - int dim2 = std::max((inj_dim_idx_ + 1) % 3, (inj_dim_idx_ + 2) % 3); + int dim1 = std::min((INJECT_DIM_IDX_ + 1) % 3, (INJECT_DIM_IDX_ + 2) % 3); + int dim2 = std::max((INJECT_DIM_IDX_ + 1) % 3, (INJECT_DIM_IDX_ + 2) % 3); for (Int3 cell_idx = ilo; cell_idx[dim1] < ihi[dim1]; cell_idx[dim1]++) { for (cell_idx[dim2] = ilo[dim2]; cell_idx[dim2] < ihi[dim2]; @@ -177,7 +177,7 @@ class Inflow for (int patch_idx = 0; patch_idx < grid_.n_patches(); patch_idx++) { const auto& patch = grid_.patches[patch_idx]; - if (patch.off[inj_dim_idx_] != 0) { + if (patch.off[INJECT_DIM_IDX_] != 0) { continue; } @@ -191,10 +191,10 @@ class Inflow SetupParticles setup_particles_; psc_particle_np np_; double (*offset_in_cell_dist_)(); - const static int inj_dim_idx_ = (!Dim::InvarX::value ? 0 - : !Dim::InvarY::value ? 1 - : !Dim::InvarZ::value ? 2 - : -1); + const static int INJECT_DIM_IDX_ = (!Dim::InvarX::value ? 0 + : !Dim::InvarY::value ? 1 + : !Dim::InvarZ::value ? 2 + : -1); }; class TestInjector From 05e2f55d01d42cec2134592e3343094c64b9654c Mon Sep 17 00:00:00 2001 From: James McClung Date: Fri, 28 Mar 2025 11:32:50 -0400 Subject: [PATCH 29/29] test_bnd_prt_inflow: make half() static --- src/libpsc/tests/test_bnd_prt_inflow.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/libpsc/tests/test_bnd_prt_inflow.cxx b/src/libpsc/tests/test_bnd_prt_inflow.cxx index 82a2d02c8c..ff16293403 100644 --- a/src/libpsc/tests/test_bnd_prt_inflow.cxx +++ b/src/libpsc/tests/test_bnd_prt_inflow.cxx @@ -205,7 +205,7 @@ class TestInjector std::vector prts; }; -double half() { return 0.5; } +static double half() { return 0.5; } TEST(TestSetupParticlesInflow, Advance) {