--- /dev/null
+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)
* 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
*/
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.
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.
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.
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
* @{
+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"
+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"
+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"
+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"
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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);
+ }
+}
--- /dev/null
+
+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
+