diff --git a/src/core/numerics/interpolator/interpolator.h b/src/core/numerics/interpolator/interpolator.h index a4ff6493c..501b29350 100644 --- a/src/core/numerics/interpolator/interpolator.h +++ b/src/core/numerics/interpolator/interpolator.h @@ -1,6 +1,8 @@ #ifndef PHARE_CORE_NUMERICS_INTERPOLATOR_INTERPOLATOR_H #define PHARE_CORE_NUMERICS_INTERPOLATOR_INTERPOLATOR_H + + #include #include @@ -360,17 +362,17 @@ namespace core - auto const& yDenStartIndex = startIndex[static_cast(densityCentering[0])][0]; - auto const& yDenWeights = weights[static_cast(densityCentering[0])][0]; + auto const& yDenStartIndex = startIndex[static_cast(densityCentering[0])][1]; + auto const& yDenWeights = weights[static_cast(densityCentering[0])][1]; - auto const& yXFluxStartIndex = startIndex[static_cast(fluxCentering[0][1])][0]; - auto const& yXFluxWeights = weights[static_cast(fluxCentering[0][1])][0]; + auto const& yXFluxStartIndex = startIndex[static_cast(fluxCentering[0][1])][1]; + auto const& yXFluxWeights = weights[static_cast(fluxCentering[0][1])][1]; - auto const& yYFluxStartIndex = startIndex[static_cast(fluxCentering[1][1])][0]; - auto const& yYFluxWeights = weights[static_cast(fluxCentering[1][1])][0]; + auto const& yYFluxStartIndex = startIndex[static_cast(fluxCentering[1][1])][1]; + auto const& yYFluxWeights = weights[static_cast(fluxCentering[1][1])][1]; - auto const& yZFluxStartIndex = startIndex[static_cast(fluxCentering[2][1])][0]; - auto const& yZFluxWeights = weights[static_cast(fluxCentering[2][1])][0]; + auto const& yZFluxStartIndex = startIndex[static_cast(fluxCentering[2][1])][1]; + auto const& yZFluxWeights = weights[static_cast(fluxCentering[2][1])][1]; @@ -441,32 +443,32 @@ namespace core - auto const& yDenStartIndex = startIndex[static_cast(densityCentering[0])][0]; - auto const& yDenWeights = weights[static_cast(densityCentering[0])][0]; + auto const& yDenStartIndex = startIndex[static_cast(densityCentering[0])][1]; + auto const& yDenWeights = weights[static_cast(densityCentering[0])][1]; - auto const& yXFluxStartIndex = startIndex[static_cast(fluxCentering[0][1])][0]; - auto const& yXFluxWeights = weights[static_cast(fluxCentering[0][1])][0]; + auto const& yXFluxStartIndex = startIndex[static_cast(fluxCentering[0][1])][1]; + auto const& yXFluxWeights = weights[static_cast(fluxCentering[0][1])][1]; - auto const& yYFluxStartIndex = startIndex[static_cast(fluxCentering[1][1])][0]; - auto const& yYFluxWeights = weights[static_cast(fluxCentering[1][1])][0]; + auto const& yYFluxStartIndex = startIndex[static_cast(fluxCentering[1][1])][1]; + auto const& yYFluxWeights = weights[static_cast(fluxCentering[1][1])][1]; - auto const& yZFluxStartIndex = startIndex[static_cast(fluxCentering[2][1])][0]; - auto const& yZFluxWeights = weights[static_cast(fluxCentering[2][1])][0]; + auto const& yZFluxStartIndex = startIndex[static_cast(fluxCentering[2][1])][1]; + auto const& yZFluxWeights = weights[static_cast(fluxCentering[2][1])][1]; - auto const& zDenStartIndex = startIndex[static_cast(densityCentering[0])][0]; - auto const& zDenWeights = weights[static_cast(densityCentering[0])][0]; + auto const& zDenStartIndex = startIndex[static_cast(densityCentering[0])][2]; + auto const& zDenWeights = weights[static_cast(densityCentering[0])][2]; - auto const& zXFluxStartIndex = startIndex[static_cast(fluxCentering[0][2])][0]; - auto const& zXFluxWeights = weights[static_cast(fluxCentering[0][2])][0]; + auto const& zXFluxStartIndex = startIndex[static_cast(fluxCentering[0][2])][2]; + auto const& zXFluxWeights = weights[static_cast(fluxCentering[0][2])][2]; - auto const& zYFluxStartIndex = startIndex[static_cast(fluxCentering[1][2])][0]; - auto const& zYFluxWeights = weights[static_cast(fluxCentering[1][2])][0]; + auto const& zYFluxStartIndex = startIndex[static_cast(fluxCentering[1][2])][2]; + auto const& zYFluxWeights = weights[static_cast(fluxCentering[1][2])][2]; - auto const& zZFluxStartIndex = startIndex[static_cast(fluxCentering[2][2])][0]; - auto const& zZFluxWeights = weights[static_cast(fluxCentering[2][2])][0]; + auto const& zZFluxStartIndex = startIndex[static_cast(fluxCentering[2][2])][2]; + auto const& zZFluxWeights = weights[static_cast(fluxCentering[2][2])][2]; auto const partRho = particle.weight * coef / cellVolume; auto const xPartFlux = particle.v[0] * particle.weight * coef / cellVolume; diff --git a/src/core/utilities/types.h b/src/core/utilities/types.h index 2e43350ed..1cab23298 100644 --- a/src/core/utilities/types.h +++ b/src/core/utilities/types.h @@ -152,6 +152,15 @@ namespace core return sized_array(arr); } + template + constexpr decltype(auto) ConstArray(Type val = 0) + { + std::array arr{}; + for (uint8_t i = 0; i < size; i++) + arr[i] = val; + return arr; + } + } // namespace core } // namespace PHARE diff --git a/tests/core/numerics/interpolator/test_main.cpp b/tests/core/numerics/interpolator/test_main.cpp index 1287d35b7..b3e89233e 100644 --- a/tests/core/numerics/interpolator/test_main.cpp +++ b/tests/core/numerics/interpolator/test_main.cpp @@ -1,3 +1,5 @@ + + #include "gmock/gmock.h" #include "gtest/gtest.h" @@ -525,7 +527,7 @@ TYPED_TEST(A3DInterpolator, canComputeAllEMfieldsAtParticle) template -class ACollectionOfParticles : public ::testing::Test +class ACollectionOfParticles_1d : public ::testing::Test { public: static constexpr uint32_t nx = 40; @@ -543,7 +545,7 @@ class ACollectionOfParticles : public ::testing::Test - ACollectionOfParticles() + ACollectionOfParticles_1d() : part{} , particles{} , rho{"field", HybridQuantity::Scalar::rho, nx} @@ -669,38 +671,77 @@ class ACollectionOfParticles : public ::testing::Test protected: Interpolator interpolator; }; +TYPED_TEST_SUITE_P(ACollectionOfParticles_1d); - - -TYPED_TEST_SUITE_P(ACollectionOfParticles); - - -TYPED_TEST_P(ACollectionOfParticles, DepositCorrectlyTheirWeight) +TYPED_TEST_P(ACollectionOfParticles_1d, DepositCorrectlyTheirWeight_1d) { - EXPECT_DOUBLE_EQ(this->rho(25u), 1.0); - EXPECT_DOUBLE_EQ(this->vx(25u), 2.0); - EXPECT_DOUBLE_EQ(this->vy(25u), -1.0); - EXPECT_DOUBLE_EQ(this->vz(25u), 1.0); + EXPECT_DOUBLE_EQ(this->rho(25), 1.0); + EXPECT_DOUBLE_EQ(this->vx(25), 2.0); + EXPECT_DOUBLE_EQ(this->vy(25), -1.0); + EXPECT_DOUBLE_EQ(this->vz(25), 1.0); } +REGISTER_TYPED_TEST_SUITE_P(ACollectionOfParticles_1d, DepositCorrectlyTheirWeight_1d); + +using MyTypes = ::testing::Types, Interpolator<1, 2>, Interpolator<1, 3>>; +INSTANTIATE_TYPED_TEST_SUITE_P(testInterpolator, ACollectionOfParticles_1d, MyTypes); +template +struct ACollectionOfParticles_2d : public ::testing::Test +{ + static constexpr size_t dim = 2; + static constexpr uint32_t nx = 15, ny = 15; + static constexpr int start = 0, end = 5; + + ACollectionOfParticles_2d() + : rho{"field", HybridQuantity::Scalar::rho, nx, ny} + , vx{"v_x", HybridQuantity::Scalar::Vx, nx, ny} + , vy{"v_y", HybridQuantity::Scalar::Vy, nx, ny} + , vz{"v_z", HybridQuantity::Scalar::Vz, nx, ny} + , v{"v", HybridQuantity::Vector::V} + { + v.setBuffer("v_x", &vx); + v.setBuffer("v_y", &vy); + v.setBuffer("v_z", &vz); -REGISTER_TYPED_TEST_SUITE_P(ACollectionOfParticles, DepositCorrectlyTheirWeight); + for (int i = start; i < end; i++) + for (int j = start; j < end; j++) + { + auto& part = particles.emplace_back(); + part.iCell = {i, j}; + part.delta = ConstArray(.5); + part.weight = .1 * layout.meshSize()[0]; + part.v[0] = +2.; + part.v[1] = -1.; + part.v[2] = +1.; + } + interpolator(std::begin(particles), std::end(particles), rho, v, layout); + } + GridLayout> layout{ + ConstArray(.1), {nx, ny}, ConstArray(0)}; -/*using GridLayoutYee2DO1 = GridLayout>; -using GridLayoutYee2DO2 = GridLayout>; -using GridLayoutYee2DO3 = GridLayout>; -using GridLayoutYee3DO1 = GridLayout>; -using GridLayoutYee3DO2 = GridLayout>; -using GridLayoutYee3DO3 = GridLayout>;*/ + ParticleArray particles; + Field, typename HybridQuantity::Scalar> rho, vx, vy, vz; + VecField, HybridQuantity> v; + Interpolator interpolator; +}; +TYPED_TEST_SUITE_P(ACollectionOfParticles_2d); +TYPED_TEST_P(ACollectionOfParticles_2d, DepositCorrectlyTheirWeight_2d) +{ + EXPECT_DOUBLE_EQ(this->rho(7, 7), 1.0); + EXPECT_DOUBLE_EQ(this->vx(7, 7), 2.0); + EXPECT_DOUBLE_EQ(this->vy(7, 7), -1.0); + EXPECT_DOUBLE_EQ(this->vz(7, 7), 1.0); +} +REGISTER_TYPED_TEST_SUITE_P(ACollectionOfParticles_2d, DepositCorrectlyTheirWeight_2d); -using MyTypes = ::testing::Types, Interpolator<1, 2>, Interpolator<1, 3>>; -INSTANTIATE_TYPED_TEST_SUITE_P(testInterpolator, ACollectionOfParticles, MyTypes); +using My2dTypes = ::testing::Types, Interpolator<2, 2>, Interpolator<2, 3>>; +INSTANTIATE_TYPED_TEST_SUITE_P(testInterpolator, ACollectionOfParticles_2d, My2dTypes); int main(int argc, char** argv)