try TC build with rebalance_coasrest always == true

pyphare/pyphare/core/box.py

@@ -76,7 +76,7 @@ class BoxIterator:
     def __init__(self, box):
         self.l = box.lower
         self.u = box.upper
-        self.c = l.copy()
+        self.c = self.l.copy()
         self.c[0] -= 1  # hax for first iteration
 
     def __iter__(self):

pyphare/pyphare/pharesee/particles.py

@@ -191,25 +191,46 @@ def as_tuples(self):
         ]
 
 
+def delta_check_or_show_icells(part1, part2, idx1, idx2, deltol):
+    try:
+        np.testing.assert_allclose(part1.deltas[idx1], part2.deltas[idx2], atol=deltol)
+    except AssertionError as e:
+        diff = part1.deltas[idx1] - part2.deltas[idx2] > deltol
+        idxs = np.where(diff > 0)
+        print("icells ", part1.iCells[idxs], part2.iCells[idxs])
+        raise e
+
+
+def v_check_or_show_icells(part1, part2, idx1, idx2, vi):
+    try:
+        np.testing.assert_allclose(part1.v[idx1, vi], part2.v[idx2, vi], atol=1e-18)
+    except AssertionError as e:
+        diff = part1.v[idx1, vi] - part2.v[idx2, vi] > 1e-18
+        idxs = np.where(diff > 0)
+        print("icells ", part1.iCells[idxs], part2.iCells[idxs])
+        raise e
+
+
 def all_assert_sorted(part1, part2):
     idx1 = _arg_sort(part1)
     idx2 = _arg_sort(part2)
 
     np.testing.assert_equal(part1.ndim, part2.ndim)
     np.testing.assert_equal(part1.size(), part2.size())
 
+    np.testing.assert_array_equal(part1.iCells[idx1], part2.iCells[idx2])
+
     deltol = (
         1e-6
         if any([part.deltas.dtype == np.float32 for part in [part1, part2]])
         else 1e-12
     )
 
-    np.testing.assert_array_equal(part1.iCells[idx1], part2.iCells[idx2])
-    np.testing.assert_allclose(part1.deltas[idx1], part2.deltas[idx2], atol=deltol)
+    delta_check_or_show_icells(part1, part2, idx1, idx2, deltol)
 
-    np.testing.assert_allclose(part1.v[idx1, 0], part2.v[idx2, 0], atol=1e-12)
-    np.testing.assert_allclose(part1.v[idx1, 1], part2.v[idx2, 1], atol=1e-12)
-    np.testing.assert_allclose(part1.v[idx1, 2], part2.v[idx2, 2], atol=1e-12)
+    v_check_or_show_icells(part1, part2, idx1, idx2, 0)
+    v_check_or_show_icells(part1, part2, idx1, idx2, 1)
+    v_check_or_show_icells(part1, part2, idx1, idx2, 2)
 
 
 def any_assert(part1, part2):

src/amr/level_initializer/hybrid_level_initializer.hpp

@@ -52,24 +52,24 @@ namespace solver
             auto& hybridModel = static_cast<HybridModel&>(model);
             auto& level       = amr_types::getLevel(*hierarchy, levelNumber);
 
-            auto& hybMessenger = dynamic_cast<HybridMessenger&>(messenger);
+            auto& hybMessenger  = dynamic_cast<HybridMessenger&>(messenger);
+            bool isRegriddingL0 = isRegridding and levelNumber == 0;
 
-            if (isRootLevel(levelNumber))
+            if (isRegridding)
             {
-                PHARE_LOG_START("hybridLevelInitializer::initialize : root level init");
-                model.initialize(level);
-                messenger.fillRootGhosts(model, level, initDataTime);
-                PHARE_LOG_STOP("hybridLevelInitializer::initialize : root level init");
+                std::cout << "regriding level " << levelNumber << "\n";
+                PHARE_LOG_START("hybridLevelInitializer::initialize : regriding block");
+                messenger.regrid(hierarchy, levelNumber, oldLevel, model, initDataTime);
+                PHARE_LOG_STOP("hybridLevelInitializer::initialize : regriding block");
             }
-
             else
             {
-                if (isRegridding)
+                if (isRootLevel(levelNumber))
                 {
-                    std::cout << "regriding level " << levelNumber << "\n";
-                    PHARE_LOG_START("hybridLevelInitializer::initialize : regriding block");
-                    messenger.regrid(hierarchy, levelNumber, oldLevel, model, initDataTime);
-                    PHARE_LOG_STOP("hybridLevelInitializer::initialize : regriding block");
+                    PHARE_LOG_START("hybridLevelInitializer::initialize : root level init");
+                    model.initialize(level);
+                    messenger.fillRootGhosts(model, level, initDataTime);
+                    PHARE_LOG_STOP("hybridLevelInitializer::initialize : root level init");
                 }
                 else
                 {
@@ -109,55 +109,59 @@ namespace solver
             // we are at a sync time across levels and that the time interpolation
             // is not needed. But is still seems to use the messenger tempoeraries like
             // NiOld etc. so prepareStep() must be called, see end of the function.
-            hybMessenger.fillIonMomentGhosts(hybridModel.state.ions, level, initDataTime);
+
+            if (!isRegriddingL0)
+                hybMessenger.fillIonMomentGhosts(hybridModel.state.ions, level, initDataTime);
 
 
             // now moments are known everywhere, compute J and E
             // via Ampere and Ohm
             // this only needs to be done for the root level
             // since otherwise initLevel has done it already
-
-            if (isRootLevel(levelNumber))
-            {
-                auto& B = hybridModel.state.electromag.B;
-                auto& J = hybridModel.state.J;
-
-                for (auto& patch : level)
-                {
-                    auto _      = hybridModel.resourcesManager->setOnPatch(*patch, B, J);
-                    auto layout = PHARE::amr::layoutFromPatch<GridLayoutT>(*patch);
-                    auto __     = core::SetLayout(&layout, ampere_);
-                    ampere_(B, J);
-
-                    hybridModel.resourcesManager->setTime(J, *patch, 0.);
-                }
-                hybMessenger.fillCurrentGhosts(J, levelNumber, 0.);
-
-                auto& electrons = hybridModel.state.electrons;
-                auto& E         = hybridModel.state.electromag.E;
-
-                for (auto& patch : level)
+            if (!isRegriddingL0)
+                if (isRootLevel(levelNumber))
                 {
-                    auto layout = PHARE::amr::layoutFromPatch<GridLayoutT>(*patch);
-                    auto _ = hybridModel.resourcesManager->setOnPatch(*patch, B, E, J, electrons);
-                    electrons.update(layout);
-                    auto& Ve = electrons.velocity();
-                    auto& Ne = electrons.density();
-                    auto& Pe = electrons.pressure();
-                    auto __  = core::SetLayout(&layout, ohm_);
-                    ohm_(Ne, Ve, Pe, B, J, E);
-                    hybridModel.resourcesManager->setTime(E, *patch, 0.);
+                    auto& B = hybridModel.state.electromag.B;
+                    auto& J = hybridModel.state.J;
+
+                    for (auto& patch : level)
+                    {
+                        auto _      = hybridModel.resourcesManager->setOnPatch(*patch, B, J);
+                        auto layout = PHARE::amr::layoutFromPatch<GridLayoutT>(*patch);
+                        auto __     = core::SetLayout(&layout, ampere_);
+                        ampere_(B, J);
+
+                        hybridModel.resourcesManager->setTime(J, *patch, 0.);
+                    }
+                    hybMessenger.fillCurrentGhosts(J, levelNumber, 0.);
+
+                    auto& electrons = hybridModel.state.electrons;
+                    auto& E         = hybridModel.state.electromag.E;
+
+                    for (auto& patch : level)
+                    {
+                        auto layout = PHARE::amr::layoutFromPatch<GridLayoutT>(*patch);
+                        auto _
+                            = hybridModel.resourcesManager->setOnPatch(*patch, B, E, J, electrons);
+                        electrons.update(layout);
+                        auto& Ve = electrons.velocity();
+                        auto& Ne = electrons.density();
+                        auto& Pe = electrons.pressure();
+                        auto __  = core::SetLayout(&layout, ohm_);
+                        ohm_(Ne, Ve, Pe, B, J, E);
+                        hybridModel.resourcesManager->setTime(E, *patch, 0.);
+                    }
+
+                    hybMessenger.fillElectricGhosts(E, levelNumber, 0.);
                 }
 
-                hybMessenger.fillElectricGhosts(E, levelNumber, 0.);
-            }
-
             // quantities have been computed on the level,like the moments and J
             // that we later in the code need to get on level ghost nodes via
             // space and TIME interpolation. We thus need to save current values
             // in "old" messenger temporaries.
             // NOTE :  this may probably be skipped for finest level since, TBC at some point
-            hybMessenger.prepareStep(hybridModel, level, initDataTime);
+            if (!isRegriddingL0)
+                hybMessenger.prepareStep(hybridModel, level, initDataTime);
         }
     };
 } // namespace solver

src/amr/messengers/hybrid_hybrid_messenger_strategy.hpp

@@ -172,23 +172,25 @@ namespace amr
          * not doing so will result in communication to/from that level being impossible
          */
         void registerLevel(std::shared_ptr<SAMRAI::hier::PatchHierarchy> const& hierarchy,
-                           int const levelNumber) override
+                           int const levelNumber,
+                           std::shared_ptr<SAMRAI::hier::PatchLevel> const& oldLevel) override
         {
             auto const level = hierarchy->getPatchLevel(levelNumber);
 
-            magSharedNodesRefiners_.registerLevel(hierarchy, level);
-            elecSharedNodesRefiners_.registerLevel(hierarchy, level);
-            currentSharedNodesRefiners_.registerLevel(hierarchy, level);
+            magSharedNodesRefiners_.registerLevel(hierarchy, level, oldLevel);
+            elecSharedNodesRefiners_.registerLevel(hierarchy, level, oldLevel);
+            currentSharedNodesRefiners_.registerLevel(hierarchy, level, oldLevel);
 
-            magPatchGhostsRefiners_.registerLevel(hierarchy, level);
-            magGhostsRefiners_.registerLevel(hierarchy, level);
-            elecGhostsRefiners_.registerLevel(hierarchy, level);
-            currentGhostsRefiners_.registerLevel(hierarchy, level);
+            magPatchGhostsRefiners_.registerLevel(hierarchy, level, oldLevel);
+            magGhostsRefiners_.registerLevel(hierarchy, level, oldLevel);
+            elecGhostsRefiners_.registerLevel(hierarchy, level, oldLevel);
+            currentGhostsRefiners_.registerLevel(hierarchy, level, oldLevel);
 
-            rhoGhostsRefiners_.registerLevel(hierarchy, level);
-            velGhostsRefiners_.registerLevel(hierarchy, level);
+            rhoGhostsRefiners_.registerLevel(hierarchy, level, oldLevel);
+            velGhostsRefiners_.registerLevel(hierarchy, level, oldLevel);
+
+            patchGhostPartRefiners_.registerLevel(hierarchy, level, oldLevel);
 
-            patchGhostPartRefiners_.registerLevel(hierarchy, level);
 
             // root level is not initialized with a schedule using coarser level data
             // so we don't create these schedules if root level
@@ -198,9 +200,9 @@ namespace amr
                 // those are for refinement
                 magneticInitRefiners_.registerLevel(hierarchy, level);
                 electricInitRefiners_.registerLevel(hierarchy, level);
-                domainParticlesRefiners_.registerLevel(hierarchy, level);
-                lvlGhostPartOldRefiners_.registerLevel(hierarchy, level);
-                lvlGhostPartNewRefiners_.registerLevel(hierarchy, level);
+                domainParticlesRefiners_.registerLevel(hierarchy, level, oldLevel);
+                lvlGhostPartOldRefiners_.registerLevel(hierarchy, level, oldLevel);
+                lvlGhostPartNewRefiners_.registerLevel(hierarchy, level, oldLevel);
 
                 // and these for coarsening
                 magnetoSynchronizers_.registerLevel(hierarchy, level);
@@ -223,9 +225,13 @@ namespace amr
         {
             auto& hybridModel = dynamic_cast<HybridModel&>(model);
             auto level        = hierarchy->getPatchLevel(levelNumber);
+
+            bool isRegriddingL0 = levelNumber == 0 and oldLevel;
+
             magneticInitRefiners_.regrid(hierarchy, levelNumber, oldLevel, initDataTime);
             electricInitRefiners_.regrid(hierarchy, levelNumber, oldLevel, initDataTime);
             domainParticlesRefiners_.regrid(hierarchy, levelNumber, oldLevel, initDataTime);
+
             patchGhostPartRefiners_.fill(levelNumber, initDataTime);
 
             // regriding will fill the new level wherever it has points that overlap
@@ -236,14 +242,18 @@ namespace amr
             // Specifically, we need all fine faces to have equal magnetic field and also
             // equal to that of the shared coarse face.
             // This means that we now need to fill ghosts and border included
-            auto& B = hybridModel.state.electromag.B;
-            auto& E = hybridModel.state.electromag.E;
-            // magSharedNodesRefiners_.fill(B, levelNumber, initDataTime);
-            magGhostsRefiners_.fill(B, levelNumber, initDataTime);
-            // elecSharedNodesRefiners_.fill(E, levelNumber, initDataTime);
-            elecGhostsRefiners_.fill(E, levelNumber, initDataTime);
 
-            fix_magnetic_divergence_(*hierarchy, levelNumber, B);
+            if (!isRegriddingL0)
+            {
+                auto& B = hybridModel.state.electromag.B;
+                auto& E = hybridModel.state.electromag.E;
+                // magSharedNodesRefiners_.fill(B, levelNumber, initDataTime);
+                magGhostsRefiners_.fill(B, levelNumber, initDataTime);
+                // elecSharedNodesRefiners_.fill(E, levelNumber, initDataTime);
+                elecGhostsRefiners_.fill(E, levelNumber, initDataTime);
+
+                fix_magnetic_divergence_(*hierarchy, levelNumber, B);
+            }
 
             // we now call only levelGhostParticlesOld.fill() and not .regrid()
             // regrid() would refine from next coarser in regions of level not overlaping
@@ -253,8 +263,14 @@ namespace amr
             // https://github.com/PHAREHUB/PHARE/issues/604 calling .fill() ensures that
             // levelGhostParticlesOld particles are filled exclusively from spliting next
             // coarser domain ones like when a new finest level is created.
-            lvlGhostPartOldRefiners_.fill(levelNumber, initDataTime);
-            copyLevelGhostOldToPushable_(*level, model);
+
+
+            if (levelNumber != rootLevelNumber)
+            {
+                lvlGhostPartOldRefiners_.fill(levelNumber, initDataTime);
+                copyLevelGhostOldToPushable_(*level, model);
+            }
+
 
             // computeIonMoments_(*level, model);
             // levelGhostNew will be refined in next firstStep

src/amr/messengers/hybrid_messenger.hpp

@@ -110,9 +110,10 @@ namespace amr
          * @brief see IMessenger::registerLevel
          */
         void registerLevel(std::shared_ptr<SAMRAI::hier::PatchHierarchy> const& hierarchy,
-                           int const levelNumber) override
+                           int const levelNumber,
+                           std::shared_ptr<SAMRAI::hier::PatchLevel> const& oldLevel) override
         {
-            strat_->registerLevel(hierarchy, levelNumber);
+            strat_->registerLevel(hierarchy, levelNumber, oldLevel);
         }
 
 

src/amr/messengers/hybrid_messenger_strategy.hpp

@@ -50,7 +50,9 @@ namespace amr
 
 
         virtual void registerLevel(std::shared_ptr<SAMRAI::hier::PatchHierarchy> const& hierarchy,
-                                   int const levelNumber)
+                                   int const levelNumber,
+                                   std::shared_ptr<SAMRAI::hier::PatchLevel> const& oldLevel
+                                   = nullptr)
             = 0;
 
         virtual void regrid(std::shared_ptr<SAMRAI::hier::PatchHierarchy> const& hierarchy,

src/amr/messengers/messenger.hpp

@@ -118,7 +118,9 @@ namespace amr
          * given hierarchy.
          */
         virtual void registerLevel(std::shared_ptr<SAMRAI::hier::PatchHierarchy> const& hierarchy,
-                                   int const levelNumber)
+                                   int const levelNumber,
+                                   std::shared_ptr<SAMRAI::hier::PatchLevel> const& oldLevel
+                                   = nullptr)
             = 0;
 
 

src/amr/messengers/mhd_hybrid_messenger_strategy.hpp

@@ -51,7 +51,9 @@ namespace amr
         }
 
         void registerLevel(std::shared_ptr<SAMRAI::hier::PatchHierarchy> const& /*hierarchy*/,
-                           int const /*levelNumber*/) override
+                           int const /*levelNumber*/,
+                           std::shared_ptr<SAMRAI::hier::PatchLevel> const& /*oldLevel*/
+                           ) override
         {
         }
 

src/amr/messengers/mhd_messenger.hpp

@@ -43,7 +43,8 @@ namespace amr
 
 
         void registerLevel(std::shared_ptr<SAMRAI::hier::PatchHierarchy> const& /*hierarchy*/,
-                           int const /*levelNumber*/) override
+                           int const /*levelNumber*/,
+                           std::shared_ptr<SAMRAI::hier::PatchLevel> const& /*oldLevel*/) override
         {
         }