The neighbour-searching functionality within Aboria uses the underlying spatial data structure (i.e. cell list, octtree, kdtree) to perform its task. You can directly access this underlying data structure as well, in order for you to write your own spatial algorithms.
Aboria uses iterators to provide a generic way to interact with the spatial
data structure, so that the code that you write is independent of the particlular
data structure that you use. An iterator is a generic C++ object that "iterates"
through a 1D container. For example, the classic iterator is the loop index
i
in the code below.
size_t N = 10; std::vector<double> v(N); for (size_t i = 0; i < N; ++i) { v[i] = i; }
The STL library abstracts away the particular type of i
,
and defines a set of iterators for each container. For example, the std::vector
class has its own iterator, which you can use as follows.
size_t index = 0; for (std::vector<double>::iterator i = v.begin(); i != v.end(); ++i) { *i = index++; }
Or in more compact notation
index = 0; for (auto i = v.begin(); i != v.end(); ++i) { *i = index++; }
The iterators in Aboria are similar to STL iterators in that they can step
through a given range of objects using the ++
operator, with some slight differences that we will describe below.
Before we can start using the data structure iterators in Aboria, we need a particle set. Lets create a random set of points in 2D:
N = 100; typedef Particles<std::tuple<>, 2> Particles_t; typedef typename Particles_t::position position; Particles_t particles(N); std::uniform_real_distribution<double> uniform(0, 1); for (size_t i = 0; i < N; ++i) { auto &gen = get<generator>(particles)[i]; get<position>(particles)[i] = vdouble2(uniform(gen), uniform(gen)); } particles.init_neighbour_search(vdouble2(0, 0), vdouble2(1, 1), vdouble2(false, false));
In order to start interacting with the spatial data structures, we need to
get its query object from the particle set. A query object is a lightweight
object that has all the information neccessary to access either the particle
set itself or the underlying spatial data structures. It was designed because
the Aboria::Particles
container
and the neighbour search classes (e.g. Aboria::NeighbourQueryBase
)
were unsuitable for copying to a gpu in order to perform calculations there,
so a simpler class, the query class, was created with the required functionality.
typedef Particles_t::query_type Query_t; Query_t query = particles.get_query();
The base class for all the query objects is Aboria::NeighbourQueryBase
,
and all the query classes for the individual data structures are derived from
this. Now that we have query
,
we can create a Aboria::NeighbourQueryBase::child_iterator
.
This is the lowest level data structure iterator, and allows you to iterate
through a set of child nodes attached to a single parent node within a tree
structure.
Note | |
---|---|
All the spatial data structures in Aboria are considered trees. For the HyperOctree and kdtree data structures, this description is obvious, but the cell list is also treated as a tree, in this case a tree with one root node having N children, where N is the total number of buckets in the cell list. |
Note | |
---|---|
Aboria tends to use the terms nodes and buckets fairly interchangably. |
You can create a child_iterator by using the Aboria::NeighbourQueryBase::get_children
function. This creates a child_iterator that loops through the children of
the root node of the tree. We will use this iterator to loop through all these
children and print out the spatial bounds of each node. In order to determine
when we have reached the end of the children, we can compare the iterator to
false
. This pattern is widely
used in Aboria, rather than specific end
iterators as used in the STL.
for (auto i = query.get_children(); i != false; ++i) { std::cout << query.get_bounds(i) << std::endl; }
Above we use the Aboria::NeighbourQueryBase::get_bounds
function to get the bounds of the child iterator. This returns a Aboria::bbox
class that contains the minimum
and maximum spatial extents of the node pointed to by i
.
Note that here we are using the default spatial data structure, a cell list
provided by Aboria::CellList
,
so the "tree" here will only have 2 levels, and the loop above will
loop through the second (i.e. non-root) level. We can also create a proper
tree structure using the hyper oct-tree data structure given by Aboria::HyperOctree
, like so:
typedef Particles<std::tuple<>, 2, std::vector, HyperOctree> ParticlesOcttree_t; ParticlesOcttree_t particles_octtree(N); for (size_t i = 0; i < N; ++i) { auto &gen = get<generator>(particles_octtree)[i]; get<position>(particles_octtree)[i] = vdouble2(uniform(gen), uniform(gen)); } particles_octtree.init_neighbour_search(vdouble2(0, 0), vdouble2(1, 1), vdouble2(false, false)); auto query_octtree = particles_octtree.get_query();
Now particles_octtree
contains
a full oct-tree, dividing the spatial domain into a hierarchical set of boxes
that make up our tree data structure. The simplest iteration we might want
to do on the tree is a depth-first iteration, which is easiest achieved by
recursion. The Aboria::NeighbourQueryBase::get_children
function can be used to get the children of a Aboria::NeighbourQueryBase::child_iterator
,
and using a C++ lambda function to provide the recursion we can implement a
depth-first iteration like so
std::cout << "recursive depth-first" << std::endl; for (auto i = query_octtree.get_children(); i != false; ++i) { std::function<void(decltype(i) &)> depth_first; depth_first = [&](const auto &parent) { std::cout << query_octtree.get_bounds(parent) << std::endl; for (auto i = query_octtree.get_children(parent); i != false; ++i) { depth_first(i); } }; depth_first(i); }
This construction might be a bit clumsy to use in practice however, so Aboria
provides a special depth-first iterator Aboria::NeighbourQueryBase::all_iterator
to allow you to write a loop equivalent to the recursive depth-first code given
above.
The Aboria::NeighbourQueryBase::get_subtree
function returns a Aboria::NeighbourQueryBase::all_iterator
that performs a depth-first iteration over the tree. Note that you can also
pass in a child_iterator to Aboria::NeighbourQueryBase::get_subtree
to iterate over the sub-tree below a particular node of the tree.
std::cout << "subtree depth-first" << std::endl; for (auto i = query_octtree.get_subtree(); i != false; ++i) { std::cout << query_octtree.get_bounds(i.get_child_iterator()) << std::endl; }
You might also want to distinguish between leaf nodes (nodes with no children)
and non-leaf nodes. You can do this with the Aboria::NeighbourQueryBase::is_leaf_node
function, which takes a reference to a node (rather than an iterator), and
can be used like so
std::cout << "subtree depth-first showing leaf nodes" << std::endl; for (auto i = query_octtree.get_subtree(); i != false; ++i) { auto ci = i.get_child_iterator(); if (query_octtree.is_leaf_node(*ci)) { std::cout << "leaf node with bounds = " << query_octtree.get_bounds(ci) << std::endl; } else { std::cout << "non-leaf node with bounds = " << query_octtree.get_bounds(ci) << std::endl; } }
Leaf nodes in the tree are the only nodes that contain particles. You can loop
through all the particles in a given leaf node using the Aboria::NeighbourQueryBase::get_bucket_particles
function, which returns an iterator. Note for non-leaf nodes, the Aboria::NeighbourQueryBase::get_bucket_particles
will return an iterator that is immediatelly false, so this loop is safe even
for non-leaf nodes.
std::cout << "subtree depth-first showing leaf nodes and particles" << std::endl; for (auto i = query_octtree.get_subtree(); i != false; ++i) { auto ci = i.get_child_iterator(); if (query_octtree.is_leaf_node(*ci)) { std::cout << "leaf node with bounds = " << query_octtree.get_bounds(ci) << std::endl; for (auto j = query_octtree.get_bucket_particles(*ci); j != false; ++j) { std::cout << "\t has particle with position" << get<position>(*j) << std::endl; } } else { std::cout << "non-leaf node with bounds = " << query_octtree.get_bounds(ci) << std::endl; } }
Aboria also provides functions to query leaf nodes, or buckets, within a certain
distance of a point, and these are used internally for the neighbour search
functionality discussed in earlier sections. You can use the Aboria::NeighbourQueryBase::get_buckets_near_point
function, which returns a Aboria::NeighbourQueryBase::query_iterator
of all the buckets with a given distance of a point. This function also takes
a template argument P
, which
refers to the p-norm distance that it uses (i.e. P=2 is the standard euclidean
distance).
Caution | |
---|---|
The distance search provided by |
const int P = 2; const vdouble2 search_point = vdouble2(0.5, 0.5); const double search_radius = 0.1; std::cout << "searching within " << search_point << " of point " << search_point << std::endl; for (auto i = query_octtree.get_buckets_near_point<P>(search_point, search_radius); i != false; ++i) { auto ci = i.get_child_iterator(); std::cout << "\t found bucket at " << query_octtree.get_bounds(ci) << std::endl; for (auto j = query_octtree.get_bucket_particles(*ci); j != false; ++j) { std::cout << "\t\t found particle at " << get<position>(*j) << std::endl; } }