Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Function to get particle cell index. #5118

Merged
merged 30 commits into from
Aug 21, 2024
Merged
Show file tree
Hide file tree
Changes from 12 commits
Commits
Show all changes
30 commits
Select commit Hold shift + click to select a range
a7d358f
Add particle cell id function
archermarx Aug 7, 2024
89216c6
Add header
archermarx Aug 7, 2024
a71fad0
Rewrite in header
archermarx Aug 7, 2024
1a8a045
Remove commented-out code
archermarx Aug 7, 2024
68459d8
Merge branch 'development' into particle_cell_id
archermarx Aug 7, 2024
ccde01d
better ignoring of unused variables
archermarx Aug 7, 2024
7b80c85
Const decls for clang-tidy
archermarx Aug 7, 2024
fcfbb45
Add overload for getParticleCellIndex using SuperParticleType
archermarx Aug 7, 2024
68f278a
Use ignore_unsued for index variables to appease clang-tidy
archermarx Aug 7, 2024
cd4cd15
Update Source/Utils/ParticleUtils.cpp
archermarx Aug 8, 2024
2eed9cd
Merge branch 'ECP-WarpX:development' into particle_cell_id
archermarx Aug 8, 2024
8d1d32f
Remove unused variable 'lo'
archermarx Aug 8, 2024
737d2cf
Update Source/Utils/ParticleUtils.H
archermarx Aug 12, 2024
002e8ce
Add use of function to ParticleReductionFunctor.cpp
archermarx Aug 12, 2024
de86103
Update Source/Utils/ParticleUtils.H
archermarx Aug 12, 2024
6c7b222
Update Source/Utils/ParticleUtils.H
archermarx Aug 12, 2024
bc28bde
Merge branch 'particle_cell_id' of github.com:archermarx/WarpX into p…
archermarx Aug 12, 2024
15b7670
Revert use of new function in ParticleUtils.cpp to get CI to pass
archermarx Aug 14, 2024
9bda9dd
Merge branch 'development' into particle_cell_id
archermarx Aug 14, 2024
b5059c4
Remove parameter
archermarx Aug 14, 2024
25111a4
fix typo
archermarx Aug 14, 2024
c177530
Use new overload of amrex::getParticleCell
archermarx Aug 14, 2024
970f8a6
remove particleutils.h header include in ParticleReductionFunctor
archermarx Aug 14, 2024
513bd14
Update ParticleUtils.H
archermarx Aug 20, 2024
8e13732
Merge branch 'ECP-WarpX:development' into particle_cell_id
archermarx Aug 20, 2024
055776b
Fix structured bindings
archermarx Aug 21, 2024
b248dc0
Remove new function entirely and rely on getParticleCell
archermarx Aug 21, 2024
4e8f688
Try using amrex::getParticleCell in ParticleUtils.cpp
archermarx Aug 21, 2024
fef22fe
Revert "Try using amrex::getParticleCell in ParticleUtils.cpp"
archermarx Aug 21, 2024
b20c5ce
remove unnecessary include
archermarx Aug 21, 2024
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
37 changes: 9 additions & 28 deletions Source/Diagnostics/ComputeDiagFunctors/TemperatureFunctor.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@
#include "Particles/MultiParticleContainer.H"
#include "Particles/WarpXParticleContainer.H"
#include "Utils/Parser/ParserUtils.H"
#include "Utils/ParticleUtils.H"
archermarx marked this conversation as resolved.
Show resolved Hide resolved
#include "WarpX.H"

#include <ablastr/coarsen/sample.H>
Expand Down Expand Up @@ -53,22 +54,10 @@ TemperatureFunctor::operator() (amrex::MultiFab& mf_dst, const int dcomp, const
amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dxi)
{
// Get position in AMReX convention to calculate corresponding index.
// Ideally this will be replaced with the AMReX NGP interpolator
// Always do x direction.
int ii = 0, jj = 0, kk = 0;
const amrex::ParticleReal x = p.pos(0);
const amrex::Real lx = (x - plo[0]) * dxi[0];
ii = static_cast<int>(amrex::Math::floor(lx));
#if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_3D) || defined(WARPX_DIM_RZ)
const amrex::ParticleReal y = p.pos(1);
const amrex::Real ly = (y - plo[1]) * dxi[1];
jj = static_cast<int>(amrex::Math::floor(ly));
#endif
#if defined(WARPX_DIM_3D)
const amrex::ParticleReal z = p.pos(2);
const amrex::Real lz = (z - plo[2]) * dxi[2];
kk = static_cast<int>(amrex::Math::floor(lz));
#endif
const auto cellIndex = ParticleUtils::getParticleCellIndex(p, plo, dxi, amrex::Dim3{0, 0, 0}).dim3();
const int ii = cellIndex.x;
const int jj = cellIndex.y;
const int kk = cellIndex.z;

const amrex::ParticleReal w = p.rdata(PIdx::w);
const amrex::ParticleReal ux = p.rdata(PIdx::ux);
Expand Down Expand Up @@ -119,18 +108,10 @@ TemperatureFunctor::operator() (amrex::MultiFab& mf_dst, const int dcomp, const
GetPosition.AsStored(ip, xp, yp, zp);

// Get position in AMReX convention to calculate corresponding index.
int ii = 0, jj = 0, kk = 0;
const amrex::Real lx = (xp - plo[0]) * dxi[0];
ii = static_cast<int>(amrex::Math::floor(lx));
#if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
const amrex::Real lz = (zp - plo[1]) * dxi[1];
jj = static_cast<int>(amrex::Math::floor(lz));
#elif defined(WARPX_DIM_3D)
const amrex::Real ly = (yp - plo[1]) * dxi[1];
jj = static_cast<int>(amrex::Math::floor(ly));
const amrex::Real lz = (zp - plo[2]) * dxi[2];
kk = static_cast<int>(amrex::Math::floor(lz));
#endif
const auto cellIndex = ParticleUtils::getParticleCellIndex(xp, yp, zp, plo, dxi, amrex::Dim3{0, 0, 0}).dim3();
const int ii = cellIndex.x;
const int jj = cellIndex.y;
const int kk = cellIndex.z;

const amrex::ParticleReal w = wp[ip];
const amrex::ParticleReal ux = uxp[ip] - out_array(ii, jj, kk, 1);
Expand Down
69 changes: 69 additions & 0 deletions Source/Utils/ParticleUtils.H
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,75 @@

namespace ParticleUtils {

/**
* \brief Determine which cell a particle lies in.
archermarx marked this conversation as resolved.
Show resolved Hide resolved
* Returns an amrex::IntVect with the i, j, k indices of the cell the particle lies within
*
* @param[in] x the particle's x position
* @param[in] y the particle's y position
* @param[in] z the particle's z position
* @param[in] plo the particle lower bounds
* @param[in] dxi the inverse cell size
* @param[in] lo
*/
AMREX_GPU_HOST_DEVICE AMREX_INLINE
amrex::IntVect getParticleCellIndex (amrex::ParticleReal xp,
[[maybe_unused]] amrex::ParticleReal yp,
[[maybe_unused]] amrex::ParticleReal zp,
amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& plo,
amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dxi,
amrex::Dim3 lo) {
int ii = 0, jj = 0, kk = 0;
amrex::ignore_unused(jj, kk);

const amrex::Real lx = (xp - plo[0]) * dxi[0] - lo.x;
ii = static_cast<int>(amrex::Math::floor(lx));
#if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_3D) || defined(WARPX_DIM_RZ)
const amrex::Real ly = (yp - plo[1]) * dxi[1] - lo.y;
jj = static_cast<int>(amrex::Math::floor(ly));
#endif
#if defined(WARPX_DIM_3D)
const amrex::Real lz = (zp - plo[2]) * dxi[2] - lo.z;
kk = static_cast<int>(amrex::Math::floor(lz));
#endif

return amrex::IntVect(AMREX_D_DECL(ii, jj, kk));
}

/**
* \brief Determine which cell a particle lies in.
archermarx marked this conversation as resolved.
Show resolved Hide resolved
* Returns an amrex::IntVect with the i, j, k indices of the cell the particle lies within
*
* @param[in] p the particle whose cell index we are finding
* @param[in] plo the particle lower bounds
* @param[in] dxi the inverse cell size
* @param[in] lo
*/
AMREX_GPU_HOST_DEVICE AMREX_INLINE
amrex::IntVect getParticleCellIndex (const WarpXParticleContainer::ParticleType& p,
amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& plo,
amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dxi,
amrex::Dim3 lo) {
return getParticleCellIndex(p.pos(0), p.pos(1), p.pos(2), plo, dxi, lo);
}

/**
* \brief Determine which cell a particle lies in.
archermarx marked this conversation as resolved.
Show resolved Hide resolved
* Returns an amrex::IntVect with the i, j, k indices of the cell the particle lies within
*
* @param[in] p the particle whose cell index we are finding
* @param[in] plo the particle lower bounds
* @param[in] dxi the inverse cell size
* @param[in] lo
*/
AMREX_GPU_HOST_DEVICE AMREX_INLINE
amrex::IntVect getParticleCellIndex (const WarpXParticleContainer::SuperParticleType& p,
amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& plo,
amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dxi,
amrex::Dim3 lo) {
return getParticleCellIndex(p.pos(0), p.pos(1), p.pos(2), plo, dxi, lo);
}

/**
* \brief Find the particles and count the particles that are in each cell. More specifically
* this function returns an amrex::DenseBins object containing an offset array and a permutation
Expand Down
8 changes: 2 additions & 6 deletions Source/Utils/ParticleUtils.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,6 @@ namespace ParticleUtils
// Extract box properties
Geometry const& geom = WarpX::GetInstance().Geom(lev);
Box const& cbx = mfi.tilebox(IntVect::TheZeroVector()); //Cell-centered box
const auto lo = lbound(cbx);
const auto dxi = geom.InvCellSizeArray();
const auto plo = geom.ProbLoArray();

Expand All @@ -57,12 +56,9 @@ namespace ParticleUtils
ParticleBins bins;
bins.build(np, ptd, cbx,
// Pass lambda function that returns the cell index
[=] AMREX_GPU_DEVICE (ParticleType const & p) noexcept -> amrex::IntVect
[=] AMREX_GPU_DEVICE (ParticleType const & p) noexcept -> IntVect
{
return IntVect{AMREX_D_DECL(
static_cast<int>((p.pos(0)-plo[0])*dxi[0] - lo.x),
static_cast<int>((p.pos(1)-plo[1])*dxi[1] - lo.y),
static_cast<int>((p.pos(2)-plo[2])*dxi[2] - lo.z))};
return getParticleCell(p, plo, dxi, cbx);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I see most (maybe all) of the examples that involve sorting the particles per cell are still failing. Looking closer, I see that the current development function here calculates the cell indices with

ii = (p.pos(0) - plo[0]) * dxi[0] - lo.x

whereas the AMReX function uses

ii = (p.pos(0) - plo[0]) * dxi[0] + cbx.smallend[0]

So the difference is that the development version calculates a cell index starting from 0 for each tile whereas the AMReX function calculates a global cell index for the entire domain. It's not clear to me why this breaks the code though since the binning should only care about which particles have the same cell index, not about what the value of the indices are.

});

return bins;
Expand Down
Loading