minimize maxwellian vector capacity somewhat

PHAREHUB · Oct 3, 2024 · 4590a2b · 4590a2b
1 parent 6fca2f7
commit 4590a2b
Show file tree

Hide file tree

Showing 4 changed files with 86 additions and 19 deletions.
diff --git a/src/core/data/ions/particle_initializers/maxwellian_particle_initializer.hpp b/src/core/data/ions/particle_initializers/maxwellian_particle_initializer.hpp
@@ -178,6 +178,20 @@ void MaxwellianParticleInitializer<ParticleArray, GridLayout>::loadParticles(
     auto randGen           = getRNG(rngSeed_);
     ParticleDeltaDistribution<double> deltaDistrib;
 
+    std::size_t const expected_size = [&]() {
+        std::size_t size = 0;
+        for (std::size_t flatCellIdx = 0; flatCellIdx < ndCellIndices.size(); ++flatCellIdx)
+            if (n[flatCellIdx] >= densityCutOff_)
+                ++size;
+        return size * nbrParticlePerCell_;
+    }();
+
+    double const over_alloc_factor = .1; // todo from dict
+    std::size_t const allocate_size
+        = expected_size
+          + (layout.AMRBox().surface_cell_count() * (nbrParticlePerCell_ * over_alloc_factor));
+    particles.reserve(allocate_size);
+
     for (std::size_t flatCellIdx = 0; flatCellIdx < ndCellIndices.size(); ++flatCellIdx)
     {
         if (n[flatCellIdx] < densityCutOff_)

diff --git a/src/core/data/particles/particle.hpp b/src/core/data/particles/particle.hpp
@@ -39,6 +39,9 @@ NO_DISCARD auto cellAsPoint(Particle const& particle)
 template<size_t dim>
 struct Particle
 {
+    std::array<std::uint8_t, 3> constexpr static sizes = {52, 64, 76};
+    static_assert(sizeof(Particle<dim>) == sizes[dim]);
+
     static_assert(dim > 0 and dim < 4, "Only dimensions 1,2,3 are supported.");
     static const size_t dimension = dim;
 
@@ -54,12 +57,11 @@ struct Particle
 
     Particle() = default;
 
-    double weight = 0;
-    double charge = 0;
-
-    std::array<int, dim> iCell    = ConstArray<int, dim>();
-    std::array<double, dim> delta = ConstArray<double, dim>();
-    std::array<double, 3> v       = ConstArray<double, 3>();
+    double weight                 = 0;                         // 8
+    double charge                 = 0;                         // 8
+    std::array<int, dim> iCell    = ConstArray<int, dim>();    // 4 * dim
+    std::array<double, dim> delta = ConstArray<double, dim>(); // 8 * dim
+    std::array<double, 3> v       = ConstArray<double, 3>();   // 8 * 3
 
     NO_DISCARD bool operator==(Particle<dim> const& that) const
     {
@@ -121,9 +123,9 @@ inline constexpr auto is_phare_particle_type
 
 template<std::size_t dim, template<std::size_t> typename ParticleA,
          template<std::size_t> typename ParticleB>
-NO_DISCARD typename std::enable_if_t<
-    is_phare_particle_type<dim, ParticleA<dim>> and is_phare_particle_type<dim, ParticleB<dim>>,
-    bool>
+NO_DISCARD typename std::enable_if_t<is_phare_particle_type<dim, ParticleA<dim>>
+                                         and is_phare_particle_type<dim, ParticleB<dim>>,
+                                     bool>
 operator==(ParticleA<dim> const& particleA, ParticleB<dim> const& particleB)
 {
     return particleA.weight == particleB.weight and //

diff --git a/src/core/utilities/box/box.hpp b/src/core/utilities/box/box.hpp
@@ -26,7 +26,7 @@ template<typename Type, std::size_t dim>
 struct Box
 {
     static const size_t dimension = dim;
-
+    using value_type              = Type;
 
     Point<Type, dim> lower;
     Point<Type, dim> upper;
@@ -130,9 +130,16 @@ struct Box
             return iterator{this, {upper[0] + 1, upper[1] + 1, upper[2] + 1}};
         }
     }
-    using value_type = Type;
 
 
+    NO_DISCARD auto surface_cell_count() const
+    {
+        // assumes box never smaller than 3 in any direction
+        auto const shape_ = shape();
+        auto const nested = shape_ - 2;
+        return core::product(shape_) - core::product(nested);
+    }
+
     NO_DISCARD constexpr static std::size_t nbrRemainBoxes()
     {
         if constexpr (dim == 1)

diff --git a/tests/core/data/maxwellian_particle_initializer/test_maxwellian_particle_initializer.cpp b/tests/core/data/maxwellian_particle_initializer/test_maxwellian_particle_initializer.cpp
@@ -9,22 +9,23 @@
 #include "core/utilities/box/box.hpp"
 #include "core/utilities/point/point.hpp"
 
-#include "gmock/gmock.h"
 #include "gtest/gtest.h"
 
 #include "tests/initializer/init_functions.hpp"
-using namespace PHARE::initializer::test_fn::func_1d; // density/etc are here
+#include "tests/core/data/gridlayout/test_gridlayout.hpp"
 
 using namespace PHARE::core;
-using namespace PHARE::initializer;
+
+namespace PHARE::initializer::test_fn::func_1d
+{
 
 
 class AMaxwellianParticleInitializer1D : public ::testing::Test
 {
 private:
-    using GridLayoutT       = GridLayout<GridLayoutImplYee<1, 1>>;
     using ParticleArrayT    = ParticleArray<1>;
     using InitFunctionArray = std::array<InitFunction<1>, 3>;
+    using GridLayoutT       = GridLayout<GridLayoutImplYee<1, 1>>;
 
 public:
     AMaxwellianParticleInitializer1D()
@@ -34,10 +35,8 @@ class AMaxwellianParticleInitializer1D : public ::testing::Test
               density, InitFunctionArray{vx, vy, vz}, InitFunctionArray{vthx, vthy, vthz}, 1.,
               nbrParticlesPerCell)}
     {
-        //
     }
 
-
     GridLayoutT layout;
     ParticleArrayT particles;
     std::uint32_t nbrParticlesPerCell{10000};
@@ -46,7 +45,6 @@ class AMaxwellianParticleInitializer1D : public ::testing::Test
 
 
 
-
 TEST_F(AMaxwellianParticleInitializer1D, loadsTheCorrectNbrOfParticles)
 {
     auto nbrCells             = layout.nbrCells();
@@ -57,7 +55,6 @@ TEST_F(AMaxwellianParticleInitializer1D, loadsTheCorrectNbrOfParticles)
 
 
 
-
 TEST_F(AMaxwellianParticleInitializer1D, loadsParticlesInTheDomain)
 {
     initializer->loadParticles(particles, layout);
@@ -73,8 +70,55 @@ TEST_F(AMaxwellianParticleInitializer1D, loadsParticlesInTheDomain)
     }
 }
 
+} // namespace PHARE::initializer::test_fn::func_1d
+
+
+namespace PHARE::initializer::test_fn::func_2d
+{
+
+class AMaxwellianParticleInitializer2D : public ::testing::Test
+{
+private:
+    using ParticleArrayT    = ParticleArray<2>;
+    using InitFunctionArray = std::array<InitFunction<2>, 3>;
+    using GridLayoutT       = GridLayout<GridLayoutImplYee<2, 1>>;
+
+
+
+public:
+    AMaxwellianParticleInitializer2D()
+        : layout{50}
+        , initializer{std::make_unique<MaxwellianParticleInitializer<ParticleArrayT, GridLayoutT>>(
+              density, InitFunctionArray{vx, vy, vz}, InitFunctionArray{vthx, vthy, vthz}, 1.,
+              nbrParticlesPerCell)}
+    {
+    }
+
+    TestGridLayout<GridLayoutT> layout;
+    ParticleArrayT particles{layout.AMRBox()};
+    std::uint32_t nbrParticlesPerCell{600};
+    std::unique_ptr<MaxwellianParticleInitializer<ParticleArrayT, GridLayoutT>> initializer;
+};
+
+TEST_F(AMaxwellianParticleInitializer2D, loadsTheCorrectNbrOfParticles)
+{
+    // vector push back allocation observations
+    // 100 ppc = 262144 - 250000 = 12144 == 12144 * 64 / 1e6 == .7MB overallocated
+    // 600 ppc = 2097152 - 1500000 = 597152 * 64 / 1e6 == 38MB overallocated
+
+    auto const expectedNbrParticles = nbrParticlesPerCell * product(layout.AMRBox().shape());
+    initializer->loadParticles(particles, layout);
+    EXPECT_EQ(expectedNbrParticles, particles.size());
+    auto outer_cell_count = std::pow(50, 2) - std::pow(48, 2);
+    EXPECT_EQ(particles.capacity(), particles.size() + (outer_cell_count * 10));
+
+    // new method
+    // 100 ppc = (1511760 - 1500000) * 64 / 1e6 == 0.12544 overallocated
+    // 600 ppc = (1511760 - 1500000) * 64 / 1e6 == 0.75264 overallocated
+}
 
 
+} // namespace PHARE::initializer::test_fn::func_2d
 
 int main(int argc, char** argv)
 {