Skip to content

Commit

Permalink
Add more comments and cleanup code
Browse files Browse the repository at this point in the history
  • Loading branch information
aprokop committed May 24, 2021
1 parent 71b9710 commit 63328cc
Show file tree
Hide file tree
Showing 3 changed files with 107 additions and 76 deletions.
3 changes: 2 additions & 1 deletion examples/dbscan/ArborX_DBSCANVerification.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -273,7 +273,8 @@ bool verifyDBSCAN(ExecutionSpace exec_space, Primitives const &primitives,

ArborX::BVH<MemorySpace> bvh(exec_space, primitives);

auto const predicates = buildPredicates(primitives, eps);
auto const predicates =
DBSCAN::PrimitivesWithRadius<Primitives>{primitives, eps};

Kokkos::View<int *, MemorySpace> indices("ArborX::DBSCAN::indices", 0);
Kokkos::View<int *, MemorySpace> offset("ArborX::DBSCAN::offset", 0);
Expand Down
144 changes: 83 additions & 61 deletions src/ArborX_DBSCAN.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -23,58 +23,12 @@
namespace ArborX
{

template <typename Primitives, typename FilterAndPermuter>
struct PrimitivesWithRadius
{
Primitives _primitives;
FilterAndPermuter _filter;
double _r;
};

struct Iota
{
size_t _n;
KOKKOS_FUNCTION unsigned int operator()(int const i) const { return i; }
KOKKOS_FUNCTION size_t extent(int) const { return _n; }
};

template <typename Primitives>
auto buildPredicates(Primitives const &v, double r)
{
using Access = AccessTraits<Primitives, PrimitivesTag>;

return PrimitivesWithRadius<Primitives, Iota>{
v, Iota{(size_t)Access::size(v)}, r};
}

template <typename Primitives, typename FilterAndPermuter>
auto buildPredicates(Primitives const &v, double r, FilterAndPermuter filter)
{
return PrimitivesWithRadius<Primitives, FilterAndPermuter>{v, filter, r};
}

template <typename Primitives, typename FilterAndPermuter>
struct AccessTraits<PrimitivesWithRadius<Primitives, FilterAndPermuter>,
PredicatesTag>
{
using PrimitivesAccess = AccessTraits<Primitives, PrimitivesTag>;

using memory_space = typename PrimitivesAccess::memory_space;
using Predicates = PrimitivesWithRadius<Primitives, FilterAndPermuter>;

static size_t size(Predicates const &w) { return w._filter.extent(0); }
static KOKKOS_FUNCTION auto get(Predicates const &w, size_t i)
{
int index = w._filter(i);
return attach(
intersects(Sphere{PrimitivesAccess::get(w._primitives, index), w._r}),
(int)index);
}
};

namespace DBSCAN
{

// All points are marked as if they were core points minpts = 2 case.
// Obviously, this is not true. However, in the algorithms it is used only for
// pairs of points within the distance eps, in which case it is correct.
struct CCSCorePoints
{
KOKKOS_FUNCTION bool operator()(int) const { return true; }
Expand Down Expand Up @@ -102,6 +56,7 @@ struct Parameters
{
// Print timers to standard output
bool _print_timers = false;
// Algorithm implementation (FDBSCAN or FDBSCAN-DenseBox)
Implementation _implementation = Implementation::FDBSCAN;

Parameters &setPrintTimers(bool print_timers)
Expand All @@ -116,6 +71,23 @@ struct Parameters
}
};

template <typename Primitives>
struct PrimitivesWithRadius
{
Primitives _primitives;
double _r;
};

template <typename Primitives, typename PermuteFilter>
struct PrimitivesWithRadiusReorderedAndFiltered
{
Primitives _primitives;
double _r;
PermuteFilter _filter;
};

// Mixed primitives consist of a set of boxes corresponding to dense cells,
// following by boxes corresponding to points in non-dense cells.
template <typename PointPrimitives, typename DenseCellOffsets,
typename CellIndices, typename Permutation>
struct MixedBoxPrimitives
Expand All @@ -130,6 +102,48 @@ struct MixedBoxPrimitives

} // namespace DBSCAN

template <typename Primitives>
struct AccessTraits<DBSCAN::PrimitivesWithRadius<Primitives>, PredicatesTag>
{
using PrimitivesAccess = AccessTraits<Primitives, PrimitivesTag>;

using memory_space = typename PrimitivesAccess::memory_space;
using Predicates = DBSCAN::PrimitivesWithRadius<Primitives>;

static size_t size(Predicates const &w)
{
return PrimitivesAccess::size(w._primitives);
}
static KOKKOS_FUNCTION auto get(Predicates const &w, size_t i)
{
return attach(
intersects(Sphere{PrimitivesAccess::get(w._primitives, i), w._r}),
(int)i);
}
};

template <typename Primitives, typename PermuteFilter>
struct AccessTraits<
DBSCAN::PrimitivesWithRadiusReorderedAndFiltered<Primitives, PermuteFilter>,
PredicatesTag>
{
using PrimitivesAccess = AccessTraits<Primitives, PrimitivesTag>;

using memory_space = typename PrimitivesAccess::memory_space;
using Predicates =
DBSCAN::PrimitivesWithRadiusReorderedAndFiltered<Primitives,
PermuteFilter>;

static size_t size(Predicates const &w) { return w._filter.extent(0); }
static KOKKOS_FUNCTION auto get(Predicates const &w, size_t i)
{
int index = w._filter(i);
return attach(
intersects(Sphere{PrimitivesAccess::get(w._primitives, index), w._r}),
(int)index);
}
};

template <typename PointPrimitives, typename MixedOffsets, typename CellIndices,
typename Permutation>
struct AccessTraits<DBSCAN::MixedBoxPrimitives<PointPrimitives, MixedOffsets,
Expand All @@ -155,11 +169,16 @@ struct AccessTraits<DBSCAN::MixedBoxPrimitives<PointPrimitives, MixedOffsets,
auto num_dense_primitives = dco.size() - 1;
if (i < num_dense_primitives)
{
// For a primitive corresponding to a dense cell, use that cell's box. It
// may not be tight around the points inside, but is cheap to compute.
auto cell_index = w._sorted_cell_indices(dco(i));
return w._grid.cellBox(cell_index);
}
else
{
// For a primitive corresponding to a point in a non-dense cell, use that
// point. But first, figure out its index, which requires some
// computations.
using Access = AccessTraits<PointPrimitives, PrimitivesTag>;

i = (i - num_dense_primitives) + w._num_points_in_dense_cells;
Expand Down Expand Up @@ -214,8 +233,6 @@ dbscan(ExecutionSpace const &exec_space, Primitives const &primitives,

if (parameters._implementation == DBSCAN::Implementation::FDBSCAN)
{
auto const predicates = buildPredicates(primitives, eps);

// Build the tree
timer_start(timer);
Kokkos::Profiling::pushRegion("ArborX::dbscan::tree_construction");
Expand All @@ -225,6 +242,8 @@ dbscan(ExecutionSpace const &exec_space, Primitives const &primitives,

timer_start(timer);
Kokkos::Profiling::pushRegion("ArborX::dbscan::clusters");
auto const predicates =
DBSCAN::PrimitivesWithRadius<Primitives>{primitives, eps};
if (is_special_case)
{
// Perform the queries and build clusters through callback
Expand Down Expand Up @@ -264,8 +283,6 @@ dbscan(ExecutionSpace const &exec_space, Primitives const &primitives,
else if (parameters._implementation ==
DBSCAN::Implementation::FDBSCAN_DenseBox)
{
auto const predicates = buildPredicates(primitives, eps);

// Find dense boxes
timer_start(timer);
Kokkos::Profiling::pushRegion("ArborX::dbscan::dense_cells");
Expand Down Expand Up @@ -313,9 +330,7 @@ dbscan(ExecutionSpace const &exec_space, Primitives const &primitives,
(100.f * num_dense_cells) / num_nonempty_cells);
printf("#dense cell points : %10d [%.2f%%]\n", num_points_in_dense_cells,
(100.f * num_points_in_dense_cells) / n);
printf("#dense primitives : %10d\n", num_dense_cells);
printf("#point primitives : %10d\n", num_points_in_sparse_cells);
printf("#total primitives : %10d\n",
printf("#mixed primitives : %10d\n",
num_dense_cells + num_points_in_sparse_cells);
}

Expand Down Expand Up @@ -346,6 +361,8 @@ dbscan(ExecutionSpace const &exec_space, Primitives const &primitives,
// Perform the queries and build clusters through callback
using CorePoints = DBSCAN::CCSCorePoints;
Kokkos::Profiling::pushRegion("ArborX::dbscan::clusters::query");
auto const predicates =
DBSCAN::PrimitivesWithRadius<Primitives>{primitives, eps};
bvh.query(
exec_space, predicates,
Details::FDBSCANDenseBoxCallback<MemorySpace, CorePoints, Primitives,
Expand All @@ -371,16 +388,19 @@ dbscan(ExecutionSpace const &exec_space, Primitives const &primitives,
num_points_in_dense_cells),
KOKKOS_LAMBDA(int i) { num_neigh(permute(i)) = INT_MAX; });
// Count neighbors for points in sparse cells
auto const sparse_predicates = buildPredicates(
primitives, eps,
Kokkos::subview(permute,
Kokkos::make_pair(num_points_in_dense_cells, n)));
auto sparse_permute = Kokkos::subview(
permute, Kokkos::make_pair(num_points_in_dense_cells, n));

auto const sparse_predicates =
DBSCAN::PrimitivesWithRadiusReorderedAndFiltered<
Primitives, decltype(sparse_permute)>{primitives, eps,
sparse_permute};
bvh.query(exec_space, sparse_predicates,
Details::CountUpToN_DenseBox<MemorySpace, Primitives,
decltype(dense_cell_offsets),
decltype(permute)>(
num_neigh, primitives, dense_cell_offsets, permute,
core_min_size, eps));
core_min_size, eps, core_min_size));
Kokkos::Profiling::popRegion();
elapsed["neigh"] = timer_seconds(timer_local);

Expand All @@ -389,6 +409,8 @@ dbscan(ExecutionSpace const &exec_space, Primitives const &primitives,
// Perform the queries and build clusters through callback
timer_start(timer_local);
Kokkos::Profiling::pushRegion("ArborX::dbscan::clusters:query");
auto const predicates =
DBSCAN::PrimitivesWithRadius<Primitives>{primitives, eps};
bvh.query(
exec_space, predicates,
Details::FDBSCANDenseBoxCallback<MemorySpace, CorePoints, Primitives,
Expand Down
36 changes: 22 additions & 14 deletions src/details/ArborX_DetailsFDBSCANDenseBox.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,9 @@ struct CartesianGrid
_ny = std::ceil((max_corner[1] - min_corner[1]) / h);
_nz = std::ceil((max_corner[2] - min_corner[2]) / h);

// Catch overflow in grid cell indices
// Catch potential overflow in grid cell indices early. This is a
// conservative check as an actual overflow may not occur, depending on
// which cells are filled.
size_t constexpr max_size_t = std::numeric_limits<size_t>::max();
ARBORX_ASSERT(_nx == 0 || _ny == 0 || _nz == 0 ||
(_ny < max_size_t / _nx && _nz < max_size_t / (_nx * _ny)));
Expand Down Expand Up @@ -87,19 +89,21 @@ struct CountUpToN_DenseBox
Permutation _permute;
int core_min_size;
float eps;
int _n;

CountUpToN_DenseBox(Kokkos::View<int *, MemorySpace> const &counts,
Primitives const &primitives,
DenseCellOffsets const &dense_cell_offsets,
Permutation const &permute, int core_min_size_in,
float eps_in)
float eps_in, int n)
: _counts(counts)
, _primitives(primitives)
, _dense_cell_offsets(dense_cell_offsets)
, _num_dense_cells(dense_cell_offsets.size() - 1)
, _permute(permute)
, core_min_size(core_min_size_in)
, eps(eps_in)
, _n(n)
{
}

Expand All @@ -125,15 +129,15 @@ struct CountUpToN_DenseBox
if (distance(query_point, Access::get(_primitives, j)) <= eps)
{
Kokkos::atomic_fetch_add(&count, 1);
if (count >= core_min_size)
if (count >= _n)
return ArborX::CallbackTreeTraversalControl::early_exit;
}
}
}
else
{
Kokkos::atomic_fetch_add(&count, 1);
if (count >= core_min_size)
if (count >= _n)
return ArborX::CallbackTreeTraversalControl::early_exit;
}

Expand Down Expand Up @@ -199,10 +203,9 @@ struct FDBSCANDenseBoxCallback
{
int j = _permute(jj);

// As soon as a pair is found, all other potentially search pairs for
// these two boxes will stop too
// TODO: can this be optimized out, i.e., not reload
// _union_find.representative(i)
// As soon as a pair is found, stop the search. If it is a case of
// merging two dense cells, this will stop all other threads from
// processing the same merge.
if (_union_find.representative(i) == _union_find.representative(j))
break;

Expand Down Expand Up @@ -355,21 +358,26 @@ int reorderDenseAndSparseCells(ExecutionSpace const &exec_space,
template <typename ExecutionSpace, typename CellIndices, typename Permutation,
typename Labels>
void unionFindWithinEachDenseCell(ExecutionSpace const &exec_space,
CellIndices sorted_cell_indices,
CellIndices sorted_dense_cell_indices,
Permutation permute, Labels labels)
{
using MemorySpace = typename Permutation::memory_space;

UnionFind<MemorySpace> union_find{labels};

auto const n = sorted_cell_indices.size();
// The algorithm relies on the fact that the cell indices array only contains
// dense cells. Thus, as long as two cell indices are the same, a) they
// belong to the same cell, and b) that cell is dense, thus they should be in
// the same cluster. If, on the other hand, the array also contained
// non-dense cells, that would not have been possible, as an additional
// computations would have to be done to figure out if the points belong to a
// dense cell, which would have required a linear scan.
auto const n = sorted_dense_cell_indices.size();
Kokkos::parallel_for("ArborX::dbscan::union_find_within_each_dense_box",
Kokkos::RangePolicy<ExecutionSpace>(exec_space, 1, n),
KOKKOS_LAMBDA(int i) {
// Each thread union-merges a pair if they belong to
// the same cell
if (sorted_cell_indices(i) ==
sorted_cell_indices(i - 1))
if (sorted_dense_cell_indices(i) ==
sorted_dense_cell_indices(i - 1))
union_find.merge(permute(i), permute(i - 1));
});
}
Expand Down

0 comments on commit 63328cc

Please sign in to comment.