From: Timo Heister <timo.heister@gmail.com>
Date: Fri, 22 Jan 2021 17:04:59 +0000 (-0500)
Subject: Introduce Mapping::is_compatible_with
X-Git-Tag: v9.3.0-rc1~550^2~2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=19fbe8cebc3266ca35b4bec0761a88257a44b277;p=dealii.git

Introduce Mapping::is_compatible_with

- Introduce Mapping::is_compatible_with()
- Implement for all Mapping classes
- Use in FEValues::reinit()
---

diff --git a/include/deal.II/fe/mapping.h b/include/deal.II/fe/mapping.h
index b6d587aa04..ede2ebe400 100644
--- a/include/deal.II/fe/mapping.h
+++ b/include/deal.II/fe/mapping.h
@@ -401,6 +401,14 @@ public:
   virtual bool
   preserves_vertex_locations() const = 0;
 
+  /**
+   * Returns is this instance of Mapping is compatible with the type of cell
+   * in @p cell_type.
+   *
+   */
+  virtual bool
+  is_compatible_with(const ReferenceCell::Type &cell_type) const = 0;
+
   /**
    * @name Mapping points between reference and real cells
    * @{
diff --git a/include/deal.II/fe/mapping_cartesian.h b/include/deal.II/fe/mapping_cartesian.h
index 5fbb61f554..543ab649dd 100644
--- a/include/deal.II/fe/mapping_cartesian.h
+++ b/include/deal.II/fe/mapping_cartesian.h
@@ -88,6 +88,9 @@ public:
   virtual bool
   preserves_vertex_locations() const override;
 
+  virtual bool
+  is_compatible_with(const ReferenceCell::Type &cell_type) const override;
+
   /**
    * @name Mapping points between reference and real cells
    * @{
diff --git a/include/deal.II/fe/mapping_fe.h b/include/deal.II/fe/mapping_fe.h
index 95c8591884..0e0bb66592 100644
--- a/include/deal.II/fe/mapping_fe.h
+++ b/include/deal.II/fe/mapping_fe.h
@@ -85,6 +85,9 @@ public:
                      &cell) const override;
 
 
+  virtual bool
+  is_compatible_with(const ReferenceCell::Type &cell_type) const override;
+
   /**
    * Always returns @p true because the default implementation of functions in
    * this class preserves vertex locations.
@@ -604,6 +607,18 @@ MappingFE<dim, spacedim>::preserves_vertex_locations() const
   return true;
 }
 
+
+template <int dim, int spacedim>
+bool
+MappingFE<dim, spacedim>::is_compatible_with(
+  const ReferenceCell::Type &cell_type) const
+{
+  if (cell_type.get_dimension() != dim)
+    return false; // TODO: or is this an error?
+
+  return true;
+}
+
 #endif // DOXYGEN
 
 /* -------------- declaration of explicit specializations ------------- */
diff --git a/include/deal.II/fe/mapping_fe_field.h b/include/deal.II/fe/mapping_fe_field.h
index 7f36f19cf9..c0bf9d8e2b 100644
--- a/include/deal.II/fe/mapping_fe_field.h
+++ b/include/deal.II/fe/mapping_fe_field.h
@@ -210,6 +210,9 @@ public:
   virtual bool
   preserves_vertex_locations() const override;
 
+  virtual bool
+  is_compatible_with(const ReferenceCell::Type &cell_type) const override;
+
   /**
    * Return the mapped vertices of a cell.
    *
diff --git a/include/deal.II/fe/mapping_manifold.h b/include/deal.II/fe/mapping_manifold.h
index 60c3f4b206..aeb7ef723b 100644
--- a/include/deal.II/fe/mapping_manifold.h
+++ b/include/deal.II/fe/mapping_manifold.h
@@ -79,6 +79,9 @@ public:
   virtual bool
   preserves_vertex_locations() const override;
 
+  virtual bool
+  is_compatible_with(const ReferenceCell::Type &cell_type) const override;
+
   /**
    * @name Mapping points between reference and real cells
    * @{
@@ -434,6 +437,23 @@ MappingManifold<dim, spacedim>::preserves_vertex_locations() const
   return true;
 }
 
+
+template <int dim, int spacedim>
+bool
+MappingManifold<dim, spacedim>::is_compatible_with(
+  const ReferenceCell::Type &cell_type) const
+{
+  if (cell_type.get_dimension() != dim)
+    return false; // TODO: or is this an error?
+
+  if (cell_type.is_hyper_cube())
+    return true;
+
+  return false;
+}
+
+
+
 #endif // DOXYGEN
 
 /* -------------- declaration of explicit specializations ------------- */
diff --git a/include/deal.II/fe/mapping_q.h b/include/deal.II/fe/mapping_q.h
index 168191bce4..54154b470d 100644
--- a/include/deal.II/fe/mapping_q.h
+++ b/include/deal.II/fe/mapping_q.h
@@ -124,6 +124,9 @@ public:
   get_bounding_box(const typename Triangulation<dim, spacedim>::cell_iterator
                      &cell) const override;
 
+  virtual bool
+  is_compatible_with(const ReferenceCell::Type &cell_type) 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.
diff --git a/include/deal.II/fe/mapping_q_generic.h b/include/deal.II/fe/mapping_q_generic.h
index 8a2d6ddb0e..d456a3ad68 100644
--- a/include/deal.II/fe/mapping_q_generic.h
+++ b/include/deal.II/fe/mapping_q_generic.h
@@ -170,6 +170,9 @@ public:
   get_bounding_box(const typename Triangulation<dim, spacedim>::cell_iterator
                      &cell) const override;
 
+  virtual bool
+  is_compatible_with(const ReferenceCell::Type &cell_type) const override;
+
   /**
    * @name Mapping points between reference and real cells
    * @{
diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc
index 4f5348d568..7000d921dd 100644
--- a/source/fe/fe_values.cc
+++ b/source/fe/fe_values.cc
@@ -4499,6 +4499,13 @@ void
 FEValues<dim, spacedim>::reinit(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell)
 {
+  // Check that mapping and reference cell type are compatible:
+  Assert(this->get_mapping().is_compatible_with(cell->reference_cell_type()),
+         ExcMessage(
+           "You are trying to call FEValues::reinit() with a cell of type " +
+           cell->reference_cell_type().to_string() +
+           " with a Mapping that is not compatible with it."));
+
   // no FE in this cell, so no assertion
   // necessary here
   this->maybe_invalidate_previous_present_cell(cell);
@@ -4528,6 +4535,13 @@ FEValues<dim, spacedim>::reinit(
            static_cast<const FiniteElementData<dim> &>(cell->get_fe()),
          (typename FEValuesBase<dim, spacedim>::ExcFEDontMatch()));
 
+  // Check that mapping and reference cell type are compatible:
+  Assert(this->get_mapping().is_compatible_with(cell->reference_cell_type()),
+         ExcMessage(
+           "You are trying to call FEValues::reinit() with a cell of type " +
+           cell->reference_cell_type().to_string() +
+           " with a Mapping that is not compatible with it."));
+
   this->maybe_invalidate_previous_present_cell(cell);
   this->check_cell_similarity(cell);
 
diff --git a/source/fe/mapping_cartesian.cc b/source/fe/mapping_cartesian.cc
index dbb65d925e..1a138613d1 100644
--- a/source/fe/mapping_cartesian.cc
+++ b/source/fe/mapping_cartesian.cc
@@ -70,6 +70,22 @@ MappingCartesian<dim, spacedim>::preserves_vertex_locations() const
 
 
 
+template <int dim, int spacedim>
+bool
+MappingCartesian<dim, spacedim>::is_compatible_with(
+  const ReferenceCell::Type &cell_type) const
+{
+  if (cell_type.get_dimension() != dim)
+    return false; // TODO: or is this an error?
+
+  if (cell_type.is_hyper_cube())
+    return true;
+
+  return false;
+}
+
+
+
 template <int dim, int spacedim>
 UpdateFlags
 MappingCartesian<dim, spacedim>::requires_update_flags(
diff --git a/source/fe/mapping_fe_field.cc b/source/fe/mapping_fe_field.cc
index 0cad06d6f4..c9815c61b3 100644
--- a/source/fe/mapping_fe_field.cc
+++ b/source/fe/mapping_fe_field.cc
@@ -342,6 +342,19 @@ MappingFEField<dim, spacedim, VectorType, void>::preserves_vertex_locations()
 
 
 
+template <int dim, int spacedim, typename VectorType>
+bool
+MappingFEField<dim, spacedim, VectorType, void>::is_compatible_with(
+  const ReferenceCell::Type &cell_type) const
+{
+  if (cell_type.get_dimension() != dim)
+    return false; // TODO: or is this an error?
+
+  return true;
+}
+
+
+
 template <int dim, int spacedim, typename VectorType>
 boost::container::small_vector<Point<spacedim>,
                                GeometryInfo<dim>::vertices_per_cell>
@@ -2221,6 +2234,7 @@ MappingFEField<dim, spacedim, VectorType, void>::get_degree() const
 }
 
 
+
 template <int dim, int spacedim, typename VectorType>
 ComponentMask
 MappingFEField<dim, spacedim, VectorType, void>::get_component_mask() const
diff --git a/source/fe/mapping_q.cc b/source/fe/mapping_q.cc
index 38f412e293..185bf102dc 100644
--- a/source/fe/mapping_q.cc
+++ b/source/fe/mapping_q.cc
@@ -559,6 +559,22 @@ MappingQ<dim, spacedim>::get_bounding_box(
 
 
 
+template <int dim, int spacedim>
+bool
+MappingQ<dim, spacedim>::is_compatible_with(
+  const ReferenceCell::Type &cell_type) const
+{
+  if (cell_type.get_dimension() != dim)
+    return false; // TODO: or is this an error?
+
+  if (cell_type.is_hyper_cube())
+    return true;
+
+  return false;
+}
+
+
+
 // explicit instantiations
 #include "mapping_q.inst"
 
diff --git a/source/fe/mapping_q_generic.cc b/source/fe/mapping_q_generic.cc
index 02dd963685..0f90658d45 100644
--- a/source/fe/mapping_q_generic.cc
+++ b/source/fe/mapping_q_generic.cc
@@ -1652,6 +1652,22 @@ MappingQGeneric<dim, spacedim>::get_bounding_box(
 
 
 
+template <int dim, int spacedim>
+bool
+MappingQGeneric<dim, spacedim>::is_compatible_with(
+  const ReferenceCell::Type &cell_type) const
+{
+  if (cell_type.get_dimension() != dim)
+    return false; // TODO: or is this an error?
+
+  if (cell_type.is_hyper_cube())
+    return true;
+
+  return false;
+}
+
+
+
 //--------------------------- Explicit instantiations -----------------------
 #include "mapping_q_generic.inst"