Simplify implementation and fix bugs

EZoni · EZoni · commit d16803f92ba2 · 2025-10-24T16:05:24.000-07:00
diff --git a/Source/Evolve/WarpXEvolve.cpp b/Source/Evolve/WarpXEvolve.cpp
@@ -504,29 +504,48 @@ WarpX::OneStep_nosub (
     ExecutePythonCallback("particlescraper");
     ExecutePythonCallback("beforedeposition");
 
-    // push particles (half position and full momentum)
-    PushParticlesandDeposit(
-        a_cur_time,
-        /*skip_current=*/true,
-        PositionPushType::FirstHalf,
-        MomentumPushType::Full
-    );
+    // with collisions placed in the middle of the position push and after the momentum push
+    if (m_collisions_split_position_push) {
+        // push particles (half position and full momentum)
+        PushParticlesandDeposit(
+            a_cur_time,
+            /*skip_current=*/true,
+            PositionPushType::FirstHalf,
+            MomentumPushType::Full
+        );
 
-    // communicate particle data
-    mypc->Redistribute();
+        // communicate particle data
+        // FIXME Copy local communication from WarpX::HandleParticlesAtBoundaries
+        mypc->Redistribute();
 
-    // perform particle collisions
-    ExecutePythonCallback("beforecollisions");
-    mypc->doCollisions(a_step, a_cur_time, a_dt);
-    ExecutePythonCallback("aftercollisions");
+        // perform particle collisions
+        ExecutePythonCallback("beforecollisions");
+        mypc->doCollisions(a_step, a_cur_time, a_dt);
+        ExecutePythonCallback("aftercollisions");
 
-    // push particles (half position)
-    PushParticlesandDeposit(
-        a_cur_time,
-        /*skip_current=*/false,
-        PositionPushType::SecondHalf,
-        MomentumPushType::None
-    );
+        // push particles (half position)
+        PushParticlesandDeposit(
+            a_cur_time,
+            /*skip_current=*/false,
+            PositionPushType::SecondHalf,
+            MomentumPushType::None
+        );
+    }
+    // with collisions placed before the position and momentum push, or without collisions
+    else {
+        // perform particle collisions
+        ExecutePythonCallback("beforecollisions");
+        mypc->doCollisions(a_step, a_cur_time, a_dt);
+        ExecutePythonCallback("aftercollisions");
+
+        // push particles (half position)
+        PushParticlesandDeposit(
+            a_cur_time,
+            /*skip_current=*/false,
+            PositionPushType::Full,
+            MomentumPushType::Full
+        );
+    }
 
     ExecutePythonCallback("afterdeposition");
 
diff --git a/Source/Particles/MultiParticleContainer.H b/Source/Particles/MultiParticleContainer.H
@@ -69,10 +69,7 @@ class MultiParticleContainer
 
 public:
 
-    MultiParticleContainer (
-        amrex::AmrCore* amr_core,
-        bool const collisions_split_position_push
-    );
+    MultiParticleContainer (amrex::AmrCore* amr_core);
 
     ~MultiParticleContainer() = default;
 
diff --git a/Source/Particles/MultiParticleContainer.cpp b/Source/Particles/MultiParticleContainer.cpp
@@ -92,11 +92,9 @@ namespace
     };
 }
 
-MultiParticleContainer::MultiParticleContainer (
-    AmrCore* amr_core,
-    bool const collisions_split_position_push
-)
+MultiParticleContainer::MultiParticleContainer (AmrCore* amr_core)
 {
+
     ReadParameters();
 
     auto const nspecies = static_cast<int>(species_names.size());
@@ -105,28 +103,13 @@ MultiParticleContainer::MultiParticleContainer (
     allcontainers.resize(nspecies + nlasers);
     for (int i = 0; i < nspecies; ++i) {
         if (species_types[i] == PCTypes::Physical) {
-            allcontainers[i] = std::make_unique<PhysicalParticleContainer>(
-                amr_core,
-                /*ispecies=*/i,
-                /*name=*/species_names[i],
-                collisions_split_position_push
-            );
+            allcontainers[i] = std::make_unique<PhysicalParticleContainer>(amr_core, i, species_names[i]);
         }
         else if (species_types[i] == PCTypes::RigidInjected) {
-            allcontainers[i] = std::make_unique<RigidInjectedParticleContainer>(
-                amr_core,
-                /*ispecies=*/i,
-                /*name=*/species_names[i],
-                collisions_split_position_push
-            );
+            allcontainers[i] = std::make_unique<RigidInjectedParticleContainer>(amr_core, i, species_names[i]);
         }
         else if (species_types[i] == PCTypes::Photon) {
-            allcontainers[i] = std::make_unique<PhotonParticleContainer>(
-                amr_core,
-                /*ispecies=*/i,
-                /*name=*/species_names[i],
-                collisions_split_position_push
-            );
+            allcontainers[i] = std::make_unique<PhotonParticleContainer>(amr_core, i, species_names[i]);
         }
         allcontainers[i]->m_deposit_on_main_grid = m_deposit_on_main_grid[i];
         allcontainers[i]->m_gather_from_main_grid = m_gather_from_main_grid[i];
diff --git a/Source/Particles/PhotonParticleContainer.H b/Source/Particles/PhotonParticleContainer.H
@@ -33,13 +33,9 @@ class PhotonParticleContainer
     : public PhysicalParticleContainer
 {
 public:
-    PhotonParticleContainer (
-        amrex::AmrCore* amr_core,
-        int ispecies,
-        const std::string& name,
-        bool const collisions_split_position_push
-    );
-
+    PhotonParticleContainer (amrex::AmrCore* amr_core,
+                                    int ispecies,
+                                    const std::string& name);
     ~PhotonParticleContainer () override = default;
 
     PhotonParticleContainer ( PhotonParticleContainer const &)             = delete;
diff --git a/Source/Particles/PhotonParticleContainer.cpp b/Source/Particles/PhotonParticleContainer.cpp
@@ -43,12 +43,9 @@
 
 using namespace amrex;
 
-PhotonParticleContainer::PhotonParticleContainer (
-    AmrCore* amr_core,
-    int ispecies,
-    const std::string& name,
-    bool const collisions_split_position_push
-) : PhysicalParticleContainer(amr_core, ispecies, name, collisions_split_position_push)
+PhotonParticleContainer::PhotonParticleContainer (AmrCore* amr_core, int ispecies,
+                                                  const std::string& name)
+    : PhysicalParticleContainer(amr_core, ispecies, name)
 {
     const ParmParse pp_species_name(species_name);
 
diff --git a/Source/Particles/PhysicalParticleContainer.H b/Source/Particles/PhysicalParticleContainer.H
@@ -49,12 +49,9 @@ class PhysicalParticleContainer
 {
 public:
 
-    PhysicalParticleContainer (
-        amrex::AmrCore* amr_core,
-        int ispecies,
-        const std::string& name,
-        bool const collisions_split_position_push
-    );
+    PhysicalParticleContainer (amrex::AmrCore* amr_core,
+                               int ispecies,
+                               const std::string& name);
 
     PhysicalParticleContainer (amrex::AmrCore* amr_core);
 
@@ -498,13 +495,6 @@ protected:
     /* Vector of user-defined parser for initializing user-defined real attributes */
     amrex::Vector< std::unique_ptr<amrex::Parser> > m_user_real_attrib_parser;
 
-    /**
-     * \brief Local copy of the WarpX class attribute that controls whether collisions
-     *        are placed in the middle of the position push and after the momentum push,
-     *        or before both the position and momentum push. The default value
-     *        is set in WarpX::ReadParameters().
-     */
-    bool m_collisions_split_position_push;
 };
 
 #endif
diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp
@@ -116,14 +116,10 @@
 
 using namespace amrex;
 
-PhysicalParticleContainer::PhysicalParticleContainer (
-    AmrCore* amr_core,
-    int ispecies,
-    const std::string& name,
-    bool const collisions_split_position_push
-) : WarpXParticleContainer(amr_core, ispecies),
-    species_name(name),
-    m_collisions_split_position_push(collisions_split_position_push)
+PhysicalParticleContainer::PhysicalParticleContainer (AmrCore* amr_core, int ispecies,
+                                                      const std::string& name)
+    : WarpXParticleContainer(amr_core, ispecies),
+      species_name(name)
 {
     BackwardCompatibility();
 
@@ -341,13 +337,11 @@ PhysicalParticleContainer::PhysicalParticleContainer (
         m_boundary_conditions.SetThermalVelocity(boundary_uth);
     }
 
-    // If the position push is split to perform collisions, add the
-    // average momentum components to deposit the current correctly
-    if (m_collisions_split_position_push) {
-        this->AddRealComp("ux_avg", /*communicate=*/0);
-        this->AddRealComp("uy_avg", /*communicate=*/0);
-        this->AddRealComp("uz_avg", /*communicate=*/0);
-    }
+    // Add the average momentum components to deposit the current correctly
+    // if the position push is split to perform collisions
+    this->AddRealComp("ux_avg", /*communicate=*/0);
+    this->AddRealComp("uy_avg", /*communicate=*/0);
+    this->AddRealComp("uz_avg", /*communicate=*/0);
 }
 
 void
@@ -477,7 +471,6 @@ PhysicalParticleContainer::Evolve (ablastr::fields::MultiFabRegister& fields,
         !do_not_deposit &&
         (position_push_type == PositionPushType::Full || position_push_type == PositionPushType::FirstHalf)
     );
-    // FIXME Check the logic here
     bool const deposit_current = (
         !skip_deposition &&
         !do_not_deposit &&
@@ -516,9 +509,9 @@ PhysicalParticleContainer::Evolve (ablastr::fields::MultiFabRegister& fields,
             // Extract particle data
             auto& attribs = pti.GetAttribs();
             auto&  wp = attribs[PIdx::w];
-            auto& uxp = (m_collisions_split_position_push) ? pti.GetAttribs("ux_avg") : attribs[PIdx::ux];
-            auto& uyp = (m_collisions_split_position_push) ? pti.GetAttribs("uy_avg") : attribs[PIdx::uy];
-            auto& uzp = (m_collisions_split_position_push) ? pti.GetAttribs("uz_avg") : attribs[PIdx::uz];
+            auto& uxp = pti.GetAttribs("ux_avg");
+            auto& uyp = pti.GetAttribs("uy_avg");
+            auto& uzp = pti.GetAttribs("uz_avg");
 
             const long np = pti.numParticles();
 
@@ -1347,9 +1340,9 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
     ParticleReal* const AMREX_RESTRICT uz = attribs[PIdx::uz].dataPtr() + offset;
 
     // Average momentum
-    amrex::ParticleReal* const AMREX_RESTRICT ux_avg = pti.GetAttribs("ux_avg").dataPtr() + offset;
-    amrex::ParticleReal* const AMREX_RESTRICT uy_avg = pti.GetAttribs("uy_avg").dataPtr() + offset;
-    amrex::ParticleReal* const AMREX_RESTRICT uz_avg = pti.GetAttribs("uz_avg").dataPtr() + offset;
+    amrex::ParticleReal* AMREX_RESTRICT ux_avg = pti.GetAttribs("ux_avg").dataPtr() + offset;
+    amrex::ParticleReal* AMREX_RESTRICT uy_avg = pti.GetAttribs("uy_avg").dataPtr() + offset;
+    amrex::ParticleReal* AMREX_RESTRICT uz_avg = pti.GetAttribs("uz_avg").dataPtr() + offset;
 
     CopyParticleAttribs copyAttribs;
     if (copy_particle_attribs) {
@@ -1491,6 +1484,7 @@ PhysicalParticleContainer::PushPX (WarpXParIter& pti,
         }
 #endif
         // Update average momentum
+        // FIXME Improve the logic here
         if (momentum_push_type != MomentumPushType::None) {
             ux_avg[ip] = ux[ip];
             uy_avg[ip] = uy[ip];
diff --git a/Source/Particles/RigidInjectedParticleContainer.H b/Source/Particles/RigidInjectedParticleContainer.H
@@ -48,13 +48,9 @@ class RigidInjectedParticleContainer
     : public PhysicalParticleContainer
 {
 public:
-    RigidInjectedParticleContainer (
-        amrex::AmrCore* amr_core,
-        int ispecies,
-        const std::string& name,
-        bool const collisions_split_position_push
-    );
-
+    RigidInjectedParticleContainer (amrex::AmrCore* amr_core,
+                                    int ispecies,
+                                    const std::string& name);
     ~RigidInjectedParticleContainer () override = default;
 
     RigidInjectedParticleContainer ( RigidInjectedParticleContainer const &)             = delete;
diff --git a/Source/Particles/RigidInjectedParticleContainer.cpp b/Source/Particles/RigidInjectedParticleContainer.cpp
@@ -56,12 +56,9 @@
 
 using namespace amrex;
 
-RigidInjectedParticleContainer::RigidInjectedParticleContainer (
-    AmrCore* amr_core,
-    int ispecies,
-    const std::string& name,
-    bool const collisions_split_position_push
-) : PhysicalParticleContainer(amr_core, ispecies, name, collisions_split_position_push)
+RigidInjectedParticleContainer::RigidInjectedParticleContainer (AmrCore* amr_core, int ispecies,
+                                                                const std::string& name)
+    : PhysicalParticleContainer(amr_core, ispecies, name)
 {
 
 #if defined(WARPX_DIM_RCYLINDER) || defined(WARPX_DIM_RSPHERE)
diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp
@@ -346,7 +346,7 @@ WarpX::WarpX ()
     t_old.resize(nlevs_max, std::numeric_limits<Real>::lowest());
     dt.resize(nlevs_max, std::numeric_limits<Real>::max());
 
-    mypc = std::make_unique<MultiParticleContainer>(this, m_collisions_split_position_push);
+    mypc = std::make_unique<MultiParticleContainer>(this);
 
     // Loop over species (particles and lasers)
     // and set current injection position per species
