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

Use reduction to compute min and max particle distances in NeighborParticles test. #3212

Merged
merged 9 commits into from
Jul 22, 2024
4 changes: 2 additions & 2 deletions Tests/Particles/NeighborParticles/Constants.H
Original file line number Diff line number Diff line change
Expand Up @@ -8,8 +8,8 @@ namespace Params
// This is designed to represent MFiX-like conditions where the grid spacing is
// roughly 2.5 times the particle diameter. In main.cpp we set grid spacing to 1
// so here we set cutoff to diameter = 1/2.5 --> cutoff = 0.2
static constexpr amrex::Real cutoff = 0.2 ;
static constexpr amrex::Real min_r = 1.e-4;
static constexpr amrex::ParticleReal cutoff = 0.2 ;
static constexpr amrex::ParticleReal min_r = 1.e-4;
}

#endif
59 changes: 30 additions & 29 deletions Tests/Particles/NeighborParticles/MDParticleContainer.cpp
Original file line number Diff line number Diff line change
@@ -1,8 +1,9 @@
#include <AMReX_Reduce.H>
#include <AMReX_SPACE.H>

#include "MDParticleContainer.H"
#include "Constants.H"

#include "CheckPair.H"
#include <AMReX_SPACE.H>

using namespace amrex;

Expand Down Expand Up @@ -148,8 +149,9 @@ std::pair<Real, Real> MDParticleContainer::minAndMaxDistance()
const int lev = 0;
auto& plev = GetParticles(lev);

Real min_d = std::numeric_limits<Real>::max();
Real max_d = std::numeric_limits<Real>::min();
ReduceOps<ReduceOpMin, ReduceOpMax> reduce_op;
ReduceData<ParticleReal, ParticleReal> reduce_data(reduce_op);
using ReduceTuple = typename decltype(reduce_data)::Type;

for(MFIter mfi = MakeMFIter(lev); mfi.isValid(); ++mfi)
{
Expand All @@ -159,41 +161,40 @@ std::pair<Real, Real> MDParticleContainer::minAndMaxDistance()

auto& ptile = plev[index];
auto& aos = ptile.GetArrayOfStructs();
const size_t np = aos.numParticles();

auto nbor_data = m_neighbor_list[lev][index].data();
ParticleType* pstruct = aos().dataPtr();

Gpu::DeviceScalar<Real> min_d_gpu(min_d);
Gpu::DeviceScalar<Real> max_d_gpu(max_d);

Real* pmin_d = min_d_gpu.dataPtr();
Real* pmax_d = max_d_gpu.dataPtr();

AMREX_FOR_1D ( np, i,
{
ParticleType& p1 = pstruct[i];
ParticleReal min_start = std::numeric_limits<ParticleReal>::max();
ParticleReal max_start = std::numeric_limits<ParticleReal>::lowest();

for (const auto& p2 : nbor_data.getNeighbors(i))
{
AMREX_D_TERM(Real dx = p1.pos(0) - p2.pos(0);,
Real dy = p1.pos(1) - p2.pos(1);,
Real dz = p1.pos(2) - p2.pos(2);)
reduce_op.eval(aos.numParticles(), reduce_data,
[=] AMREX_GPU_DEVICE (int i) -> ReduceTuple
{
ParticleType& p1 = pstruct[i];

Real r2 = AMREX_D_TERM(dx*dx, + dy*dy, + dz*dz);
r2 = amrex::max(r2, Params::min_r*Params::min_r);
Real r = sqrt(r2);
ParticleReal min_d = min_start;
ParticleReal max_d = max_start;

Gpu::Atomic::Min(pmin_d, r);
Gpu::Atomic::Max(pmax_d, r);
}
});
for (const auto& p2 : nbor_data.getNeighbors(i))
{
AMREX_D_TERM(ParticleReal dx = p1.pos(0) - p2.pos(0);,
ParticleReal dy = p1.pos(1) - p2.pos(1);,
ParticleReal dz = p1.pos(2) - p2.pos(2);)

Gpu::Device::streamSynchronize();
ParticleReal r2 = AMREX_D_TERM(dx*dx, + dy*dy, + dz*dz);
r2 = amrex::max(r2, Params::min_r*Params::min_r);
auto r = ParticleReal(std::sqrt(r2));

min_d = std::min(min_d, min_d_gpu.dataValue());
max_d = std::max(max_d, max_d_gpu.dataValue());
min_d = std::min(min_d, r);
max_d = std::max(max_d, r);
}
return {min_d, max_d};
});
}

ParticleReal min_d = amrex::get<0>(reduce_data.value(reduce_op));
ParticleReal max_d = amrex::get<1>(reduce_data.value(reduce_op));
ParallelDescriptor::ReduceRealMin(min_d, ParallelDescriptor::IOProcessorNumber());
ParallelDescriptor::ReduceRealMax(max_d, ParallelDescriptor::IOProcessorNumber());

Expand Down
Loading