Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
25 commits
Select commit Hold shift + click to select a range
230a1a7
Test runtime only particle container
AlexanderSinn Apr 22, 2025
00455fa
fix ab5
AlexanderSinn Apr 22, 2025
fe8c4ec
clean GetStructOfArrays
AlexanderSinn Apr 22, 2025
2e1f95c
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Apr 23, 2025
933de27
Merge branch 'development' into Test_runtime_only_particle_container
AlexanderSinn May 17, 2025
25bd871
fix spin data
AlexanderSinn May 17, 2025
51ceb00
merge dev
AlexanderSinn Jun 15, 2025
3c3e2c6
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Jul 26, 2025
6b1410b
fix define
AlexanderSinn Jul 26, 2025
12b2b4f
change name
AlexanderSinn Jul 26, 2025
b16a6a1
use amrex ReorderParticles
AlexanderSinn Jul 26, 2025
33bd8eb
Merge branch 'development' into Test_runtime_only_particle_container
AlexanderSinn Aug 13, 2025
b8b7ed3
Merge branch 'development' into Test_runtime_only_particle_container
AlexanderSinn Oct 17, 2025
e73b304
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Nov 10, 2025
6ddd246
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Nov 25, 2025
4455a2e
Merge branch 'development' into Test_runtime_only_particle_container
AlexanderSinn Dec 2, 2025
e83c0c9
Merge branch 'development' into Test_runtime_only_particle_container
AlexanderSinn Jan 6, 2026
5aac640
Merge branch 'development' into Test_runtime_only_particle_container
AlexanderSinn Mar 20, 2026
e395fc6
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Mar 27, 2026
a0a00f2
Use AMReX development
AlexanderSinn Apr 3, 2026
373e3d1
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Apr 3, 2026
af4059b
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Apr 14, 2026
ef95927
fix memory leak
AlexanderSinn Apr 20, 2026
42cabac
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn Apr 27, 2026
03fd5ba
Merge branch 'Hi-PACE:development' into Test_runtime_only_particle_co…
AlexanderSinn May 21, 2026
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 6 additions & 6 deletions src/diagnostics/OpenPMDWriter.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -381,7 +381,7 @@ OpenPMDWriter::CopyBeams (MultiBeam& beams, const amrex::Vector< std::string > b

if (np != 0) {
// copy data from GPU to IO buffer
auto& soa = beam.getBeamSlice(WhichBeamSlice::This).GetStructOfArrays();
auto& slice = beam.getBeamSlice(WhichBeamSlice::This);

for (std::size_t idx=0; idx<m_uint64_beam_data[ibeam].size(); idx++) {
const auto old_size = m_uint64_beam_data[ibeam][idx].size();
Expand All @@ -391,13 +391,13 @@ OpenPMDWriter::CopyBeams (MultiBeam& beams, const amrex::Vector< std::string > b
);
}
amrex::Gpu::copyAsync(amrex::Gpu::deviceToHost,
soa.GetIdCPUData().begin(),
soa.GetIdCPUData().begin() + np,
slice.GetIdCPUData().begin(),
slice.GetIdCPUData().begin() + np,
m_uint64_beam_data[ibeam][idx].data() + m_offset[ibeam]);
}

AMREX_ALWAYS_ASSERT_WITH_MESSAGE(
int(m_real_beam_data[ibeam].size()) == soa.NumRealComps(),
int(m_real_beam_data[ibeam].size()) == slice.NumRealComps(),
"List of real names in openPMD Writer class does not match the beam");

for (std::size_t idx=0; idx<m_real_beam_data[ibeam].size(); idx++) {
Expand All @@ -408,8 +408,8 @@ OpenPMDWriter::CopyBeams (MultiBeam& beams, const amrex::Vector< std::string > b
);
}
amrex::Gpu::copyAsync(amrex::Gpu::deviceToHost,
soa.GetRealData(idx).begin(),
soa.GetRealData(idx).begin() + np,
slice.GetRealData(idx).begin(),
slice.GetRealData(idx).begin() + np,
m_real_beam_data[ibeam][idx].data() + m_offset[ibeam]);
}
}
Expand Down
28 changes: 5 additions & 23 deletions src/particles/beam/BeamParticleContainer.H
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,8 @@ struct BeamIdx
w, // weight
ux, uy, uz, // momentum
real_nattribs_in_buffer,
real_nattribs=real_nattribs_in_buffer
real_nattribs=real_nattribs_in_buffer,
sx=real_nattribs, sy, sz
};
enum {
// no extra components stored in MultiBuffer, besides 64bit idcpu
Expand All @@ -45,26 +46,7 @@ struct WhichBeamSlice {
enum beam_slice : int { Next=0, This, N };
};

using BeamTile = amrex::ParticleTile<
amrex::SoAParticle<
BeamIdx::real_nattribs,
BeamIdx::int_nattribs
>,
BeamIdx::real_nattribs,
BeamIdx::int_nattribs,
amrex::PolymorphicArenaAllocator
>;

using BeamTileInit = amrex::ParticleTile<
amrex::SoAParticle<
BeamIdx::real_nattribs_in_buffer,
BeamIdx::int_nattribs_in_buffer
>,
BeamIdx::real_nattribs_in_buffer,
BeamIdx::int_nattribs_in_buffer,
// use PolymorphicArenaAllocator to either use Pinned or Device memory at runtime
amrex::PolymorphicArenaAllocator
>;
using BeamTile = amrex::ParticleTileRT<amrex::Real, int>;

/** \brief Container for particles of 1 beam species. */
class BeamParticleContainer
Expand Down Expand Up @@ -185,7 +167,7 @@ public:

void resize (int which_slice, int num_particles, int num_slipped_particles);

BeamTileInit& getBeamInitSlice () {
BeamTile& getBeamInitSlice () {
return m_init_slice;
}

Expand Down Expand Up @@ -215,7 +197,7 @@ private:
std::array<BeamTile, WhichBeamSlice::N> m_slices {};
std::array<int, WhichBeamSlice::N> m_num_particles_without_slipped {};
std::array<int, WhichBeamSlice::N> m_num_particles_with_slipped {};
BeamTileInit m_init_slice {};
BeamTile m_init_slice {};
BoxSorter m_init_sorter;
uint64_t m_total_num_particles = 0;
public:
Expand Down
96 changes: 29 additions & 67 deletions src/particles/beam/BeamParticleContainer.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -109,11 +109,22 @@ BeamParticleContainer::ReadParameters ()
queryWithParserAlt(pp, "spin_anom", m_spin_anom, pp_alt);
}

getBeamInitSlice().define(m_do_spin_tracking ? 3 : 0, 0, nullptr, nullptr,
m_initialize_on_cpu ? amrex::The_Pinned_Arena() : amrex::The_Arena());
getBeamInitSlice().define(
BeamIdx::real_nattribs_in_buffer + (m_do_spin_tracking ? 3 : 0),
BeamIdx::int_nattribs_in_buffer,
nullptr,
nullptr,
m_initialize_on_cpu ? amrex::The_Pinned_Arena() : amrex::The_Arena()
);

for (auto& beam_tile : m_slices) {
beam_tile.define(m_do_spin_tracking ? 3 : 0, 0, nullptr, nullptr, amrex::The_Arena());
beam_tile.define(
BeamIdx::real_nattribs + (m_do_spin_tracking ? 3 : 0),
BeamIdx::int_nattribs,
nullptr,
nullptr,
amrex::The_Arena()
);
}
}

Expand Down Expand Up @@ -271,7 +282,7 @@ BeamParticleContainer::InitData (const amrex::Geometry& geom)
m_total_num_particles = getBeamInitSlice().size();
if (Hipace::HeadRank()) {
m_init_sorter.sortParticlesByBox(
getBeamInitSlice().GetStructOfArrays().GetRealData(BeamIdx::z).dataPtr(),
getBeamInitSlice().GetRealData(BeamIdx::z).dataPtr(),
getBeamInitSlice().size(), m_initialize_on_cpu, geom);
}
#else
Expand Down Expand Up @@ -333,10 +344,10 @@ void BeamParticleContainer::TagByLevel (const int current_N_level,
{
HIPACE_PROFILE("BeamParticleContainer::TagByLevel()");

auto& soa = getBeamSlice(which_slice).GetStructOfArrays();
const amrex::Real * const pos_x = soa.GetRealData(BeamIdx::x).data();
const amrex::Real * const pos_y = soa.GetRealData(BeamIdx::y).data();
int * const p_mr_level = soa.GetIntData(BeamIdx::mr_level).data();
auto& slice = getBeamSlice(which_slice);
const amrex::Real * const pos_x = slice.GetRealData(BeamIdx::x).data();
const amrex::Real * const pos_y = slice.GetRealData(BeamIdx::y).data();
int * const p_mr_level = slice.GetIntData(BeamIdx::mr_level).data();

const int lev1_idx = std::min(1, current_N_level-1);
const int lev2_idx = std::min(2, current_N_level-1);
Expand Down Expand Up @@ -395,9 +406,9 @@ BeamParticleContainer::initializeSlice (int slice, int which_slice) {
ptd.rdata(BeamIdx::uy)[ip] = ptd_init.rdata(BeamIdx::uy)[idx_src];
ptd.rdata(BeamIdx::uz)[ip] = ptd_init.rdata(BeamIdx::uz)[idx_src];
if (do_spin_tracking) {
ptd.m_runtime_rdata[0][ip] = ptd_init.m_runtime_rdata[0][idx_src];
ptd.m_runtime_rdata[1][ip] = ptd_init.m_runtime_rdata[1][idx_src];
ptd.m_runtime_rdata[2][ip] = ptd_init.m_runtime_rdata[2][idx_src];
ptd.rdata(BeamIdx::sx)[ip] = ptd_init.rdata(BeamIdx::sx)[idx_src];
ptd.rdata(BeamIdx::sy)[ip] = ptd_init.rdata(BeamIdx::sy)[idx_src];
ptd.rdata(BeamIdx::sz)[ip] = ptd_init.rdata(BeamIdx::sz)[idx_src];
}
ptd.idcpu(ip) = ptd_init.idcpu(idx_src);
ptd.idata(BeamIdx::nsubcycles)[ip] = 0;
Expand All @@ -413,9 +424,9 @@ BeamParticleContainer::initializeSlice (int slice, int which_slice) {
const amrex::RealVect initial_spin_norm = m_initial_spin / m_initial_spin.vectorLength();
amrex::ParallelFor(getNumParticles(which_slice),
[=] AMREX_GPU_DEVICE (const int ip) {
ptd.m_runtime_rdata[0][ip] = initial_spin_norm[0];
ptd.m_runtime_rdata[1][ip] = initial_spin_norm[1];
ptd.m_runtime_rdata[2][ip] = initial_spin_norm[2];
ptd.rdata(BeamIdx::sx)[ip] = initial_spin_norm[0];
ptd.rdata(BeamIdx::sy)[ip] = initial_spin_norm[1];
ptd.rdata(BeamIdx::sz)[ip] = initial_spin_norm[2];
}
);
}
Expand All @@ -442,60 +453,11 @@ BeamParticleContainer::ReorderParticles (int beam_slice, int step, amrex::Geomet
HIPACE_PROFILE("BeamParticleContainer::ReorderParticles()");

const int np = getNumParticles(beam_slice);
const int np_total = getNumParticlesIncludingSlipped(beam_slice);
auto& ptile = getBeamSlice(beam_slice);
amrex::Gpu::DeviceVector<unsigned int> perm;
amrex::PermutationForDeposition<unsigned int>(perm, np, ptile, slice_geom.Domain(),
slice_geom, m_reorder_idx_type);
const unsigned int* permutations = perm.dataPtr();
auto& soa = ptile.GetStructOfArrays();

{
typename BeamTile::SoA::IdCPU tmp_idcpu;
tmp_idcpu.setArena(amrex::The_Arena());
tmp_idcpu.resize(np_total);
auto src = soa.GetIdCPUData().data();
uint64_t* dst = tmp_idcpu.data();
amrex::ParallelFor(np_total,
[=] AMREX_GPU_DEVICE (int i) {
dst[i] = i < np ? src[permutations[i]] : src[i];
});

amrex::Gpu::streamSynchronize();
soa.GetIdCPUData().swap(tmp_idcpu);
}

{ // Create a scope for the temporary vector below
BeamTile::RealVector tmp_real;
tmp_real.setArena(amrex::The_Arena());
tmp_real.resize(np_total);
for (int comp = 0; comp < soa.NumRealComps(); ++comp) {
auto src = soa.GetRealData(comp).data();
amrex::ParticleReal* dst = tmp_real.data();
amrex::ParallelFor(np_total,
[=] AMREX_GPU_DEVICE (int i) {
dst[i] = i < np ? src[permutations[i]] : src[i];
});

amrex::Gpu::streamSynchronize();
soa.GetRealData(comp).swap(tmp_real);
}
}

BeamTile::IntVector tmp_int;
tmp_int.setArena(amrex::The_Arena());
tmp_int.resize(np_total);
for (int comp = 0; comp < soa.NumIntComps(); ++comp) {
auto src = soa.GetIntData(comp).data();
int* dst = tmp_int.data();
amrex::ParallelFor(np_total,
[=] AMREX_GPU_DEVICE (int i) {
dst[i] = i < np ? src[permutations[i]] : src[i];
});

amrex::Gpu::streamSynchronize();
soa.GetIntData(comp).swap(tmp_int);
}
amrex::ReorderParticles(ptile, perm.dataPtr());
}
}

Expand Down Expand Up @@ -590,9 +552,9 @@ BeamParticleContainer::InSituComputeDiags (int islice)
{
const amrex::Real x = ptd.pos(0, ip);
const amrex::Real y = ptd.pos(1, ip);
const amrex::Real sx = ptd.m_runtime_rdata[0][ip];
const amrex::Real sy = ptd.m_runtime_rdata[1][ip];
const amrex::Real sz = ptd.m_runtime_rdata[2][ip];
const amrex::Real sx = ptd.rdata(BeamIdx::sx)[ip];
const amrex::Real sy = ptd.rdata(BeamIdx::sy)[ip];
const amrex::Real sz = ptd.rdata(BeamIdx::sz)[ip];
const amrex::Real w = ptd.rdata(BeamIdx::w)[ip];

if (!ptd.id(ip).is_valid() || x*x + y*y > insitu_radius_sq) {
Expand Down
8 changes: 4 additions & 4 deletions src/particles/beam/BeamParticleContainerInit.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,7 @@ namespace
*/
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
void AddOneBeamParticle (
const BeamTileInit::ParticleTileDataType& ptd,
const BeamTile::ParticleTileDataType& ptd,
const amrex::Real& x, const amrex::Real& y, const amrex::Real& z,
const amrex::Real& ux, const amrex::Real& uy, const amrex::Real& uz,
const amrex::Real& sx, const amrex::Real& sy, const amrex::Real& sz,
Expand All @@ -60,9 +60,9 @@ namespace
ptd.rdata(BeamIdx::uy)[ip] = uyp;
ptd.rdata(BeamIdx::uz)[ip] = uz;
if (do_spin) {
ptd.m_runtime_rdata[0][ip] = sx;
ptd.m_runtime_rdata[1][ip] = sy;
ptd.m_runtime_rdata[2][ip] = sz;
ptd.rdata(BeamIdx::sx)[ip] = sx;
ptd.rdata(BeamIdx::sy)[ip] = sy;
ptd.rdata(BeamIdx::sz)[ip] = sz;
}
ptd.rdata(BeamIdx::w )[ip] = std::abs(weight);

Expand Down
57 changes: 28 additions & 29 deletions src/particles/collisions/CoulombCollision.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -88,12 +88,12 @@ CoulombCollision::doPlasmaPlasmaCoulombCollision (
for (PlasmaParticleIterator pti(species1); pti.isValid(); ++pti) {

// Get particles SoA data
auto& soa1 = pti.GetStructOfArrays();
amrex::Real* const ux1 = soa1.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy1 = soa1.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi1 = soa1.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w1 = soa1.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev1 = soa1.GetIntData(PlasmaIdx::ion_lev).data();
auto& ptile1 = pti.GetParticleTile();
amrex::Real* const ux1 = ptile1.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy1 = ptile1.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi1 = ptile1.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w1 = ptile1.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev1 = ptile1.GetIntData(PlasmaIdx::ion_lev).data();
PlasmaBins::index_type * const indices1 = bins1.permutationPtr();
PlasmaBins::index_type const * const offsets1 = bins1.offsetsPtr();
amrex::Real q1 = species1.GetCharge();
Expand Down Expand Up @@ -154,12 +154,12 @@ CoulombCollision::doPlasmaPlasmaCoulombCollision (
for (PlasmaParticleIterator pti(species1); pti.isValid(); ++pti) {

// Get particles SoA data for species 1
auto& soa1 = pti.GetStructOfArrays();
amrex::Real* const ux1 = soa1.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy1 = soa1.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi1 = soa1.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w1 = soa1.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev1 = soa1.GetIntData(PlasmaIdx::ion_lev).data();
auto& ptile1 = pti.GetParticleTile();
amrex::Real* const ux1 = ptile1.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy1 = ptile1.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi1 = ptile1.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w1 = ptile1.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev1 = ptile1.GetIntData(PlasmaIdx::ion_lev).data();
PlasmaBins::index_type * const indices1 = bins1.permutationPtr();
PlasmaBins::index_type const * const offsets1 = bins1.offsetsPtr();
amrex::Real q1 = species1.GetCharge();
Expand All @@ -168,12 +168,11 @@ CoulombCollision::doPlasmaPlasmaCoulombCollision (

// Get particles SoA data for species 2
auto& ptile2 = species2.ParticlesAt(lev, pti.index(), pti.LocalTileIndex());
auto& soa2 = ptile2.GetStructOfArrays();
amrex::Real* const ux2 = soa2.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy2 = soa2.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi2= soa2.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w2 = soa2.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev2 = soa2.GetIntData(PlasmaIdx::ion_lev).data();
amrex::Real* const ux2 = ptile2.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy2 = ptile2.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi2= ptile2.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w2 = ptile2.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev2 = ptile2.GetIntData(PlasmaIdx::ion_lev).data();
PlasmaBins::index_type * const indices2 = bins2.permutationPtr();
PlasmaBins::index_type const * const offsets2 = bins2.offsetsPtr();
amrex::Real q2 = species2.GetCharge();
Expand Down Expand Up @@ -265,11 +264,11 @@ CoulombCollision::doBeamPlasmaCoulombCollision (
for (PlasmaParticleIterator pti(species2); pti.isValid(); ++pti) {

// // Get particles SoA data for species 1
auto& soa1 = species1.getBeamSlice(WhichBeamSlice::This).GetStructOfArrays();
amrex::Real* const ux1 = soa1.GetRealData(BeamIdx::ux).data();
amrex::Real* const uy1 = soa1.GetRealData(BeamIdx::uy).data();
amrex::Real* const psi1 = soa1.GetRealData(BeamIdx::uz).data();
const amrex::Real* const w1 = soa1.GetRealData(BeamIdx::w).data();
auto& ptile1 = species1.getBeamSlice(WhichBeamSlice::This);
amrex::Real* const ux1 = ptile1.GetRealData(BeamIdx::ux).data();
amrex::Real* const uy1 = ptile1.GetRealData(BeamIdx::uy).data();
amrex::Real* const psi1 = ptile1.GetRealData(BeamIdx::uz).data();
const amrex::Real* const w1 = ptile1.GetRealData(BeamIdx::w).data();
BeamBins::index_type * const indices1 = bins1.permutationPtr();
BeamBins::index_type const * const offsets1 = bins1.offsetsPtr();
amrex::Real q1 = species1.GetCharge();
Expand All @@ -278,12 +277,12 @@ CoulombCollision::doBeamPlasmaCoulombCollision (

// Get particles SoA data for species 2
//auto& ptile2 = species2.ParticlesAt(lev, pti.index(), pti.LocalTileIndex());
auto& soa2 = pti.GetStructOfArrays();
amrex::Real* const ux2 = soa2.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy2 = soa2.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi2= soa2.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w2 = soa2.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev2 = soa2.GetIntData(PlasmaIdx::ion_lev).data();
auto& ptile2 = pti.GetParticleTile();
amrex::Real* const ux2 = ptile2.GetRealData(PlasmaIdx::ux_half_step).data();
amrex::Real* const uy2 = ptile2.GetRealData(PlasmaIdx::uy_half_step).data();
amrex::Real* const psi2= ptile2.GetRealData(PlasmaIdx::psi_half_step).data();
const amrex::Real* const w2 = ptile2.GetRealData(PlasmaIdx::w).data();
const int* const ion_lev2 = ptile2.GetIntData(PlasmaIdx::ion_lev).data();
PlasmaBins::index_type * const indices2 = bins2.permutationPtr();
PlasmaBins::index_type const * const offsets2 = bins2.offsetsPtr();
amrex::Real q2 = species2.GetCharge();
Expand Down
Loading
Loading