]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add mapping argument to find_closest_vertex_of_cell(). Closes #7922.
authorLuca Heltai <luca.heltai@sissa.it>
Sun, 10 May 2020 00:21:05 +0000 (02:21 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Sun, 10 May 2020 08:19:24 +0000 (10:19 +0200)
doc/news/changes/minor/20200510LucaHeltai [new file with mode: 0644]
include/deal.II/grid/grid_tools.h
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/grid_tools_find_closest_vertex_02.cc [new file with mode: 0644]
tests/grid/grid_tools_find_closest_vertex_02.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20200510LucaHeltai b/doc/news/changes/minor/20200510LucaHeltai
new file mode 100644 (file)
index 0000000..d4b9fc7
--- /dev/null
@@ -0,0 +1,4 @@
+Fixed: GridTools::find_closest_vertex_of_cell() now takes an optional mapping argument, to take into
+account of all those cases where a mapping may modify the location of the cell vertices.
+<br>
+(Luca Heltai, 2020/05/10)
index b3935cd15abf6f60dc7f0a203ae999f890c25b79..9a8016766c4e0853c9e256e9a9d55db5ab100d33 100644 (file)
@@ -1824,15 +1824,20 @@ namespace GridTools
 
   /**
    * Return the local vertex index of cell @p cell that is closest to
-   * the given location @p position.
+   * the given location @p position. The location of the vertices is extracted
+   * from the (optional) @p mapping argument, to guarantee that the correct
+   * answer is returned when the underlying mapping modifies the position of the
+   * vertices.
    *
-   * @author Rene Gassmoeller, Luca Heltai, 2017.
+   * @author Rene Gassmoeller, Luca Heltai, 2017, 2020.
    */
   template <int dim, int spacedim>
   unsigned int
   find_closest_vertex_of_cell(
     const typename Triangulation<dim, spacedim>::active_cell_iterator &cell,
-    const Point<spacedim> &position);
+    const Point<spacedim> &                                            position,
+    const Mapping<dim, spacedim> &                                     mapping =
+      StaticMappingQ1<dim, spacedim>::mapping);
 
   /**
    * Compute a globally unique index for each vertex and hanging node
index 685c975582751c8e96dc6752e7da374fc4d3e5b6..6a78bf705e28f29581dcd653e09b5969db95c58b 100644 (file)
@@ -1604,9 +1604,12 @@ namespace GridTools
         // invalid cell, then query for the closest global vertex
         if (current_cell.state() == IteratorState::valid)
           {
+            const auto cell_vertices = mapping.get_vertices(current_cell);
             const unsigned int closest_vertex =
-              find_closest_vertex_of_cell<dim, spacedim>(current_cell, p);
-            vertex_to_point      = p - current_cell->vertex(closest_vertex);
+              find_closest_vertex_of_cell<dim, spacedim>(current_cell,
+                                                         p,
+                                                         mapping);
+            vertex_to_point      = p - cell_vertices[closest_vertex];
             closest_vertex_index = current_cell->vertex_index(closest_vertex);
           }
         else
@@ -1736,15 +1739,16 @@ namespace GridTools
   unsigned int
   find_closest_vertex_of_cell(
     const typename Triangulation<dim, spacedim>::active_cell_iterator &cell,
-    const Point<spacedim> &                                            position)
+    const Point<spacedim> &                                            position,
+    const Mapping<dim, spacedim> &                                     mapping)
   {
-    double       minimum_distance = position.distance_square(cell->vertex(0));
+    const auto   vertices         = mapping.get_vertices(cell);
+    double       minimum_distance = position.distance_square(vertices[0]);
     unsigned int closest_vertex   = 0;
 
     for (unsigned int v = 1; v < GeometryInfo<dim>::vertices_per_cell; ++v)
       {
-        const double vertex_distance =
-          position.distance_square(cell->vertex(v));
+        const double vertex_distance = position.distance_square(vertices[v]);
         if (vertex_distance < minimum_distance)
           {
             closest_vertex   = v;
index c08037543632c9074a30b0155aeeb9ee8bf21c75..77b083198df159ba400f0e1b25570100450bf873 100644 (file)
@@ -68,7 +68,8 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
         const typename Triangulation<
           deal_II_dimension,
           deal_II_space_dimension>::active_cell_iterator &,
-        const Point<deal_II_space_dimension> &);
+        const Point<deal_II_space_dimension> &,
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &);
 
       template std::map<unsigned int, types::global_vertex_index>
       compute_local_to_global_vertex_index_map(
diff --git a/tests/grid/grid_tools_find_closest_vertex_02.cc b/tests/grid/grid_tools_find_closest_vertex_02.cc
new file mode 100644 (file)
index 0000000..28ed479
--- /dev/null
@@ -0,0 +1,78 @@
+// ---------------------------------------------------------------------
+//
+// 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 find_closest_vertex_of_cell when a mapping is present
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/mapping_q_eulerian.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+
+template <int dim, int spacedim>
+void
+test(const Point<spacedim> &p, double displacement)
+{
+  deallog << "dim: " << dim << ", spacedim: " << spacedim
+          << " with displacement " << displacement << std::endl;
+
+  Triangulation<dim, spacedim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global();
+
+  // Vector fe
+  FESystem<dim, spacedim>   fe{FE_Q<dim, spacedim>(1), spacedim};
+  DoFHandler<dim, spacedim> dh(tria);
+  dh.distribute_dofs(fe);
+
+  // Displacement vector.
+  Vector<double> disp(dh.n_dofs());
+  for (unsigned int i = 0; i < dh.n_dofs(); ++i)
+    disp[i] = displacement;
+
+  // mapping
+  MappingQEulerian<dim, Vector<double>, spacedim> mapping(1, dh, disp);
+
+  for (const auto cell : tria.active_cell_iterators())
+    {
+      const auto i =
+        GridTools::find_closest_vertex_of_cell<dim, spacedim>(cell, p, mapping);
+      const auto &v = mapping.get_vertices(cell);
+
+      deallog << "Closest vertex of cell " << cell << " to " << p
+              << " is cell->vertex(" << i << ") :" << v[i] << std::endl;
+    }
+};
+
+
+int
+main()
+{
+  initlog();
+  test<2, 2>(Point<2>(.45, .45), 0);
+  test<2, 2>(Point<2>(.45, .45), 0.5);
+
+  test<2, 3>(Point<3>(.9, .9, .1), 0.0);
+  test<2, 3>(Point<3>(.9, .9, .1), 1.0);
+  return 0;
+}
diff --git a/tests/grid/grid_tools_find_closest_vertex_02.output b/tests/grid/grid_tools_find_closest_vertex_02.output
new file mode 100644 (file)
index 0000000..3839578
--- /dev/null
@@ -0,0 +1,21 @@
+
+DEAL::dim: 2, spacedim: 2 with displacement 0.00000
+DEAL::Closest vertex of cell 1.0 to 0.450000 0.450000 is cell->vertex(3) :0.500000 0.500000
+DEAL::Closest vertex of cell 1.1 to 0.450000 0.450000 is cell->vertex(2) :0.500000 0.500000
+DEAL::Closest vertex of cell 1.2 to 0.450000 0.450000 is cell->vertex(1) :0.500000 0.500000
+DEAL::Closest vertex of cell 1.3 to 0.450000 0.450000 is cell->vertex(0) :0.500000 0.500000
+DEAL::dim: 2, spacedim: 2 with displacement 0.500000
+DEAL::Closest vertex of cell 1.0 to 0.450000 0.450000 is cell->vertex(0) :0.500000 0.500000
+DEAL::Closest vertex of cell 1.1 to 0.450000 0.450000 is cell->vertex(0) :1.00000 0.500000
+DEAL::Closest vertex of cell 1.2 to 0.450000 0.450000 is cell->vertex(0) :0.500000 1.00000
+DEAL::Closest vertex of cell 1.3 to 0.450000 0.450000 is cell->vertex(0) :1.00000 1.00000
+DEAL::dim: 2, spacedim: 3 with displacement 0.00000
+DEAL::Closest vertex of cell 1.0 to 0.900000 0.900000 0.100000 is cell->vertex(3) :0.500000 0.500000 0.00000
+DEAL::Closest vertex of cell 1.1 to 0.900000 0.900000 0.100000 is cell->vertex(3) :1.00000 0.500000 0.00000
+DEAL::Closest vertex of cell 1.2 to 0.900000 0.900000 0.100000 is cell->vertex(3) :0.500000 1.00000 0.00000
+DEAL::Closest vertex of cell 1.3 to 0.900000 0.900000 0.100000 is cell->vertex(3) :1.00000 1.00000 0.00000
+DEAL::dim: 2, spacedim: 3 with displacement 1.00000
+DEAL::Closest vertex of cell 1.0 to 0.900000 0.900000 0.100000 is cell->vertex(0) :1.00000 1.00000 1.00000
+DEAL::Closest vertex of cell 1.1 to 0.900000 0.900000 0.100000 is cell->vertex(0) :1.50000 1.00000 1.00000
+DEAL::Closest vertex of cell 1.2 to 0.900000 0.900000 0.100000 is cell->vertex(0) :1.00000 1.50000 1.00000
+DEAL::Closest vertex of cell 1.3 to 0.900000 0.900000 0.100000 is cell->vertex(0) :1.50000 1.50000 1.00000

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.