Update Particle Container to Pure SoA

Transition particle containers to pure SoA layouts.

Source/Diagnostics/ComputeDiagFunctors/BackTransformParticleFunctor.H

@@ -143,19 +143,19 @@ struct LorentzTransformParticles
         const amrex::ParticleReal uzp = uz_old_p * weight_old
                                       + uz_new_p * weight_new;
 #if defined (WARPX_DIM_3D)
-        dst.m_aos[i_dst].pos(0) = xp;
-        dst.m_aos[i_dst].pos(1) = yp;
-        dst.m_aos[i_dst].pos(2) = zp;
+        dst.m_rdata[PIdx::x][i_dst] = xp;
+        dst.m_rdata[PIdx::y][i_dst] = yp;
+        dst.m_rdata[PIdx::z][i_dst] = zp;
 #elif defined (WARPX_DIM_RZ)
-        dst.m_aos[i_dst].pos(0) = std::sqrt(xp*xp + yp*yp);
-        dst.m_aos[i_dst].pos(1) = zp;
+        dst.m_rdata[PIdx::x][i_dst] = std::sqrt(xp*xp + yp*yp);
+        dst.m_rdata[PIdx::z][i_dst] = zp;
         dst.m_rdata[PIdx::theta][i_dst] = std::atan2(yp, xp);
 #elif defined (WARPX_DIM_XZ)
-        dst.m_aos[i_dst].pos(0) = xp;
-        dst.m_aos[i_dst].pos(1) = zp;
+        dst.m_rdata[PIdx::x][i_dst] = xp;
+        dst.m_rdata[PIdx::z][i_dst] = zp;
         amrex::ignore_unused(yp);
 #elif defined (WARPX_DIM_1D)
-        dst.m_aos[i_dst].pos(0) = zp;
+        dst.m_rdata[PIdx::z][i_dst] = zp;
         amrex::ignore_unused(xp, yp);
 #else
         amrex::ignore_unused(xp, yp, zp);

Source/Diagnostics/FlushFormats/FlushFormatCheckpoint.cpp

@@ -200,8 +200,9 @@ FlushFormatCheckpoint::CheckpointParticles (
         auto runtime_inames = pc->getParticleRuntimeiComps();
         for (auto const& x : runtime_inames) { int_names[x.second+0] = x.first; }
 
-        pc->Checkpoint(dir, part_diag.getSpeciesName(), true,
-                       real_names, int_names);
+        // TODO
+        //pc->Checkpoint(dir, part_diag.getSpeciesName(), true,
+        //               real_names, int_names);
     }
 }
 

Source/Diagnostics/FlushFormats/FlushFormatPlotfile.cpp

@@ -406,10 +406,11 @@ FlushFormatPlotfile::WriteParticles(const std::string& dir,
         }
         // real_names contains a list of all particle attributes.
         // real_flags & int_flags are 1 or 0, whether quantity is dumped or not.
-        tmp.WritePlotFile(
-            dir, part_diag.getSpeciesName(),
-            real_flags, int_flags,
-            real_names, int_names);
+        // TODO
+        //tmp.WritePlotFile(
+        //    dir, part_diag.getSpeciesName(),
+        //    real_flags, int_flags,
+        //    real_names, int_names);
     }
 }
 

Source/Diagnostics/ParticleIO.cpp

@@ -206,10 +206,10 @@ MultiParticleContainer::Restart (const std::string& dir)
             }
         }
 
-        pc->Restart(dir, species_names.at(i));
+        //pc->Restart(dir, species_names.at(i)); // TODO
     }
     for (unsigned i = species_names.size(); i < species_names.size()+lasers_names.size(); ++i) {
-        allcontainers.at(i)->Restart(dir, lasers_names.at(i-species_names.size()));
+        //allcontainers.at(i)->Restart(dir, lasers_names.at(i-species_names.size())); // TODO
     }
 }
 

Source/Diagnostics/ReducedDiags/FieldProbeParticleContainer.H

@@ -24,7 +24,8 @@ struct FieldProbePIdx
 {
     enum
     {
-        Ex = 0, Ey, Ez,
+        x, y, z,
+        Ex, Ey, Ez,
         Bx, By, Bz,
         S, //!< the Poynting vector
         nattribs
@@ -37,16 +38,21 @@ struct FieldProbePIdx
  * nattribs tells the particle container to allot 7 SOA values.
  */
 class FieldProbeParticleContainer
-    : public amrex::ParticleContainer<0, 0, FieldProbePIdx::nattribs>
+    : public amrex::ParticleContainerPureSoA<FieldProbePIdx::nattribs, 0>
 {
 public:
+    static constexpr int NStructReal = 0;
+    static constexpr int NStructInt = 0;
+    static constexpr int NReal = FieldProbePIdx::nattribs;
+    static constexpr int NInt = 0;
+
     FieldProbeParticleContainer (amrex::AmrCore* amr_core);
     virtual ~FieldProbeParticleContainer() {}
 
     //! amrex iterator for our number of attributes
-    using iterator = amrex::ParIter<0, 0, FieldProbePIdx::nattribs, 0>;
+    using iterator = amrex::ParIterSoA<FieldProbePIdx::nattribs, 0>;
     //! amrex iterator for our number of attributes (read-only)
-    using const_iterator = amrex::ParConstIter<0, 0, FieldProbePIdx::nattribs, 0>;
+    using const_iterator = amrex::ParConstIterSoA<FieldProbePIdx::nattribs, 0>;
 
     //! similar to WarpXParticleContainer::AddNParticles but does not include u(x,y,z)
     void AddNParticles (int lev, amrex::Vector<amrex::ParticleReal> const & x, amrex::Vector<amrex::ParticleReal> const & y, amrex::Vector<amrex::ParticleReal> const & z);

Source/Diagnostics/ReducedDiags/FieldProbeParticleContainer.cpp

@@ -60,7 +60,7 @@
 using namespace amrex;
 
 FieldProbeParticleContainer::FieldProbeParticleContainer (AmrCore* amr_core)
-    : ParticleContainer<0, 0, FieldProbePIdx::nattribs>(amr_core->GetParGDB())
+    : ParticleContainerPureSoA<FieldProbePIdx::nattribs, 0>(amr_core->GetParGDB())
 {
     SetParticleSize();
 }
@@ -90,37 +90,23 @@ FieldProbeParticleContainer::AddNParticles (int lev,
      * is then coppied to the permament tile which is stored on the particle
      * (particle_tile).
      */
+    using PinnedTile = typename ContainerLike<amrex::PinnedArenaAllocator>::ParticleTileType;
 
-    using PinnedTile = ParticleTile<amrex::Particle<NStructReal, NStructInt>,
-                                    NArrayReal, NArrayInt,
-                                    amrex::PinnedArenaAllocator>;
     PinnedTile pinned_tile;
     pinned_tile.define(NumRuntimeRealComps(), NumRuntimeIntComps());
 
     for (int i = 0; i < np; i++)
     {
-        ParticleType p;
-        p.id() = ParticleType::NextID();
-        p.cpu() = ParallelDescriptor::MyProc();
-#if defined(WARPX_DIM_3D)
-        p.pos(0) = x[i];
-        p.pos(1) = y[i];
-        p.pos(2) = z[i];
-#elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
-        amrex::ignore_unused(y);
-        p.pos(0) = x[i];
-        p.pos(1) = z[i];
-#elif defined(WARPX_DIM_1D_Z)
-        amrex::ignore_unused(x, y);
-        p.pos(0) = z[i];
-#endif
-        // write position, cpu id, and particle id to particle
-        pinned_tile.push_back(p);
+        pinned_tile.push_back_int(0, ParticleType::NextID());
+        pinned_tile.push_back_int(1, amrex::ParallelDescriptor::MyProc());
     }
 
     // write Real attributes (SoA) to particle initialized zero
     DefineAndReturnParticleTile(0, 0, 0);
 
+    pinned_tile.push_back_real(FieldProbePIdx::x, x);
+    pinned_tile.push_back_real(FieldProbePIdx::y, y);
+    pinned_tile.push_back_real(FieldProbePIdx::z, z);
     pinned_tile.push_back_real(FieldProbePIdx::Ex, np, 0.0);
     pinned_tile.push_back_real(FieldProbePIdx::Ey, np, 0.0);
     pinned_tile.push_back_real(FieldProbePIdx::Ez, np, 0.0);

Source/Diagnostics/WarpXOpenPMD.H

@@ -62,7 +62,7 @@ class WarpXParticleCounter
 {
 public:
   using ParticleContainer = typename WarpXParticleContainer::ContainerLike<amrex::PinnedArenaAllocator>;
-  using ParticleIter = typename amrex::ParIter<0, 0, PIdx::nattribs, 0, amrex::PinnedArenaAllocator>;
+  using ParticleIter = typename amrex::ParIterSoA<PIdx::nattribs, 0, amrex::PinnedArenaAllocator>;
 
   WarpXParticleCounter (ParticleContainer* pc);
   unsigned long GetTotalNumParticles () {return m_Total;}
@@ -98,7 +98,7 @@ class WarpXOpenPMDPlot
 {
 public:
   using ParticleContainer = typename WarpXParticleContainer::ContainerLike<amrex::PinnedArenaAllocator>;
-  using ParticleIter = typename amrex::ParConstIter<0, 0, PIdx::nattribs, 0, amrex::PinnedArenaAllocator>;
+  using ParticleIter = typename amrex::ParConstIterSoA<PIdx::nattribs, 0, amrex::PinnedArenaAllocator>;
 
   /** Initialize openPMD I/O routines
    *

Source/Diagnostics/WarpXOpenPMD.cpp

@@ -939,31 +939,7 @@ WarpXOpenPMDPlot::SaveRealProperty (ParticleIter& pti,
 {
   auto const numParticleOnTile = pti.numParticles();
   uint64_t const numParticleOnTile64 = static_cast<uint64_t>( numParticleOnTile );
-  auto const& aos = pti.GetArrayOfStructs();  // size =  numParticlesOnTile
   auto const& soa = pti.GetStructOfArrays();
-  // first we concatinate the AoS into contiguous arrays
-  {
-    // note: WarpX does not yet use extra AoS Real attributes
-    for( auto idx=0; idx<ParticleIter::ContainerType::NStructReal; idx++ ) {  // lgtm [cpp/constant-comparison]
-      if( write_real_comp[idx] ) {
-          // handle scalar and non-scalar records by name
-          const auto [record_name, component_name] = detail::name2openPMD(real_comp_names[idx]);
-          auto currRecord = currSpecies[record_name];
-          auto currRecordComp = currRecord[component_name];
-
-          std::shared_ptr< amrex::ParticleReal > d(
-              new amrex::ParticleReal[numParticleOnTile],
-              [](amrex::ParticleReal const *p){ delete[] p; }
-          );
-
-          for( auto kk=0; kk<numParticleOnTile; kk++ )
-               d.get()[kk] = aos[kk].rdata(idx);
-
-          currRecordComp.storeChunk(d,
-               {offset}, {numParticleOnTile64});
-      }
-    }
-  }
 
   auto const getComponentRecord = [&currSpecies](std::string const comp_name) {
     // handle scalar and non-scalar records by name
@@ -975,9 +951,8 @@ WarpXOpenPMDPlot::SaveRealProperty (ParticleIter& pti,
   {
     auto const real_counter = std::min(write_real_comp.size(), real_comp_names.size());
     for (auto idx=0; idx<real_counter; idx++) {
-      auto ii = ParticleIter::ContainerType::NStructReal + idx;  // jump over extra AoS names
-      if (write_real_comp[ii]) {
-        getComponentRecord(real_comp_names[ii]).storeChunkRaw(
+      if (write_real_comp[idx]) {
+        getComponentRecord(real_comp_names[idx]).storeChunkRaw(
           soa.GetRealData(idx).data(), {offset}, {numParticleOnTile64});
       }
     }

Source/EmbeddedBoundary/ParticleBoundaryProcess.H

@@ -25,12 +25,12 @@ struct NoOp {
 struct Absorb {
     template <typename PData>
     AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
-    void operator() (const PData& ptd, int i,
+    void operator() (PData& ptd, int i,
                      const amrex::RealVect& /*pos*/, const amrex::RealVect& /*normal*/,
                      amrex::RandomEngine const& /*engine*/) const noexcept
     {
-        auto& p = ptd.m_aos[i];
-        p.id() = -p.id();
+        //ptd.id(i) = -ptd.id(i);  // FIXME: in AMReX - why does that not work?
+        ptd.m_idata[0][i] = -ptd.id(i);
     }
 };
 }

Source/EmbeddedBoundary/ParticleScraper.H

@@ -167,13 +167,12 @@ scrapeParticles (PC& pc, const amrex::Vector<const amrex::MultiFab*>& distance_t
             auto& tile = pti.GetParticleTile();
             auto ptd = tile.getParticleTileData();
             const auto np = tile.numParticles();
-            amrex::Particle<0,0> * const particles = tile.GetArrayOfStructs()().data();
             auto phi = (*distance_to_eb[lev])[pti].array();  // signed distance function
             amrex::ParallelForRNG( np,
             [=] AMREX_GPU_DEVICE (const int ip, amrex::RandomEngine const& engine) noexcept
             {
                 // skip particles that are already flagged for removal
-                if (particles[ip].id() < 0) return;
+                if (ptd.id(ip) < 0) return;
 
                 amrex::ParticleReal xp, yp, zp;
                 getPosition(ip, xp, yp, zp);

Source/Particles/Collision/BinaryCollision/ParticleCreationFunc.H

@@ -204,14 +204,21 @@ public:
 
                 // If the colliding particle weight decreases to zero, remove particle by
                 // setting its id to -1
-                constexpr amrex::Long minus_one_long = -1;
                 if (w1[p_pair_indices_1[i]] <= amrex::ParticleReal(0.))
                 {
-                    particle_ptr_1[p_pair_indices_1[i]].atomicSetID(minus_one_long);
+                    //particle_ptr_1[p_pair_indices_1[i]].atomicSetID(minus_one_long); // TODO
+                    // or?
+                    //soa_1.m_idata[PIdxInt::id] = -1;
+                    // or?
+                    particle_ptr_1->id() = -1;
                 }
                 if (w2[p_pair_indices_2[i]] <= amrex::ParticleReal(0.))
                 {
-                    particle_ptr_2[p_pair_indices_2[i]].atomicSetID(minus_one_long);
+                    //particle_ptr_2[p_pair_indices_2[i]].atomicSetID(minus_one_long); // TODO
+                    // or?
+                    //soa_2.m_idata[PIdxInt::id] = -1;
+                    // or?
+                    particle_ptr_2->id() = -1;
                 }
 
                 // Initialize the product particles' momentum, using a function depending on the

Source/Particles/ElementaryProcess/QEDPairGeneration.H

@@ -160,7 +160,7 @@ public:
             p_ux, p_uy, p_uz,
             engine);
 
-        src.m_aos[i_src].id() = -1; //destroy photon after pair generation
+        src.m_idata[0][i_src] = -1;  // destroy photon after pair generation
     }
 
 private:

Source/Particles/LaserParticleContainer.H

@@ -50,11 +50,10 @@ public:
      * \brief Method to initialize runtime attributes. Does nothing for LaserParticleContainer.
      */
     virtual void DefaultInitializeRuntimeAttributes (
-                        amrex::ParticleTile<amrex::Particle<NStructReal, NStructInt>,
-                                            NArrayReal, NArrayInt, amrex::PinnedArenaAllocator>& /*pinned_tile*/,
-                        const int /*n_external_attr_real*/,
-                        const int /*n_external_attr_int*/,
-                        const amrex::RandomEngine& /*engine*/) override final {}
+        typename ContainerLike<amrex::PinnedArenaAllocator>::ParticleTileType& /*pinned_tile*/,
+        const int /*n_external_attr_real*/,
+        const int /*n_external_attr_int*/,
+        const amrex::RandomEngine& /*engine*/) override final {}
 
     virtual void ReadHeader (std::istream& is) final;
 

Source/Particles/MultiParticleContainer.cpp

@@ -1400,6 +1400,7 @@ MultiParticleContainer::doQEDSchwinger ()
 
         const auto Transform = SchwingerTransformFunc{m_qed_schwinger_y_size, PIdx::w};
 
+
         const auto num_added = filterCreateTransformFromFAB<1>( dst_ele_tile,
                                dst_pos_tile, box, fieldsEB, np_ele_dst,
                                np_pos_dst,Filter, CreateEle, CreatePos,
@@ -1551,6 +1552,7 @@ void MultiParticleContainer::doQedBreitWheeler (int lev,
 
             const auto np_dst_ele = dst_ele_tile.numParticles();
             const auto np_dst_pos = dst_pos_tile.numParticles();
+
             const auto num_added = filterCopyTransformParticles<1>(
                                                       dst_ele_tile, dst_pos_tile,
                                                       src_tile, np_dst_ele, np_dst_pos,