]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make the class MappingQ equivalent to MappingQGeneric
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 9 Jun 2021 18:08:07 +0000 (20:08 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 5 Jul 2021 10:08:26 +0000 (12:08 +0200)
include/deal.II/fe/mapping_c1.h
include/deal.II/fe/mapping_q.h
include/deal.II/fe/mapping_q_eulerian.h
source/fe/mapping_c1.cc
source/fe/mapping_q.cc
source/fe/mapping_q_eulerian.cc

index b5f98850ddf339dd5867f39da6e2803cad9edaee..621d4547d1e30244de21fc78670c3cf5aa56c3dd 100644 (file)
@@ -28,14 +28,14 @@ DEAL_II_NAMESPACE_OPEN
 
 /**
  * Mapping class that uses C1 (continuously differentiable) cubic mappings of
- * the boundary. This class is built atop of MappingQ by simply determining
- * the interpolation points for a cubic mapping of the boundary differently:
- * MappingQ chooses them such that they interpolate the boundary, while this
- * class chooses them such that the discretized boundary is globally
- * continuously differentiable.
+ * the boundary. This class is built atop of MappingQGeneric by simply
+ * determining the interpolation points for a cubic mapping of the boundary
+ * differently: MappingQ chooses them such that they interpolate the boundary,
+ * while this class chooses them such that the discretized boundary is
+ * globally continuously differentiable.
  */
 template <int dim, int spacedim = dim>
-class MappingC1 : public MappingQ<dim, spacedim>
+class MappingC1 : public MappingQGeneric<dim, spacedim>
 {
 public:
   /**
@@ -51,48 +51,33 @@ public:
   virtual std::unique_ptr<Mapping<dim, spacedim>>
   clone() const override;
 
-protected:
   /**
-   * A class derived from MappingQGeneric that provides the generic mapping
-   * with support points on boundary objects so that the corresponding Q3
-   * mapping ends up being C1.
+   * For <tt>dim=2,3</tt>. Append the support points of all shape functions
+   * located on bounding lines to the vector @p a. Points located on the
+   * line but on vertices are not included.
+   *
+   * This function chooses the respective points not such that they are
+   * interpolating the boundary (as does the base class), but rather such
+   * that the resulting cubic mapping is a continuous one.
    */
-  class MappingC1Generic : public MappingQGeneric<dim, spacedim>
-  {
-  public:
-    /**
-     * Constructor.
-     */
-    MappingC1Generic();
-
-    /**
-     * For <tt>dim=2,3</tt>. Append the support points of all shape functions
-     * located on bounding lines to the vector @p a. Points located on the
-     * line but on vertices are not included.
-     *
-     * This function chooses the respective points not such that they are
-     * interpolating the boundary (as does the base class), but rather such
-     * that the resulting cubic mapping is a continuous one.
-     */
-    virtual void
-    add_line_support_points(
-      const typename Triangulation<dim>::cell_iterator &cell,
-      std::vector<Point<dim>> &                         a) const override;
-
-    /**
-     * For <tt>dim=3</tt>. Append the support points of all shape functions
-     * located on bounding faces (quads in 3d) to the vector @p a. Points
-     * located on the line but on vertices are not included.
-     *
-     * This function chooses the respective points not such that they are
-     * interpolating the boundary (as does the base class), but rather such
-     * that the resulting cubic mapping is a continuous one.
-     */
-    virtual void
-    add_quad_support_points(
-      const typename Triangulation<dim>::cell_iterator &cell,
-      std::vector<Point<dim>> &                         a) const override;
-  };
+  virtual void
+  add_line_support_points(
+    const typename Triangulation<dim>::cell_iterator &cell,
+    std::vector<Point<dim>> &                         a) const override;
+
+  /**
+   * For <tt>dim=3</tt>. Append the support points of all shape functions
+   * located on bounding faces (quads in 3d) to the vector @p a. Points
+   * located on the line but on vertices are not included.
+   *
+   * This function chooses the respective points not such that they are
+   * interpolating the boundary (as does the base class), but rather such
+   * that the resulting cubic mapping is a continuous one.
+   */
+  virtual void
+  add_quad_support_points(
+    const typename Triangulation<dim>::cell_iterator &cell,
+    std::vector<Point<dim>> &                         a) const override;
 };
 
 /*@}*/
@@ -103,25 +88,22 @@ protected:
 
 template <>
 void
-MappingC1<1>::MappingC1Generic::add_line_support_points(
-  const Triangulation<1>::cell_iterator &,
-  std::vector<Point<1>> &) const;
+MappingC1<1>::add_line_support_points(const Triangulation<1>::cell_iterator &,
+                                      std::vector<Point<1>> &) const;
 template <>
 void
-MappingC1<2>::MappingC1Generic::add_line_support_points(
+MappingC1<2>::add_line_support_points(
   const Triangulation<2>::cell_iterator &cell,
   std::vector<Point<2>> &                a) const;
 
 template <>
 void
-MappingC1<1>::MappingC1Generic::add_quad_support_points(
-  const Triangulation<1>::cell_iterator &,
-  std::vector<Point<1>> &) const;
+MappingC1<1>::add_quad_support_points(const Triangulation<1>::cell_iterator &,
+                                      std::vector<Point<1>> &) const;
 template <>
 void
-MappingC1<2>::MappingC1Generic::add_quad_support_points(
-  const Triangulation<2>::cell_iterator &,
-  std::vector<Point<2>> &) const;
+MappingC1<2>::add_quad_support_points(const Triangulation<2>::cell_iterator &,
+                                      std::vector<Point<2>> &) const;
 
 
 #endif // DOXYGEN
index abe2f7b43b174973339e9a8f747c29d560b64f07..5c8262aaa472740447c4e80dd350c1609120d32b 100644 (file)
@@ -27,354 +27,32 @@ DEAL_II_NAMESPACE_OPEN
 /*@{*/
 
 /**
- * A class that implements a polynomial mapping $Q_p$ of degree $p$ on cells
- * at the boundary of the domain (or, if requested in the constructor, for all
- * cells) and linear mappings for interior cells.
- *
- * The class is in fact poorly named since (unless explicitly specified during
- * the construction of the object, see below), it does not actually use
- * mappings of degree $p$ <i>everywhere</i>, but only on cells at the
- * boundary. This is in contrast to the MappingQGeneric class which indeed
- * does use a polynomial mapping $Q_p$ of degree $p$ everywhere. The point of
- * the current class is that in many situations, curved domains are only
- * provided with information about how exactly edges at the boundary are
- * shaped, but we do not know anything about internal edges. Thus, in the
- * absence of other information, we can only assume that internal edges are
- * straight lines, and in that case internal cells may as well be treated is
- * bilinear quadrilaterals or trilinear hexahedra. (An example of how such
- * meshes look is shown in step-1 already, but it is also discussed in the
- * "Results" section of step-6.) Because bi-/trilinear mappings are
- * significantly cheaper to compute than higher order mappings, it is
- * advantageous in such situations to use the higher order mapping only on
- * cells at the boundary of the domain. This class implements exactly this
- * behavior.
- *
- * There are a number of special cases worth considering:
- * - If you want to use a higher order mapping for all cells, you can
- * achieve this by setting the second argument to the constructor to true.
- * This only makes sense if you can actually provide information about how
- * interior edges and faces of the mesh should be curved. This is typically
- * done by associating a Manifold with interior cells and edges. A simple
- * example of this is discussed in the "Results" section of step-6; a full
- * discussion of manifolds is provided in step-53.
- * - If you pass true as the second argument to this class, then it
- * is in fact completely equivalent to generating a MappingQGeneric object
- * right away.
- * - This class is also entirely equivalent to MappingQGeneric if the
- * polynomial degree provided is one. This is because in that case, no
- * distinction between the mapping used on cells in the interior and on the
- * boundary of the domain can be made.
- * - If you are working on meshes embedded in higher space dimensions,
- * i.e., if dim!=spacedim, then every cell is considered to be at the boundary
- * of the domain and consequently a higher order mapping is used for all
- * cells; again this class is then equivalent to using MappingQGeneric right
- * away.
- *
- * <h4>Behavior along curved boundaries and with different manifolds</h4>
- *
- * For the behavior of the mapping and convergence rates in case of mixing
- * different manifolds, please consult the respective section of
- * MappingQGeneric.
+ * A class that implements a polynomial mapping $Q_p$ of degree $p$ on all
+ * cells. This class is completely equivalent to the MappingQGeneric class.
  */
 template <int dim, int spacedim = dim>
-class MappingQ : public Mapping<dim, spacedim>
+class MappingQ : public MappingQGeneric<dim, spacedim>
 {
 public:
   /**
    * Constructor.  @p polynomial_degree denotes the polynomial degree of the
    * polynomials that are used to map cells boundary.
-   *
-   * The second argument determines whether the higher order mapping should
-   * also be used on interior cells. If its value is <code>false</code> (the
-   * default), then a lower order mapping is used in the interior. This is
-   * sufficient for most cases where higher order mappings are only used to
-   * better approximate the boundary. In that case, cells bounded by straight
-   * lines are acceptable in the interior. However, there are cases where one
-   * would also like to use a higher order mapping in the interior. The
-   * MappingQEulerian class is one such case.
-   *
-   * The value of @p use_mapping_q_on_all_cells is ignored if @p dim is not
-   * equal to @p spacedim, i.e., if we are considering meshes on surfaces
-   * embedded into higher dimensional spaces.
-   */
-  MappingQ(const unsigned int polynomial_degree,
-           const bool         use_mapping_q_on_all_cells = false);
-
-  /**
-   * Copy constructor.
-   */
-  MappingQ(const MappingQ<dim, spacedim> &mapping);
-
-  /**
-   * Return the degree of the mapping, i.e. the value which was passed to the
-   * constructor.
-   */
-  unsigned int
-  get_degree() const;
-
-  /**
-   * 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;
-
-  virtual bool
-  is_compatible_with(const ReferenceCell &reference_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 Point<spacedim>
-  transform_unit_to_real_cell(
-    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-    const Point<dim> &p) const override;
-
-  /**
-   * Transform the point @p p on the real cell to the point @p p_unit on the
-   * unit cell @p cell and returns @p p_unit.
-   *
-   * Uses Newton iteration and the @p transform_unit_to_real_cell function.
-   *
-   * In the codimension one case, this function returns the normal projection
-   * of the real point @p p on the curve or surface identified by the @p cell.
-   *
-   * @note Polynomial mappings from the reference (unit) cell coordinates to
-   * the coordinate system of a real cell are not always invertible if the
-   * point for which the inverse mapping is to be computed lies outside the
-   * cell's boundaries.  In such cases, the current function may fail to
-   * compute a point on the reference cell whose image under the mapping
-   * equals the given point @p p.  If this is the case then this function
-   * throws an exception of type Mapping::ExcTransformationFailed .  Whether
-   * the given point @p p lies outside the cell can therefore be determined by
-   * checking whether the return reference coordinates lie inside of outside
-   * the reference cell (e.g., using GeometryInfo::is_inside_unit_cell) or
-   * whether the exception mentioned above has been thrown.
-   */
-  virtual Point<dim>
-  transform_real_to_unit_cell(
-    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-    const Point<spacedim> &p) const override;
-
-  // for documentation, see the Mapping base class
-  virtual void
-  transform(const ArrayView<const Tensor<1, dim>> &                  input,
-            const MappingKind                                        kind,
-            const typename Mapping<dim, spacedim>::InternalDataBase &internal,
-            const ArrayView<Tensor<1, spacedim>> &output) const override;
-
-  // for documentation, see the Mapping base class
-  virtual void
-  transform(const ArrayView<const DerivativeForm<1, dim, spacedim>> &input,
-            const MappingKind                                        kind,
-            const typename Mapping<dim, spacedim>::InternalDataBase &internal,
-            const ArrayView<Tensor<2, spacedim>> &output) const override;
-
-  // for documentation, see the Mapping base class
-  virtual void
-  transform(const ArrayView<const Tensor<2, dim>> &                  input,
-            const MappingKind                                        kind,
-            const typename Mapping<dim, spacedim>::InternalDataBase &internal,
-            const ArrayView<Tensor<2, spacedim>> &output) const override;
-
-  // for documentation, see the Mapping base class
-  virtual void
-  transform(const ArrayView<const DerivativeForm<2, dim, spacedim>> &input,
-            const MappingKind                                        kind,
-            const typename Mapping<dim, spacedim>::InternalDataBase &internal,
-            const ArrayView<Tensor<3, spacedim>> &output) const override;
-
-  // for documentation, see the Mapping base class
-  virtual void
-  transform(const ArrayView<const Tensor<3, dim>> &                  input,
-            const MappingKind                                        kind,
-            const typename Mapping<dim, spacedim>::InternalDataBase &internal,
-            const ArrayView<Tensor<3, spacedim>> &output) const override;
-
-  /**
-   * Return a pointer to a copy of the present object. The caller of this copy
-   * then assumes ownership of it.
-   */
-
-  virtual std::unique_ptr<Mapping<dim, spacedim>>
-  clone() const override;
-
-
-  /**
-   * @name Interface with FEValues
-   * @{
    */
+  MappingQ(const unsigned int polynomial_degree);
 
-protected:
   /**
-   * Storage for internal data of this mapping. See Mapping::InternalDataBase
-   * for an extensive description.
-   *
-   * This includes data that is computed once when the object is created (in
-   * get_data()) as well as data the class wants to store from between the
-   * call to fill_fe_values(), fill_fe_face_values(), or
-   * fill_fe_subface_values() until possible later calls from the finite
-   * element to functions such as transform(). The latter class of member
-   * variables are marked as 'mutable'.
-   *
-   * The current class uses essentially the same fields for storage as the
-   * MappingQGeneric class. Consequently, it inherits from
-   * MappingQGeneric::InternalData, rather than from
-   * Mapping::InternalDataBase. The principal difference to
-   * MappingQGeneric::InternalData is that MappingQ switches between $Q_1$ and
-   * $Q_p$ mappings depending on the cell we are on, so the internal data
-   * object needs to also store a pointer to an InternalData object that
-   * pertains to a $Q_1$ mapping.
+   * The second argument is here for backward compatibility with previous
+   * versions of deal.II, but it does not have any effect on the workings of
+   * this class.
    */
-  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
-  {
-  public:
-    /**
-     * Constructor.
-     */
-    InternalData();
-
-
-    /**
-     * Return an estimate (in bytes) for the memory consumption of this object.
-     */
-    virtual std::size_t
-    memory_consumption() const override;
-
-    /**
-     * Flag that is set by the <tt>fill_fe_[[sub]face]_values</tt> function.
-     *
-     * If this flag is @p true we are on an interior cell and the @p
-     * mapping_q1_data is used.
-     */
-    mutable bool use_mapping_q1_on_current_cell;
-
-    /**
-     * A pointer to a structure to store the information for the pure $Q_1$
-     * mapping that is, by default, used on all interior cells.
-     */
-    std::unique_ptr<typename MappingQGeneric<dim, spacedim>::InternalData>
-      mapping_q1_data;
-
-    /**
-     * A pointer to a structure to store the information for the full $Q_p$
-     * mapping that is, by default, used on all boundary cells.
-     */
-    std::unique_ptr<typename MappingQGeneric<dim, spacedim>::InternalData>
-      mapping_qp_data;
-  };
-
-protected:
-  // documentation can be found in Mapping::requires_update_flags()
-  virtual UpdateFlags
-  requires_update_flags(const UpdateFlags update_flags) const override;
-
-  // documentation can be found in Mapping::get_data()
-  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
-  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
-
-  using Mapping<dim, spacedim>::get_face_data;
-
-  // documentation can be found in Mapping::get_face_data()
-  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
-  get_face_data(const UpdateFlags               flags,
-                const hp::QCollection<dim - 1> &quadrature) const override;
-
-  // documentation can be found in Mapping::get_subface_data()
-  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
-  get_subface_data(const UpdateFlags          flags,
-                   const Quadrature<dim - 1> &quadrature) const override;
-
-  // documentation can be found in Mapping::fill_fe_values()
-  virtual CellSimilarity::Similarity
-  fill_fe_values(
-    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-    const CellSimilarity::Similarity                            cell_similarity,
-    const Quadrature<dim> &                                     quadrature,
-    const typename Mapping<dim, spacedim>::InternalDataBase &   internal_data,
-    internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
-      &output_data) const override;
-
-  using Mapping<dim, spacedim>::fill_fe_face_values;
-
-  // documentation can be found in Mapping::fill_fe_face_values()
-  virtual void
-  fill_fe_face_values(
-    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-    const unsigned int                                          face_no,
-    const hp::QCollection<dim - 1> &                            quadrature,
-    const typename Mapping<dim, spacedim>::InternalDataBase &   internal_data,
-    internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
-      &output_data) const override;
-
-  // documentation can be found in Mapping::fill_fe_subface_values()
-  virtual void
-  fill_fe_subface_values(
-    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-    const unsigned int                                          face_no,
-    const unsigned int                                          subface_no,
-    const Quadrature<dim - 1> &                                 quadrature,
-    const typename Mapping<dim, spacedim>::InternalDataBase &   internal_data,
-    internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
-      &output_data) const override;
-
-  /**
-   * @}
-   */
-
-protected:
-  /**
-   * The polynomial degree of the cells to be used on all cells at the
-   * boundary of the domain, or everywhere if so specified.
-   */
-  const unsigned int polynomial_degree;
-
-  /**
-   * If this flag is set @p true then @p MappingQ is used on all cells, not
-   * only on boundary cells.
-   */
-  const bool use_mapping_q_on_all_cells;
-
-  /**
-   * Pointer to a Q1 mapping. This mapping is used on interior cells unless
-   * use_mapping_q_on_all_cells was set in the call to the constructor. The
-   * mapping is also used on any cell in the transform_real_to_unit_cell() to
-   * compute a cheap initial guess for the position of the point before we
-   * employ the more expensive Newton iteration using the full mapping.
-   *
-   * @note MappingQEulerian resets this pointer to an object of type
-   * MappingQ1Eulerian to ensure that the Q1 mapping also knows about the
-   * proper shifts and transformations of the Eulerian displacements. This
-   * also means that we really need to store our own Q1 mapping here, rather
-   * than simply resorting to StaticMappingQ1::mapping.
-   *
-   * @note If the polynomial degree used for the current object is one, then
-   * the qp_mapping and q1_mapping variables point to the same underlying
-   * object.
-   */
-  std::shared_ptr<const MappingQGeneric<dim, spacedim>> q1_mapping;
+  DEAL_II_DEPRECATED_EARLY
+  MappingQ(const unsigned int polynomial_degree,
+           const bool         use_mapping_q_on_all_cells);
 
   /**
-   * Pointer to a Q_p mapping. This mapping is used on boundary cells unless
-   * use_mapping_q_on_all_cells was set in the call to the constructor (in
-   * which case it is used for all cells).
-   *
-   * @note MappingQEulerian and MappingC1 reset this pointer to an object of
-   * their own implementation to ensure that the Q_p mapping also knows about
-   * the proper shifts and transformations of the Eulerian displacements
-   * (Eulerian case) and proper choice of support points (C1 case).
-   *
-   * @note If the polynomial degree used for the current object is one, then
-   * the qp_mapping and q1_mapping variables point to the same underlying
-   * object.
+   * Copy constructor.
    */
-  std::shared_ptr<const MappingQGeneric<dim, spacedim>> qp_mapping;
+  MappingQ(const MappingQ<dim, spacedim> &mapping);
 };
 
 /*@}*/
index 1e19e037d2a656a7aa958715b1f1b61681013a00..bd1a4686e61d8f661c0e77bbd6e4a146b0597e70 100644 (file)
@@ -90,7 +90,7 @@ class Vector;
  * of the vector can be specified as template parameter <tt>VectorType</tt>.
  */
 template <int dim, typename VectorType = Vector<double>, int spacedim = dim>
-class MappingQEulerian : public MappingQ<dim, spacedim>
+class MappingQEulerian : public MappingQGeneric<dim, spacedim>
 {
 public:
   /**
@@ -141,9 +141,10 @@ public:
   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;
+  virtual std::vector<Point<spacedim>>
+  compute_mapping_support_points(
+    const typename Triangulation<dim, spacedim>::cell_iterator &cell)
+    const override;
 
   /**
    * Exception which is thrown when the mapping is being evaluated at
@@ -182,7 +183,6 @@ protected:
                MappingQEulerian<dim, VectorType, spacedim>>
     euler_dof_handler;
 
-
 private:
   /**
    * Multigrid level at which the mapping is to be used.
@@ -190,89 +190,37 @@ private:
   const unsigned int level;
 
   /**
-   * A class derived from MappingQGeneric that provides the generic mapping
-   * with support points on boundary objects so that the corresponding Q3
-   * mapping ends up being C1.
+   * Special quadrature rule used to define the support points in the
+   * reference configuration.
    */
-  class MappingQEulerianGeneric : public MappingQGeneric<dim, spacedim>
+  class SupportQuadrature : public Quadrature<dim>
   {
   public:
     /**
-     * Constructor.
-     */
-    MappingQEulerianGeneric(
-      const unsigned int                                 degree,
-      const MappingQEulerian<dim, VectorType, spacedim> &mapping_q_eulerian);
-
-    /**
-     * Return the mapped vertices of the cell. For the current class, this
-     * function does not use the support points from the geometry of the
-     * current cell but instead evaluates an externally given displacement
-     * field in addition to the geometry of the cell.
+     * Constructor, with an argument defining the desired polynomial degree.
      */
-    virtual boost::container::small_vector<Point<spacedim>,
-                                           GeometryInfo<dim>::vertices_per_cell>
-    get_vertices(const typename Triangulation<dim, spacedim>::cell_iterator
-                   &cell) const override;
-
-    /**
-     * Compute the positions of the support points in the current
-     * configuration. See the documentation of
-     * MappingQGeneric::compute_mapping_support_points() for more information.
-     */
-    virtual std::vector<Point<spacedim>>
-    compute_mapping_support_points(
-      const typename Triangulation<dim, spacedim>::cell_iterator &cell)
-      const override;
-
-    /**
-     * Always return @p false because MappingQEulerianGeneric does not in general
-     * preserve vertex locations (unless the translation vector happens to
-     * provide for zero displacements at vertex locations).
-     */
-    virtual bool
-    preserves_vertex_locations() const override;
-
-  private:
-    /**
-     * Reference to the surrounding object off of which we live.
-     */
-    const MappingQEulerian<dim, VectorType, spacedim> &mapping_q_eulerian;
-
-
-    /**
-     * Special quadrature rule used to define the support points in the
-     * reference configuration.
-     */
-    class SupportQuadrature : public Quadrature<dim>
-    {
-    public:
-      /**
-       * Constructor, with an argument defining the desired polynomial degree.
-       */
-      SupportQuadrature(const unsigned int map_degree);
-    };
+    SupportQuadrature(const unsigned int map_degree);
+  };
 
-    /**
-     * A member variable holding the quadrature points in the right order.
-     */
-    const SupportQuadrature support_quadrature;
+  /**
+   * A member variable holding the quadrature points in the right order.
+   */
+  const SupportQuadrature support_quadrature;
 
-    /**
-     * FEValues object used to query the given finite element field at the
-     * support points in the reference configuration.
-     *
-     * The variable is marked as mutable since we have to call
-     * FEValues::reinit from compute_mapping_support_points, a function that
-     * is 'const'.
-     */
-    mutable FEValues<dim, spacedim> fe_values;
+  /**
+   * FEValues object used to query the given finite element field at the
+   * support points in the reference configuration.
+   *
+   * The variable is marked as mutable since we have to call
+   * FEValues::reinit from compute_mapping_support_points, a function that
+   * is 'const'.
+   */
+  mutable FEValues<dim, spacedim> fe_values;
 
-    /**
-     * A variable to guard access to the fe_values variable.
-     */
-    mutable Threads::Mutex fe_values_mutex;
-  };
+  /**
+   * A variable to guard access to the fe_values variable.
+   */
+  mutable Threads::Mutex fe_values_mutex;
 };
 
 /*@}*/
index de998fc4d732959e964e5b52c5c3f7db27ea271d..d5165e3887835e7da046ae1eb6fd70f6dc975a41 100644 (file)
@@ -27,35 +27,19 @@ DEAL_II_NAMESPACE_OPEN
 
 
 
-template <int dim, int spacedim>
-MappingC1<dim, spacedim>::MappingC1Generic::MappingC1Generic()
-  : MappingQGeneric<dim, spacedim>(3)
-{}
-
-
-
 template <int dim, int spacedim>
 MappingC1<dim, spacedim>::MappingC1()
-  : MappingQ<dim, spacedim>(3)
+  : MappingQGeneric<dim, spacedim>(3)
 {
   Assert(dim > 1, ExcImpossibleInDim(dim));
-
-  // replace the mapping_qp objects of the base class by something
-  // that knows about generating data points based on the geometry
-  //
-  // we only need to replace the Qp mapping because that's the one that's
-  // used on boundary cells where it matters
-  this->qp_mapping =
-    std::make_shared<MappingC1<dim, spacedim>::MappingC1Generic>();
 }
 
 
 
 template <>
 void
-MappingC1<1>::MappingC1Generic::add_line_support_points(
-  const Triangulation<1>::cell_iterator &,
-  std::vector<Point<1>> &) const
+MappingC1<1>::add_line_support_points(const Triangulation<1>::cell_iterator &,
+                                      std::vector<Point<1>> &) const
 {
   const unsigned int dim = 1;
   (void)dim;
@@ -66,7 +50,7 @@ MappingC1<1>::MappingC1Generic::add_line_support_points(
 
 template <>
 void
-MappingC1<2>::MappingC1Generic::add_line_support_points(
+MappingC1<2>::add_line_support_points(
   const Triangulation<2>::cell_iterator &cell,
   std::vector<Point<2>> &                a) const
 {
@@ -172,7 +156,7 @@ MappingC1<2>::MappingC1Generic::add_line_support_points(
 
 template <int dim, int spacedim>
 void
-MappingC1<dim, spacedim>::MappingC1Generic::add_line_support_points(
+MappingC1<dim, spacedim>::add_line_support_points(
   const typename Triangulation<dim>::cell_iterator &,
   std::vector<Point<dim>> &) const
 {
@@ -183,9 +167,8 @@ MappingC1<dim, spacedim>::MappingC1Generic::add_line_support_points(
 
 template <>
 void
-MappingC1<1>::MappingC1Generic::add_quad_support_points(
-  const Triangulation<1>::cell_iterator &,
-  std::vector<Point<1>> &) const
+MappingC1<1>::add_quad_support_points(const Triangulation<1>::cell_iterator &,
+                                      std::vector<Point<1>> &) const
 {
   const unsigned int dim = 1;
   (void)dim;
@@ -196,9 +179,8 @@ MappingC1<1>::MappingC1Generic::add_quad_support_points(
 
 template <>
 void
-MappingC1<2>::MappingC1Generic::add_quad_support_points(
-  const Triangulation<2>::cell_iterator &,
-  std::vector<Point<2>> &) const
+MappingC1<2>::add_quad_support_points(const Triangulation<2>::cell_iterator &,
+                                      std::vector<Point<2>> &) const
 {
   const unsigned int dim = 2;
   (void)dim;
@@ -209,7 +191,7 @@ MappingC1<2>::MappingC1Generic::add_quad_support_points(
 
 template <int dim, int spacedim>
 void
-MappingC1<dim, spacedim>::MappingC1Generic::add_quad_support_points(
+MappingC1<dim, spacedim>::add_quad_support_points(
   const typename Triangulation<dim>::cell_iterator &,
   std::vector<Point<dim>> &) const
 {
index 06da719bb17d1af258b7f8cbba964517d4389ba2..aaf181dd8f92c56873b76bca5adee54c5d556031 100644 (file)
@@ -38,542 +38,23 @@ DEAL_II_NAMESPACE_OPEN
 
 
 template <int dim, int spacedim>
-MappingQ<dim, spacedim>::InternalData::InternalData()
-  : use_mapping_q1_on_current_cell(false)
+MappingQ<dim, spacedim>::MappingQ(const unsigned int degree)
+  : MappingQGeneric<dim, spacedim>(degree)
 {}
 
 
 
 template <int dim, int spacedim>
-std::size_t
-MappingQ<dim, spacedim>::InternalData::memory_consumption() const
-{
-  return (
-    Mapping<dim, spacedim>::InternalDataBase::memory_consumption() +
-    MemoryConsumption::memory_consumption(use_mapping_q1_on_current_cell) +
-    MemoryConsumption::memory_consumption(mapping_q1_data) +
-    MemoryConsumption::memory_consumption(mapping_qp_data));
-}
-
-
-
-template <int dim, int spacedim>
-MappingQ<dim, spacedim>::MappingQ(const unsigned int degree,
-                                  const bool         use_mapping_q_on_all_cells)
-  : polynomial_degree(degree)
-  ,
-
-  // see whether we want to use *this* mapping objects on *all* cells,
-  // or defer to an explicit Q1 mapping on interior cells. if
-  // degree==1, then we are already that Q1 mapping, so we don't need
-  // it; if dim!=spacedim, there is also no need for anything because
-  // we're most likely on a curved manifold
-  use_mapping_q_on_all_cells(degree == 1 || use_mapping_q_on_all_cells ||
-                             (dim != spacedim))
-  ,
-  // create a Q1 mapping for use on interior cells (if necessary)
-  // or to create a good initial guess in transform_real_to_unit_cell()
-  q1_mapping(std::make_shared<MappingQGeneric<dim, spacedim>>(1))
-  ,
-
-  // create a Q_p mapping; if p=1, simply share the Q_1 mapping already
-  // created via the shared_ptr objects
-  qp_mapping(this->polynomial_degree > 1 ?
-               std::make_shared<const MappingQGeneric<dim, spacedim>>(degree) :
-               q1_mapping)
+MappingQ<dim, spacedim>::MappingQ(const unsigned int degree, const bool)
+  : MappingQGeneric<dim, spacedim>(degree)
 {}
 
 
 
 template <int dim, int spacedim>
 MappingQ<dim, spacedim>::MappingQ(const MappingQ<dim, spacedim> &mapping)
-  : polynomial_degree(mapping.polynomial_degree)
-  , use_mapping_q_on_all_cells(mapping.use_mapping_q_on_all_cells)
-{
-  // Note that we really do have to use clone() here, since mapping.q1_mapping
-  // may be MappingQ1Eulerian and mapping.qp_mapping may be MappingQEulerian.
-  std::shared_ptr<const Mapping<dim, spacedim>> other_q1_map =
-    mapping.q1_mapping->clone();
-  q1_mapping = std::dynamic_pointer_cast<const MappingQGeneric<dim, spacedim>>(
-    other_q1_map);
-  Assert(q1_mapping != nullptr, ExcInternalError());
-  Assert(q1_mapping->get_degree() == 1, ExcInternalError());
-
-  // Same as the other constructor: if possible reuse the Q1 mapping
-  if (this->polynomial_degree == 1)
-    {
-      qp_mapping = q1_mapping;
-    }
-  else
-    {
-      std::shared_ptr<const Mapping<dim, spacedim>> other_qp_map =
-        mapping.qp_mapping->clone();
-      qp_mapping =
-        std::dynamic_pointer_cast<const MappingQGeneric<dim, spacedim>>(
-          other_qp_map);
-      Assert(qp_mapping != nullptr, ExcInternalError());
-    }
-}
-
-
-
-template <int dim, int spacedim>
-unsigned int
-MappingQ<dim, spacedim>::get_degree() const
-{
-  return polynomial_degree;
-}
-
-
-
-template <int dim, int spacedim>
-inline bool
-MappingQ<dim, spacedim>::preserves_vertex_locations() const
-{
-  return true;
-}
-
-
-
-template <int dim, int spacedim>
-UpdateFlags
-MappingQ<dim, spacedim>::requires_update_flags(const UpdateFlags in) const
-{
-  return (q1_mapping->requires_update_flags(in) |
-          qp_mapping->requires_update_flags(in));
-}
-
-
-
-template <int dim, int spacedim>
-std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
-MappingQ<dim, spacedim>::get_data(const UpdateFlags      update_flags,
-                                  const Quadrature<dim> &quadrature) const
-{
-  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
-    std::make_unique<InternalData>();
-  auto &data = dynamic_cast<InternalData &>(*data_ptr);
-
-  // build the Q1 and Qp internal data objects in parallel
-  Threads::Task<
-    std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>>
-    do_get_data = Threads::new_task(&MappingQGeneric<dim, spacedim>::get_data,
-                                    *qp_mapping,
-                                    update_flags,
-                                    quadrature);
-
-  if (!use_mapping_q_on_all_cells)
-    data.mapping_q1_data = Utilities::dynamic_unique_cast<
-      typename MappingQGeneric<dim, spacedim>::InternalData>(
-      std::move(q1_mapping->get_data(update_flags, quadrature)));
-
-  // wait for the task above to finish and use returned value
-  data.mapping_qp_data = Utilities::dynamic_unique_cast<
-    typename MappingQGeneric<dim, spacedim>::InternalData>(
-    std::move(do_get_data.return_value()));
-  return data_ptr;
-}
-
-
-
-template <int dim, int spacedim>
-std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
-MappingQ<dim, spacedim>::get_face_data(
-  const UpdateFlags               update_flags,
-  const hp::QCollection<dim - 1> &quadrature) const
-{
-  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
-    std::make_unique<InternalData>();
-  auto &data = dynamic_cast<InternalData &>(*data_ptr);
-
-  std::unique_ptr<typename MappingQGeneric<dim, spacedim>::InternalDataBase> (
-    MappingQGeneric<dim, spacedim>::*mapping_get_face_data)(
-    const UpdateFlags, const hp::QCollection<dim - 1> &) const =
-    &MappingQGeneric<dim, spacedim>::get_face_data;
-
-  // build the Q1 and Qp internal data objects in parallel
-  Threads::Task<
-    std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>>
-    do_get_data = Threads::new_task(mapping_get_face_data,
-                                    *qp_mapping,
-                                    update_flags,
-                                    quadrature);
-
-  if (!use_mapping_q_on_all_cells)
-    data.mapping_q1_data = Utilities::dynamic_unique_cast<
-      typename MappingQGeneric<dim, spacedim>::InternalData>(
-      std::move(q1_mapping->get_face_data(update_flags, quadrature)));
-
-  // wait for the task above to finish and use returned value
-  data.mapping_qp_data = Utilities::dynamic_unique_cast<
-    typename MappingQGeneric<dim, spacedim>::InternalData>(
-    std::move(do_get_data.return_value()));
-  return data_ptr;
-}
-
-
-
-template <int dim, int spacedim>
-std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
-MappingQ<dim, spacedim>::get_subface_data(
-  const UpdateFlags          update_flags,
-  const Quadrature<dim - 1> &quadrature) const
-{
-  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
-    std::make_unique<InternalData>();
-  auto &data = dynamic_cast<InternalData &>(*data_ptr);
-
-  // build the Q1 and Qp internal data objects in parallel
-  Threads::Task<
-    std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>>
-    do_get_data =
-      Threads::new_task(&MappingQGeneric<dim, spacedim>::get_subface_data,
-                        *qp_mapping,
-                        update_flags,
-                        quadrature);
-
-  if (!use_mapping_q_on_all_cells)
-    data.mapping_q1_data = Utilities::dynamic_unique_cast<
-      typename MappingQGeneric<dim, spacedim>::InternalData>(
-      std::move(q1_mapping->get_subface_data(update_flags, quadrature)));
-
-  // wait for the task above to finish and use returned value
-  data.mapping_qp_data = Utilities::dynamic_unique_cast<
-    typename MappingQGeneric<dim, spacedim>::InternalData>(
-    std::move(do_get_data.return_value()));
-  return data_ptr;
-}
-
-
-// Note that the CellSimilarity flag is modifiable, since MappingQ can need to
-// recalculate data even when cells are similar.
-template <int dim, int spacedim>
-CellSimilarity::Similarity
-MappingQ<dim, spacedim>::fill_fe_values(
-  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-  const CellSimilarity::Similarity                            cell_similarity,
-  const Quadrature<dim> &                                     quadrature,
-  const typename Mapping<dim, spacedim>::InternalDataBase &   internal_data,
-  internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
-    &output_data) const
-{
-  // convert data object to internal data for this class. fails with an
-  // exception if that is not possible
-  Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
-         ExcInternalError());
-  const InternalData &data = static_cast<const InternalData &>(internal_data);
-
-  // check whether this cell needs the full mapping or can be treated by a
-  // reduced Q1 mapping, e.g. if the cell is in the interior of the domain
-  data.use_mapping_q1_on_current_cell =
-    !(use_mapping_q_on_all_cells || cell->has_boundary_lines());
-
-
-  // call the base class. we need to ensure that the flag indicating whether
-  // we can use some similarity has to be modified - for a general MappingQ,
-  // the data needs to be recomputed anyway since then the mapping changes the
-  // data. this needs to be known also for later operations, so modify the
-  // variable here. this also affects the calculation of the next cell -- if
-  // we use Q1 data on the next cell, the data will still be invalid.
-  const CellSimilarity::Similarity updated_cell_similarity =
-    ((data.use_mapping_q1_on_current_cell == false) &&
-         (this->polynomial_degree > 1) ?
-       CellSimilarity::invalid_next_cell :
-       cell_similarity);
-
-  // depending on the results above, decide whether the Q1 mapping or
-  // the Qp mapping needs to handle this cell
-  if (data.use_mapping_q1_on_current_cell)
-    q1_mapping->fill_fe_values(cell,
-                               updated_cell_similarity,
-                               quadrature,
-                               *data.mapping_q1_data,
-                               output_data);
-  else
-    qp_mapping->fill_fe_values(cell,
-                               updated_cell_similarity,
-                               quadrature,
-                               *data.mapping_qp_data,
-                               output_data);
-
-  return updated_cell_similarity;
-}
-
-
-
-template <int dim, int spacedim>
-void
-MappingQ<dim, spacedim>::fill_fe_face_values(
-  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-  const unsigned int                                          face_no,
-  const hp::QCollection<dim - 1> &                            quadrature,
-  const typename Mapping<dim, spacedim>::InternalDataBase &   internal_data,
-  internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
-    &output_data) const
-{
-  // convert data object to internal data for this class. fails with an
-  // exception if that is not possible
-  Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
-         ExcInternalError());
-  const InternalData &data = static_cast<const InternalData &>(internal_data);
-
-  // check whether this cell needs the full mapping or can be treated by a
-  // reduced Q1 mapping, e.g. if the cell is entirely in the interior of the
-  // domain. note that it is not sufficient to ask whether the present _face_
-  // is in the interior, as the mapping on the face depends on the mapping of
-  // the cell, which in turn depends on the fact whether _any_ of the faces of
-  // this cell is at the boundary, not only the present face
-  data.use_mapping_q1_on_current_cell =
-    !(use_mapping_q_on_all_cells || cell->has_boundary_lines());
-
-  // depending on the results above, decide whether the Q1 mapping or
-  // the Qp mapping needs to handle this cell
-  if (data.use_mapping_q1_on_current_cell)
-    q1_mapping->fill_fe_face_values(
-      cell, face_no, quadrature, *data.mapping_q1_data, output_data);
-  else
-    qp_mapping->fill_fe_face_values(
-      cell, face_no, quadrature, *data.mapping_qp_data, output_data);
-}
-
-
-template <int dim, int spacedim>
-void
-MappingQ<dim, spacedim>::fill_fe_subface_values(
-  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-  const unsigned int                                          face_no,
-  const unsigned int                                          subface_no,
-  const Quadrature<dim - 1> &                                 quadrature,
-  const typename Mapping<dim, spacedim>::InternalDataBase &   internal_data,
-  internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
-    &output_data) const
-{
-  // convert data object to internal data for this class. fails with an
-  // exception if that is not possible
-  Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
-         ExcInternalError());
-  const InternalData &data = static_cast<const InternalData &>(internal_data);
-
-  // check whether this cell needs the full mapping or can be treated by a
-  // reduced Q1 mapping, e.g. if the cell is entirely in the interior of the
-  // domain. note that it is not sufficient to ask whether the present _face_
-  // is in the interior, as the mapping on the face depends on the mapping of
-  // the cell, which in turn depends on the fact whether _any_ of the faces of
-  // this cell is at the boundary, not only the present face
-  data.use_mapping_q1_on_current_cell =
-    !(use_mapping_q_on_all_cells || cell->has_boundary_lines());
-
-  // depending on the results above, decide whether the Q1 mapping or
-  // the Qp mapping needs to handle this cell
-  if (data.use_mapping_q1_on_current_cell)
-    q1_mapping->fill_fe_subface_values(cell,
-                                       face_no,
-                                       subface_no,
-                                       quadrature,
-                                       *data.mapping_q1_data,
-                                       output_data);
-  else
-    qp_mapping->fill_fe_subface_values(cell,
-                                       face_no,
-                                       subface_no,
-                                       quadrature,
-                                       *data.mapping_qp_data,
-                                       output_data);
-}
-
-
-
-template <int dim, int spacedim>
-void
-MappingQ<dim, spacedim>::transform(
-  const ArrayView<const Tensor<1, dim>> &                  input,
-  const MappingKind                                        mapping_kind,
-  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-  const ArrayView<Tensor<1, spacedim>> &                   output) const
-{
-  AssertDimension(input.size(), output.size());
-
-  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
-  Assert(data != nullptr, ExcInternalError());
-
-  // check whether we should in fact work on the Q1 portion of it
-  if (data->use_mapping_q1_on_current_cell)
-    q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
-  else
-    // otherwise use the full mapping
-    qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
-}
-
-
-
-template <int dim, int spacedim>
-void
-MappingQ<dim, spacedim>::transform(
-  const ArrayView<const DerivativeForm<1, dim, spacedim>> &input,
-  const MappingKind                                        mapping_kind,
-  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-  const ArrayView<Tensor<2, spacedim>> &                   output) const
-{
-  AssertDimension(input.size(), output.size());
-  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
-            &mapping_data) != nullptr),
-         ExcInternalError());
-  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
-
-  // check whether we should in fact work on the Q1 portion of it
-  if (data->use_mapping_q1_on_current_cell)
-    q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
-  else
-    // otherwise use the full mapping
-    qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
-}
-
-
-
-template <int dim, int spacedim>
-void
-MappingQ<dim, spacedim>::transform(
-  const ArrayView<const Tensor<2, dim>> &                  input,
-  const MappingKind                                        mapping_kind,
-  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-  const ArrayView<Tensor<2, spacedim>> &                   output) const
-{
-  AssertDimension(input.size(), output.size());
-  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
-            &mapping_data) != nullptr),
-         ExcInternalError());
-  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
-
-  // check whether we should in fact work on the Q1 portion of it
-  if (data->use_mapping_q1_on_current_cell)
-    q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
-  else
-    // otherwise use the full mapping
-    qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
-}
-
-
-
-template <int dim, int spacedim>
-void
-MappingQ<dim, spacedim>::transform(
-  const ArrayView<const DerivativeForm<2, dim, spacedim>> &input,
-  const MappingKind                                        mapping_kind,
-  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-  const ArrayView<Tensor<3, spacedim>> &                   output) const
-{
-  AssertDimension(input.size(), output.size());
-  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
-            &mapping_data) != nullptr),
-         ExcInternalError());
-  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
-
-  // check whether we should in fact work on the Q1 portion of it
-  if (data->use_mapping_q1_on_current_cell)
-    q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
-  else
-    // otherwise use the full mapping
-    qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
-}
-
-
-
-template <int dim, int spacedim>
-void
-MappingQ<dim, spacedim>::transform(
-  const ArrayView<const Tensor<3, dim>> &                  input,
-  const MappingKind                                        mapping_kind,
-  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-  const ArrayView<Tensor<3, spacedim>> &                   output) const
-{
-  AssertDimension(input.size(), output.size());
-  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
-            &mapping_data) != nullptr),
-         ExcInternalError());
-  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
-
-  // check whether we should in fact work on the Q1 portion of it
-  if (data->use_mapping_q1_on_current_cell)
-    q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
-  else
-    // otherwise use the full mapping
-    qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
-}
-
-
-
-template <int dim, int spacedim>
-Point<spacedim>
-MappingQ<dim, spacedim>::transform_unit_to_real_cell(
-  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-  const Point<dim> &                                          p) const
-{
-  // first see, whether we want to use a linear or a higher order
-  // mapping, then either use our own facilities or that of the Q1
-  // mapping we store
-  if (use_mapping_q_on_all_cells || cell->has_boundary_lines())
-    return qp_mapping->transform_unit_to_real_cell(cell, p);
-  else
-    return q1_mapping->transform_unit_to_real_cell(cell, p);
-}
-
-
-
-template <int dim, int spacedim>
-Point<dim>
-MappingQ<dim, spacedim>::transform_real_to_unit_cell(
-  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
-  const Point<spacedim> &                                     p) const
-{
-  if (cell->has_boundary_lines() || use_mapping_q_on_all_cells ||
-      (dim != spacedim))
-    return qp_mapping->transform_real_to_unit_cell(cell, p);
-  else
-    return q1_mapping->transform_real_to_unit_cell(cell, p);
-}
-
-
-
-template <int dim, int spacedim>
-std::unique_ptr<Mapping<dim, spacedim>>
-MappingQ<dim, spacedim>::clone() const
-{
-  return std::make_unique<MappingQ<dim, spacedim>>(
-    this->polynomial_degree, this->use_mapping_q_on_all_cells);
-}
-
-
-
-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));
-}
-
-
-
-template <int dim, int spacedim>
-bool
-MappingQ<dim, spacedim>::is_compatible_with(
-  const ReferenceCell &reference_cell) const
-{
-  Assert(dim == reference_cell.get_dimension(),
-         ExcMessage("The dimension of your mapping (" +
-                    Utilities::to_string(dim) +
-                    ") and the reference cell cell_type (" +
-                    Utilities::to_string(reference_cell.get_dimension()) +
-                    " ) do not agree."));
-
-  return reference_cell.is_hyper_cube();
-}
-
+  : MappingQGeneric<dim, spacedim>(mapping)
+{}
 
 
 // explicit instantiations
index 98675809de22f96b53f7a43ace1462156c46d1ee..91021526780a25dc7e73a7b7b477bb6ae8cdd20d 100644 (file)
@@ -44,20 +44,6 @@ DEAL_II_NAMESPACE_OPEN
 
 // .... MAPPING Q EULERIAN CONSTRUCTOR
 
-template <int dim, class VectorType, int spacedim>
-MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
-  MappingQEulerianGeneric(
-    const unsigned int                                 degree,
-    const MappingQEulerian<dim, VectorType, spacedim> &mapping_q_eulerian)
-  : MappingQGeneric<dim, spacedim>(degree)
-  , mapping_q_eulerian(mapping_q_eulerian)
-  , support_quadrature(degree)
-  , fe_values(mapping_q_eulerian.euler_dof_handler->get_fe(),
-              support_quadrature,
-              update_values | update_quadrature_points)
-{}
-
-
 
 template <int dim, class VectorType, int spacedim>
 MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerian(
@@ -65,21 +51,15 @@ MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerian(
   const DoFHandler<dim, spacedim> &euler_dof_handler,
   const VectorType &               euler_vector,
   const unsigned int               level)
-  : MappingQ<dim, spacedim>(degree, true)
+  : MappingQGeneric<dim, spacedim>(degree)
   , euler_vector(&euler_vector)
   , euler_dof_handler(&euler_dof_handler)
   , level(level)
-{
-  // reset the q1 mapping we use for interior cells (and previously
-  // set by the MappingQ constructor) to a MappingQ1Eulerian with the
-  // current vector
-  this->q1_mapping =
-    std::make_shared<MappingQ1Eulerian<dim, VectorType, spacedim>>(
-      euler_dof_handler, euler_vector);
-
-  // also reset the qp mapping pointer with our own class
-  this->qp_mapping = std::make_shared<MappingQEulerianGeneric>(degree, *this);
-}
+  , support_quadrature(degree)
+  , fe_values(euler_dof_handler.get_fe(),
+              support_quadrature,
+              update_values | update_quadrature_points)
+{}
 
 
 
@@ -96,8 +76,8 @@ MappingQEulerian<dim, VectorType, spacedim>::clone() const
 // .... SUPPORT QUADRATURE CONSTRUCTOR
 
 template <int dim, class VectorType, int spacedim>
-MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
-  SupportQuadrature::SupportQuadrature(const unsigned int map_degree)
+MappingQEulerian<dim, VectorType, spacedim>::SupportQuadrature::
+  SupportQuadrature(const unsigned int map_degree)
   : Quadrature<dim>(Utilities::fixed_power<dim>(map_degree + 1))
 {
   // first we determine the support points on the unit cell in lexicographic
@@ -127,9 +107,7 @@ MappingQEulerian<dim, VectorType, spacedim>::get_vertices(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
 {
   // get the vertices as the first 2^dim mapping support points
-  const std::vector<Point<spacedim>> a =
-    dynamic_cast<const MappingQEulerianGeneric &>(*this->qp_mapping)
-      .compute_mapping_support_points(cell);
+  const std::vector<Point<spacedim>> a = compute_mapping_support_points(cell);
 
   boost::container::small_vector<Point<spacedim>,
                                  GeometryInfo<dim>::vertices_per_cell>
@@ -141,43 +119,16 @@ MappingQEulerian<dim, VectorType, spacedim>::get_vertices(
 
 
 
-template <int dim, class VectorType, int spacedim>
-boost::container::small_vector<Point<spacedim>,
-                               GeometryInfo<dim>::vertices_per_cell>
-MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
-  get_vertices(
-    const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
-{
-  return mapping_q_eulerian.get_vertices(cell);
-}
-
-
-
-template <int dim, class VectorType, int spacedim>
-bool
-MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
-  preserves_vertex_locations() const
-{
-  return false;
-}
-
-
-
 template <int dim, class VectorType, int spacedim>
 std::vector<Point<spacedim>>
-MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
-  compute_mapping_support_points(
-    const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
+MappingQEulerian<dim, VectorType, spacedim>::compute_mapping_support_points(
+  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
 {
-  const bool mg_vector =
-    mapping_q_eulerian.level != numbers::invalid_unsigned_int;
+  const bool mg_vector = level != numbers::invalid_unsigned_int;
 
   const types::global_dof_index n_dofs =
-    mg_vector ?
-      mapping_q_eulerian.euler_dof_handler->n_dofs(mapping_q_eulerian.level) :
-      mapping_q_eulerian.euler_dof_handler->n_dofs();
-  const types::global_dof_index vector_size =
-    mapping_q_eulerian.euler_vector->size();
+    mg_vector ? euler_dof_handler->n_dofs(level) : euler_dof_handler->n_dofs();
+  const types::global_dof_index vector_size = euler_vector->size();
 
   (void)n_dofs;
   (void)vector_size;
@@ -186,8 +137,8 @@ MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
 
   // we then transform our tria iterator into a dof iterator so we can access
   // data not associated with triangulations
-  typename DoFHandler<dim, spacedim>::cell_iterator dof_cell(
-    *cell, mapping_q_eulerian.euler_dof_handler);
+  typename DoFHandler<dim, spacedim>::cell_iterator dof_cell(*cell,
+                                                             euler_dof_handler);
 
   Assert(mg_vector || dof_cell->is_active() == true, ExcInactiveCell());
 
@@ -203,8 +154,7 @@ MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
   // that. This implies that the user should order components appropriately,
   // or create a separate dof handler for the displacements.
   const unsigned int n_support_pts = support_quadrature.size();
-  const unsigned int n_components =
-    mapping_q_eulerian.euler_dof_handler->get_fe(0).n_components();
+  const unsigned int n_components = euler_dof_handler->get_fe(0).n_components();
 
   Assert(n_components >= spacedim,
          ExcDimensionMismatch(n_components, spacedim));
@@ -213,7 +163,7 @@ MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
     n_support_pts, Vector<typename VectorType::value_type>(n_components));
 
   std::vector<types::global_dof_index> dof_indices(
-    mapping_q_eulerian.euler_dof_handler->get_fe(0).n_dofs_per_cell());
+    euler_dof_handler->get_fe(0).n_dofs_per_cell());
   // fill shift vector for each support point using an fe_values object. make
   // sure that the fe_values variable isn't used simultaneously from different
   // threads
@@ -222,13 +172,10 @@ MappingQEulerian<dim, VectorType, spacedim>::MappingQEulerianGeneric::
   if (mg_vector)
     {
       dof_cell->get_mg_dof_indices(dof_indices);
-      fe_values.get_function_values(*mapping_q_eulerian.euler_vector,
-                                    dof_indices,
-                                    shift_vector);
+      fe_values.get_function_values(*euler_vector, dof_indices, shift_vector);
     }
   else
-    fe_values.get_function_values(*mapping_q_eulerian.euler_vector,
-                                  shift_vector);
+    fe_values.get_function_values(*euler_vector, shift_vector);
 
   // and finally compute the positions of the support points in the deformed
   // configuration.
@@ -258,11 +205,12 @@ MappingQEulerian<dim, VectorType, spacedim>::fill_fe_values(
   // call the function of the base class, but ignoring
   // any potentially detected cell similarity between
   // the current and the previous cell
-  MappingQ<dim, spacedim>::fill_fe_values(cell,
-                                          CellSimilarity::invalid_next_cell,
-                                          quadrature,
-                                          internal_data,
-                                          output_data);
+  MappingQGeneric<dim, spacedim>::fill_fe_values(
+    cell,
+    CellSimilarity::invalid_next_cell,
+    quadrature,
+    internal_data,
+    output_data);
   // also return the updated flag since any detected
   // similarity wasn't based on the mapped field, but
   // the original vertices which are meaningless
@@ -270,19 +218,6 @@ 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"
 

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.