]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Return the box of support points instead of vertices. 11570/head
authorLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Jan 2021 21:19:16 +0000 (22:19 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Jan 2021 21:19:16 +0000 (22:19 +0100)
12 files changed:
doc/news/changes/minor/20210116LucaHeltai [new file with mode: 0644]
include/deal.II/fe/mapping.h
include/deal.II/fe/mapping_fe.h
include/deal.II/fe/mapping_q.h
include/deal.II/fe/mapping_q_eulerian.h
include/deal.II/fe/mapping_q_generic.h
source/fe/mapping_fe.cc
source/fe/mapping_q.cc
source/fe/mapping_q_eulerian.cc
source/fe/mapping_q_generic.cc
tests/mappings/mapping_get_bounding_box_01.cc [new file with mode: 0644]
tests/mappings/mapping_get_bounding_box_01.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20210116LucaHeltai b/doc/news/changes/minor/20210116LucaHeltai
new file mode 100644 (file)
index 0000000..5f54d90
--- /dev/null
@@ -0,0 +1,4 @@
+New: Mapping::get_bounding_box() now returns the bounding box of the support points, when the mapping
+is of type MappingQ, MappingQGeneric, and MappingQEulerian.
+<br>
+(Luca Heltai, 2021/01/16)
index f04ee0890a0f9242a31254aa7ff5dc3c6effbec1..b6d587aa04be71e7393f2091b67bf66982f34882 100644 (file)
@@ -375,11 +375,11 @@ public:
    * displacements or choose completely different locations, e.g.,
    * MappingQEulerian, MappingQ1Eulerian, or MappingFEField.
    *
-   * This function returns the bounding box containing all the vertices of the
-   * cell, as returned by the get_vertices() method. Beware of the fact that
-   * for higher order mappings this bounding box is only an approximation of the
-   * true bounding box, since it does not take into account curved faces, and it
-   * may be smaller than the true bounding box.
+   * For linear mappings, this function returns the bounding box containing all
+   * the vertices of the cell, as returned by the get_vertices() method. For
+   * higher order mappings defined through support points, the bounding box is
+   * only guaranteed to contain all the support points, and it is, in general,
+   * only an approximation of the true bounding box, which may be larger.
    *
    * @param[in] cell The cell for which you want to compute the bounding box
    */
index 0865491ee0c5db86344e60643f47a9306395c9b7..24f4ad439cdcbacc94be88c44ece2ad4366900e1 100644 (file)
@@ -79,6 +79,12 @@ public:
   unsigned int
   get_degree() const;
 
+  // for documentation, see the Mapping base class
+  virtual BoundingBox<spacedim>
+  get_bounding_box(const typename Triangulation<dim, spacedim>::cell_iterator
+                     &cell) const override;
+
+
   /**
    * Always returns @p true because the default implementation of functions in
    * this class preserves vertex locations.
index 47ad2bf5d80de1604cbb07139a5bc7f3d4bfd207..168191bce434cd95701f1f60adebe30a00fe23e5 100644 (file)
@@ -119,6 +119,11 @@ public:
   virtual bool
   preserves_vertex_locations() const override;
 
+  // for documentation, see the Mapping base class
+  virtual BoundingBox<spacedim>
+  get_bounding_box(const typename Triangulation<dim, spacedim>::cell_iterator
+                     &cell) const override;
+
   /**
    * Transform the point @p p on the unit cell to the point @p p_real on the
    * real cell @p cell and returns @p p_real.
index 10bd13b13f888084c9586be54c86bc880bce5423..bf4d09041751ebb78cee4bf66f52de413dbdfeea 100644 (file)
@@ -140,6 +140,11 @@ public:
   virtual bool
   preserves_vertex_locations() const override;
 
+  // for documentation, see the Mapping base class
+  virtual BoundingBox<spacedim>
+  get_bounding_box(const typename Triangulation<dim, spacedim>::cell_iterator
+                     &cell) const override;
+
   /**
    * Exception which is thrown when the mapping is being evaluated at
    * non-active cell.
index a3000b187be9c308964b0e32d5d5137944cb7354..8a2d6ddb0e3df4e3f33ada353f511d59215079dc 100644 (file)
@@ -165,6 +165,11 @@ public:
   virtual bool
   preserves_vertex_locations() const override;
 
+  // for documentation, see the Mapping base class
+  virtual BoundingBox<spacedim>
+  get_bounding_box(const typename Triangulation<dim, spacedim>::cell_iterator
+                     &cell) const override;
+
   /**
    * @name Mapping points between reference and real cells
    * @{
index 3f7148dbca07893aad0e98d37c88592c8311fd86..1aea54c729fd6852c53059fb59dfa43408be3b42 100644 (file)
@@ -2157,6 +2157,16 @@ MappingFE<dim, spacedim>::compute_mapping_support_points(
 
 
 
+template <int dim, int spacedim>
+BoundingBox<spacedim>
+MappingFE<dim, spacedim>::get_bounding_box(
+  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
+{
+  return BoundingBox<spacedim>(this->compute_mapping_support_points(cell));
+}
+
+
+
 //--------------------------- Explicit instantiations -----------------------
 #include "mapping_fe.inst"
 
index 69276944a313854daabc12817a7824d6dff3b8fe..38f412e293a3bf80ac90accfe8c9f049f96cab63 100644 (file)
@@ -544,6 +544,21 @@ MappingQ<dim, spacedim>::clone() const
 
 
 
+template <int dim, int spacedim>
+BoundingBox<spacedim>
+MappingQ<dim, spacedim>::get_bounding_box(
+  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
+{
+  if (cell->has_boundary_lines() || use_mapping_q_on_all_cells ||
+      (dim != spacedim))
+    return BoundingBox<spacedim>(
+      qp_mapping->compute_mapping_support_points(cell));
+  else
+    return BoundingBox<spacedim>(q1_mapping->get_vertices(cell));
+}
+
+
+
 // explicit instantiations
 #include "mapping_q.inst"
 
index d2990e8cb9f1d128d8898697482e4e0d74e35fa5..dcd0d57c4f89444daeef86c7425c9b7cb785bbd8 100644 (file)
@@ -271,6 +271,18 @@ MappingQEulerian<dim, VectorType, spacedim>::fill_fe_values(
 
 
 
+template <int dim, class VectorType, int spacedim>
+BoundingBox<spacedim>
+MappingQEulerian<dim, VectorType, spacedim>::get_bounding_box(
+  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
+{
+  return BoundingBox<spacedim>(
+    dynamic_cast<const MappingQEulerianGeneric &>(*this->qp_mapping)
+      .compute_mapping_support_points(cell));
+}
+
+
+
 // explicit instantiations
 #include "mapping_q_eulerian.inst"
 
index 4559fb6c9c26c58d06a4e5e5ec3c72d301d9e6ae..73c167865428b73e73203faced6e8897675ed4dc 100644 (file)
@@ -1643,6 +1643,16 @@ MappingQGeneric<dim, spacedim>::compute_mapping_support_points(
 
 
 
+template <int dim, int spacedim>
+BoundingBox<spacedim>
+MappingQGeneric<dim, spacedim>::get_bounding_box(
+  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
+{
+  return BoundingBox<spacedim>(this->compute_mapping_support_points(cell));
+}
+
+
+
 //--------------------------- Explicit instantiations -----------------------
 #include "mapping_q_generic.inst"
 
diff --git a/tests/mappings/mapping_get_bounding_box_01.cc b/tests/mappings/mapping_get_bounding_box_01.cc
new file mode 100644 (file)
index 0000000..6c495f3
--- /dev/null
@@ -0,0 +1,80 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2006 - 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 bounding boxes are created correctly with high order mappings
+
+#include <deal.II/base/bounding_box_data_out.h>
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/fe/mapping_q.h>
+#include <deal.II/fe/mapping_q_generic.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools_cache.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include <boost/core/demangle.hpp>
+
+#include "../tests.h"
+
+
+template <int dim, int spacedim>
+void
+test_bounding_boxes(const Mapping<dim, spacedim> &mapping,
+                    const unsigned int            degree)
+{
+  deallog << "Testing " << boost::core::demangle(typeid(mapping).name()) << "("
+          << degree << ")" << std::endl;
+
+  Triangulation<dim, spacedim> triangulation;
+  GridGenerator::hyper_ball(triangulation);
+  GridTools::Cache<dim, spacedim> cache(triangulation, mapping);
+
+  const auto boxes = cache.get_cell_bounding_boxes_rtree();
+
+  {
+    std::string fname = "boxes_" + std::to_string(dim) + "_" +
+                        std::to_string(spacedim) + "_" +
+                        std::to_string(degree) + ".vtk";
+    std::ofstream           ofile(fname);
+    BoundingBoxDataOut<dim> data_out;
+    DataOutBase::VtkFlags   flags;
+    flags.print_date_and_time = false;
+    data_out.set_flags(flags);
+    data_out.build_patches(boxes);
+    data_out.write_vtk(ofile);
+    cat_file(fname.c_str());
+  }
+}
+
+
+int
+main()
+{
+  initlog();
+  {
+    unsigned int degree = 2;
+    MappingQ<2>  mapping(degree);
+    test_bounding_boxes(mapping, degree);
+  }
+  {
+    unsigned int       degree = 2;
+    MappingQGeneric<2> mapping(degree);
+    test_bounding_boxes(mapping, degree);
+  }
+}
diff --git a/tests/mappings/mapping_get_bounding_box_01.output b/tests/mappings/mapping_get_bounding_box_01.output
new file mode 100644 (file)
index 0000000..6ceb716
--- /dev/null
@@ -0,0 +1,79 @@
+
+DEAL::Testing dealii::MappingQ<2, 2>(2)
+# vtk DataFile Version 3.0
+#This file was generated by the deal.II library.
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 20 double
+-0.707107 -1 0
+0.707107 -1 0
+-0.707107 -0.292893 0
+0.707107 -0.292893 0
+-1 -0.707107 0
+-0.292893 -0.707107 0
+-1 0.707107 0
+-0.292893 0.707107 0
+-0.292893 -0.292893 0
+0.292893 -0.292893 0
+-0.292893 0.292893 0
+0.292893 0.292893 0
+0.292893 -0.707107 0
+1 -0.707107 0
+0.292893 0.707107 0
+1 0.707107 0
+-0.707107 0.292893 0
+0.707107 0.292893 0
+-0.707107 1 0
+0.707107 1 0
+
+CELLS 5 25
+4      0       1       3       2
+4      4       5       7       6
+4      8       9       11      10
+4      12      13      15      14
+4      16      17      19      18
+
+CELL_TYPES 5
+ 9 9 9 9 9
+POINT_DATA 20
+
+DEAL::Testing dealii::MappingQGeneric<2, 2>(2)
+# vtk DataFile Version 3.0
+#This file was generated by the deal.II library.
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 20 double
+-0.707107 -1 0
+0.707107 -1 0
+-0.707107 -0.292893 0
+0.707107 -0.292893 0
+-1 -0.707107 0
+-0.292893 -0.707107 0
+-1 0.707107 0
+-0.292893 0.707107 0
+-0.292893 -0.292893 0
+0.292893 -0.292893 0
+-0.292893 0.292893 0
+0.292893 0.292893 0
+0.292893 -0.707107 0
+1 -0.707107 0
+0.292893 0.707107 0
+1 0.707107 0
+-0.707107 0.292893 0
+0.707107 0.292893 0
+-0.707107 1 0
+0.707107 1 0
+
+CELLS 5 25
+4      0       1       3       2
+4      4       5       7       6
+4      8       9       11      10
+4      12      13      15      14
+4      16      17      19      18
+
+CELL_TYPES 5
+ 9 9 9 9 9
+POINT_DATA 20
+

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.