minimize maxwellian vector capacity somewhat

src/amr/data/field/refine/electric_field_refiner.hpp

@@ -4,13 +4,15 @@
 
 #include "core/def/phare_mpi.hpp"
 
-#include <SAMRAI/hier/Box.h>
-
+#include "amr/amr_constants.hpp"
 #include "amr/resources_manager/amr_utils.hpp"
+
 #include "core/utilities/constants.hpp"
 #include "core/data/grid/gridlayoutdefs.hpp"
 #include "core/utilities/point/point.hpp"
 
+#include <SAMRAI/hier/Box.h>
+
 #include <cstddef>
 
 namespace PHARE::amr
@@ -43,8 +45,8 @@ class ElectricFieldRefiner
     {
         TBOX_ASSERT(coarseField.physicalQuantity() == fineField.physicalQuantity());
 
-        auto const locFineIdx          = AMRToLocal(fineIndex, fineBox_);
-        auto constexpr refinementRatio = 2;
+        auto const locFineIdx = AMRToLocal(fineIndex, fineBox_);
+
         auto coarseIdx{fineIndex};
         for (auto& idx : coarseIdx)
             idx = idx / refinementRatio;

src/amr/data/particles/refine/particles_data_split.hpp

@@ -2,30 +2,30 @@
 #define PHARE_PARTICLES_DATA_SPLIT_HPP
 
 
+#include "core/def.hpp"
+#include "phare_core.hpp"
 #include "core/def/phare_mpi.hpp"
+#include "core/utilities/constants.hpp"
 
-#include "core/def.hpp"
-#include "amr/data/particles/particles_data.hpp"
-#include "amr/resources_manager/amr_utils.hpp"
 #include "split.hpp"
-#include "core/utilities/constants.hpp"
-#include "phare_core.hpp"
 #include "amr/amr_constants.hpp"
+#include "amr/utilities/box/amr_box.hpp"
+#include "amr/resources_manager/amr_utils.hpp"
+#include "amr/data/particles/particles_data.hpp"
 
 #include <SAMRAI/geom/CartesianPatchGeometry.h>
 #include <SAMRAI/hier/Box.h>
 #include <SAMRAI/hier/RefineOperator.h>
 #include <SAMRAI/pdat/CellOverlap.h>
 
-#include <functional>
 
 
 namespace PHARE
 {
 namespace amr
 {
-    enum class ParticlesDataSplitType {
-        coarseBoundary,
+    enum class ParticlesDataSplitType : std::uint8_t {
+        coarseBoundary = 0,
         interior,
         coarseBoundaryOld,
         coarseBoundaryNew
@@ -53,6 +53,12 @@ namespace amr
     template<typename ParticleArray, ParticlesDataSplitType splitType, typename Splitter>
     class ParticlesRefineOperator : public SAMRAI::hier::RefineOperator
     {
+        using Particle_t = typename ParticleArray::value_type;
+        static constexpr bool putParticlesInCoarseBoundary
+            = splitType == ParticlesDataSplitType::coarseBoundary
+              || splitType == ParticlesDataSplitType::coarseBoundaryOld
+              || splitType == ParticlesDataSplitType::coarseBoundaryNew;
+
     public:
         static constexpr auto dim           = Splitter::dimension;
         static constexpr auto interpOrder   = Splitter::interp_order;
@@ -99,7 +105,7 @@ namespace amr
                 = std::dynamic_pointer_cast<ParticlesData<ParticleArray>>(
                     source.getPatchData(sourceComponent));
 
-            // Finnaly we need the cartesion geometry of both patch.
+            // Finaly we need the cartesion geometry of both patch.
             auto patchGeomDestination
                 = std::dynamic_pointer_cast<SAMRAI::geom::CartesianPatchGeometry>(
                     destination.getPatchGeometry());
@@ -133,9 +139,15 @@ namespace amr
             }
         }
 
+        template<typename Fn>
+        static void _reserve(ParticleArray& particles, Fn&& counter)
+        { // only reserve for regular refinement of domain to start with
+            if constexpr (!putParticlesInCoarseBoundary)
+                particles.reserve(counter() * nbRefinedPart);
+        }
 
         /** @brief given two ParticlesData (destination and source),
-         * an overlap , a ratio and the geometry of both patches, perform the
+         * an overlap, a ratio and the geometry of both patches, perform the
          * splitting of coarse particles onto the destination patch
          */
         void refine_(ParticlesData<ParticleArray>& destParticlesData,
@@ -166,43 +178,48 @@ namespace amr
             // same for destinationGhostBox and destinationDomainBox the later will allow to get an
             // index relative to the interior
 
+            std::array const particlesArrays{&srcInteriorParticles, &srcGhostParticles};
+
+            auto const count_expected = [&]() {
+                std::size_t incoming_estimate = 0;
+                for (auto const& destinationBox : destBoxes)
+                {
+                    auto const coarseSplitBox
+                        = coarsen(phare_box_from<dim>(getSplitBox(destinationBox)));
+                    for (auto const& sourceParticlesArray : particlesArrays)
+                        incoming_estimate += sourceParticlesArray->nbr_particles_in(coarseSplitBox);
+                }
+                return incoming_estimate;
+            };
+
+            _reserve(destDomainParticles, count_expected);
+
             Splitter split;
 
             // The PatchLevelFillPattern had compute boxes that correspond to the expected filling.
             // In case of a coarseBoundary it will most likely give multiple boxes
             // in case of interior, this will be just one box usually
             for (auto const& destinationBox : destBoxes)
             {
-                std::array particlesArrays{&srcInteriorParticles, &srcGhostParticles};
-                auto splitBox = getSplitBox(destinationBox);
+                auto const splitBox = getSplitBox(destinationBox);
 
-                auto isInDest = [&destinationBox](auto const& particle) //
+                auto const isInDest = [&destinationBox](auto const& particle) //
                 { return isInBox(destinationBox, particle); };
 
-
                 for (auto const& sourceParticlesArray : particlesArrays)
                 {
                     for (auto const& particle : *sourceParticlesArray)
                     {
-                        std::array<typename ParticleArray::value_type, nbRefinedPart>
-                            refinedParticles;
-                        auto particleRefinedPos = toFineGrid<interpOrder>(particle);
+                        auto const particleRefinedPos = toFineGrid<interpOrder>(particle);
 
                         if (isInBox(splitBox, particleRefinedPos))
                         {
+                            std::array<Particle_t, nbRefinedPart> refinedParticles;
                             split(particleRefinedPos, refinedParticles);
 
-
                             // we need to know in which of interior or levelGhostParticlesXXXX
                             // arrays we must put particles
 
-                            bool constexpr putParticlesInCoarseBoundary
-                                = splitType == ParticlesDataSplitType::coarseBoundary
-                                  || splitType == ParticlesDataSplitType::coarseBoundaryOld
-                                  || splitType == ParticlesDataSplitType::coarseBoundaryNew;
-
-
-
                             if constexpr (putParticlesInCoarseBoundary)
                             {
                                 if constexpr (splitType == ParticlesDataSplitType::coarseBoundary)
@@ -243,9 +260,9 @@ namespace amr
                                              std::back_inserter(destDomainParticles), isInDest);
                             }
                         } // end is candidate for split
-                    }     // end loop on particles
-                }         // end loop on source particle arrays
-            }             // loop on destination box
+                    } // end loop on particles
+                } // end loop on source particle arrays
+            } // loop on destination box
         }
 
 

src/amr/utilities/box/amr_box.hpp

@@ -2,13 +2,13 @@
 #define PHARE_AMR_UTILITIES_BOX_BOX_HPP
 
 
+#include "core/def.hpp"
 #include "core/def/phare_mpi.hpp"
+#include "core/utilities/box/box.hpp"
 
+#include "amr/amr_constants.hpp"
 
 #include "SAMRAI/hier/Box.h"
-#include "core/utilities/box/box.hpp"
-#include "core/def.hpp"
-
 
 namespace PHARE::amr
 {
@@ -85,6 +85,22 @@ struct Box : public PHARE::core::Box<Type, dim>
     }
 };
 
+template<std::size_t dim>
+auto refine(core::Box<int, dim> box)
+{
+    for (std::uint8_t di = 0; di < dim; ++di)
+        box.lower[di] *= refinementRatio, box.upper[di] *= refinementRatio;
+    return box;
+}
+template<std::size_t dim>
+auto coarsen(core::Box<int, dim> box)
+{
+    for (std::uint8_t di = 0; di < dim; ++di)
+        box.lower[di] /= refinementRatio, box.upper[di] /= refinementRatio;
+
+    return box;
+}
+
 } // namespace PHARE::amr
 
 #endif

src/core/data/ions/particle_initializers/maxwellian_particle_initializer.hpp

@@ -44,7 +44,7 @@ class MaxwellianParticleInitializer : public ParticleInitializer<ParticleArray,
         std::uint32_t const& nbrParticlesPerCell, std::optional<std::size_t> seed = {},
         Basis const basis                                 = Basis::Cartesian,
         std::array<InputFunction, 3> const& magneticField = {nullptr, nullptr, nullptr},
-        double densityCutOff                              = 1e-5)
+        double const densityCutOff = 1e-5, double const over_alloc_factor = .1)
         : density_{density}
         , bulkVelocity_{bulkVelocity}
         , thermalVelocity_{thermalVelocity}
@@ -54,6 +54,7 @@ class MaxwellianParticleInitializer : public ParticleInitializer<ParticleArray,
         , nbrParticlePerCell_{nbrParticlesPerCell}
         , basis_{basis}
         , rngSeed_{seed}
+        , over_alloc_factor_{over_alloc_factor}
     {
     }
 
@@ -92,6 +93,7 @@ class MaxwellianParticleInitializer : public ParticleInitializer<ParticleArray,
     std::uint32_t nbrParticlePerCell_;
     Basis basis_;
     std::optional<std::size_t> rngSeed_;
+    double over_alloc_factor_ = .1;
 };
 
 
@@ -178,6 +180,17 @@ void MaxwellianParticleInitializer<ParticleArray, GridLayout>::loadParticles(
     auto randGen           = getRNG(rngSeed_);
     ParticleDeltaDistribution<double> deltaDistrib;
 
+    auto const expected_size = std::accumulate(n, n + ndCellIndices.size(), std::size_t{0},
+                                               [&](auto const& sum, auto const& density_value) {
+                                                   return (density_value > densityCutOff_)
+                                                              ? sum + nbrParticlePerCell_
+                                                              : sum;
+                                               });
+
+    auto const incoming_estimate
+        = (layout.AMRBox().surface_cell_count() * (nbrParticlePerCell_ * over_alloc_factor_));
+    particles.reserve(expected_size + incoming_estimate);
+
     for (std::size_t flatCellIdx = 0; flatCellIdx < ndCellIndices.size(); ++flatCellIdx)
     {
         if (n[flatCellIdx] < densityCutOff_)

src/core/data/particles/particle.hpp

@@ -39,6 +39,7 @@ 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(dim > 0 and dim < 4, "Only dimensions 1,2,3 are supported.");
     static const size_t dimension = dim;
 
@@ -54,12 +55,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 +121,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 //

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)

tests/amr/data/particles/refine/CMakeLists.txt

@@ -29,6 +29,4 @@ target_link_libraries(${PROJECT_NAME} PRIVATE
   ${GTEST_LIBS})
 
 
-add_phare_test(${PROJECT_NAME} ${CMAKE_CURRENT_BINARY_DIR})
-
-
+add_no_mpi_phare_test(${PROJECT_NAME} ${CMAKE_CURRENT_BINARY_DIR})