From: Luca Heltai Date: Sat, 23 May 2020 22:41:36 +0000 (+0200) Subject: pack_rtree_of_indices(). X-Git-Tag: v9.3.0-rc1~1560^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9ebcb7efc8b4fcfdaa00074e08e19e2c183f81b6;p=dealii.git pack_rtree_of_indices(). --- diff --git a/doc/news/changes/minor/20200524LucaHeltai b/doc/news/changes/minor/20200524LucaHeltai new file mode 100644 index 0000000000..650c8cfe0c --- /dev/null +++ b/doc/news/changes/minor/20200524LucaHeltai @@ -0,0 +1,4 @@ +New: pack_rtree_of_indices() and IndexableGetterFromIndices allow to construct a +RTree object that stores indices to existing containers of indexable objects. +
+(Luca Heltai, 2020/05/24) diff --git a/include/deal.II/numerics/rtree.h b/include/deal.II/numerics/rtree.h index 26c424ab97..d1762a9408 100644 --- a/include/deal.II/numerics/rtree.h +++ b/include/deal.II/numerics/rtree.h @@ -19,6 +19,7 @@ #include #include +#include #include #include @@ -127,25 +128,38 @@ DEAL_II_NAMESPACE_OPEN * // Returns the 3 closest points to the Segment defined above. * @endcode * - * @author Luca Heltai, 2018. + * In general, a tree does not need to store the actual objects, as long as it + * knows how to access a const reference to an indexable type. This is + * achieved by passing the optional template argument @p IndexableGetter, that + * extracts a const reference to one of the possible indexable types given an + * object of type @p LeafType. As an example, one may store points in a container, + * and only create a tree of the indices within the container, using the + * IndexableGetterFromIndices class defined below, and the function + * pack_rtree_of_indices(). + * + * @author Luca Heltai, 2018, 2020. */ template > -using RTree = boost::geometry::index::rtree; + typename IndexType = boost::geometry::index::linear<16>, + typename IndexableGetter = + boost::geometry::index::indexable> +using RTree = + boost::geometry::index::rtree; /** * Construct the correct RTree object by passing an iterator range. * * Notice that the order of the parameters is the opposite with respect to the * RTree class, since we can automatically infer the @p LeafType from the - * arguments, and we only need to specify the @p IndexType if the default is not - * adequate. + * arguments. * * @author Luca Heltai, 2018. */ template , - typename LeafTypeIterator> -RTree + typename LeafTypeIterator, + typename IndexableGetter = boost::geometry::index::indexable< + typename LeafTypeIterator::value_type>> +RTree pack_rtree(const LeafTypeIterator &begin, const LeafTypeIterator &end); /** @@ -161,10 +175,121 @@ pack_rtree(const LeafTypeIterator &begin, const LeafTypeIterator &end); * @author Luca Heltai, 2018. */ template , - typename ContainerType> -RTree + typename ContainerType, + typename IndexableGetter = boost::geometry::index::indexable< + typename ContainerType::value_type>> +RTree pack_rtree(const ContainerType &container); +/** + * A class that may be used as an @p IndexableGetter template argument for an + * RTree that stores indices to entries in a @p Container type. + * + * This class may be used as a proxy to extract an indexable type from + * compatible containers. For example: + * @code + * std::vector, double> > point_temperature = fill(); + * IndexableGetterFromIndices + * getter(point_temperature); + * + * const Point &p = getter(i); // returns point_temperature[i].first; + * @endcode + * + * This class is used by the pack_rtree_of_indices() function to construct an + * RTree where the leaves are indices pointing to the entries of the container + * passed to this class. + * + * @author Luca Heltai, 2020. + */ +template +class IndexableGetterFromIndices +{ +public: + /** + * An alias for the boost type that is used to extract a Point, Segment, or + * BoundingBox from compatible types (pairs, tuples, etc.). + */ + using IndexableGetter = + typename boost::geometry::index::indexable; + + /** + * An alias to the actual geometrical type. + */ + using result_type = typename IndexableGetter::result_type; + + /** + * An alias to the index type. + */ + using size_t = typename Container::size_type; + + /** + * Constructor. Store a const reference to a container. + */ + explicit IndexableGetterFromIndices(Container const &c) + : container(c) + {} + + /** + * Implements the @p IndexableGetter requirements of the rtree class. + */ + result_type + operator()(size_t i) const + { + return getter(container[i]); + } + +private: + /** + * A const reference to the container. + */ + const Container &container; + + /** + * An instantiation of the getter that allows easy translation from the + * container @p value_type and the actual indexable type. + */ + IndexableGetter getter; +}; + +/** + * Construct a RTree object that stores the indices of an existing container of + * indexable types. The only requirement on the container is that it supports + * operator[] for any index between 0 and the size of the container (i.e., a + * std::vector, or an std::array will do, however an std::map won't). + * + * Differently from the object created by the pack_rtree() function, in this + * case we don't store the actual geometrical types, but just their indices: + * + * @code + * namespace bgi = boost::geometry::index; + * std::vector> some_points = fill(); + * auto tree = pack_rtree(points); // the tree contains a copy of the points + * auto index_tree = pack_rtree_of_indices(points); // the tree contains only + * // the indices of the + * // points + * BoundingBox box = build_a_box(); + * + * for(const auto &p: tree | bgi::adaptors::queried(bgi::intersects(box))) + * std::cout << "Point p: " << p << std::endl; + * + * for(const auto &i: index_tree | bgi::adaptors::queried(bgi::intersects(box))) + * std::cout << "Point p: " << some_points[i] << std::endl; + * @endcode + * + * The leaves stored in the tree are the indices of the corresponding entry in + * the container. A reference to the external container is stored internally, + * but keep in mind that if you change the container, you should rebuild the + * tree. + * + * @author Luca Heltai, 2020. + */ +template , + typename ContainerType> +RTree> +pack_rtree_of_indices(const ContainerType &container); + /** * Helper structure that allows one to extract a level from an RTree as a vector * of BoundingBox objects. @@ -319,20 +444,46 @@ extract_rtree_level(const Rtree &tree, const unsigned int level); // Inline and template functions #ifndef DOXYGEN -template -RTree + +template +RTree pack_rtree(const LeafTypeIterator &begin, const LeafTypeIterator &end) { - return RTree(begin, end); + return RTree(begin, end); } -template -RTree +template +RTree pack_rtree(const ContainerType &container) { - return pack_rtree(container.begin(), container.end()); + return pack_rtree( + container.begin(), container.end()); +} + + + +template +RTree> +pack_rtree_of_indices(const ContainerType &container) +{ + std_cxx20::ranges::iota_view + indices(0, container.size()); + return RTree>( + indices.begin(), + indices.end(), + IndexType(), + IndexableGetterFromIndices(container)); } diff --git a/tests/boost/rtree_03.cc b/tests/boost/rtree_03.cc new file mode 100644 index 0000000000..94b5ec1bec --- /dev/null +++ b/tests/boost/rtree_03.cc @@ -0,0 +1,57 @@ +/* --------------------------------------------------------------------- + * + * Copyright (C) 2020 by the deal.II authors + * + * This file is part of the deal.II library. + * + * The deal.II library is free software; you can use it, redistribute + * it, and/or modify it under the terms of the GNU Lesser General + * Public License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * The full text of the license can be found in the file LICENSE.md at + * the top level directory of deal.II. + * + * --------------------------------------------------------------------- + */ + +// Check that we can construct a boost rtree based on indices. + +#include +#include + +#include +#include +#include + +#include + +#include + +#include "../tests.h" + +namespace bgi = boost::geometry::index; + +int +main() +{ + initlog(); + const unsigned int N = 30; + std::vector> points(N); + for (auto &p : points) + p = random_point<2>(); + + auto tree = pack_rtree(points); + auto tree2 = pack_rtree_of_indices(points); + + Point<2> p0(0, 0); + Point<2> p1(.4, .7111); + Segment<2> segment(p0, p1); + + for (const auto &p : tree | bgi::adaptors::queried(bgi::nearest(segment, 3))) + deallog << "P: " << p << std::endl; + + deallog << "Same with other tree?" << std::endl; + + for (const auto &i : tree2 | bgi::adaptors::queried(bgi::nearest(segment, 3))) + deallog << "P[" << i << "]: " << points[i] << std::endl; +} diff --git a/tests/boost/rtree_03.output b/tests/boost/rtree_03.output new file mode 100644 index 0000000000..369ce83c5f --- /dev/null +++ b/tests/boost/rtree_03.output @@ -0,0 +1,8 @@ + +DEAL::P: 0.129790 0.108809 +DEAL::P: 0.277775 0.553970 +DEAL::P: 0.0641713 0.0200230 +DEAL::Same with other tree? +DEAL::P[13]: 0.129790 0.108809 +DEAL::P[4]: 0.277775 0.553970 +DEAL::P[28]: 0.0641713 0.0200230