From 3036aa59762a27ee562036e29b78dc97a21053e0 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Tue, 19 Sep 2017 14:49:28 +0200 Subject: [PATCH] GridTools::find_closest_vertex(map, point) --- doc/news/changes/minor/20170919LucaHeltai | 6 ++ include/deal.II/grid/grid_tools.h | 13 ++++ source/grid/grid_tools.cc | 15 +++++ source/grid/grid_tools.inst.in | 4 ++ .../grid_tools_extract_used_vertices_02.cc | 62 +++++++++++++++++++ ...grid_tools_extract_used_vertices_02.output | 23 +++++++ 6 files changed, 123 insertions(+) create mode 100644 doc/news/changes/minor/20170919LucaHeltai create mode 100644 tests/grid/grid_tools_extract_used_vertices_02.cc create mode 100644 tests/grid/grid_tools_extract_used_vertices_02.output diff --git a/doc/news/changes/minor/20170919LucaHeltai b/doc/news/changes/minor/20170919LucaHeltai new file mode 100644 index 0000000000..e406b0a7fe --- /dev/null +++ b/doc/news/changes/minor/20170919LucaHeltai @@ -0,0 +1,6 @@ +New: A GridTools::find_closest_vertex() has been added +that takes a map of index-vertex position and a target +point, and returns the index of the closest point in +the map. +
+(Luca Heltai, 2017/09/19) diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index 4585e16b59..f6fdc81b45 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -593,6 +593,19 @@ namespace GridTools const Triangulation &container, const Mapping &mapping = StaticMappingQ1::mapping); + /** + * Find and return the number of the closest vertex to a given point in the + * map of vertices passed as the first argument. + * + * @param vertices A map of index->vertex, as returned by + * extract_used_vertices() + * @param p The target point + * @return The index of the vertex that is closest to the target point `p` + */ + template + unsigned int + find_closest_vertex (const std::map> &vertices, + const Point &p); /** * Find and return the number of the used vertex (or marked vertex) in a diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 1229339a78..cac93c628a 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -4943,6 +4943,21 @@ next_cell: return result; } + template + unsigned int + find_closest_vertex(const std::map > &vertices, + const Point &p) + { + auto id_and_v = std::min_element(vertices.begin(), vertices.end(), + [&](const std::pair> &p1, + const std::pair> &p2) -> bool + { + return p1.second.distance(p) < p2.second.distance(p); + } + ); + return id_and_v->first; + } + } /* namespace GridTools */ diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index 2fc9db2041..3938a71118 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -143,6 +143,10 @@ for (deal_II_space_dimension : SPACE_DIMENSIONS) find_active_cell_around_point (const Mapping &, const parallel::distributed::Triangulation &, const Point &); + + template unsigned int + GridTools::find_closest_vertex(const std::map >& vertices, + const Point& p); } diff --git a/tests/grid/grid_tools_extract_used_vertices_02.cc b/tests/grid/grid_tools_extract_used_vertices_02.cc new file mode 100644 index 0000000000..50f3374133 --- /dev/null +++ b/tests/grid/grid_tools_extract_used_vertices_02.cc @@ -0,0 +1,62 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2001 - 2015 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +// Check extract used_vertices and find_closest_vertex + +#include "../tests.h" +#include +#include +#include +#include + + +template +void test (const Point &p) +{ + deallog << "dim: " << dim << ", spacedim: " + << spacedim << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube(tria); + + tria.refine_global(1); + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + + for (auto cell = tria.begin_active(2); cell != tria.end(2); ++cell) + cell->set_coarsen_flag(); + + tria.execute_coarsening_and_refinement(); + + auto m = GridTools::extract_used_vertices(tria); + + for (auto &e: m) + deallog << "Vertex: " << e.first << ": " << e.second << std::endl; + + auto i = GridTools::find_closest_vertex(m,p); + deallog << "Closest vertex to " << p + << ", v["<< i << "] :" + << m[i] << std::endl; +}; + + +int main () +{ + initlog(); + test<2,2> (Point<2>(.2,.2)); + test<2,2> (Point<2>(.6,.9)); + return 0; +} + diff --git a/tests/grid/grid_tools_extract_used_vertices_02.output b/tests/grid/grid_tools_extract_used_vertices_02.output new file mode 100644 index 0000000000..de11874f0a --- /dev/null +++ b/tests/grid/grid_tools_extract_used_vertices_02.output @@ -0,0 +1,23 @@ + +DEAL::dim: 2, spacedim: 2 +DEAL::Vertex: 0: 0.00000 0.00000 +DEAL::Vertex: 1: 1.00000 0.00000 +DEAL::Vertex: 2: 0.00000 1.00000 +DEAL::Vertex: 3: 1.00000 1.00000 +DEAL::Vertex: 4: 0.500000 0.00000 +DEAL::Vertex: 5: 0.00000 0.500000 +DEAL::Vertex: 6: 1.00000 0.500000 +DEAL::Vertex: 7: 0.500000 1.00000 +DEAL::Vertex: 8: 0.500000 0.500000 +DEAL::Closest vertex to 0.200000 0.200000, v[0] :0.00000 0.00000 +DEAL::dim: 2, spacedim: 2 +DEAL::Vertex: 0: 0.00000 0.00000 +DEAL::Vertex: 1: 1.00000 0.00000 +DEAL::Vertex: 2: 0.00000 1.00000 +DEAL::Vertex: 3: 1.00000 1.00000 +DEAL::Vertex: 4: 0.500000 0.00000 +DEAL::Vertex: 5: 0.00000 0.500000 +DEAL::Vertex: 6: 1.00000 0.500000 +DEAL::Vertex: 7: 0.500000 1.00000 +DEAL::Vertex: 8: 0.500000 0.500000 +DEAL::Closest vertex to 0.600000 0.900000, v[7] :0.500000 1.00000 -- 2.39.5