From: Luca Heltai Date: Mon, 14 Aug 2017 00:32:10 +0000 (-0600) Subject: Addressed comments, renamed classes, externalized implementation X-Git-Tag: v9.0.0-rc1~1262^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9c5f6d3c80d25b420cb41cc7b2671c07d82a54f6;p=dealii.git Addressed comments, renamed classes, externalized implementation --- diff --git a/cmake/configure/configure_nanoflann.cmake b/cmake/configure/configure_nanoflann.cmake index c6bc02b693..6585109417 100644 --- a/cmake/configure/configure_nanoflann.cmake +++ b/cmake/configure/configure_nanoflann.cmake @@ -1,6 +1,6 @@ ## --------------------------------------------------------------------- ## -## Copyright (C) 2012 - 2014 by the deal.II authors +## Copyright (C) 2017 by the deal.II authors ## ## This file is part of the deal.II library. ## diff --git a/cmake/modules/FindNANOFLANN.cmake b/cmake/modules/FindNANOFLANN.cmake index b34d80a11d..79bb99c089 100644 --- a/cmake/modules/FindNANOFLANN.cmake +++ b/cmake/modules/FindNANOFLANN.cmake @@ -1,6 +1,6 @@ ## --------------------------------------------------------------------- ## -## Copyright (C) 2014 - 2015 by the deal.II authors +## Copyright (C) 2017 by the deal.II authors ## ## This file is part of the deal.II library. ## diff --git a/doc/news/changes/major/20170813LucaHeltai b/doc/news/changes/major/20170813LucaHeltai index a1ccfab047..ccccf9c534 100644 --- a/doc/news/changes/major/20170813LucaHeltai +++ b/doc/news/changes/major/20170813LucaHeltai @@ -1,7 +1,7 @@ New: A new KDTreeDistance class has been added that interfaces -to the nanoflann library (https://github.com/jlblancoc/nanoflann). +the nanoflann library (https://github.com/jlblancoc/nanoflann). This class can be used to extract nearest neighbour information -on collection of points, query for the closes points to a target +on collection of points, query for the closest points to a target point or all points contained within a given distance.
(Luca Heltai, 2017/08/13) diff --git a/include/deal.II/numerics/kdtree_distance.h b/include/deal.II/numerics/kdtree.h similarity index 65% rename from include/deal.II/numerics/kdtree_distance.h rename to include/deal.II/numerics/kdtree.h index 08e048aa10..9e9556375a 100644 --- a/include/deal.II/numerics/kdtree_distance.h +++ b/include/deal.II/numerics/kdtree.h @@ -13,8 +13,8 @@ // // --------------------------------------------------------------------- -#ifndef _dealii__numerics_kdtree_distance_h -#define _dealii__numerics_kdtree_distance_h +#ifndef _dealii__numerics_kdtree_h +#define _dealii__numerics_kdtree_h #include @@ -33,15 +33,35 @@ DEAL_II_NAMESPACE_OPEN * A wrapper for the nanoflann library, used to compute the distance from a * collection of points, and to efficiently return nearest neighbors to a * target point. This function uses nanoflann to efficiently partition the - * space in a tree. The cost of each query is then roughly of order log(n), - * where n is the number of points stored in this class. + * space in a k-dimensional tree. The cost of each query is then roughly of + * order log(n), where n is the number of points stored in this class. * - * The wrapper provides methods that give access to some of thefunctionalities + * The wrapper provides methods that give access to some of the functionalities * of the nanoflann library, like searching the n nearest neighbors, or - * searching the points that fall within a raidius of a target point. + * searching the points that fall within a radius of a target point. + * + * @quotation + * From wikipedia (https://en.wikipedia.org/wiki/K-d_tree): + * + * A k-d tree is a binary tree in which every node is a k-dimensional point. + * Every non-leaf node can be thought of as implicitly generating a splitting + * hyperplane that divides the space into two parts, known as half-spaces. + * Points to the left of this hyperplane are represented by the left subtree of + * that node and points right of the hyperplane are represented by the right + * subtree. The hyperplane direction is chosen in the following way: every node + * in the tree is associated with one of the k-dimensions, with the hyperplane + * perpendicular to that dimension's axis. So, for example, if for a particular + * split the "x" axis is chosen, all points in the subtree with a smaller "x" + * value than the node will appear in the left subtree and all points with + * larger "x" value will be in the right subtree. In such a case, the + * hyperplane would be set by the x-value of the point, and its normal would be + * the unit x-axis. + * @endquotation + * + * @author Luca Heltai, 2017. */ template -class KDTreeDistance +class KDTree { public: /** @@ -56,7 +76,7 @@ public: * a vector of points will result in an exception. Only a reference * to the points is stored, so you should make sure that the life of * the the vector you pass is longer than the life of this class, or - * you'll get undefinite behaviour. + * you'll get undefined behaviour. * * If you update the vector of points in someway, remember to call * again the set_points() method. The tree and the index are @@ -65,8 +85,8 @@ public: * your points, and do not call again set_points(), your results * will likely be wrong. */ - KDTreeDistance(const unsigned int &max_leaf_size=10, - const std::vector > &pts=std::vector >()); + KDTree(const unsigned int &max_leaf_size=10, + const std::vector > &pts=std::vector >()); /** @@ -86,45 +106,32 @@ public: * Reference to the vector of points from which we want to compute * the distance. */ - const std::vector > &points; //!< A const ref to the data set origin + const std::vector > &points; /** - * The constrcutor needs the data set source. + * The constructor needs the vector of points from which we want to build + * the tree. */ - PointCloudAdaptor(const std::vector > &_points) : points(_points) { } + PointCloudAdaptor(const std::vector > &_points); /** * Return number of points in the data set (required by nanoflann). */ - inline size_t kdtree_get_point_count() const - { - return points.size(); - } + size_t kdtree_get_point_count() const; /** * Return the L2 distance between points */ - inline coord_t kdtree_distance(const coord_t *p1, const size_t idx_p2,size_t size) const - { - AssertDimension(size, dim); - coord_t res=0.0; - for (size_t d=0; d , - PointCloudAdaptor, dim, unsigned int> KDTree; + PointCloudAdaptor, dim, unsigned int> NanoFlannKDTree; /** @@ -188,7 +195,6 @@ public: * @param[in] point: the target point * @param[in] radius: the radius of the ball * @param[in] sorted: sort the output results in ascending order with respect to distances - * @param[out] matches: * * @return vector of indices and distances of the matching points */ @@ -197,11 +203,8 @@ public: bool sorted=false) const; /** - * Fill two vectors with the indices and distances of the closest - * points to the given target point. The vectors are filled with - * indices and distances until there is space in them. You should - * resize them to the number of closest points you wish to get. An - * assertion is thrown if the vectors do not have the same size. + * Fill a vectors with the indices and distances of the closest @p n_points + * points to the given target point. * * @param[in] target: the target point * @param[in] n_points: the number of requested points @@ -227,15 +230,16 @@ private: /** * The actual kdtree. */ - std::unique_ptr kdtree; + std::unique_ptr kdtree; }; //------------ inline functions ------------- +#ifndef DOXYGEN template inline -unsigned int KDTreeDistance::size() const +unsigned int KDTree::size() const { if (adaptor) return adaptor->points.size(); @@ -243,23 +247,66 @@ unsigned int KDTreeDistance::size() const return 0; }; + + template inline const Point & -KDTreeDistance::operator[](unsigned int i) const +KDTree::operator[](unsigned int i) const { AssertIndexRange(i, size()); return adaptor->points[i]; } + +template +KDTree::PointCloudAdaptor::PointCloudAdaptor(const std::vector > &_points) + : points(_points) +{} + + + +template +inline size_t +KDTree::PointCloudAdaptor::kdtree_get_point_count() const +{ + return points.size(); +} + + + +template +inline double +KDTree::PointCloudAdaptor::kdtree_get_pt(const size_t idx, int d) const +{ + AssertIndexRange(d,dim); + return points[idx][d]; +} + + + template template inline bool -KDTreeDistance::PointCloudAdaptor::kdtree_get_bbox(BBOX &) const +KDTree::PointCloudAdaptor::kdtree_get_bbox(BBOX &) const { return false; } + + +template +inline double +KDTree::PointCloudAdaptor::kdtree_distance(const double *p1, const size_t idx_p2,size_t size) const +{ + AssertDimension(size, dim); + double res=0.0; + for (size_t d=0; d +#include #ifdef DEAL_II_WITH_NANOFLANN @@ -8,8 +8,8 @@ DEAL_II_NAMESPACE_OPEN template -KDTreeDistance::KDTreeDistance(const unsigned int &max_leaf_size, - const std::vector > &pts) +KDTree::KDTree(const unsigned int &max_leaf_size, + const std::vector > &pts) : max_leaf_size(max_leaf_size) { if (pts.size() > 0) @@ -18,7 +18,7 @@ KDTreeDistance::KDTreeDistance(const unsigned int &max_leaf_size, template -std::vector > KDTreeDistance::get_points_within_ball(const Point ¢er, +std::vector > KDTree::get_points_within_ball(const Point ¢er, const double &radius, bool sorted) const { @@ -36,7 +36,7 @@ std::vector > KDTreeDistance::get_points_wi } template -std::vector > KDTreeDistance::get_closest_points(const Point &target, +std::vector > KDTree::get_closest_points(const Point &target, const unsigned int n_points) const { Assert(adaptor, ExcNotInitialized()); @@ -52,18 +52,18 @@ std::vector > KDTreeDistance::get_closest_p } template -void KDTreeDistance::set_points(const std::vector > &pts) +void KDTree::set_points(const std::vector > &pts) { Assert(pts.size() > 0, ExcMessage("Expecting a non zero set of points.")); adaptor = std_cxx14::make_unique(pts); - kdtree = std_cxx14::make_unique(dim, *adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(max_leaf_size)); + kdtree = std_cxx14::make_unique(dim, *adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(max_leaf_size)); kdtree->buildIndex(); } -template class KDTreeDistance<1>; -template class KDTreeDistance<2>; -template class KDTreeDistance<3>; +template class KDTree<1>; +template class KDTree<2>; +template class KDTree<3>; DEAL_II_NAMESPACE_CLOSE diff --git a/tests/numerics/kdtree_distance_01.cc b/tests/numerics/kdtree_01.cc similarity index 62% rename from tests/numerics/kdtree_distance_01.cc rename to tests/numerics/kdtree_01.cc index 2c0ae731e6..9879dd38a6 100644 --- a/tests/numerics/kdtree_distance_01.cc +++ b/tests/numerics/kdtree_01.cc @@ -1,23 +1,23 @@ -//----------------------------------------------------------- +// --------------------------------------------------------------------- // -// Copyright (C) 2015 by the deal2lkit authors +// Copyright (C) 2017 by the deal.II authors // -// This file is part of the deal2lkit library. +// This file is part of the deal.II library. // -// The deal2lkit 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 at -// the top level of the deal2lkit distribution. +// 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 at +// the top level of the deal2lkit distribution. // -//----------------------------------------------------------- +//--------------------------------------------------------------------- // Create a list of points, and compute the minimum distance from some other points // to this set, using a kdtree library #include "../tests.h" -#include +#include int main () diff --git a/tests/numerics/kdtree_distance_01.output b/tests/numerics/kdtree_01.output similarity index 100% rename from tests/numerics/kdtree_distance_01.output rename to tests/numerics/kdtree_01.output diff --git a/tests/numerics/kdtree_distance_02.cc b/tests/numerics/kdtree_02.cc similarity index 63% rename from tests/numerics/kdtree_distance_02.cc rename to tests/numerics/kdtree_02.cc index 0706c20825..ff423b7f73 100644 --- a/tests/numerics/kdtree_distance_02.cc +++ b/tests/numerics/kdtree_02.cc @@ -1,22 +1,22 @@ -//----------------------------------------------------------- +// --------------------------------------------------------------------- // -// Copyright (C) 2015 by the deal2lkit authors +// Copyright (C) 2017 by the deal.II authors // -// This file is part of the deal2lkit library. +// This file is part of the deal.II library. // -// The deal2lkit 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 at -// the top level of the deal2lkit distribution. +// 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 at +// the top level of the deal2lkit distribution. // -//----------------------------------------------------------- +//--------------------------------------------------------------------- // Create a list of points, and compute the closest points to a given one. #include "../tests.h" -#include +#include int main () { diff --git a/tests/numerics/kdtree_distance_02.output b/tests/numerics/kdtree_02.output similarity index 100% rename from tests/numerics/kdtree_distance_02.output rename to tests/numerics/kdtree_02.output diff --git a/tests/numerics/kdtree_distance_03.cc b/tests/numerics/kdtree_03.cc similarity index 63% rename from tests/numerics/kdtree_distance_03.cc rename to tests/numerics/kdtree_03.cc index df87f83f2e..b7fefa610d 100644 --- a/tests/numerics/kdtree_distance_03.cc +++ b/tests/numerics/kdtree_03.cc @@ -1,22 +1,22 @@ -//----------------------------------------------------------- +// --------------------------------------------------------------------- // -// Copyright (C) 2015 by the deal2lkit authors +// Copyright (C) 2017 by the deal.II authors // -// This file is part of the deal2lkit library. +// This file is part of the deal.II library. // -// The deal2lkit 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 at -// the top level of the deal2lkit distribution. +// 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 at +// the top level of the deal2lkit distribution. // -//----------------------------------------------------------- +//--------------------------------------------------------------------- // Create a list of points, and the ones within a specified radius of the target points #include "../tests.h" -#include +#include int main () { diff --git a/tests/numerics/kdtree_distance_03.output b/tests/numerics/kdtree_03.output similarity index 100% rename from tests/numerics/kdtree_distance_03.output rename to tests/numerics/kdtree_03.output