]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Intoduce ReferenceCell in QProjector 10720/head
authorPeter Munch <peterrmuench@gmail.com>
Sat, 18 Jul 2020 20:20:46 +0000 (22:20 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Tue, 21 Jul 2020 08:54:28 +0000 (10:54 +0200)
30 files changed:
include/deal.II/base/qprojector.h
include/deal.II/fe/fe_poly.templates.h
include/deal.II/fe/fe_poly_face.h
include/deal.II/fe/fe_tools.templates.h
include/deal.II/grid/reference_cell.h
source/base/CMakeLists.txt
source/base/qprojector.cc [new file with mode: 0644]
source/base/quadrature.cc
source/fe/fe.cc
source/fe/fe_abf.cc
source/fe/fe_bdm.cc
source/fe/fe_bernstein.cc
source/fe/fe_nedelec.cc
source/fe/fe_nedelec_sz.cc
source/fe/fe_p1nc.cc
source/fe/fe_poly_tensor.cc
source/fe/fe_q_base.cc
source/fe/fe_q_bubbles.cc
source/fe/fe_q_dg0.cc
source/fe/fe_raviart_thomas.cc
source/fe/fe_raviart_thomas_nodal.cc
source/fe/fe_rt_bubbles.cc
source/fe/fe_system.cc
source/fe/fe_trace.cc
source/fe/mapping_cartesian.cc
source/fe/mapping_fe_field.cc
source/fe/mapping_manifold.cc
source/fe/mapping_q1.cc
source/fe/mapping_q_generic.cc
source/grid/grid_out.cc

index af3818566eff6d8382219ff104d83b2b0e6c584c..7e231242f7c649e54b90d5753097ba7d827d0ff2 100644 (file)
@@ -22,6 +22,8 @@
 #include <deal.II/base/geometry_info.h>
 #include <deal.II/base/quadrature.h>
 
+#include <deal.II/grid/reference_cell.h>
+
 DEAL_II_NAMESPACE_OPEN
 
 
@@ -86,20 +88,49 @@ public:
    * Compute the quadrature points on the cell if the given quadrature formula
    * is used on face <tt>face_no</tt>. For further details, see the general
    * doc for this class.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static void
+  DEAL_II_DEPRECATED static void
   project_to_face(const SubQuadrature &    quadrature,
                   const unsigned int       face_no,
                   std::vector<Point<dim>> &q_points);
 
+  /**
+   * Compute the quadrature points on the cell if the given quadrature formula
+   * is used on face <tt>face_no</tt>. For further details, see the general
+   * doc for this class.
+   */
+  static void
+  project_to_face(const ReferenceCell::Type reference_cell_type,
+                  const SubQuadrature &     quadrature,
+                  const unsigned int        face_no,
+                  std::vector<Point<dim>> & q_points);
+
   /**
    * Compute the cell quadrature formula corresponding to using
    * <tt>quadrature</tt> on face <tt>face_no</tt>. For further details, see
    * the general doc for this class.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static Quadrature<dim>
+  DEAL_II_DEPRECATED static Quadrature<dim>
   project_to_face(const SubQuadrature &quadrature, const unsigned int face_no);
 
+  /**
+   * Compute the cell quadrature formula corresponding to using
+   * <tt>quadrature</tt> on face <tt>face_no</tt>. For further details, see
+   * the general doc for this class.
+   */
+  static Quadrature<dim>
+  project_to_face(const ReferenceCell::Type reference_cell_type,
+                  const SubQuadrature &     quadrature,
+                  const unsigned int        face_no);
+
   /**
    * Compute the quadrature points on the cell if the given quadrature formula
    * is used on face <tt>face_no</tt>, subface number <tt>subface_no</tt>
@@ -108,8 +139,12 @@ public:
    *
    * @note Only the points are transformed. The quadrature weights are the
    * same as those of the original rule.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static void
+  DEAL_II_DEPRECATED static void
   project_to_subface(const SubQuadrature &          quadrature,
                      const unsigned int             face_no,
                      const unsigned int             subface_no,
@@ -117,6 +152,24 @@ public:
                      const RefinementCase<dim - 1> &ref_case =
                        RefinementCase<dim - 1>::isotropic_refinement);
 
+  /**
+   * Compute the quadrature points on the cell if the given quadrature formula
+   * is used on face <tt>face_no</tt>, subface number <tt>subface_no</tt>
+   * corresponding to RefineCase::Type <tt>ref_case</tt>. The last argument is
+   * only used in 3D.
+   *
+   * @note Only the points are transformed. The quadrature weights are the
+   * same as those of the original rule.
+   */
+  static void
+  project_to_subface(const ReferenceCell::Type      reference_cell_type,
+                     const SubQuadrature &          quadrature,
+                     const unsigned int             face_no,
+                     const unsigned int             subface_no,
+                     std::vector<Point<dim>> &      q_points,
+                     const RefinementCase<dim - 1> &ref_case =
+                       RefinementCase<dim - 1>::isotropic_refinement);
+
   /**
    * Compute the cell quadrature formula corresponding to using
    * <tt>quadrature</tt> on subface <tt>subface_no</tt> of face
@@ -125,14 +178,39 @@ public:
    *
    * @note Only the points are transformed. The quadrature weights are the
    * same as those of the original rule.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static Quadrature<dim>
+  DEAL_II_DEPRECATED static Quadrature<dim>
   project_to_subface(const SubQuadrature &          quadrature,
                      const unsigned int             face_no,
                      const unsigned int             subface_no,
                      const RefinementCase<dim - 1> &ref_case =
                        RefinementCase<dim - 1>::isotropic_refinement);
 
+  /**
+   * Compute the cell quadrature formula corresponding to using
+   * <tt>quadrature</tt> on subface <tt>subface_no</tt> of face
+   * <tt>face_no</tt> with RefinementCase<dim-1> <tt>ref_case</tt>. The last
+   * argument is only used in 3D.
+   *
+   * @note Only the points are transformed. The quadrature weights are the
+   * same as those of the original rule.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
+   */
+  static Quadrature<dim>
+  project_to_subface(const ReferenceCell::Type      reference_cell_type,
+                     const SubQuadrature &          quadrature,
+                     const unsigned int             face_no,
+                     const unsigned int             subface_no,
+                     const RefinementCase<dim - 1> &ref_case =
+                       RefinementCase<dim - 1>::isotropic_refinement);
+
   /**
    * Take a face quadrature formula and generate a cell quadrature formula
    * from it where the quadrature points of the given argument are projected
@@ -148,10 +226,34 @@ public:
    *
    * @note In 3D, this function produces eight sets of quadrature points for
    * each face, in order to cope possibly different orientations of the mesh.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static Quadrature<dim>
+  DEAL_II_DEPRECATED static Quadrature<dim>
   project_to_all_faces(const SubQuadrature &quadrature);
 
+  /**
+   * Take a face quadrature formula and generate a cell quadrature formula
+   * from it where the quadrature points of the given argument are projected
+   * on all faces.
+   *
+   * The weights of the new rule are replications of the original weights.
+   * Thus, the sum of the weights is not one, but the number of faces, which
+   * is the surface of the reference cell.
+   *
+   * This in particular allows us to extract a subset of points corresponding
+   * to a single face and use it as a quadrature on this face, as is done in
+   * FEFaceValues.
+   *
+   * @note In 3D, this function produces eight sets of quadrature points for
+   * each face, in order to cope possibly different orientations of the mesh.
+   */
+  static Quadrature<dim>
+  project_to_all_faces(const ReferenceCell::Type reference_cell_type,
+                       const SubQuadrature &     quadrature);
+
   /**
    * Take a face quadrature formula and generate a cell quadrature formula
    * from it where the quadrature points of the given argument are projected
@@ -164,10 +266,31 @@ public:
    * This in particular allows us to extract a subset of points corresponding
    * to a single subface and use it as a quadrature on this face, as is done
    * in FESubfaceValues.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static Quadrature<dim>
+  DEAL_II_DEPRECATED static Quadrature<dim>
   project_to_all_subfaces(const SubQuadrature &quadrature);
 
+  /**
+   * Take a face quadrature formula and generate a cell quadrature formula
+   * from it where the quadrature points of the given argument are projected
+   * on all subfaces.
+   *
+   * Like in project_to_all_faces(), the weights of the new rule sum up to the
+   * number of faces (not subfaces), which is the surface of the reference
+   * cell.
+   *
+   * This in particular allows us to extract a subset of points corresponding
+   * to a single subface and use it as a quadrature on this face, as is done
+   * in FESubfaceValues.
+   */
+  static Quadrature<dim>
+  project_to_all_subfaces(const ReferenceCell::Type reference_cell_type,
+                          const SubQuadrature &     quadrature);
+
   /**
    * Project a given quadrature formula to a child of a cell. You may want to
    * use this function in case you want to extend an integral only over the
@@ -177,11 +300,30 @@ public:
    * As integration using this quadrature formula now only extends over a
    * fraction of the cell, the weights of the resulting object are divided by
    * GeometryInfo<dim>::children_per_cell.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static Quadrature<dim>
+  DEAL_II_DEPRECATED static Quadrature<dim>
   project_to_child(const Quadrature<dim> &quadrature,
                    const unsigned int     child_no);
 
+  /**
+   * Project a given quadrature formula to a child of a cell. You may want to
+   * use this function in case you want to extend an integral only over the
+   * area which a potential child would occupy. The child numbering is the
+   * same as the children would be numbered upon refinement of the cell.
+   *
+   * As integration using this quadrature formula now only extends over a
+   * fraction of the cell, the weights of the resulting object are divided by
+   * GeometryInfo<dim>::children_per_cell.
+   */
+  static Quadrature<dim>
+  project_to_child(const ReferenceCell::Type reference_cell_type,
+                   const Quadrature<dim> &   quadrature,
+                   const unsigned int        child_no);
+
   /**
    * Project a quadrature rule to all children of a cell. Similarly to
    * project_to_all_subfaces(), this function replicates the formula generated
@@ -190,19 +332,50 @@ public:
    *
    * The child numbering is the same as the children would be numbered upon
    * refinement of the cell.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static Quadrature<dim>
+  DEAL_II_DEPRECATED static Quadrature<dim>
   project_to_all_children(const Quadrature<dim> &quadrature);
 
+  /**
+   * Project a quadrature rule to all children of a cell. Similarly to
+   * project_to_all_subfaces(), this function replicates the formula generated
+   * by project_to_child() for all children, such that the weights sum up to
+   * one, the volume of the total cell again.
+   *
+   * The child numbering is the same as the children would be numbered upon
+   * refinement of the cell.
+   */
+  static Quadrature<dim>
+  project_to_all_children(const ReferenceCell::Type reference_cell_type,
+                          const Quadrature<dim> &   quadrature);
+
   /**
    * Project the one dimensional rule <tt>quadrature</tt> to the straight line
    * connecting the points <tt>p1</tt> and <tt>p2</tt>.
+   *
+   * @deprecated This function makes an implicit assumption that the cell is
+   *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+   *   this function that takes the reference cell type instead.
    */
-  static Quadrature<dim>
+  DEAL_II_DEPRECATED static Quadrature<dim>
   project_to_line(const Quadrature<1> &quadrature,
                   const Point<dim> &   p1,
                   const Point<dim> &   p2);
 
+  /**
+   * Project the one dimensional rule <tt>quadrature</tt> to the straight line
+   * connecting the points <tt>p1</tt> and <tt>p2</tt>.
+   */
+  static Quadrature<dim>
+  project_to_line(const ReferenceCell::Type reference_cell_type,
+                  const Quadrature<1> &     quadrature,
+                  const Point<dim> &        p1,
+                  const Point<dim> &        p2);
+
   /**
    * Since the project_to_all_faces() and project_to_all_subfaces() functions
    * chain together the quadrature points and weights of all projections of a
@@ -241,14 +414,36 @@ public:
      * The last argument denotes the number of quadrature points the lower-
      * dimensional face quadrature formula (the one that has been projected
      * onto the faces) has.
+     *
+     * @deprecated This function makes an implicit assumption that the cell is
+     *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+     *   this function that takes the reference cell type instead.
      */
-    static DataSetDescriptor
+    DEAL_II_DEPRECATED static DataSetDescriptor
     face(const unsigned int face_no,
          const bool         face_orientation,
          const bool         face_flip,
          const bool         face_rotation,
          const unsigned int n_quadrature_points);
 
+    /**
+     * Static function to generate an offset object for a given face of a cell
+     * with the given face orientation, flip and rotation. This function of
+     * course is only allowed if <tt>dim>=2</tt>, and the face orientation,
+     * flip and rotation are ignored if the space dimension equals 2.
+     *
+     * The last argument denotes the number of quadrature points the lower-
+     * dimensional face quadrature formula (the one that has been projected
+     * onto the faces) has.
+     */
+    static DataSetDescriptor
+    face(const ReferenceCell::Type reference_cell_type,
+         const unsigned int        face_no,
+         const bool                face_orientation,
+         const bool                face_flip,
+         const bool                face_rotation,
+         const unsigned int        n_quadrature_points);
+
     /**
      * Static function to generate an offset object for a given subface of a
      * cell with the given face orientation, flip and rotation. This function
@@ -260,8 +455,12 @@ public:
      * projected onto the faces) has.
      *
      * Through the last argument anisotropic refinement can be respected.
+     *
+     * @deprecated This function makes an implicit assumption that the cell is
+     *   a line (1D), a quad (2D), or a hex (3D). Use the other version of
+     *   this function that takes the reference cell type instead.
      */
-    static DataSetDescriptor
+    DEAL_II_DEPRECATED static DataSetDescriptor
     subface(const unsigned int               face_no,
             const unsigned int               subface_no,
             const bool                       face_orientation,
@@ -271,6 +470,29 @@ public:
             const internal::SubfaceCase<dim> ref_case =
               internal::SubfaceCase<dim>::case_isotropic);
 
+    /**
+     * Static function to generate an offset object for a given subface of a
+     * cell with the given face orientation, flip and rotation. This function
+     * of course is only allowed if <tt>dim>=2</tt>, and the face orientation,
+     * flip and rotation are ignored if the space dimension equals 2.
+     *
+     * The last but one argument denotes the number of quadrature points the
+     * lower-dimensional face quadrature formula (the one that has been
+     * projected onto the faces) has.
+     *
+     * Through the last argument anisotropic refinement can be respected.
+     */
+    static DataSetDescriptor
+    subface(const ReferenceCell::Type        reference_cell_type,
+            const unsigned int               face_no,
+            const unsigned int               subface_no,
+            const bool                       face_orientation,
+            const bool                       face_flip,
+            const bool                       face_rotation,
+            const unsigned int               n_quadrature_points,
+            const internal::SubfaceCase<dim> ref_case =
+              internal::SubfaceCase<dim>::case_isotropic);
+
     /**
      * Conversion operator to an integer denoting the offset of the first
      * element of this dataset in the set of quadrature formulas all projected
@@ -365,18 +587,41 @@ QProjector<1>::project_to_face(const Quadrature<0> &,
                                std::vector<Point<1>> &);
 template <>
 void
+QProjector<1>::project_to_face(const ReferenceCell::Type reference_cell_type,
+                               const Quadrature<0> &,
+                               const unsigned int,
+                               std::vector<Point<1>> &);
+template <>
+void
 QProjector<2>::project_to_face(const Quadrature<1> &  quadrature,
                                const unsigned int     face_no,
                                std::vector<Point<2>> &q_points);
 template <>
 void
+QProjector<2>::project_to_face(const ReferenceCell::Type reference_cell_type,
+                               const Quadrature<1> &     quadrature,
+                               const unsigned int        face_no,
+                               std::vector<Point<2>> &   q_points);
+template <>
+void
 QProjector<3>::project_to_face(const Quadrature<2> &  quadrature,
                                const unsigned int     face_no,
                                std::vector<Point<3>> &q_points);
+template <>
+void
+QProjector<3>::project_to_face(const ReferenceCell::Type reference_cell_type,
+                               const Quadrature<2> &     quadrature,
+                               const unsigned int        face_no,
+                               std::vector<Point<3>> &   q_points);
 
 template <>
 Quadrature<1>
 QProjector<1>::project_to_all_faces(const Quadrature<0> &quadrature);
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_faces(
+  const ReferenceCell::Type reference_cell_type,
+  const Quadrature<0> &     quadrature);
 
 
 template <>
@@ -388,6 +633,14 @@ QProjector<1>::project_to_subface(const Quadrature<0> &,
                                   const RefinementCase<0> &);
 template <>
 void
+QProjector<1>::project_to_subface(const ReferenceCell::Type reference_cell_type,
+                                  const Quadrature<0> &,
+                                  const unsigned int,
+                                  const unsigned int,
+                                  std::vector<Point<1>> &,
+                                  const RefinementCase<0> &);
+template <>
+void
 QProjector<2>::project_to_subface(const Quadrature<1> &  quadrature,
                                   const unsigned int     face_no,
                                   const unsigned int     subface_no,
@@ -395,15 +648,36 @@ QProjector<2>::project_to_subface(const Quadrature<1> &  quadrature,
                                   const RefinementCase<1> &);
 template <>
 void
+QProjector<2>::project_to_subface(const ReferenceCell::Type reference_cell_type,
+                                  const Quadrature<1> &     quadrature,
+                                  const unsigned int        face_no,
+                                  const unsigned int        subface_no,
+                                  std::vector<Point<2>> &   q_points,
+                                  const RefinementCase<1> &);
+template <>
+void
 QProjector<3>::project_to_subface(const Quadrature<2> &    quadrature,
                                   const unsigned int       face_no,
                                   const unsigned int       subface_no,
                                   std::vector<Point<3>> &  q_points,
                                   const RefinementCase<2> &face_ref_case);
+template <>
+void
+QProjector<3>::project_to_subface(const ReferenceCell::Type reference_cell_type,
+                                  const Quadrature<2> &     quadrature,
+                                  const unsigned int        face_no,
+                                  const unsigned int        subface_no,
+                                  std::vector<Point<3>> &   q_points,
+                                  const RefinementCase<2> & face_ref_case);
 
 template <>
 Quadrature<1>
 QProjector<1>::project_to_all_subfaces(const Quadrature<0> &quadrature);
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_subfaces(
+  const ReferenceCell::Type reference_cell_type,
+  const Quadrature<0> &     quadrature);
 
 
 #endif // DOXYGEN
index c4bcb2c100581ae8abc692f0531cba089ec3a614..cf88f802a73cb8db4f5a7d6c23d23fb47dfb6044 100644 (file)
@@ -358,8 +358,9 @@ FE_Poly<dim, spacedim>::fill_fe_face_values(
   // to take (all data sets for all
   // faces are stored contiguously)
 
-  const typename QProjector<dim>::DataSetDescriptor offset =
-    QProjector<dim>::DataSetDescriptor::face(face_no,
+  const auto offset =
+    QProjector<dim>::DataSetDescriptor::face(this->reference_cell_type(),
+                                             face_no,
                                              cell->face_orientation(face_no),
                                              cell->face_flip(face_no),
                                              cell->face_rotation(face_no),
@@ -450,8 +451,9 @@ FE_Poly<dim, spacedim>::fill_fe_subface_values(
   // to take (all data sets for all
   // sub-faces are stored contiguously)
 
-  const typename QProjector<dim>::DataSetDescriptor offset =
-    QProjector<dim>::DataSetDescriptor::subface(face_no,
+  const auto offset =
+    QProjector<dim>::DataSetDescriptor::subface(this->reference_cell_type(),
+                                                face_no,
                                                 sub_no,
                                                 cell->face_orientation(face_no),
                                                 cell->face_flip(face_no),
index d2b3834c33e43c1ce581edc69bf25049774a7be8..cd2be76ff14fe91759186ed4a1eca6ddc6aaae13 100644 (file)
@@ -167,7 +167,7 @@ protected:
     return get_face_data(update_flags,
                          mapping,
                          QProjector<dim - 1>::project_to_all_children(
-                           quadrature),
+                           ReferenceCell::get_hypercube(dim - 1), quadrature),
                          output_data);
   }
 
index 02f472e41d55ae98bf29ff4af63fa7c67ef73cc0..34c7d92789867daab53bbec0e2aa2ebc338bbc39 100644 (file)
@@ -2035,7 +2035,9 @@ namespace FETools
     // hating the anisotropic implementation
     QGauss<dim - 1>       q_gauss(degree + 1);
     const Quadrature<dim> q_fine =
-      QProjector<dim>::project_to_face(q_gauss, face_fine);
+      QProjector<dim>::project_to_face(fe.reference_cell_type(),
+                                       q_gauss,
+                                       face_fine);
     const unsigned int nq = q_fine.size();
 
     FEValues<dim> fine(mapping,
@@ -2084,10 +2086,8 @@ namespace FETools
          cell_number < GeometryInfo<dim>::max_children_per_face;
          ++cell_number)
       {
-        const Quadrature<dim> q_coarse =
-          QProjector<dim>::project_to_subface(q_gauss,
-                                              face_coarse,
-                                              cell_number);
+        const Quadrature<dim> q_coarse = QProjector<dim>::project_to_subface(
+          fe.reference_cell_type(), q_gauss, face_coarse, cell_number);
         FEValues<dim> coarse(mapping, fe, q_coarse, update_values);
 
         typename Triangulation<dim, spacedim>::active_cell_iterator fine_cell =
index 8832704e63af3542c8453213ed2afef34baea433..d678ce0f61935120f2db36fabd48905bb9e13ed0 100644 (file)
@@ -45,6 +45,52 @@ namespace ReferenceCell
     Invalid = static_cast<std::uint8_t>(-1)
   };
 
+  /**
+   * Return the correct simplex reference cell type for the given dimension
+   * @p dim.
+   */
+  inline Type
+  get_simplex(const unsigned int dim)
+  {
+    switch (dim)
+      {
+        case 0:
+          return Type::Vertex;
+        case 1:
+          return Type::Line;
+        case 2:
+          return Type::Tri;
+        case 3:
+          return Type::Tet;
+        default:
+          Assert(false, ExcNotImplemented());
+          return Type::Invalid;
+      }
+  }
+
+  /**
+   * Return the correct hypercube reference cell type for the given dimension
+   * @p dim.
+   */
+  inline Type
+  get_hypercube(const unsigned int dim)
+  {
+    switch (dim)
+      {
+        case 0:
+          return Type::Vertex;
+        case 1:
+          return Type::Line;
+        case 2:
+          return Type::Quad;
+        case 3:
+          return Type::Hex;
+        default:
+          Assert(false, ExcNotImplemented());
+          return Type::Invalid;
+      }
+  }
+
   namespace internal
   {
     /**
index a046beecc8305c31bcb2f5073feebc3db7575b0d..5a56310fe77bbb68daed1faaf194b82dec6893d6 100644 (file)
@@ -76,6 +76,7 @@ SET(_unity_include_src
   polynomials_raviart_thomas.cc
   polynomials_rt_bubbles.cc
   process_grid.cc
+  qprojector.cc
   quadrature.cc
   quadrature_lib.cc
   quadrature_selector.cc
diff --git a/source/base/qprojector.cc b/source/base/qprojector.cc
new file mode 100644 (file)
index 0000000..275a935
--- /dev/null
@@ -0,0 +1,1509 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/geometry_info.h>
+#include <deal.II/base/qprojector.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+
+template <int dim>
+Quadrature<2>
+QProjector<dim>::reflect(const Quadrature<2> &q)
+{
+  std::vector<Point<2>> q_points(q.size());
+  std::vector<double>   weights(q.size());
+  for (unsigned int i = 0; i < q.size(); ++i)
+    {
+      q_points[i][0] = q.point(i)[1];
+      q_points[i][1] = q.point(i)[0];
+
+      weights[i] = q.weight(i);
+    }
+
+  return Quadrature<2>(q_points, weights);
+}
+
+
+template <int dim>
+Quadrature<2>
+QProjector<dim>::rotate(const Quadrature<2> &q, const unsigned int n_times)
+{
+  std::vector<Point<2>> q_points(q.size());
+  std::vector<double>   weights(q.size());
+  for (unsigned int i = 0; i < q.size(); ++i)
+    {
+      switch (n_times % 4)
+        {
+          case 0:
+            // 0 degree
+            q_points[i][0] = q.point(i)[0];
+            q_points[i][1] = q.point(i)[1];
+            break;
+          case 1:
+            // 90 degree counterclockwise
+            q_points[i][0] = 1.0 - q.point(i)[1];
+            q_points[i][1] = q.point(i)[0];
+            break;
+          case 2:
+            // 180 degree counterclockwise
+            q_points[i][0] = 1.0 - q.point(i)[0];
+            q_points[i][1] = 1.0 - q.point(i)[1];
+            break;
+          case 3:
+            // 270 degree counterclockwise
+            q_points[i][0] = q.point(i)[1];
+            q_points[i][1] = 1.0 - q.point(i)[0];
+            break;
+        }
+
+      weights[i] = q.weight(i);
+    }
+
+  return Quadrature<2>(q_points, weights);
+}
+
+
+template <>
+void
+QProjector<1>::project_to_face(const Quadrature<0> &  quadrature,
+                               const unsigned int     face_no,
+                               std::vector<Point<1>> &q_points)
+{
+  project_to_face(ReferenceCell::Type::Line, quadrature, face_no, q_points);
+}
+
+
+
+template <>
+void
+QProjector<1>::project_to_face(const ReferenceCell::Type reference_cell_type,
+                               const Quadrature<0> &,
+                               const unsigned int     face_no,
+                               std::vector<Point<1>> &q_points)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Line, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 1;
+  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
+  AssertDimension(q_points.size(), 1);
+
+  q_points[0] = Point<dim>(static_cast<double>(face_no));
+}
+
+
+
+template <>
+void
+QProjector<2>::project_to_face(const Quadrature<1> &  quadrature,
+                               const unsigned int     face_no,
+                               std::vector<Point<2>> &q_points)
+{
+  project_to_face(ReferenceCell::Type::Quad, quadrature, face_no, q_points);
+}
+
+
+
+template <>
+void
+QProjector<2>::project_to_face(const ReferenceCell::Type reference_cell_type,
+                               const Quadrature<1> &     quadrature,
+                               const unsigned int        face_no,
+                               std::vector<Point<2>> &   q_points)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Quad, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 2;
+  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
+  Assert(q_points.size() == quadrature.size(),
+         ExcDimensionMismatch(q_points.size(), quadrature.size()));
+
+  for (unsigned int p = 0; p < quadrature.size(); ++p)
+    switch (face_no)
+      {
+        case 0:
+          q_points[p] = Point<dim>(0, quadrature.point(p)(0));
+          break;
+        case 1:
+          q_points[p] = Point<dim>(1, quadrature.point(p)(0));
+          break;
+        case 2:
+          q_points[p] = Point<dim>(quadrature.point(p)(0), 0);
+          break;
+        case 3:
+          q_points[p] = Point<dim>(quadrature.point(p)(0), 1);
+          break;
+        default:
+          Assert(false, ExcInternalError());
+      }
+}
+
+
+
+template <>
+void
+QProjector<3>::project_to_face(const Quadrature<2> &  quadrature,
+                               const unsigned int     face_no,
+                               std::vector<Point<3>> &q_points)
+{
+  project_to_face(ReferenceCell::Type::Hex, quadrature, face_no, q_points);
+}
+
+
+
+template <>
+void
+QProjector<3>::project_to_face(const ReferenceCell::Type reference_cell_type,
+                               const Quadrature<2> &     quadrature,
+                               const unsigned int        face_no,
+                               std::vector<Point<3>> &   q_points)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Hex, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 3;
+  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
+  Assert(q_points.size() == quadrature.size(),
+         ExcDimensionMismatch(q_points.size(), quadrature.size()));
+
+  for (unsigned int p = 0; p < quadrature.size(); ++p)
+    switch (face_no)
+      {
+        case 0:
+          q_points[p] =
+            Point<dim>(0, quadrature.point(p)(0), quadrature.point(p)(1));
+          break;
+        case 1:
+          q_points[p] =
+            Point<dim>(1, quadrature.point(p)(0), quadrature.point(p)(1));
+          break;
+        case 2:
+          q_points[p] =
+            Point<dim>(quadrature.point(p)(1), 0, quadrature.point(p)(0));
+          break;
+        case 3:
+          q_points[p] =
+            Point<dim>(quadrature.point(p)(1), 1, quadrature.point(p)(0));
+          break;
+        case 4:
+          q_points[p] =
+            Point<dim>(quadrature.point(p)(0), quadrature.point(p)(1), 0);
+          break;
+        case 5:
+          q_points[p] =
+            Point<dim>(quadrature.point(p)(0), quadrature.point(p)(1), 1);
+          break;
+
+        default:
+          Assert(false, ExcInternalError());
+      }
+}
+
+
+
+template <>
+void
+QProjector<1>::project_to_subface(const Quadrature<0> &    quadrature,
+                                  const unsigned int       face_no,
+                                  const unsigned int       subface_no,
+                                  std::vector<Point<1>> &  q_points,
+                                  const RefinementCase<0> &ref_case)
+{
+  project_to_subface(ReferenceCell::Type::Line,
+                     quadrature,
+                     face_no,
+                     subface_no,
+                     q_points,
+                     ref_case);
+}
+
+
+
+template <>
+void
+QProjector<1>::project_to_subface(const ReferenceCell::Type reference_cell_type,
+                                  const Quadrature<0> &,
+                                  const unsigned int face_no,
+                                  const unsigned int,
+                                  std::vector<Point<1>> &q_points,
+                                  const RefinementCase<0> &)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Line, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 1;
+  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
+  AssertDimension(q_points.size(), 1);
+
+  q_points[0] = Point<dim>(static_cast<double>(face_no));
+}
+
+
+
+template <>
+void
+QProjector<2>::project_to_subface(const Quadrature<1> &    quadrature,
+                                  const unsigned int       face_no,
+                                  const unsigned int       subface_no,
+                                  std::vector<Point<2>> &  q_points,
+                                  const RefinementCase<1> &ref_case)
+{
+  project_to_subface(ReferenceCell::Type::Quad,
+                     quadrature,
+                     face_no,
+                     subface_no,
+                     q_points,
+                     ref_case);
+}
+
+
+
+template <>
+void
+QProjector<2>::project_to_subface(const ReferenceCell::Type reference_cell_type,
+                                  const Quadrature<1> &     quadrature,
+                                  const unsigned int        face_no,
+                                  const unsigned int        subface_no,
+                                  std::vector<Point<2>> &   q_points,
+                                  const RefinementCase<1> &)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Quad, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 2;
+  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
+  AssertIndexRange(subface_no, GeometryInfo<dim>::max_children_per_face);
+
+  Assert(q_points.size() == quadrature.size(),
+         ExcDimensionMismatch(q_points.size(), quadrature.size()));
+
+  for (unsigned int p = 0; p < quadrature.size(); ++p)
+    switch (face_no)
+      {
+        case 0:
+          switch (subface_no)
+            {
+              case 0:
+                q_points[p] = Point<dim>(0, quadrature.point(p)(0) / 2);
+                break;
+              case 1:
+                q_points[p] = Point<dim>(0, quadrature.point(p)(0) / 2 + 0.5);
+                break;
+              default:
+                Assert(false, ExcInternalError());
+            }
+          break;
+        case 1:
+          switch (subface_no)
+            {
+              case 0:
+                q_points[p] = Point<dim>(1, quadrature.point(p)(0) / 2);
+                break;
+              case 1:
+                q_points[p] = Point<dim>(1, quadrature.point(p)(0) / 2 + 0.5);
+                break;
+              default:
+                Assert(false, ExcInternalError());
+            }
+          break;
+        case 2:
+          switch (subface_no)
+            {
+              case 0:
+                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2, 0);
+                break;
+              case 1:
+                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2 + 0.5, 0);
+                break;
+              default:
+                Assert(false, ExcInternalError());
+            }
+          break;
+        case 3:
+          switch (subface_no)
+            {
+              case 0:
+                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2, 1);
+                break;
+              case 1:
+                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2 + 0.5, 1);
+                break;
+              default:
+                Assert(false, ExcInternalError());
+            }
+          break;
+
+        default:
+          Assert(false, ExcInternalError());
+      }
+}
+
+
+
+template <>
+void
+QProjector<3>::project_to_subface(const Quadrature<2> &    quadrature,
+                                  const unsigned int       face_no,
+                                  const unsigned int       subface_no,
+                                  std::vector<Point<3>> &  q_points,
+                                  const RefinementCase<2> &ref_case)
+{
+  project_to_subface(ReferenceCell::Type::Hex,
+                     quadrature,
+                     face_no,
+                     subface_no,
+                     q_points,
+                     ref_case);
+}
+
+
+
+template <>
+void
+QProjector<3>::project_to_subface(const ReferenceCell::Type reference_cell_type,
+                                  const Quadrature<2> &     quadrature,
+                                  const unsigned int        face_no,
+                                  const unsigned int        subface_no,
+                                  std::vector<Point<3>> &   q_points,
+                                  const RefinementCase<2> & ref_case)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Hex, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 3;
+  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
+  AssertIndexRange(subface_no, GeometryInfo<dim>::max_children_per_face);
+  Assert(q_points.size() == quadrature.size(),
+         ExcDimensionMismatch(q_points.size(), quadrature.size()));
+
+  // one coordinate is at a const value. for
+  // faces 0, 2 and 4 this value is 0.0, for
+  // faces 1, 3 and 5 it is 1.0
+  double const_value = face_no % 2;
+  // local 2d coordinates are xi and eta,
+  // global 3d coordinates are x, y and
+  // z. those have to be mapped. the following
+  // indices tell, which global coordinate
+  // (0->x, 1->y, 2->z) corresponds to which
+  // local one
+  unsigned int xi_index    = numbers::invalid_unsigned_int,
+               eta_index   = numbers::invalid_unsigned_int,
+               const_index = face_no / 2;
+  // the xi and eta values have to be scaled
+  // (by factor 0.5 or factor 1.0) depending on
+  // the refinement case and translated (by 0.0
+  // or 0.5) depending on the refinement case
+  // and subface_no.
+  double xi_scale = 1.0, eta_scale = 1.0, xi_translation = 0.0,
+         eta_translation = 0.0;
+  // set the index mapping between local and
+  // global coordinates
+  switch (face_no / 2)
+    {
+      case 0:
+        xi_index  = 1;
+        eta_index = 2;
+        break;
+      case 1:
+        xi_index  = 2;
+        eta_index = 0;
+        break;
+      case 2:
+        xi_index  = 0;
+        eta_index = 1;
+        break;
+    }
+  // set the scale and translation parameter
+  // for individual subfaces
+  switch (ref_case)
+    {
+      case RefinementCase<dim - 1>::cut_x:
+        xi_scale       = 0.5;
+        xi_translation = subface_no % 2 * 0.5;
+        break;
+      case RefinementCase<dim - 1>::cut_y:
+        eta_scale       = 0.5;
+        eta_translation = subface_no % 2 * 0.5;
+        break;
+      case RefinementCase<dim - 1>::cut_xy:
+        xi_scale        = 0.5;
+        eta_scale       = 0.5;
+        xi_translation  = int(subface_no % 2) * 0.5;
+        eta_translation = int(subface_no / 2) * 0.5;
+        break;
+      default:
+        Assert(false, ExcInternalError());
+        break;
+    }
+  // finally, compute the scaled, translated,
+  // projected quadrature points
+  for (unsigned int p = 0; p < quadrature.size(); ++p)
+    {
+      q_points[p][xi_index] =
+        xi_scale * quadrature.point(p)(0) + xi_translation;
+      q_points[p][eta_index] =
+        eta_scale * quadrature.point(p)(1) + eta_translation;
+      q_points[p][const_index] = const_value;
+    }
+}
+
+
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_faces(const Quadrature<0> &quadrature)
+{
+  return project_to_all_faces(ReferenceCell::Type::Line, quadrature);
+}
+
+
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_faces(
+  const ReferenceCell::Type reference_cell_type,
+  const Quadrature<0> &     quadrature)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Line, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 1;
+
+  const unsigned int n_points = 1, n_faces = GeometryInfo<dim>::faces_per_cell;
+
+  // first fix quadrature points
+  std::vector<Point<dim>> q_points;
+  q_points.reserve(n_points * n_faces);
+  std::vector<Point<dim>> help(n_points);
+
+
+  // project to each face and append
+  // results
+  for (unsigned int face = 0; face < n_faces; ++face)
+    {
+      project_to_face(quadrature, face, help);
+      std::copy(help.begin(), help.end(), std::back_inserter(q_points));
+    }
+
+  // next copy over weights
+  std::vector<double> weights;
+  weights.reserve(n_points * n_faces);
+  for (unsigned int face = 0; face < n_faces; ++face)
+    std::copy(quadrature.get_weights().begin(),
+              quadrature.get_weights().end(),
+              std::back_inserter(weights));
+
+  Assert(q_points.size() == n_points * n_faces, ExcInternalError());
+  Assert(weights.size() == n_points * n_faces, ExcInternalError());
+
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <>
+Quadrature<2>
+QProjector<2>::project_to_all_faces(
+  const ReferenceCell::Type reference_cell_type,
+  const SubQuadrature &     quadrature)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Quad, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 2;
+
+  const unsigned int n_points = quadrature.size(),
+                     n_faces  = GeometryInfo<dim>::faces_per_cell;
+
+  // first fix quadrature points
+  std::vector<Point<dim>> q_points;
+  q_points.reserve(n_points * n_faces);
+  std::vector<Point<dim>> help(n_points);
+
+  // project to each face and append
+  // results
+  for (unsigned int face = 0; face < n_faces; ++face)
+    {
+      project_to_face(quadrature, face, help);
+      std::copy(help.begin(), help.end(), std::back_inserter(q_points));
+    }
+
+  // next copy over weights
+  std::vector<double> weights;
+  weights.reserve(n_points * n_faces);
+  for (unsigned int face = 0; face < n_faces; ++face)
+    std::copy(quadrature.get_weights().begin(),
+              quadrature.get_weights().end(),
+              std::back_inserter(weights));
+
+  Assert(q_points.size() == n_points * n_faces, ExcInternalError());
+  Assert(weights.size() == n_points * n_faces, ExcInternalError());
+
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <>
+Quadrature<2>
+QProjector<2>::project_to_all_faces(const SubQuadrature &quadrature)
+{
+  return project_to_all_faces(ReferenceCell::Type::Quad, quadrature);
+}
+
+
+
+template <>
+Quadrature<3>
+QProjector<3>::project_to_all_faces(
+  const ReferenceCell::Type reference_cell_type,
+  const SubQuadrature &     quadrature)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Hex, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 3;
+
+  SubQuadrature q_reflected = reflect(quadrature);
+  SubQuadrature q[8]        = {quadrature,
+                        rotate(quadrature, 1),
+                        rotate(quadrature, 2),
+                        rotate(quadrature, 3),
+                        q_reflected,
+                        rotate(q_reflected, 3),
+                        rotate(q_reflected, 2),
+                        rotate(q_reflected, 1)};
+
+
+
+  const unsigned int n_points = quadrature.size(),
+                     n_faces  = GeometryInfo<dim>::faces_per_cell;
+
+  // first fix quadrature points
+  std::vector<Point<dim>> q_points;
+  q_points.reserve(n_points * n_faces * 8);
+  std::vector<Point<dim>> help(n_points);
+
+  std::vector<double> weights;
+  weights.reserve(n_points * n_faces * 8);
+
+  // do the following for all possible
+  // mutations of a face (mutation==0
+  // corresponds to a face with standard
+  // orientation, no flip and no rotation)
+  for (const auto &mutation : q)
+    {
+      // project to each face and append
+      // results
+      for (unsigned int face = 0; face < n_faces; ++face)
+        {
+          project_to_face(mutation, face, help);
+          std::copy(help.begin(), help.end(), std::back_inserter(q_points));
+        }
+
+      // next copy over weights
+      for (unsigned int face = 0; face < n_faces; ++face)
+        std::copy(mutation.get_weights().begin(),
+                  mutation.get_weights().end(),
+                  std::back_inserter(weights));
+    }
+
+
+  Assert(q_points.size() == n_points * n_faces * 8, ExcInternalError());
+  Assert(weights.size() == n_points * n_faces * 8, ExcInternalError());
+
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <>
+Quadrature<3>
+QProjector<3>::project_to_all_faces(const SubQuadrature &quadrature)
+{
+  return project_to_all_faces(ReferenceCell::Type::Hex, quadrature);
+}
+
+
+
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_subfaces(const Quadrature<0> &quadrature)
+{
+  return project_to_all_subfaces(ReferenceCell::Type::Line, quadrature);
+}
+
+
+
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_subfaces(
+  const ReferenceCell::Type reference_cell_type,
+  const Quadrature<0> &     quadrature)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Line, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 1;
+
+  const unsigned int n_points = 1, n_faces = GeometryInfo<dim>::faces_per_cell,
+                     subfaces_per_face =
+                       GeometryInfo<dim>::max_children_per_face;
+
+  // first fix quadrature points
+  std::vector<Point<dim>> q_points;
+  q_points.reserve(n_points * n_faces * subfaces_per_face);
+  std::vector<Point<dim>> help(n_points);
+
+  // project to each face and copy
+  // results
+  for (unsigned int face = 0; face < n_faces; ++face)
+    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
+      {
+        project_to_subface(quadrature, face, subface, help);
+        std::copy(help.begin(), help.end(), std::back_inserter(q_points));
+      }
+
+  // next copy over weights
+  std::vector<double> weights;
+  weights.reserve(n_points * n_faces * subfaces_per_face);
+  for (unsigned int face = 0; face < n_faces; ++face)
+    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
+      std::copy(quadrature.get_weights().begin(),
+                quadrature.get_weights().end(),
+                std::back_inserter(weights));
+
+  Assert(q_points.size() == n_points * n_faces * subfaces_per_face,
+         ExcInternalError());
+  Assert(weights.size() == n_points * n_faces * subfaces_per_face,
+         ExcInternalError());
+
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <>
+Quadrature<2>
+QProjector<2>::project_to_all_subfaces(
+  const ReferenceCell::Type reference_cell_type,
+  const SubQuadrature &     quadrature)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Quad, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim = 2;
+
+  const unsigned int n_points = quadrature.size(),
+                     n_faces  = GeometryInfo<dim>::faces_per_cell,
+                     subfaces_per_face =
+                       GeometryInfo<dim>::max_children_per_face;
+
+  // first fix quadrature points
+  std::vector<Point<dim>> q_points;
+  q_points.reserve(n_points * n_faces * subfaces_per_face);
+  std::vector<Point<dim>> help(n_points);
+
+  // project to each face and copy
+  // results
+  for (unsigned int face = 0; face < n_faces; ++face)
+    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
+      {
+        project_to_subface(quadrature, face, subface, help);
+        std::copy(help.begin(), help.end(), std::back_inserter(q_points));
+      }
+
+  // next copy over weights
+  std::vector<double> weights;
+  weights.reserve(n_points * n_faces * subfaces_per_face);
+  for (unsigned int face = 0; face < n_faces; ++face)
+    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
+      std::copy(quadrature.get_weights().begin(),
+                quadrature.get_weights().end(),
+                std::back_inserter(weights));
+
+  Assert(q_points.size() == n_points * n_faces * subfaces_per_face,
+         ExcInternalError());
+  Assert(weights.size() == n_points * n_faces * subfaces_per_face,
+         ExcInternalError());
+
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <>
+Quadrature<2>
+QProjector<2>::project_to_all_subfaces(const SubQuadrature &quadrature)
+{
+  return project_to_all_subfaces(ReferenceCell::Type::Quad, quadrature);
+}
+
+
+
+template <>
+Quadrature<3>
+QProjector<3>::project_to_all_subfaces(
+  const ReferenceCell::Type reference_cell_type,
+  const SubQuadrature &     quadrature)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Hex, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int dim         = 3;
+  SubQuadrature      q_reflected = reflect(quadrature);
+  SubQuadrature      q[8]        = {quadrature,
+                        rotate(quadrature, 1),
+                        rotate(quadrature, 2),
+                        rotate(quadrature, 3),
+                        q_reflected,
+                        rotate(q_reflected, 3),
+                        rotate(q_reflected, 2),
+                        rotate(q_reflected, 1)};
+
+  const unsigned int n_points = quadrature.size(),
+                     n_faces  = GeometryInfo<dim>::faces_per_cell,
+                     total_subfaces_per_face = 2 + 2 + 4;
+
+  // first fix quadrature points
+  std::vector<Point<dim>> q_points;
+  q_points.reserve(n_points * n_faces * total_subfaces_per_face * 8);
+  std::vector<Point<dim>> help(n_points);
+
+  std::vector<double> weights;
+  weights.reserve(n_points * n_faces * total_subfaces_per_face * 8);
+
+  // do the following for all possible
+  // mutations of a face (mutation==0
+  // corresponds to a face with standard
+  // orientation, no flip and no rotation)
+  for (const auto &mutation : q)
+    {
+      // project to each face and copy
+      // results
+      for (unsigned int face = 0; face < n_faces; ++face)
+        for (unsigned int ref_case = RefinementCase<dim - 1>::cut_xy;
+             ref_case >= RefinementCase<dim - 1>::cut_x;
+             --ref_case)
+          for (unsigned int subface = 0;
+               subface < GeometryInfo<dim - 1>::n_children(
+                           RefinementCase<dim - 1>(ref_case));
+               ++subface)
+            {
+              project_to_subface(mutation,
+                                 face,
+                                 subface,
+                                 help,
+                                 RefinementCase<dim - 1>(ref_case));
+              std::copy(help.begin(), help.end(), std::back_inserter(q_points));
+            }
+
+      // next copy over weights
+      for (unsigned int face = 0; face < n_faces; ++face)
+        for (unsigned int ref_case = RefinementCase<dim - 1>::cut_xy;
+             ref_case >= RefinementCase<dim - 1>::cut_x;
+             --ref_case)
+          for (unsigned int subface = 0;
+               subface < GeometryInfo<dim - 1>::n_children(
+                           RefinementCase<dim - 1>(ref_case));
+               ++subface)
+            std::copy(mutation.get_weights().begin(),
+                      mutation.get_weights().end(),
+                      std::back_inserter(weights));
+    }
+
+  Assert(q_points.size() == n_points * n_faces * total_subfaces_per_face * 8,
+         ExcInternalError());
+  Assert(weights.size() == n_points * n_faces * total_subfaces_per_face * 8,
+         ExcInternalError());
+
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <>
+Quadrature<3>
+QProjector<3>::project_to_all_subfaces(const SubQuadrature &quadrature)
+{
+  return project_to_all_subfaces(ReferenceCell::Type::Hex, quadrature);
+}
+
+
+
+// This function is not used in the library
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_child(const Quadrature<dim> &quadrature,
+                                  const unsigned int     child_no)
+{
+  return project_to_child(ReferenceCell::get_hypercube(dim),
+                          quadrature,
+                          child_no);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_child(const ReferenceCell::Type reference_cell_type,
+                                  const Quadrature<dim> &   quadrature,
+                                  const unsigned int        child_no)
+{
+  Assert(reference_cell_type == ReferenceCell::get_hypercube(dim),
+         ExcNotImplemented());
+  (void)reference_cell_type;
+
+  AssertIndexRange(child_no, GeometryInfo<dim>::max_children_per_cell);
+
+  const unsigned int n_q_points = quadrature.size();
+
+  std::vector<Point<dim>> q_points(n_q_points);
+  for (unsigned int i = 0; i < n_q_points; ++i)
+    q_points[i] =
+      GeometryInfo<dim>::child_to_cell_coordinates(quadrature.point(i),
+                                                   child_no);
+
+  // for the weights, things are
+  // equally simple: copy them and
+  // scale them
+  std::vector<double> weights = quadrature.get_weights();
+  for (unsigned int i = 0; i < n_q_points; ++i)
+    weights[i] *= (1. / GeometryInfo<dim>::max_children_per_cell);
+
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_all_children(const Quadrature<dim> &quadrature)
+{
+  return project_to_all_children(ReferenceCell::get_hypercube(dim), quadrature);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_all_children(
+  const ReferenceCell::Type reference_cell_type,
+  const Quadrature<dim> &   quadrature)
+{
+  Assert(reference_cell_type == ReferenceCell::get_hypercube(dim),
+         ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int n_points   = quadrature.size(),
+                     n_children = GeometryInfo<dim>::max_children_per_cell;
+
+  std::vector<Point<dim>> q_points(n_points * n_children);
+  std::vector<double>     weights(n_points * n_children);
+
+  // project to each child and copy
+  // results
+  for (unsigned int child = 0; child < n_children; ++child)
+    {
+      Quadrature<dim> help = project_to_child(quadrature, child);
+      for (unsigned int i = 0; i < n_points; ++i)
+        {
+          q_points[child * n_points + i] = help.point(i);
+          weights[child * n_points + i]  = help.weight(i);
+        }
+    }
+  return Quadrature<dim>(q_points, weights);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_line(const Quadrature<1> &quadrature,
+                                 const Point<dim> &   p1,
+                                 const Point<dim> &   p2)
+{
+  return project_to_line(ReferenceCell::get_hypercube(dim), quadrature, p1, p2);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_line(const ReferenceCell::Type reference_cell_type,
+                                 const Quadrature<1> &     quadrature,
+                                 const Point<dim> &        p1,
+                                 const Point<dim> &        p2)
+{
+  Assert(reference_cell_type == ReferenceCell::get_hypercube(dim),
+         ExcNotImplemented());
+  (void)reference_cell_type;
+
+  const unsigned int      n = quadrature.size();
+  std::vector<Point<dim>> points(n);
+  std::vector<double>     weights(n);
+  const double            length = p1.distance(p2);
+
+  for (unsigned int k = 0; k < n; ++k)
+    {
+      const double alpha = quadrature.point(k)(0);
+      points[k]          = alpha * p2;
+      points[k] += (1. - alpha) * p1;
+      weights[k] = length * quadrature.weight(k);
+    }
+  return Quadrature<dim>(points, weights);
+}
+
+
+
+template <int dim>
+typename QProjector<dim>::DataSetDescriptor
+QProjector<dim>::DataSetDescriptor::face(const unsigned int face_no,
+                                         const bool         face_orientation,
+                                         const bool         face_flip,
+                                         const bool         face_rotation,
+                                         const unsigned int n_quadrature_points)
+{
+  return face(ReferenceCell::get_hypercube(dim),
+              face_no,
+              face_orientation,
+              face_flip,
+              face_rotation,
+              n_quadrature_points);
+}
+
+
+
+template <int dim>
+typename QProjector<dim>::DataSetDescriptor
+QProjector<dim>::DataSetDescriptor::face(
+  const ReferenceCell::Type reference_cell_type,
+  const unsigned int        face_no,
+  const bool                face_orientation,
+  const bool                face_flip,
+  const bool                face_rotation,
+  const unsigned int        n_quadrature_points)
+{
+  Assert(reference_cell_type == ReferenceCell::get_hypercube(dim),
+         ExcNotImplemented());
+  (void)reference_cell_type;
+
+  Assert(face_no < GeometryInfo<dim>::faces_per_cell, ExcInternalError());
+
+  switch (dim)
+    {
+      case 1:
+      case 2:
+        return face_no * n_quadrature_points;
+
+
+      case 3:
+        {
+          // in 3d, we have to account for faces that
+          // have non-standard face orientation, flip
+          // and rotation. thus, we have to store
+          // _eight_ data sets per face or subface
+
+          // set up a table with the according offsets
+          // for non-standard orientation, first index:
+          // face_orientation (standard true=1), second
+          // index: face_flip (standard false=0), third
+          // index: face_rotation (standard false=0)
+          //
+          // note, that normally we should use the
+          // obvious offsets 0,1,2,3,4,5,6,7. However,
+          // prior to the changes enabling flipped and
+          // rotated faces, in many places of the
+          // library the convention was used, that the
+          // first dataset with offset 0 corresponds to
+          // a face in standard orientation. therefore
+          // we use the offsets 4,5,6,7,0,1,2,3 here to
+          // stick to that (implicit) convention
+          static const unsigned int offset[2][2][2] = {
+            {{4 * GeometryInfo<dim>::faces_per_cell,
+              5 * GeometryInfo<dim>::
+                    faces_per_cell}, // face_orientation=false; face_flip=false;
+                                     // face_rotation=false and true
+             {6 * GeometryInfo<dim>::faces_per_cell,
+              7 * GeometryInfo<dim>::
+                    faces_per_cell}}, // face_orientation=false; face_flip=true;
+                                      // face_rotation=false and true
+            {{0 * GeometryInfo<dim>::faces_per_cell,
+              1 * GeometryInfo<dim>::
+                    faces_per_cell}, // face_orientation=true;  face_flip=false;
+                                     // face_rotation=false and true
+             {2 * GeometryInfo<dim>::faces_per_cell,
+              3 * GeometryInfo<dim>::
+                    faces_per_cell}}}; // face_orientation=true; face_flip=true;
+                                       // face_rotation=false and true
+
+          return (
+            (face_no + offset[face_orientation][face_flip][face_rotation]) *
+            n_quadrature_points);
+        }
+
+      default:
+        Assert(false, ExcInternalError());
+    }
+  return numbers::invalid_unsigned_int;
+}
+
+
+
+template <>
+QProjector<1>::DataSetDescriptor
+QProjector<1>::DataSetDescriptor::subface(
+  const ReferenceCell::Type reference_cell_type,
+  const unsigned int        face_no,
+  const unsigned int        subface_no,
+  const bool,
+  const bool,
+  const bool,
+  const unsigned int n_quadrature_points,
+  const internal::SubfaceCase<1>)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Line, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  Assert(face_no < GeometryInfo<1>::faces_per_cell, ExcInternalError());
+  Assert(subface_no < GeometryInfo<1>::max_children_per_face,
+         ExcInternalError());
+
+  return ((face_no * GeometryInfo<1>::max_children_per_face + subface_no) *
+          n_quadrature_points);
+}
+
+
+
+template <>
+QProjector<1>::DataSetDescriptor
+QProjector<1>::DataSetDescriptor::subface(
+  const unsigned int             face_no,
+  const unsigned int             subface_no,
+  const bool                     face_orientation,
+  const bool                     face_flip,
+  const bool                     face_rotation,
+  const unsigned int             n_quadrature_points,
+  const internal::SubfaceCase<1> ref_case)
+{
+  return subface(ReferenceCell::Type::Line,
+                 face_no,
+                 subface_no,
+                 face_orientation,
+                 face_flip,
+                 face_rotation,
+                 n_quadrature_points,
+                 ref_case);
+}
+
+
+
+template <>
+QProjector<2>::DataSetDescriptor
+QProjector<2>::DataSetDescriptor::subface(
+  const ReferenceCell::Type reference_cell_type,
+  const unsigned int        face_no,
+  const unsigned int        subface_no,
+  const bool,
+  const bool,
+  const bool,
+  const unsigned int n_quadrature_points,
+  const internal::SubfaceCase<2>)
+{
+  Assert(reference_cell_type == ReferenceCell::Type::Quad, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  Assert(face_no < GeometryInfo<2>::faces_per_cell, ExcInternalError());
+  Assert(subface_no < GeometryInfo<2>::max_children_per_face,
+         ExcInternalError());
+
+  return ((face_no * GeometryInfo<2>::max_children_per_face + subface_no) *
+          n_quadrature_points);
+}
+
+
+
+template <>
+QProjector<2>::DataSetDescriptor
+QProjector<2>::DataSetDescriptor::subface(
+  const unsigned int             face_no,
+  const unsigned int             subface_no,
+  const bool                     face_orientation,
+  const bool                     face_flip,
+  const bool                     face_rotation,
+  const unsigned int             n_quadrature_points,
+  const internal::SubfaceCase<2> ref_case)
+{
+  return subface(ReferenceCell::Type::Quad,
+                 face_no,
+                 subface_no,
+                 face_orientation,
+                 face_flip,
+                 face_rotation,
+                 n_quadrature_points,
+                 ref_case);
+}
+
+
+template <>
+QProjector<3>::DataSetDescriptor
+QProjector<3>::DataSetDescriptor::subface(
+  const ReferenceCell::Type      reference_cell_type,
+  const unsigned int             face_no,
+  const unsigned int             subface_no,
+  const bool                     face_orientation,
+  const bool                     face_flip,
+  const bool                     face_rotation,
+  const unsigned int             n_quadrature_points,
+  const internal::SubfaceCase<3> ref_case)
+{
+  const unsigned int dim = 3;
+
+  Assert(reference_cell_type == ReferenceCell::Type::Hex, ExcNotImplemented());
+  (void)reference_cell_type;
+
+  Assert(face_no < GeometryInfo<dim>::faces_per_cell, ExcInternalError());
+  Assert(subface_no < GeometryInfo<dim>::max_children_per_face,
+         ExcInternalError());
+
+  // As the quadrature points created by
+  // QProjector are on subfaces in their
+  // "standard location" we have to use a
+  // permutation of the equivalent subface
+  // number in order to respect face
+  // orientation, flip and rotation. The
+  // information we need here is exactly the
+  // same as the
+  // GeometryInfo<3>::child_cell_on_face info
+  // for the bottom face (face 4) of a hex, as
+  // on this the RefineCase of the cell matches
+  // that of the face and the subfaces are
+  // numbered in the same way as the child
+  // cells.
+
+  // in 3d, we have to account for faces that
+  // have non-standard face orientation, flip
+  // and rotation. thus, we have to store
+  // _eight_ data sets per face or subface
+  // already for the isotropic
+  // case. Additionally, we have three
+  // different refinement cases, resulting in
+  // <tt>4 + 2 + 2 = 8</tt> different subfaces
+  // for each face.
+  const unsigned int total_subfaces_per_face = 8;
+
+  // set up a table with the according offsets
+  // for non-standard orientation, first index:
+  // face_orientation (standard true=1), second
+  // index: face_flip (standard false=0), third
+  // index: face_rotation (standard false=0)
+  //
+  // note, that normally we should use the
+  // obvious offsets 0,1,2,3,4,5,6,7. However,
+  // prior to the changes enabling flipped and
+  // rotated faces, in many places of the
+  // library the convention was used, that the
+  // first dataset with offset 0 corresponds to
+  // a face in standard orientation. therefore
+  // we use the offsets 4,5,6,7,0,1,2,3 here to
+  // stick to that (implicit) convention
+  static const unsigned int orientation_offset[2][2][2] = {
+    {// face_orientation=false; face_flip=false; face_rotation=false and true
+     {4 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
+      5 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face},
+     // face_orientation=false; face_flip=true;  face_rotation=false and true
+     {6 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
+      7 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face}},
+    {// face_orientation=true;  face_flip=false; face_rotation=false and true
+     {0 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
+      1 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face},
+     // face_orientation=true;  face_flip=true;  face_rotation=false and true
+     {2 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
+      3 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face}}};
+
+  // set up a table with the offsets for a
+  // given refinement case respecting the
+  // corresponding number of subfaces. the
+  // index corresponds to (RefineCase::Type - 1)
+
+  // note, that normally we should use the
+  // obvious offsets 0,2,6. However, prior to
+  // the implementation of anisotropic
+  // refinement, in many places of the library
+  // the convention was used, that the first
+  // dataset with offset 0 corresponds to a
+  // standard (isotropic) face
+  // refinement. therefore we use the offsets
+  // 6,4,0 here to stick to that (implicit)
+  // convention
+  static const unsigned int ref_case_offset[3] = {
+    6, // cut_x
+    4, // cut_y
+    0  // cut_xy
+  };
+
+
+  // for each subface of a given FaceRefineCase
+  // there is a corresponding equivalent
+  // subface number of one of the "standard"
+  // RefineCases (cut_x, cut_y, cut_xy). Map
+  // the given values to those equivalent
+  // ones.
+
+  // first, define an invalid number
+  static const unsigned int e = numbers::invalid_unsigned_int;
+
+  static const RefinementCase<dim - 1>
+    equivalent_refine_case[internal::SubfaceCase<dim>::case_isotropic + 1]
+                          [GeometryInfo<3>::max_children_per_face] = {
+                            // case_none. there should be only
+                            // invalid values here. However, as
+                            // this function is also called (in
+                            // tests) for cells which have no
+                            // refined faces, use isotropic
+                            // refinement instead
+                            {RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy},
+                            // case_x
+                            {RefinementCase<dim - 1>::cut_x,
+                             RefinementCase<dim - 1>::cut_x,
+                             RefinementCase<dim - 1>::no_refinement,
+                             RefinementCase<dim - 1>::no_refinement},
+                            // case_x1y
+                            {RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_x,
+                             RefinementCase<dim - 1>::no_refinement},
+                            // case_x2y
+                            {RefinementCase<dim - 1>::cut_x,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::no_refinement},
+                            // case_x1y2y
+                            {RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy},
+                            // case_y
+                            {RefinementCase<dim - 1>::cut_y,
+                             RefinementCase<dim - 1>::cut_y,
+                             RefinementCase<dim - 1>::no_refinement,
+                             RefinementCase<dim - 1>::no_refinement},
+                            // case_y1x
+                            {RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_y,
+                             RefinementCase<dim - 1>::no_refinement},
+                            // case_y2x
+                            {RefinementCase<dim - 1>::cut_y,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::no_refinement},
+                            // case_y1x2x
+                            {RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy},
+                            // case_xy (case_isotropic)
+                            {RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy,
+                             RefinementCase<dim - 1>::cut_xy}};
+
+  static const unsigned int
+    equivalent_subface_number[internal::SubfaceCase<dim>::case_isotropic + 1]
+                             [GeometryInfo<3>::max_children_per_face] = {
+                               // case_none, see above
+                               {0, 1, 2, 3},
+                               // case_x
+                               {0, 1, e, e},
+                               // case_x1y
+                               {0, 2, 1, e},
+                               // case_x2y
+                               {0, 1, 3, e},
+                               // case_x1y2y
+                               {0, 2, 1, 3},
+                               // case_y
+                               {0, 1, e, e},
+                               // case_y1x
+                               {0, 1, 1, e},
+                               // case_y2x
+                               {0, 2, 3, e},
+                               // case_y1x2x
+                               {0, 1, 2, 3},
+                               // case_xy (case_isotropic)
+                               {0, 1, 2, 3}};
+
+  // If face-orientation or face_rotation are
+  // non-standard, cut_x and cut_y have to be
+  // exchanged.
+  static const RefinementCase<dim - 1> ref_case_permutation[4] = {
+    RefinementCase<dim - 1>::no_refinement,
+    RefinementCase<dim - 1>::cut_y,
+    RefinementCase<dim - 1>::cut_x,
+    RefinementCase<dim - 1>::cut_xy};
+
+  // set a corresponding (equivalent)
+  // RefineCase and subface number
+  const RefinementCase<dim - 1> equ_ref_case =
+    equivalent_refine_case[ref_case][subface_no];
+  const unsigned int equ_subface_no =
+    equivalent_subface_number[ref_case][subface_no];
+  // make sure, that we got a valid subface and RefineCase
+  Assert(equ_ref_case != RefinementCase<dim>::no_refinement,
+         ExcInternalError());
+  Assert(equ_subface_no != e, ExcInternalError());
+  // now, finally respect non-standard faces
+  const RefinementCase<dim - 1> final_ref_case =
+    (face_orientation == face_rotation ? ref_case_permutation[equ_ref_case] :
+                                         equ_ref_case);
+
+  // what we have now is the number of
+  // the subface in the natural
+  // orientation of the *face*. what we
+  // need to know is the number of the
+  // subface concerning the standard face
+  // orientation as seen from the *cell*.
+
+  // this mapping is not trivial, but we
+  // have done exactly this stuff in the
+  // child_cell_on_face function. in
+  // order to reduce the amount of code
+  // as well as to make maintaining the
+  // functionality easier we want to
+  // reuse that information. So we note
+  // that on the bottom face (face 4) of
+  // a hex cell the local x and y
+  // coordinates of the face and the cell
+  // coincide, thus also the refinement
+  // case of the face corresponds to the
+  // refinement case of the cell
+  // (ignoring cell refinement along the
+  // z direction). Using this knowledge
+  // we can (ab)use the
+  // child_cell_on_face function to do
+  // exactly the transformation we are in
+  // need of now
+  const unsigned int final_subface_no =
+    GeometryInfo<dim>::child_cell_on_face(RefinementCase<dim>(final_ref_case),
+                                          4,
+                                          equ_subface_no,
+                                          face_orientation,
+                                          face_flip,
+                                          face_rotation,
+                                          equ_ref_case);
+
+  return (((face_no * total_subfaces_per_face +
+            ref_case_offset[final_ref_case - 1] + final_subface_no) +
+           orientation_offset[face_orientation][face_flip][face_rotation]) *
+          n_quadrature_points);
+}
+
+
+template <>
+QProjector<3>::DataSetDescriptor
+QProjector<3>::DataSetDescriptor::subface(
+  const unsigned int             face_no,
+  const unsigned int             subface_no,
+  const bool                     face_orientation,
+  const bool                     face_flip,
+  const bool                     face_rotation,
+  const unsigned int             n_quadrature_points,
+  const internal::SubfaceCase<3> ref_case)
+{
+  return subface(ReferenceCell::Type::Hex,
+                 face_no,
+                 subface_no,
+                 face_orientation,
+                 face_flip,
+                 face_rotation,
+                 n_quadrature_points,
+                 ref_case);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_face(const SubQuadrature &quadrature,
+                                 const unsigned int   face_no)
+{
+  return project_to_face(ReferenceCell::get_hypercube(dim),
+                         quadrature,
+                         face_no);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_face(const ReferenceCell::Type reference_cell_type,
+                                 const SubQuadrature &     quadrature,
+                                 const unsigned int        face_no)
+{
+  Assert(reference_cell_type == ReferenceCell::get_hypercube(dim),
+         ExcNotImplemented());
+  (void)reference_cell_type;
+
+  std::vector<Point<dim>> points(quadrature.size());
+  project_to_face(quadrature, face_no, points);
+  return Quadrature<dim>(points, quadrature.get_weights());
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_subface(const SubQuadrature &          quadrature,
+                                    const unsigned int             face_no,
+                                    const unsigned int             subface_no,
+                                    const RefinementCase<dim - 1> &ref_case)
+{
+  return project_to_subface(ReferenceCell::get_hypercube(dim),
+                            quadrature,
+                            face_no,
+                            subface_no,
+                            ref_case);
+}
+
+
+
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_subface(
+  const ReferenceCell::Type      reference_cell_type,
+  const SubQuadrature &          quadrature,
+  const unsigned int             face_no,
+  const unsigned int             subface_no,
+  const RefinementCase<dim - 1> &ref_case)
+{
+  Assert(reference_cell_type == ReferenceCell::get_hypercube(dim),
+         ExcNotImplemented());
+  (void)reference_cell_type;
+
+  std::vector<Point<dim>> points(quadrature.size());
+  project_to_subface(quadrature, face_no, subface_no, points, ref_case);
+  return Quadrature<dim>(points, quadrature.get_weights());
+}
+
+
+// explicit instantiations; note: we need them all for all dimensions
+template class QProjector<1>;
+template class QProjector<2>;
+template class QProjector<3>;
+
+DEAL_II_NAMESPACE_CLOSE
index 25829085e2f6263c63ffd503bfe3f15e7e6dff2f..9e79b4f98795aefbff384481842cfd930b1ce2a5 100644 (file)
@@ -424,1105 +424,6 @@ QAnisotropic<3>::QAnisotropic(const Quadrature<1> &qx,
 
 
 
-//---------------------------------------------------------------------------
-
-
-
-template <int dim>
-Quadrature<2>
-QProjector<dim>::reflect(const Quadrature<2> &q)
-{
-  std::vector<Point<2>> q_points(q.size());
-  std::vector<double>   weights(q.size());
-  for (unsigned int i = 0; i < q.size(); ++i)
-    {
-      q_points[i][0] = q.point(i)[1];
-      q_points[i][1] = q.point(i)[0];
-
-      weights[i] = q.weight(i);
-    }
-
-  return Quadrature<2>(q_points, weights);
-}
-
-
-template <int dim>
-Quadrature<2>
-QProjector<dim>::rotate(const Quadrature<2> &q, const unsigned int n_times)
-{
-  std::vector<Point<2>> q_points(q.size());
-  std::vector<double>   weights(q.size());
-  for (unsigned int i = 0; i < q.size(); ++i)
-    {
-      switch (n_times % 4)
-        {
-          case 0:
-            // 0 degree
-            q_points[i][0] = q.point(i)[0];
-            q_points[i][1] = q.point(i)[1];
-            break;
-          case 1:
-            // 90 degree counterclockwise
-            q_points[i][0] = 1.0 - q.point(i)[1];
-            q_points[i][1] = q.point(i)[0];
-            break;
-          case 2:
-            // 180 degree counterclockwise
-            q_points[i][0] = 1.0 - q.point(i)[0];
-            q_points[i][1] = 1.0 - q.point(i)[1];
-            break;
-          case 3:
-            // 270 degree counterclockwise
-            q_points[i][0] = q.point(i)[1];
-            q_points[i][1] = 1.0 - q.point(i)[0];
-            break;
-        }
-
-      weights[i] = q.weight(i);
-    }
-
-  return Quadrature<2>(q_points, weights);
-}
-
-
-template <>
-void
-QProjector<1>::project_to_face(const Quadrature<0> &,
-                               const unsigned int     face_no,
-                               std::vector<Point<1>> &q_points)
-{
-  const unsigned int dim = 1;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  AssertDimension(q_points.size(), 1);
-
-  q_points[0] = Point<dim>(static_cast<double>(face_no));
-}
-
-
-
-template <>
-void
-QProjector<2>::project_to_face(const Quadrature<1> &  quadrature,
-                               const unsigned int     face_no,
-                               std::vector<Point<2>> &q_points)
-{
-  const unsigned int dim = 2;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  Assert(q_points.size() == quadrature.size(),
-         ExcDimensionMismatch(q_points.size(), quadrature.size()));
-
-  for (unsigned int p = 0; p < quadrature.size(); ++p)
-    switch (face_no)
-      {
-        case 0:
-          q_points[p] = Point<dim>(0, quadrature.point(p)(0));
-          break;
-        case 1:
-          q_points[p] = Point<dim>(1, quadrature.point(p)(0));
-          break;
-        case 2:
-          q_points[p] = Point<dim>(quadrature.point(p)(0), 0);
-          break;
-        case 3:
-          q_points[p] = Point<dim>(quadrature.point(p)(0), 1);
-          break;
-        default:
-          Assert(false, ExcInternalError());
-      }
-}
-
-
-
-template <>
-void
-QProjector<3>::project_to_face(const Quadrature<2> &  quadrature,
-                               const unsigned int     face_no,
-                               std::vector<Point<3>> &q_points)
-{
-  const unsigned int dim = 3;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  Assert(q_points.size() == quadrature.size(),
-         ExcDimensionMismatch(q_points.size(), quadrature.size()));
-
-  for (unsigned int p = 0; p < quadrature.size(); ++p)
-    switch (face_no)
-      {
-        case 0:
-          q_points[p] =
-            Point<dim>(0, quadrature.point(p)(0), quadrature.point(p)(1));
-          break;
-        case 1:
-          q_points[p] =
-            Point<dim>(1, quadrature.point(p)(0), quadrature.point(p)(1));
-          break;
-        case 2:
-          q_points[p] =
-            Point<dim>(quadrature.point(p)(1), 0, quadrature.point(p)(0));
-          break;
-        case 3:
-          q_points[p] =
-            Point<dim>(quadrature.point(p)(1), 1, quadrature.point(p)(0));
-          break;
-        case 4:
-          q_points[p] =
-            Point<dim>(quadrature.point(p)(0), quadrature.point(p)(1), 0);
-          break;
-        case 5:
-          q_points[p] =
-            Point<dim>(quadrature.point(p)(0), quadrature.point(p)(1), 1);
-          break;
-
-        default:
-          Assert(false, ExcInternalError());
-      }
-}
-
-
-
-template <>
-void
-QProjector<1>::project_to_subface(const Quadrature<0> &,
-                                  const unsigned int face_no,
-                                  const unsigned int,
-                                  std::vector<Point<1>> &q_points,
-                                  const RefinementCase<0> &)
-{
-  const unsigned int dim = 1;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  AssertDimension(q_points.size(), 1);
-
-  q_points[0] = Point<dim>(static_cast<double>(face_no));
-}
-
-
-
-template <>
-void
-QProjector<2>::project_to_subface(const Quadrature<1> &  quadrature,
-                                  const unsigned int     face_no,
-                                  const unsigned int     subface_no,
-                                  std::vector<Point<2>> &q_points,
-                                  const RefinementCase<1> &)
-{
-  const unsigned int dim = 2;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  AssertIndexRange(subface_no, GeometryInfo<dim>::max_children_per_face);
-
-  Assert(q_points.size() == quadrature.size(),
-         ExcDimensionMismatch(q_points.size(), quadrature.size()));
-
-  for (unsigned int p = 0; p < quadrature.size(); ++p)
-    switch (face_no)
-      {
-        case 0:
-          switch (subface_no)
-            {
-              case 0:
-                q_points[p] = Point<dim>(0, quadrature.point(p)(0) / 2);
-                break;
-              case 1:
-                q_points[p] = Point<dim>(0, quadrature.point(p)(0) / 2 + 0.5);
-                break;
-              default:
-                Assert(false, ExcInternalError());
-            }
-          break;
-        case 1:
-          switch (subface_no)
-            {
-              case 0:
-                q_points[p] = Point<dim>(1, quadrature.point(p)(0) / 2);
-                break;
-              case 1:
-                q_points[p] = Point<dim>(1, quadrature.point(p)(0) / 2 + 0.5);
-                break;
-              default:
-                Assert(false, ExcInternalError());
-            }
-          break;
-        case 2:
-          switch (subface_no)
-            {
-              case 0:
-                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2, 0);
-                break;
-              case 1:
-                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2 + 0.5, 0);
-                break;
-              default:
-                Assert(false, ExcInternalError());
-            }
-          break;
-        case 3:
-          switch (subface_no)
-            {
-              case 0:
-                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2, 1);
-                break;
-              case 1:
-                q_points[p] = Point<dim>(quadrature.point(p)(0) / 2 + 0.5, 1);
-                break;
-              default:
-                Assert(false, ExcInternalError());
-            }
-          break;
-
-        default:
-          Assert(false, ExcInternalError());
-      }
-}
-
-
-
-template <>
-void
-QProjector<3>::project_to_subface(const Quadrature<2> &    quadrature,
-                                  const unsigned int       face_no,
-                                  const unsigned int       subface_no,
-                                  std::vector<Point<3>> &  q_points,
-                                  const RefinementCase<2> &ref_case)
-{
-  const unsigned int dim = 3;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  AssertIndexRange(subface_no, GeometryInfo<dim>::max_children_per_face);
-  Assert(q_points.size() == quadrature.size(),
-         ExcDimensionMismatch(q_points.size(), quadrature.size()));
-
-  // one coordinate is at a const value. for
-  // faces 0, 2 and 4 this value is 0.0, for
-  // faces 1, 3 and 5 it is 1.0
-  double const_value = face_no % 2;
-  // local 2d coordinates are xi and eta,
-  // global 3d coordinates are x, y and
-  // z. those have to be mapped. the following
-  // indices tell, which global coordinate
-  // (0->x, 1->y, 2->z) corresponds to which
-  // local one
-  unsigned int xi_index    = numbers::invalid_unsigned_int,
-               eta_index   = numbers::invalid_unsigned_int,
-               const_index = face_no / 2;
-  // the xi and eta values have to be scaled
-  // (by factor 0.5 or factor 1.0) depending on
-  // the refinement case and translated (by 0.0
-  // or 0.5) depending on the refinement case
-  // and subface_no.
-  double xi_scale = 1.0, eta_scale = 1.0, xi_translation = 0.0,
-         eta_translation = 0.0;
-  // set the index mapping between local and
-  // global coordinates
-  switch (face_no / 2)
-    {
-      case 0:
-        xi_index  = 1;
-        eta_index = 2;
-        break;
-      case 1:
-        xi_index  = 2;
-        eta_index = 0;
-        break;
-      case 2:
-        xi_index  = 0;
-        eta_index = 1;
-        break;
-    }
-  // set the scale and translation parameter
-  // for individual subfaces
-  switch (ref_case)
-    {
-      case RefinementCase<dim - 1>::cut_x:
-        xi_scale       = 0.5;
-        xi_translation = subface_no % 2 * 0.5;
-        break;
-      case RefinementCase<dim - 1>::cut_y:
-        eta_scale       = 0.5;
-        eta_translation = subface_no % 2 * 0.5;
-        break;
-      case RefinementCase<dim - 1>::cut_xy:
-        xi_scale        = 0.5;
-        eta_scale       = 0.5;
-        xi_translation  = int(subface_no % 2) * 0.5;
-        eta_translation = int(subface_no / 2) * 0.5;
-        break;
-      default:
-        Assert(false, ExcInternalError());
-        break;
-    }
-  // finally, compute the scaled, translated,
-  // projected quadrature points
-  for (unsigned int p = 0; p < quadrature.size(); ++p)
-    {
-      q_points[p][xi_index] =
-        xi_scale * quadrature.point(p)(0) + xi_translation;
-      q_points[p][eta_index] =
-        eta_scale * quadrature.point(p)(1) + eta_translation;
-      q_points[p][const_index] = const_value;
-    }
-}
-
-
-template <>
-Quadrature<1>
-QProjector<1>::project_to_all_faces(const Quadrature<0> &quadrature)
-{
-  const unsigned int dim = 1;
-
-  const unsigned int n_points = 1, n_faces = GeometryInfo<dim>::faces_per_cell;
-
-  // first fix quadrature points
-  std::vector<Point<dim>> q_points;
-  q_points.reserve(n_points * n_faces);
-  std::vector<Point<dim>> help(n_points);
-
-
-  // project to each face and append
-  // results
-  for (unsigned int face = 0; face < n_faces; ++face)
-    {
-      project_to_face(quadrature, face, help);
-      std::copy(help.begin(), help.end(), std::back_inserter(q_points));
-    }
-
-  // next copy over weights
-  std::vector<double> weights;
-  weights.reserve(n_points * n_faces);
-  for (unsigned int face = 0; face < n_faces; ++face)
-    std::copy(quadrature.get_weights().begin(),
-              quadrature.get_weights().end(),
-              std::back_inserter(weights));
-
-  Assert(q_points.size() == n_points * n_faces, ExcInternalError());
-  Assert(weights.size() == n_points * n_faces, ExcInternalError());
-
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-
-template <>
-Quadrature<2>
-QProjector<2>::project_to_all_faces(const SubQuadrature &quadrature)
-{
-  const unsigned int dim = 2;
-
-  const unsigned int n_points = quadrature.size(),
-                     n_faces  = GeometryInfo<dim>::faces_per_cell;
-
-  // first fix quadrature points
-  std::vector<Point<dim>> q_points;
-  q_points.reserve(n_points * n_faces);
-  std::vector<Point<dim>> help(n_points);
-
-  // project to each face and append
-  // results
-  for (unsigned int face = 0; face < n_faces; ++face)
-    {
-      project_to_face(quadrature, face, help);
-      std::copy(help.begin(), help.end(), std::back_inserter(q_points));
-    }
-
-  // next copy over weights
-  std::vector<double> weights;
-  weights.reserve(n_points * n_faces);
-  for (unsigned int face = 0; face < n_faces; ++face)
-    std::copy(quadrature.get_weights().begin(),
-              quadrature.get_weights().end(),
-              std::back_inserter(weights));
-
-  Assert(q_points.size() == n_points * n_faces, ExcInternalError());
-  Assert(weights.size() == n_points * n_faces, ExcInternalError());
-
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-
-template <>
-Quadrature<3>
-QProjector<3>::project_to_all_faces(const SubQuadrature &quadrature)
-{
-  const unsigned int dim = 3;
-
-  SubQuadrature q_reflected = reflect(quadrature);
-  SubQuadrature q[8]        = {quadrature,
-                        rotate(quadrature, 1),
-                        rotate(quadrature, 2),
-                        rotate(quadrature, 3),
-                        q_reflected,
-                        rotate(q_reflected, 3),
-                        rotate(q_reflected, 2),
-                        rotate(q_reflected, 1)};
-
-
-
-  const unsigned int n_points = quadrature.size(),
-                     n_faces  = GeometryInfo<dim>::faces_per_cell;
-
-  // first fix quadrature points
-  std::vector<Point<dim>> q_points;
-  q_points.reserve(n_points * n_faces * 8);
-  std::vector<Point<dim>> help(n_points);
-
-  std::vector<double> weights;
-  weights.reserve(n_points * n_faces * 8);
-
-  // do the following for all possible
-  // mutations of a face (mutation==0
-  // corresponds to a face with standard
-  // orientation, no flip and no rotation)
-  for (const auto &mutation : q)
-    {
-      // project to each face and append
-      // results
-      for (unsigned int face = 0; face < n_faces; ++face)
-        {
-          project_to_face(mutation, face, help);
-          std::copy(help.begin(), help.end(), std::back_inserter(q_points));
-        }
-
-      // next copy over weights
-      for (unsigned int face = 0; face < n_faces; ++face)
-        std::copy(mutation.get_weights().begin(),
-                  mutation.get_weights().end(),
-                  std::back_inserter(weights));
-    }
-
-
-  Assert(q_points.size() == n_points * n_faces * 8, ExcInternalError());
-  Assert(weights.size() == n_points * n_faces * 8, ExcInternalError());
-
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-
-template <>
-Quadrature<1>
-QProjector<1>::project_to_all_subfaces(const Quadrature<0> &quadrature)
-{
-  const unsigned int dim = 1;
-
-  const unsigned int n_points = 1, n_faces = GeometryInfo<dim>::faces_per_cell,
-                     subfaces_per_face =
-                       GeometryInfo<dim>::max_children_per_face;
-
-  // first fix quadrature points
-  std::vector<Point<dim>> q_points;
-  q_points.reserve(n_points * n_faces * subfaces_per_face);
-  std::vector<Point<dim>> help(n_points);
-
-  // project to each face and copy
-  // results
-  for (unsigned int face = 0; face < n_faces; ++face)
-    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
-      {
-        project_to_subface(quadrature, face, subface, help);
-        std::copy(help.begin(), help.end(), std::back_inserter(q_points));
-      }
-
-  // next copy over weights
-  std::vector<double> weights;
-  weights.reserve(n_points * n_faces * subfaces_per_face);
-  for (unsigned int face = 0; face < n_faces; ++face)
-    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
-      std::copy(quadrature.get_weights().begin(),
-                quadrature.get_weights().end(),
-                std::back_inserter(weights));
-
-  Assert(q_points.size() == n_points * n_faces * subfaces_per_face,
-         ExcInternalError());
-  Assert(weights.size() == n_points * n_faces * subfaces_per_face,
-         ExcInternalError());
-
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-
-template <>
-Quadrature<2>
-QProjector<2>::project_to_all_subfaces(const SubQuadrature &quadrature)
-{
-  const unsigned int dim = 2;
-
-  const unsigned int n_points = quadrature.size(),
-                     n_faces  = GeometryInfo<dim>::faces_per_cell,
-                     subfaces_per_face =
-                       GeometryInfo<dim>::max_children_per_face;
-
-  // first fix quadrature points
-  std::vector<Point<dim>> q_points;
-  q_points.reserve(n_points * n_faces * subfaces_per_face);
-  std::vector<Point<dim>> help(n_points);
-
-  // project to each face and copy
-  // results
-  for (unsigned int face = 0; face < n_faces; ++face)
-    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
-      {
-        project_to_subface(quadrature, face, subface, help);
-        std::copy(help.begin(), help.end(), std::back_inserter(q_points));
-      }
-
-  // next copy over weights
-  std::vector<double> weights;
-  weights.reserve(n_points * n_faces * subfaces_per_face);
-  for (unsigned int face = 0; face < n_faces; ++face)
-    for (unsigned int subface = 0; subface < subfaces_per_face; ++subface)
-      std::copy(quadrature.get_weights().begin(),
-                quadrature.get_weights().end(),
-                std::back_inserter(weights));
-
-  Assert(q_points.size() == n_points * n_faces * subfaces_per_face,
-         ExcInternalError());
-  Assert(weights.size() == n_points * n_faces * subfaces_per_face,
-         ExcInternalError());
-
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-
-template <>
-Quadrature<3>
-QProjector<3>::project_to_all_subfaces(const SubQuadrature &quadrature)
-{
-  const unsigned int dim         = 3;
-  SubQuadrature      q_reflected = reflect(quadrature);
-  SubQuadrature      q[8]        = {quadrature,
-                        rotate(quadrature, 1),
-                        rotate(quadrature, 2),
-                        rotate(quadrature, 3),
-                        q_reflected,
-                        rotate(q_reflected, 3),
-                        rotate(q_reflected, 2),
-                        rotate(q_reflected, 1)};
-
-  const unsigned int n_points = quadrature.size(),
-                     n_faces  = GeometryInfo<dim>::faces_per_cell,
-                     total_subfaces_per_face = 2 + 2 + 4;
-
-  // first fix quadrature points
-  std::vector<Point<dim>> q_points;
-  q_points.reserve(n_points * n_faces * total_subfaces_per_face * 8);
-  std::vector<Point<dim>> help(n_points);
-
-  std::vector<double> weights;
-  weights.reserve(n_points * n_faces * total_subfaces_per_face * 8);
-
-  // do the following for all possible
-  // mutations of a face (mutation==0
-  // corresponds to a face with standard
-  // orientation, no flip and no rotation)
-  for (const auto &mutation : q)
-    {
-      // project to each face and copy
-      // results
-      for (unsigned int face = 0; face < n_faces; ++face)
-        for (unsigned int ref_case = RefinementCase<dim - 1>::cut_xy;
-             ref_case >= RefinementCase<dim - 1>::cut_x;
-             --ref_case)
-          for (unsigned int subface = 0;
-               subface < GeometryInfo<dim - 1>::n_children(
-                           RefinementCase<dim - 1>(ref_case));
-               ++subface)
-            {
-              project_to_subface(mutation,
-                                 face,
-                                 subface,
-                                 help,
-                                 RefinementCase<dim - 1>(ref_case));
-              std::copy(help.begin(), help.end(), std::back_inserter(q_points));
-            }
-
-      // next copy over weights
-      for (unsigned int face = 0; face < n_faces; ++face)
-        for (unsigned int ref_case = RefinementCase<dim - 1>::cut_xy;
-             ref_case >= RefinementCase<dim - 1>::cut_x;
-             --ref_case)
-          for (unsigned int subface = 0;
-               subface < GeometryInfo<dim - 1>::n_children(
-                           RefinementCase<dim - 1>(ref_case));
-               ++subface)
-            std::copy(mutation.get_weights().begin(),
-                      mutation.get_weights().end(),
-                      std::back_inserter(weights));
-    }
-
-  Assert(q_points.size() == n_points * n_faces * total_subfaces_per_face * 8,
-         ExcInternalError());
-  Assert(weights.size() == n_points * n_faces * total_subfaces_per_face * 8,
-         ExcInternalError());
-
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-
-// This function is not used in the library
-template <int dim>
-Quadrature<dim>
-QProjector<dim>::project_to_child(const Quadrature<dim> &quadrature,
-                                  const unsigned int     child_no)
-{
-  AssertIndexRange(child_no, GeometryInfo<dim>::max_children_per_cell);
-
-  const unsigned int n_q_points = quadrature.size();
-
-  std::vector<Point<dim>> q_points(n_q_points);
-  for (unsigned int i = 0; i < n_q_points; ++i)
-    q_points[i] =
-      GeometryInfo<dim>::child_to_cell_coordinates(quadrature.point(i),
-                                                   child_no);
-
-  // for the weights, things are
-  // equally simple: copy them and
-  // scale them
-  std::vector<double> weights = quadrature.get_weights();
-  for (unsigned int i = 0; i < n_q_points; ++i)
-    weights[i] *= (1. / GeometryInfo<dim>::max_children_per_cell);
-
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-template <int dim>
-Quadrature<dim>
-QProjector<dim>::project_to_all_children(const Quadrature<dim> &quadrature)
-{
-  const unsigned int n_points   = quadrature.size(),
-                     n_children = GeometryInfo<dim>::max_children_per_cell;
-
-  std::vector<Point<dim>> q_points(n_points * n_children);
-  std::vector<double>     weights(n_points * n_children);
-
-  // project to each child and copy
-  // results
-  for (unsigned int child = 0; child < n_children; ++child)
-    {
-      Quadrature<dim> help = project_to_child(quadrature, child);
-      for (unsigned int i = 0; i < n_points; ++i)
-        {
-          q_points[child * n_points + i] = help.point(i);
-          weights[child * n_points + i]  = help.weight(i);
-        }
-    }
-  return Quadrature<dim>(q_points, weights);
-}
-
-
-
-template <int dim>
-Quadrature<dim>
-QProjector<dim>::project_to_line(const Quadrature<1> &quadrature,
-                                 const Point<dim> &   p1,
-                                 const Point<dim> &   p2)
-{
-  const unsigned int      n = quadrature.size();
-  std::vector<Point<dim>> points(n);
-  std::vector<double>     weights(n);
-  const double            length = p1.distance(p2);
-
-  for (unsigned int k = 0; k < n; ++k)
-    {
-      const double alpha = quadrature.point(k)(0);
-      points[k]          = alpha * p2;
-      points[k] += (1. - alpha) * p1;
-      weights[k] = length * quadrature.weight(k);
-    }
-  return Quadrature<dim>(points, weights);
-}
-
-
-
-template <int dim>
-typename QProjector<dim>::DataSetDescriptor
-QProjector<dim>::DataSetDescriptor::face(const unsigned int face_no,
-                                         const bool         face_orientation,
-                                         const bool         face_flip,
-                                         const bool         face_rotation,
-                                         const unsigned int n_quadrature_points)
-{
-  Assert(face_no < GeometryInfo<dim>::faces_per_cell, ExcInternalError());
-
-  switch (dim)
-    {
-      case 1:
-      case 2:
-        return face_no * n_quadrature_points;
-
-
-      case 3:
-        {
-          // in 3d, we have to account for faces that
-          // have non-standard face orientation, flip
-          // and rotation. thus, we have to store
-          // _eight_ data sets per face or subface
-
-          // set up a table with the according offsets
-          // for non-standard orientation, first index:
-          // face_orientation (standard true=1), second
-          // index: face_flip (standard false=0), third
-          // index: face_rotation (standard false=0)
-          //
-          // note, that normally we should use the
-          // obvious offsets 0,1,2,3,4,5,6,7. However,
-          // prior to the changes enabling flipped and
-          // rotated faces, in many places of the
-          // library the convention was used, that the
-          // first dataset with offset 0 corresponds to
-          // a face in standard orientation. therefore
-          // we use the offsets 4,5,6,7,0,1,2,3 here to
-          // stick to that (implicit) convention
-          static const unsigned int offset[2][2][2] = {
-            {{4 * GeometryInfo<dim>::faces_per_cell,
-              5 * GeometryInfo<dim>::
-                    faces_per_cell}, // face_orientation=false; face_flip=false;
-                                     // face_rotation=false and true
-             {6 * GeometryInfo<dim>::faces_per_cell,
-              7 * GeometryInfo<dim>::
-                    faces_per_cell}}, // face_orientation=false; face_flip=true;
-                                      // face_rotation=false and true
-            {{0 * GeometryInfo<dim>::faces_per_cell,
-              1 * GeometryInfo<dim>::
-                    faces_per_cell}, // face_orientation=true;  face_flip=false;
-                                     // face_rotation=false and true
-             {2 * GeometryInfo<dim>::faces_per_cell,
-              3 * GeometryInfo<dim>::
-                    faces_per_cell}}}; // face_orientation=true; face_flip=true;
-                                       // face_rotation=false and true
-
-          return (
-            (face_no + offset[face_orientation][face_flip][face_rotation]) *
-            n_quadrature_points);
-        }
-
-      default:
-        Assert(false, ExcInternalError());
-    }
-  return numbers::invalid_unsigned_int;
-}
-
-
-
-template <>
-QProjector<1>::DataSetDescriptor
-QProjector<1>::DataSetDescriptor::subface(
-  const unsigned int face_no,
-  const unsigned int subface_no,
-  const bool,
-  const bool,
-  const bool,
-  const unsigned int n_quadrature_points,
-  const internal::SubfaceCase<1>)
-{
-  Assert(face_no < GeometryInfo<1>::faces_per_cell, ExcInternalError());
-  Assert(subface_no < GeometryInfo<1>::max_children_per_face,
-         ExcInternalError());
-
-  return ((face_no * GeometryInfo<1>::max_children_per_face + subface_no) *
-          n_quadrature_points);
-}
-
-
-
-template <>
-QProjector<2>::DataSetDescriptor
-QProjector<2>::DataSetDescriptor::subface(
-  const unsigned int face_no,
-  const unsigned int subface_no,
-  const bool,
-  const bool,
-  const bool,
-  const unsigned int n_quadrature_points,
-  const internal::SubfaceCase<2>)
-{
-  Assert(face_no < GeometryInfo<2>::faces_per_cell, ExcInternalError());
-  Assert(subface_no < GeometryInfo<2>::max_children_per_face,
-         ExcInternalError());
-
-  return ((face_no * GeometryInfo<2>::max_children_per_face + subface_no) *
-          n_quadrature_points);
-}
-
-
-template <>
-QProjector<3>::DataSetDescriptor
-QProjector<3>::DataSetDescriptor::subface(
-  const unsigned int             face_no,
-  const unsigned int             subface_no,
-  const bool                     face_orientation,
-  const bool                     face_flip,
-  const bool                     face_rotation,
-  const unsigned int             n_quadrature_points,
-  const internal::SubfaceCase<3> ref_case)
-{
-  const unsigned int dim = 3;
-
-  Assert(face_no < GeometryInfo<dim>::faces_per_cell, ExcInternalError());
-  Assert(subface_no < GeometryInfo<dim>::max_children_per_face,
-         ExcInternalError());
-
-  // As the quadrature points created by
-  // QProjector are on subfaces in their
-  // "standard location" we have to use a
-  // permutation of the equivalent subface
-  // number in order to respect face
-  // orientation, flip and rotation. The
-  // information we need here is exactly the
-  // same as the
-  // GeometryInfo<3>::child_cell_on_face info
-  // for the bottom face (face 4) of a hex, as
-  // on this the RefineCase of the cell matches
-  // that of the face and the subfaces are
-  // numbered in the same way as the child
-  // cells.
-
-  // in 3d, we have to account for faces that
-  // have non-standard face orientation, flip
-  // and rotation. thus, we have to store
-  // _eight_ data sets per face or subface
-  // already for the isotropic
-  // case. Additionally, we have three
-  // different refinement cases, resulting in
-  // <tt>4 + 2 + 2 = 8</tt> different subfaces
-  // for each face.
-  const unsigned int total_subfaces_per_face = 8;
-
-  // set up a table with the according offsets
-  // for non-standard orientation, first index:
-  // face_orientation (standard true=1), second
-  // index: face_flip (standard false=0), third
-  // index: face_rotation (standard false=0)
-  //
-  // note, that normally we should use the
-  // obvious offsets 0,1,2,3,4,5,6,7. However,
-  // prior to the changes enabling flipped and
-  // rotated faces, in many places of the
-  // library the convention was used, that the
-  // first dataset with offset 0 corresponds to
-  // a face in standard orientation. therefore
-  // we use the offsets 4,5,6,7,0,1,2,3 here to
-  // stick to that (implicit) convention
-  static const unsigned int orientation_offset[2][2][2] = {
-    {// face_orientation=false; face_flip=false; face_rotation=false and true
-     {4 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
-      5 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face},
-     // face_orientation=false; face_flip=true;  face_rotation=false and true
-     {6 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
-      7 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face}},
-    {// face_orientation=true;  face_flip=false; face_rotation=false and true
-     {0 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
-      1 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face},
-     // face_orientation=true;  face_flip=true;  face_rotation=false and true
-     {2 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face,
-      3 * GeometryInfo<dim>::faces_per_cell * total_subfaces_per_face}}};
-
-  // set up a table with the offsets for a
-  // given refinement case respecting the
-  // corresponding number of subfaces. the
-  // index corresponds to (RefineCase::Type - 1)
-
-  // note, that normally we should use the
-  // obvious offsets 0,2,6. However, prior to
-  // the implementation of anisotropic
-  // refinement, in many places of the library
-  // the convention was used, that the first
-  // dataset with offset 0 corresponds to a
-  // standard (isotropic) face
-  // refinement. therefore we use the offsets
-  // 6,4,0 here to stick to that (implicit)
-  // convention
-  static const unsigned int ref_case_offset[3] = {
-    6, // cut_x
-    4, // cut_y
-    0  // cut_xy
-  };
-
-
-  // for each subface of a given FaceRefineCase
-  // there is a corresponding equivalent
-  // subface number of one of the "standard"
-  // RefineCases (cut_x, cut_y, cut_xy). Map
-  // the given values to those equivalent
-  // ones.
-
-  // first, define an invalid number
-  static const unsigned int e = numbers::invalid_unsigned_int;
-
-  static const RefinementCase<dim - 1>
-    equivalent_refine_case[internal::SubfaceCase<dim>::case_isotropic + 1]
-                          [GeometryInfo<3>::max_children_per_face] = {
-                            // case_none. there should be only
-                            // invalid values here. However, as
-                            // this function is also called (in
-                            // tests) for cells which have no
-                            // refined faces, use isotropic
-                            // refinement instead
-                            {RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy},
-                            // case_x
-                            {RefinementCase<dim - 1>::cut_x,
-                             RefinementCase<dim - 1>::cut_x,
-                             RefinementCase<dim - 1>::no_refinement,
-                             RefinementCase<dim - 1>::no_refinement},
-                            // case_x1y
-                            {RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_x,
-                             RefinementCase<dim - 1>::no_refinement},
-                            // case_x2y
-                            {RefinementCase<dim - 1>::cut_x,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::no_refinement},
-                            // case_x1y2y
-                            {RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy},
-                            // case_y
-                            {RefinementCase<dim - 1>::cut_y,
-                             RefinementCase<dim - 1>::cut_y,
-                             RefinementCase<dim - 1>::no_refinement,
-                             RefinementCase<dim - 1>::no_refinement},
-                            // case_y1x
-                            {RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_y,
-                             RefinementCase<dim - 1>::no_refinement},
-                            // case_y2x
-                            {RefinementCase<dim - 1>::cut_y,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::no_refinement},
-                            // case_y1x2x
-                            {RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy},
-                            // case_xy (case_isotropic)
-                            {RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy,
-                             RefinementCase<dim - 1>::cut_xy}};
-
-  static const unsigned int
-    equivalent_subface_number[internal::SubfaceCase<dim>::case_isotropic + 1]
-                             [GeometryInfo<3>::max_children_per_face] = {
-                               // case_none, see above
-                               {0, 1, 2, 3},
-                               // case_x
-                               {0, 1, e, e},
-                               // case_x1y
-                               {0, 2, 1, e},
-                               // case_x2y
-                               {0, 1, 3, e},
-                               // case_x1y2y
-                               {0, 2, 1, 3},
-                               // case_y
-                               {0, 1, e, e},
-                               // case_y1x
-                               {0, 1, 1, e},
-                               // case_y2x
-                               {0, 2, 3, e},
-                               // case_y1x2x
-                               {0, 1, 2, 3},
-                               // case_xy (case_isotropic)
-                               {0, 1, 2, 3}};
-
-  // If face-orientation or face_rotation are
-  // non-standard, cut_x and cut_y have to be
-  // exchanged.
-  static const RefinementCase<dim - 1> ref_case_permutation[4] = {
-    RefinementCase<dim - 1>::no_refinement,
-    RefinementCase<dim - 1>::cut_y,
-    RefinementCase<dim - 1>::cut_x,
-    RefinementCase<dim - 1>::cut_xy};
-
-  // set a corresponding (equivalent)
-  // RefineCase and subface number
-  const RefinementCase<dim - 1> equ_ref_case =
-    equivalent_refine_case[ref_case][subface_no];
-  const unsigned int equ_subface_no =
-    equivalent_subface_number[ref_case][subface_no];
-  // make sure, that we got a valid subface and RefineCase
-  Assert(equ_ref_case != RefinementCase<dim>::no_refinement,
-         ExcInternalError());
-  Assert(equ_subface_no != e, ExcInternalError());
-  // now, finally respect non-standard faces
-  const RefinementCase<dim - 1> final_ref_case =
-    (face_orientation == face_rotation ? ref_case_permutation[equ_ref_case] :
-                                         equ_ref_case);
-
-  // what we have now is the number of
-  // the subface in the natural
-  // orientation of the *face*. what we
-  // need to know is the number of the
-  // subface concerning the standard face
-  // orientation as seen from the *cell*.
-
-  // this mapping is not trivial, but we
-  // have done exactly this stuff in the
-  // child_cell_on_face function. in
-  // order to reduce the amount of code
-  // as well as to make maintaining the
-  // functionality easier we want to
-  // reuse that information. So we note
-  // that on the bottom face (face 4) of
-  // a hex cell the local x and y
-  // coordinates of the face and the cell
-  // coincide, thus also the refinement
-  // case of the face corresponds to the
-  // refinement case of the cell
-  // (ignoring cell refinement along the
-  // z direction). Using this knowledge
-  // we can (ab)use the
-  // child_cell_on_face function to do
-  // exactly the transformation we are in
-  // need of now
-  const unsigned int final_subface_no =
-    GeometryInfo<dim>::child_cell_on_face(RefinementCase<dim>(final_ref_case),
-                                          4,
-                                          equ_subface_no,
-                                          face_orientation,
-                                          face_flip,
-                                          face_rotation,
-                                          equ_ref_case);
-
-  return (((face_no * total_subfaces_per_face +
-            ref_case_offset[final_ref_case - 1] + final_subface_no) +
-           orientation_offset[face_orientation][face_flip][face_rotation]) *
-          n_quadrature_points);
-}
-
-
-
-template <int dim>
-Quadrature<dim>
-QProjector<dim>::project_to_face(const SubQuadrature &quadrature,
-                                 const unsigned int   face_no)
-{
-  std::vector<Point<dim>> points(quadrature.size());
-  project_to_face(quadrature, face_no, points);
-  return Quadrature<dim>(points, quadrature.get_weights());
-}
-
-
-template <int dim>
-Quadrature<dim>
-QProjector<dim>::project_to_subface(const SubQuadrature &          quadrature,
-                                    const unsigned int             face_no,
-                                    const unsigned int             subface_no,
-                                    const RefinementCase<dim - 1> &ref_case)
-{
-  std::vector<Point<dim>> points(quadrature.size());
-  project_to_subface(quadrature, face_no, subface_no, points, ref_case);
-  return Quadrature<dim>(points, quadrature.get_weights());
-}
-
-
 // ------------------------------------------------------------ //
 
 namespace internal
@@ -1792,8 +693,5 @@ template class QAnisotropic<3>;
 template class QIterated<1>;
 template class QIterated<2>;
 template class QIterated<3>;
-template class QProjector<1>;
-template class QProjector<2>;
-template class QProjector<3>;
 
 DEAL_II_NAMESPACE_CLOSE
index f95944c77af5fc01984f6091bdb480adf315b17f..abf401520489fcef44517b34f8dc76c65974abff 100644 (file)
@@ -1218,7 +1218,8 @@ FiniteElement<dim, spacedim>::get_face_data(
 {
   return get_data(flags,
                   mapping,
-                  QProjector<dim>::project_to_all_faces(quadrature),
+                  QProjector<dim>::project_to_all_faces(
+                    this->reference_cell_type(), quadrature),
                   output_data);
 }
 
@@ -1236,7 +1237,8 @@ FiniteElement<dim, spacedim>::get_subface_data(
 {
   return get_data(flags,
                   mapping,
-                  QProjector<dim>::project_to_all_subfaces(quadrature),
+                  QProjector<dim>::project_to_all_subfaces(
+                    this->reference_cell_type(), quadrature),
                   output_data);
 }
 
index 432f8918062d85861bf75ee07dc07cb8ce18cfdc..83732630e5bbe8b863ed0878d6ec159bad317777 100644 (file)
@@ -207,7 +207,8 @@ FE_ABF<dim>::initialize_support_points(const unsigned int deg)
         }
 
       Quadrature<dim> faces =
-        QProjector<dim>::project_to_all_faces(face_points);
+        QProjector<dim>::project_to_all_faces(this->reference_cell_type(),
+                                              face_points);
       for (; current < GeometryInfo<dim>::faces_per_cell * n_face_points;
            ++current)
         {
@@ -339,7 +340,10 @@ FE_ABF<dim>::initialize_restriction()
       // child cell are evaluated
       // in the quadrature points
       // of a full face.
-      Quadrature<dim> q_face = QProjector<dim>::project_to_face(q_base, face);
+      Quadrature<dim> q_face =
+        QProjector<dim>::project_to_face(this->reference_cell_type(),
+                                         q_base,
+                                         face);
       // Store shape values, since the
       // evaluation suffers if not
       // ordered by point
@@ -357,8 +361,8 @@ FE_ABF<dim>::initialize_restriction()
           // the coarse face are
           // evaluated on the subface
           // only.
-          Quadrature<dim> q_sub =
-            QProjector<dim>::project_to_subface(q_base, face, sub);
+          Quadrature<dim> q_sub = QProjector<dim>::project_to_subface(
+            this->reference_cell_type(), q_base, face, sub);
           const unsigned int child = GeometryInfo<dim>::child_cell_on_face(
             RefinementCase<dim>::isotropic_refinement, face, sub);
 
@@ -434,7 +438,10 @@ FE_ABF<dim>::initialize_restriction()
   for (unsigned int child = 0; child < GeometryInfo<dim>::max_children_per_cell;
        ++child)
     {
-      Quadrature<dim> q_sub = QProjector<dim>::project_to_child(q_cell, child);
+      Quadrature<dim> q_sub =
+        QProjector<dim>::project_to_child(this->reference_cell_type(),
+                                          q_cell,
+                                          child);
 
       for (unsigned int k = 0; k < q_sub.size(); ++k)
         for (unsigned int i_child = 0; i_child < this->n_dofs_per_cell();
@@ -592,7 +599,12 @@ FE_ABF<dim>::convert_generalized_support_point_values_to_dof_values(
           // TODO: Check what the face_orientation, face_flip and face_rotation
           // have to be in 3D
           unsigned int k = QProjector<dim>::DataSetDescriptor::face(
-            face, false, false, false, n_face_points);
+            this->reference_cell_type(),
+            face,
+            false,
+            false,
+            false,
+            n_face_points);
           for (unsigned int i = 0; i < boundary_weights_abf.size(1); ++i)
             nodal_values[start_abf_dofs + i] +=
               n_orient * boundary_weights_abf(k + fp, i) *
index 264203cfb329a4cce723775da7ca33a3e701958a..a85ede4162464f433e1cc80a04d136dcb00f8d46 100644 (file)
@@ -349,13 +349,19 @@ FE_BDM<dim>::initialize_support_points(const unsigned int deg)
 
   this->generalized_support_points.resize(npoints);
 
-  Quadrature<dim> faces = QProjector<dim>::project_to_all_faces(face_points);
+  Quadrature<dim> faces =
+    QProjector<dim>::project_to_all_faces(this->reference_cell_type(),
+                                          face_points);
   for (unsigned int k = 0;
        k < face_points.size() * GeometryInfo<dim>::faces_per_cell;
        ++k)
-    this->generalized_support_points[k] =
-      faces.point(k + QProjector<dim>::DataSetDescriptor::face(
-                        0, true, false, false, this->n_dofs_per_face()));
+    this->generalized_support_points[k] = faces.point(
+      k + QProjector<dim>::DataSetDescriptor::face(this->reference_cell_type(),
+                                                   0,
+                                                   true,
+                                                   false,
+                                                   false,
+                                                   this->n_dofs_per_face()));
 
   // Currently, for backward compatibility, we do not use moments, but
   // point values on faces in 2D. In 3D, this is impossible, since the
index f8078b51730cf7ef3d7a284f253f6fc2f361b7d0..587301ecbbd87447bf915c32f08d5b1df43371fa 100644 (file)
@@ -152,8 +152,13 @@ FE_Bernstein<dim, spacedim>::get_subface_interpolation_matrix(
       // be done for the face orientation flag in 3D.
       const Quadrature<dim> subface_quadrature =
         subface == numbers::invalid_unsigned_int ?
-          QProjector<dim>::project_to_face(quad_face_support, 0) :
-          QProjector<dim>::project_to_subface(quad_face_support, 0, subface);
+          QProjector<dim>::project_to_face(this->reference_cell_type(),
+                                           quad_face_support,
+                                           0) :
+          QProjector<dim>::project_to_subface(this->reference_cell_type(),
+                                              quad_face_support,
+                                              0,
+                                              subface);
 
       for (unsigned int i = 0; i < source_fe->n_dofs_per_face(); ++i)
         {
index a6b9717f791356b76276a882fbced7167f211d12..9645e6a0c9748cd102b257621dec6ff4baddbf8e 100644 (file)
@@ -270,7 +270,8 @@ FE_Nedelec<2>::initialize_support_points(const unsigned int order)
   const unsigned int    n_boundary_points =
     GeometryInfo<dim>::lines_per_cell * n_edge_points;
   const Quadrature<dim> edge_quadrature =
-    QProjector<dim>::project_to_all_faces(reference_edge_quadrature);
+    QProjector<dim>::project_to_all_faces(this->reference_cell_type(),
+                                          reference_edge_quadrature);
 
   this->generalized_face_support_points.resize(n_edge_points);
 
@@ -297,7 +298,12 @@ FE_Nedelec<2>::initialize_support_points(const unsigned int order)
                ++line)
             this->generalized_support_points[line * n_edge_points + q_point] =
               edge_quadrature.point(QProjector<dim>::DataSetDescriptor::face(
-                                      line, true, false, false, n_edge_points) +
+                                      this->reference_cell_type(),
+                                      line,
+                                      true,
+                                      false,
+                                      false,
+                                      n_edge_points) +
                                     q_point);
 
           for (unsigned int i = 0; i < order; ++i)
@@ -323,7 +329,12 @@ FE_Nedelec<2>::initialize_support_points(const unsigned int order)
         for (unsigned int q_point = 0; q_point < n_edge_points; ++q_point)
           this->generalized_support_points[line * n_edge_points + q_point] =
             edge_quadrature.point(QProjector<dim>::DataSetDescriptor::face(
-                                    line, true, false, false, n_edge_points) +
+                                    this->reference_cell_type(),
+                                    line,
+                                    true,
+                                    false,
+                                    false,
+                                    n_edge_points) +
                                   q_point);
     }
 }
@@ -350,7 +361,9 @@ FE_Nedelec<3>::initialize_support_points(const unsigned int order)
   const QGauss<1>            reference_edge_quadrature(order + 1);
   const unsigned int         n_edge_points = reference_edge_quadrature.size();
   const Quadrature<dim - 1> &edge_quadrature =
-    QProjector<dim - 1>::project_to_all_faces(reference_edge_quadrature);
+    QProjector<dim - 1>::project_to_all_faces(ReferenceCell::get_hypercube(dim -
+                                                                           1),
+                                              reference_edge_quadrature);
 
   if (order > 0)
     {
@@ -382,7 +395,12 @@ FE_Nedelec<3>::initialize_support_points(const unsigned int order)
                                                   q_point] =
               edge_quadrature.point(
                 QProjector<dim - 1>::DataSetDescriptor::face(
-                  line, true, false, false, n_edge_points) +
+                  ReferenceCell::get_hypercube(dim - 1),
+                  line,
+                  true,
+                  false,
+                  false,
+                  n_edge_points) +
                 q_point);
 
           for (unsigned int i = 0; i < 2; ++i)
@@ -440,7 +458,8 @@ FE_Nedelec<3>::initialize_support_points(const unsigned int order)
         }
 
       const Quadrature<dim> &face_quadrature =
-        QProjector<dim>::project_to_all_faces(reference_face_quadrature);
+        QProjector<dim>::project_to_all_faces(this->reference_cell_type(),
+                                              reference_face_quadrature);
 
       for (const unsigned int face : GeometryInfo<dim>::face_indices())
         for (unsigned int q_point = 0; q_point < n_face_points; ++q_point)
@@ -449,7 +468,12 @@ FE_Nedelec<3>::initialize_support_points(const unsigned int order)
                                              GeometryInfo<dim>::lines_per_cell *
                                                n_edge_points] =
               face_quadrature.point(QProjector<dim>::DataSetDescriptor::face(
-                                      face, true, false, false, n_face_points) +
+                                      this->reference_cell_type(),
+                                      face,
+                                      true,
+                                      false,
+                                      false,
+                                      n_face_points) +
                                     q_point);
           }
 
@@ -474,7 +498,12 @@ FE_Nedelec<3>::initialize_support_points(const unsigned int order)
                                                   q_point] =
               edge_quadrature.point(
                 QProjector<dim - 1>::DataSetDescriptor::face(
-                  line, true, false, false, n_edge_points) +
+                  ReferenceCell::get_hypercube(dim - 1),
+                  line,
+                  true,
+                  false,
+                  false,
+                  n_edge_points) +
                 q_point);
 
           for (unsigned int i = 0; i < 2; ++i)
index ebf33c332b4200e3a149dd4f670d7bdc8baea554..7a9d9bb646cdd3dbc7cd1fb4ce323807708599da 100644 (file)
@@ -2056,19 +2056,22 @@ FE_NedelecSZ<dim, spacedim>::fill_fe_face_values(
   // This will fill in the missing items in the InternalData
   // (fe_internal/fe_data) which was not filled in by get_data.
   fill_edge_values(cell,
-                   QProjector<dim>::project_to_all_faces(quadrature),
+                   QProjector<dim>::project_to_all_faces(
+                     this->reference_cell_type(), quadrature),
                    fe_data);
   if (dim == 3 && this->degree > 1)
     {
       fill_face_values(cell,
-                       QProjector<dim>::project_to_all_faces(quadrature),
+                       QProjector<dim>::project_to_all_faces(
+                         this->reference_cell_type(), quadrature),
                        fe_data);
     }
 
   const UpdateFlags  flags(fe_data.update_each);
   const unsigned int n_q_points = quadrature.size();
-  const typename QProjector<dim>::DataSetDescriptor offset =
-    QProjector<dim>::DataSetDescriptor::face(face_no,
+  const auto         offset =
+    QProjector<dim>::DataSetDescriptor::face(this->reference_cell_type(),
+                                             face_no,
                                              cell->face_orientation(face_no),
                                              cell->face_flip(face_no),
                                              cell->face_rotation(face_no),
index e66652bebc0d2e1b4429d38dac00bd4b43c4bb8e..e86950d1c7da899e42744f347166f973a7e2b1af 100644 (file)
@@ -262,7 +262,9 @@ FE_P1NC::fill_fe_face_values(
 
   // compute on the face
   const Quadrature<2> quadrature_on_face =
-    QProjector<2>::project_to_face(quadrature, face_no);
+    QProjector<2>::project_to_face(this->reference_cell_type(),
+                                   quadrature,
+                                   face_no);
 
   if (flags & update_values)
     for (unsigned int i = 0; i < quadrature_on_face.size(); ++i)
@@ -306,8 +308,8 @@ FE_P1NC::fill_fe_subface_values(
     get_linear_shape_coefficients(cell);
 
   // compute on the subface
-  const Quadrature<2> quadrature_on_subface =
-    QProjector<2>::project_to_subface(quadrature, face_no, sub_no);
+  const Quadrature<2> quadrature_on_subface = QProjector<2>::project_to_subface(
+    this->reference_cell_type(), quadrature, face_no, sub_no);
 
   if (flags & update_values)
     for (unsigned int i = 0; i < quadrature_on_subface.size(); ++i)
index 0f7abca45d54bda54d24eec8142fcf9da40bfcd9..d6f1794bf37713a32c4d4e535594aa94e5bdb999 100644 (file)
@@ -969,8 +969,9 @@ FE_PolyTensor<dim, spacedim>::fill_fe_face_values(
   // to take (all data sets for all
   // faces are stored contiguously)
 
-  const typename QProjector<dim>::DataSetDescriptor offset =
-    QProjector<dim>::DataSetDescriptor::face(face_no,
+  const auto offset =
+    QProjector<dim>::DataSetDescriptor::face(this->reference_cell_type(),
+                                             face_no,
                                              cell->face_orientation(face_no),
                                              cell->face_flip(face_no),
                                              cell->face_rotation(face_no),
@@ -1595,8 +1596,9 @@ FE_PolyTensor<dim, spacedim>::fill_fe_subface_values(
   // offset determines which data set
   // to take (all data sets for all
   // sub-faces are stored contiguously)
-  const typename QProjector<dim>::DataSetDescriptor offset =
-    QProjector<dim>::DataSetDescriptor::subface(face_no,
+  const auto offset =
+    QProjector<dim>::DataSetDescriptor::subface(this->reference_cell_type(),
+                                                face_no,
                                                 sub_no,
                                                 cell->face_orientation(face_no),
                                                 cell->face_flip(face_no),
index 03e69ae8a08d3320d1dde7a13b3c3411934921a6..d6664d6ea9f792ae33008ecace5c3eb7492b6834 100644 (file)
@@ -260,19 +260,23 @@ struct FE_Q_Base<PolynomialType, xdim, xspacedim>::Implementation
         // Add nodes of lines interior in the "mother-face"
 
         // line 5: use line 9
-        QProjector<dim - 1>::project_to_subface(qline, 0, 0, p_line);
+        QProjector<dim - 1>::project_to_subface(
+          ReferenceCell::get_hypercube(dim - 1), qline, 0, 0, p_line);
         for (unsigned int i = 0; i < n; ++i)
           constraint_points.push_back(p_line[i] + Point<dim - 1>(0.5, 0));
         // line 6: use line 10
-        QProjector<dim - 1>::project_to_subface(qline, 0, 1, p_line);
+        QProjector<dim - 1>::project_to_subface(
+          ReferenceCell::get_hypercube(dim - 1), qline, 0, 1, p_line);
         for (unsigned int i = 0; i < n; ++i)
           constraint_points.push_back(p_line[i] + Point<dim - 1>(0.5, 0));
         // line 7: use line 13
-        QProjector<dim - 1>::project_to_subface(qline, 2, 0, p_line);
+        QProjector<dim - 1>::project_to_subface(
+          ReferenceCell::get_hypercube(dim - 1), qline, 2, 0, p_line);
         for (unsigned int i = 0; i < n; ++i)
           constraint_points.push_back(p_line[i] + Point<dim - 1>(0, 0.5));
         // line 8: use line 14
-        QProjector<dim - 1>::project_to_subface(qline, 2, 1, p_line);
+        QProjector<dim - 1>::project_to_subface(
+          ReferenceCell::get_hypercube(dim - 1), qline, 2, 1, p_line);
         for (unsigned int i = 0; i < n; ++i)
           constraint_points.push_back(p_line[i] + Point<dim - 1>(0, 0.5));
 
@@ -284,10 +288,12 @@ struct FE_Q_Base<PolynomialType, xdim, xspacedim>::Implementation
                subface < GeometryInfo<dim - 1>::max_children_per_face;
                ++subface)
             {
-              QProjector<dim - 1>::project_to_subface(qline,
-                                                      face,
-                                                      subface,
-                                                      p_line);
+              QProjector<dim - 1>::project_to_subface(
+                ReferenceCell::get_hypercube(dim - 1),
+                qline,
+                face,
+                subface,
+                p_line);
               constraint_points.insert(constraint_points.end(),
                                        p_line.begin(),
                                        p_line.end());
@@ -658,8 +664,13 @@ FE_Q_Base<PolynomialType, dim, spacedim>::get_subface_interpolation_matrix(
       // be done for the face orientation flag in 3D.
       const Quadrature<dim> subface_quadrature =
         subface == numbers::invalid_unsigned_int ?
-          QProjector<dim>::project_to_face(quad_face_support, 0) :
-          QProjector<dim>::project_to_subface(quad_face_support, 0, subface);
+          QProjector<dim>::project_to_face(this->reference_cell_type(),
+                                           quad_face_support,
+                                           0) :
+          QProjector<dim>::project_to_subface(this->reference_cell_type(),
+                                              quad_face_support,
+                                              0,
+                                              subface);
       for (unsigned int i = 0; i < source_fe->n_dofs_per_face(); ++i)
         {
           const Point<dim> &p = subface_quadrature.point(i);
index dc38e9d687001a31ab87314b456089a91f0cc965..ce1b950f1478bd90578702d6c29c797eaf9affb8 100644 (file)
@@ -14,7 +14,6 @@
 // ---------------------------------------------------------------------
 
 
-#include <deal.II/base/qprojector.h>
 #include <deal.II/base/quadrature.h>
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/template_constraints.h>
index 3f1b209c35d1f2641bdb0f00df2a01438c74b38d..e85acf4dc51ff03d0556ef0994374b686b8d24ca 100644 (file)
@@ -14,7 +14,6 @@
 // ---------------------------------------------------------------------
 
 
-#include <deal.II/base/qprojector.h>
 #include <deal.II/base/quadrature.h>
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/template_constraints.h>
index 19e05820615f7ce2623e166e19ba5cc495f286a5..310806e5bc4f4d1150d8b7c41091f28fa35afb19 100644 (file)
@@ -191,7 +191,8 @@ FE_RaviartThomas<dim>::initialize_support_points(const unsigned int deg)
         }
 
       Quadrature<dim> faces =
-        QProjector<dim>::project_to_all_faces(face_points);
+        QProjector<dim>::project_to_all_faces(this->reference_cell_type(),
+                                              face_points);
       for (; current < GeometryInfo<dim>::faces_per_cell * n_face_points;
            ++current)
         {
@@ -199,7 +200,12 @@ FE_RaviartThomas<dim>::initialize_support_points(const unsigned int deg)
           // into the vector
           this->generalized_support_points[current] =
             faces.point(current + QProjector<dim>::DataSetDescriptor::face(
-                                    0, true, false, false, n_face_points));
+                                    this->reference_cell_type(),
+                                    0,
+                                    true,
+                                    false,
+                                    false,
+                                    n_face_points));
         }
     }
 
@@ -276,7 +282,10 @@ FE_RaviartThomas<dim>::initialize_restriction()
       // child cell are evaluated
       // in the quadrature points
       // of a full face.
-      Quadrature<dim> q_face = QProjector<dim>::project_to_face(q_base, face);
+      Quadrature<dim> q_face =
+        QProjector<dim>::project_to_face(this->reference_cell_type(),
+                                         q_base,
+                                         face);
       // Store shape values, since the
       // evaluation suffers if not
       // ordered by point
@@ -294,8 +303,8 @@ FE_RaviartThomas<dim>::initialize_restriction()
           // the coarse face are
           // evaluated on the subface
           // only.
-          Quadrature<dim> q_sub =
-            QProjector<dim>::project_to_subface(q_base, face, sub);
+          Quadrature<dim> q_sub = QProjector<dim>::project_to_subface(
+            this->reference_cell_type(), q_base, face, sub);
           const unsigned int child = GeometryInfo<dim>::child_cell_on_face(
             RefinementCase<dim>::isotropic_refinement, face, sub);
 
@@ -372,7 +381,10 @@ FE_RaviartThomas<dim>::initialize_restriction()
   for (unsigned int child = 0; child < GeometryInfo<dim>::max_children_per_cell;
        ++child)
     {
-      Quadrature<dim> q_sub = QProjector<dim>::project_to_child(q_cell, child);
+      Quadrature<dim> q_sub =
+        QProjector<dim>::project_to_child(this->reference_cell_type(),
+                                          q_cell,
+                                          child);
 
       for (unsigned int k = 0; k < q_sub.size(); ++k)
         for (unsigned int i_child = 0; i_child < this->n_dofs_per_cell();
index 19b150ef5d518a8735b99bb56e256b54681264d5..0ddeb835023c2b94bbf98ee87ccf1bda0f88082d 100644 (file)
@@ -166,13 +166,19 @@ FE_RaviartThomasNodal<dim>::initialize_support_points(const unsigned int deg)
       for (unsigned int k = 0; k < this->n_dofs_per_face(); ++k)
         this->generalized_face_support_points[k] = face_points.point(k);
       Quadrature<dim> faces =
-        QProjector<dim>::project_to_all_faces(face_points);
+        QProjector<dim>::project_to_all_faces(this->reference_cell_type(),
+                                              face_points);
       for (unsigned int k = 0;
            k < this->n_dofs_per_face() * GeometryInfo<dim>::faces_per_cell;
            ++k)
-        this->generalized_support_points[k] =
-          faces.point(k + QProjector<dim>::DataSetDescriptor::face(
-                            0, true, false, false, this->n_dofs_per_face()));
+        this->generalized_support_points[k] = faces.point(
+          k +
+          QProjector<dim>::DataSetDescriptor::face(this->reference_cell_type(),
+                                                   0,
+                                                   true,
+                                                   false,
+                                                   false,
+                                                   this->n_dofs_per_face()));
 
       current = this->n_dofs_per_face() * GeometryInfo<dim>::faces_per_cell;
     }
@@ -621,7 +627,9 @@ FE_RaviartThomasNodal<dim>::get_face_interpolation_matrix(
   // matrix by simply taking the
   // value at the support points.
   const Quadrature<dim> face_projection =
-    QProjector<dim>::project_to_face(quad_face_support, 0);
+    QProjector<dim>::project_to_face(this->reference_cell_type(),
+                                     quad_face_support,
+                                     0);
 
   for (unsigned int i = 0; i < source_fe.n_dofs_per_face(); ++i)
     {
@@ -727,7 +735,10 @@ FE_RaviartThomasNodal<dim>::get_subface_interpolation_matrix(
   // value at the support points.
 
   const Quadrature<dim> subface_projection =
-    QProjector<dim>::project_to_subface(quad_face_support, 0, subface);
+    QProjector<dim>::project_to_subface(this->reference_cell_type(),
+                                        quad_face_support,
+                                        0,
+                                        subface);
 
   for (unsigned int i = 0; i < source_fe.n_dofs_per_face(); ++i)
     {
index aacdccfcfd5b80f734b2d51bdce4cd2ea49383fa..b7c66bdf63223bc451d2e24671381772846a0a9f 100644 (file)
@@ -147,13 +147,19 @@ FE_RT_Bubbles<dim>::initialize_support_points(const unsigned int deg)
       for (unsigned int k = 0; k < this->n_dofs_per_face(); ++k)
         this->generalized_face_support_points[k] = face_points.point(k);
       Quadrature<dim> faces =
-        QProjector<dim>::project_to_all_faces(face_points);
+        QProjector<dim>::project_to_all_faces(this->reference_cell_type(),
+                                              face_points);
       for (unsigned int k = 0;
            k < this->n_dofs_per_face() * GeometryInfo<dim>::faces_per_cell;
            ++k)
-        this->generalized_support_points[k] =
-          faces.point(k + QProjector<dim>::DataSetDescriptor::face(
-                            0, true, false, false, this->n_dofs_per_face()));
+        this->generalized_support_points[k] = faces.point(
+          k +
+          QProjector<dim>::DataSetDescriptor::face(this->reference_cell_type(),
+                                                   0,
+                                                   true,
+                                                   false,
+                                                   false,
+                                                   this->n_dofs_per_face()));
 
       current = this->n_dofs_per_face() * GeometryInfo<dim>::faces_per_cell;
     }
index 2c169eb84e2dbf9e43bc67279ded9d16b183d584..ef668c01715de2b46bb473b093a17ffbeea5ff17 100644 (file)
@@ -14,7 +14,6 @@
 // ---------------------------------------------------------------------
 
 #include <deal.II/base/memory_consumption.h>
-#include <deal.II/base/qprojector.h>
 #include <deal.II/base/quadrature.h>
 
 #include <deal.II/dofs/dof_accessor.h>
index 97cb68e8b832d0db2a13b3cc49b11161656ffabf..577c3ae778b20019b6ca42492f06789497fdc900 100644 (file)
@@ -15,7 +15,6 @@
 
 #include <deal.II/base/config.h>
 
-#include <deal.II/base/qprojector.h>
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/tensor_product_polynomials.h>
 
index c6e5918588302c6aed8651e868ad617e463d40c8..2affce2eb522c287529cce74e97be1e8a5630d86 100644 (file)
@@ -116,7 +116,8 @@ MappingCartesian<dim, spacedim>::get_face_data(
 {
   std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
     std::make_unique<InternalData>(
-      QProjector<dim>::project_to_all_faces(quadrature));
+      QProjector<dim>::project_to_all_faces(ReferenceCell::get_hypercube(dim),
+                                            quadrature));
   auto &data = dynamic_cast<InternalData &>(*data_ptr);
 
   // verify that we have computed the transitive hull of the required
@@ -140,8 +141,8 @@ MappingCartesian<dim, spacedim>::get_subface_data(
   const Quadrature<dim - 1> &quadrature) const
 {
   std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
-    std::make_unique<InternalData>(
-      QProjector<dim>::project_to_all_subfaces(quadrature));
+    std::make_unique<InternalData>(QProjector<dim>::project_to_all_subfaces(
+      ReferenceCell::get_hypercube(dim), quadrature));
   auto &data = dynamic_cast<InternalData &>(*data_ptr);
 
   // verify that we have computed the transitive hull of the required
@@ -203,8 +204,7 @@ MappingCartesian<dim, spacedim>::maybe_update_cell_quadrature_points(
 {
   if (data.update_each & update_quadrature_points)
     {
-      const typename QProjector<dim>::DataSetDescriptor offset =
-        QProjector<dim>::DataSetDescriptor::cell();
+      const auto offset = QProjector<dim>::DataSetDescriptor::cell();
 
       transform_quadrature_points(cell, data, offset, quadrature_points);
     }
@@ -224,13 +224,13 @@ MappingCartesian<dim, spacedim>::maybe_update_face_quadrature_points(
 
   if (data.update_each & update_quadrature_points)
     {
-      const typename QProjector<dim>::DataSetDescriptor offset =
-        QProjector<dim>::DataSetDescriptor::face(face_no,
-                                                 cell->face_orientation(
-                                                   face_no),
-                                                 cell->face_flip(face_no),
-                                                 cell->face_rotation(face_no),
-                                                 quadrature_points.size());
+      const auto offset = QProjector<dim>::DataSetDescriptor::face(
+        ReferenceCell::get_hypercube(dim),
+        face_no,
+        cell->face_orientation(face_no),
+        cell->face_flip(face_no),
+        cell->face_rotation(face_no),
+        quadrature_points.size());
 
 
       transform_quadrature_points(cell, data, offset, quadrature_points);
@@ -257,15 +257,15 @@ MappingCartesian<dim, spacedim>::maybe_update_subface_quadrature_points(
 
   if (data.update_each & update_quadrature_points)
     {
-      const typename QProjector<dim>::DataSetDescriptor offset =
-        QProjector<dim>::DataSetDescriptor::subface(
-          face_no,
-          sub_no,
-          cell->face_orientation(face_no),
-          cell->face_flip(face_no),
-          cell->face_rotation(face_no),
-          quadrature_points.size(),
-          cell->subface_case(face_no));
+      const auto offset = QProjector<dim>::DataSetDescriptor::subface(
+        ReferenceCell::get_hypercube(dim),
+        face_no,
+        sub_no,
+        cell->face_orientation(face_no),
+        cell->face_flip(face_no),
+        cell->face_rotation(face_no),
+        quadrature_points.size(),
+        cell->subface_case(face_no));
 
       transform_quadrature_points(cell, data, offset, quadrature_points);
     }
index 62824c3f647fc2e8b7bb7ecf0d94728aea8fa3e4..163719f12e6d3f35a984db0319319506f3aa4a76 100644 (file)
@@ -628,7 +628,9 @@ MappingFEField<dim, spacedim, VectorType, DoFHandlerType>::get_face_data(
   std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
     std::make_unique<InternalData>(euler_dof_handler->get_fe(), fe_mask);
   auto &                data = dynamic_cast<InternalData &>(*data_ptr);
-  const Quadrature<dim> q(QProjector<dim>::project_to_all_faces(quadrature));
+  const Quadrature<dim> q(
+    QProjector<dim>::project_to_all_faces(ReferenceCell::get_hypercube(dim),
+                                          quadrature));
   this->compute_face_data(update_flags, q, quadrature.size(), data);
 
   return data_ptr;
@@ -644,7 +646,9 @@ MappingFEField<dim, spacedim, VectorType, DoFHandlerType>::get_subface_data(
   std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
     std::make_unique<InternalData>(euler_dof_handler->get_fe(), fe_mask);
   auto &                data = dynamic_cast<InternalData &>(*data_ptr);
-  const Quadrature<dim> q(QProjector<dim>::project_to_all_subfaces(quadrature));
+  const Quadrature<dim> q(
+    QProjector<dim>::project_to_all_subfaces(ReferenceCell::get_hypercube(dim),
+                                             quadrature));
   this->compute_face_data(update_flags, q, quadrature.size(), data);
 
   return data_ptr;
@@ -1763,7 +1767,9 @@ MappingFEField<dim, spacedim, VectorType, DoFHandlerType>::fill_fe_face_values(
       cell,
       face_no,
       numbers::invalid_unsigned_int,
-      QProjector<dim>::DataSetDescriptor::face(face_no,
+      QProjector<dim>::DataSetDescriptor::face(ReferenceCell::get_hypercube(
+                                                 dim),
+                                               face_no,
                                                cell->face_orientation(face_no),
                                                cell->face_flip(face_no),
                                                cell->face_rotation(face_no),
@@ -1803,14 +1809,15 @@ MappingFEField<dim, spacedim, VectorType, DoFHandlerType>::
       cell,
       face_no,
       numbers::invalid_unsigned_int,
-      QProjector<dim>::DataSetDescriptor::subface(face_no,
-                                                  subface_no,
-                                                  cell->face_orientation(
-                                                    face_no),
-                                                  cell->face_flip(face_no),
-                                                  cell->face_rotation(face_no),
-                                                  quadrature.size(),
-                                                  cell->subface_case(face_no)),
+      QProjector<dim>::DataSetDescriptor::subface(
+        ReferenceCell::get_hypercube(dim),
+        face_no,
+        subface_no,
+        cell->face_orientation(face_no),
+        cell->face_flip(face_no),
+        cell->face_rotation(face_no),
+        quadrature.size(),
+        cell->subface_case(face_no)),
       quadrature,
       data,
       euler_dof_handler->get_fe(),
index 2a22a445b8e265c310d0b550f757502c2921eb40..b8dc5845733fe84273b3eb022209c31482c0f09e 100644 (file)
@@ -289,7 +289,8 @@ MappingManifold<dim, spacedim>::get_face_data(
     std::make_unique<InternalData>();
   auto &data = dynamic_cast<InternalData &>(*data_ptr);
   data.initialize_face(this->requires_update_flags(update_flags),
-                       QProjector<dim>::project_to_all_faces(quadrature),
+                       QProjector<dim>::project_to_all_faces(
+                         ReferenceCell::get_hypercube(dim), quadrature),
                        quadrature.size());
 
   return data_ptr;
@@ -307,7 +308,8 @@ MappingManifold<dim, spacedim>::get_subface_data(
     std::make_unique<InternalData>();
   auto &data = dynamic_cast<InternalData &>(*data_ptr);
   data.initialize_face(this->requires_update_flags(update_flags),
-                       QProjector<dim>::project_to_all_subfaces(quadrature),
+                       QProjector<dim>::project_to_all_subfaces(
+                         ReferenceCell::get_hypercube(dim), quadrature),
                        quadrature.size());
 
   return data_ptr;
@@ -1221,7 +1223,8 @@ MappingManifold<dim, spacedim>::fill_fe_face_values(
     cell,
     face_no,
     numbers::invalid_unsigned_int,
-    QProjector<dim>::DataSetDescriptor::face(face_no,
+    QProjector<dim>::DataSetDescriptor::face(ReferenceCell::get_hypercube(dim),
+                                             face_no,
                                              cell->face_orientation(face_no),
                                              cell->face_flip(face_no),
                                              cell->face_rotation(face_no),
@@ -1254,7 +1257,9 @@ MappingManifold<dim, spacedim>::fill_fe_subface_values(
     cell,
     face_no,
     subface_no,
-    QProjector<dim>::DataSetDescriptor::subface(face_no,
+    QProjector<dim>::DataSetDescriptor::subface(ReferenceCell::get_hypercube(
+                                                  dim),
+                                                face_no,
                                                 subface_no,
                                                 cell->face_orientation(face_no),
                                                 cell->face_flip(face_no),
index e812cb30c86efedbbd731b90836850c3e6d89950..090d9a23b5cc68dd5fefcee79d128c06fa060a0b 100644 (file)
@@ -16,7 +16,6 @@
 
 #include <deal.II/base/derivative_form.h>
 #include <deal.II/base/memory_consumption.h>
-#include <deal.II/base/qprojector.h>
 #include <deal.II/base/quadrature.h>
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/tensor_product_polynomials.h>
index 014403e94bce153cef758e87e2b8ac0a07a91250..49a0eea8e4455cbc5136a3edc9eb9c0d1df7aada 100644 (file)
@@ -2695,7 +2695,8 @@ MappingQGeneric<dim, spacedim>::get_face_data(
     std::make_unique<InternalData>(polynomial_degree);
   auto &data = dynamic_cast<InternalData &>(*data_ptr);
   data.initialize_face(this->requires_update_flags(update_flags),
-                       QProjector<dim>::project_to_all_faces(quadrature),
+                       QProjector<dim>::project_to_all_faces(
+                         ReferenceCell::get_hypercube(dim), quadrature),
                        quadrature.size());
 
   return data_ptr;
@@ -2713,7 +2714,8 @@ MappingQGeneric<dim, spacedim>::get_subface_data(
     std::make_unique<InternalData>(polynomial_degree);
   auto &data = dynamic_cast<InternalData &>(*data_ptr);
   data.initialize_face(this->requires_update_flags(update_flags),
-                       QProjector<dim>::project_to_all_subfaces(quadrature),
+                       QProjector<dim>::project_to_all_subfaces(
+                         ReferenceCell::get_hypercube(dim), quadrature),
                        quadrature.size());
 
   return data_ptr;
@@ -3235,7 +3237,8 @@ MappingQGeneric<dim, spacedim>::fill_fe_face_values(
     cell,
     face_no,
     numbers::invalid_unsigned_int,
-    QProjector<dim>::DataSetDescriptor::face(face_no,
+    QProjector<dim>::DataSetDescriptor::face(ReferenceCell::get_hypercube(dim),
+                                             face_no,
                                              cell->face_orientation(face_no),
                                              cell->face_flip(face_no),
                                              cell->face_rotation(face_no),
@@ -3281,7 +3284,9 @@ MappingQGeneric<dim, spacedim>::fill_fe_subface_values(
     cell,
     face_no,
     subface_no,
-    QProjector<dim>::DataSetDescriptor::subface(face_no,
+    QProjector<dim>::DataSetDescriptor::subface(ReferenceCell::get_hypercube(
+                                                  dim),
+                                                face_no,
                                                 subface_no,
                                                 cell->face_orientation(face_no),
                                                 cell->face_flip(face_no),
index 330655c47af1fd21a8213942a192df7293b23205..2462f1ab7d33011e7b1d30f671111d1f467eaa58 100644 (file)
@@ -4031,7 +4031,9 @@ namespace internal
           std::vector<double> dummy_weights(n_points, 1. / n_points);
           Quadrature<dim - 1> quadrature(boundary_points, dummy_weights);
 
-          q_projector = QProjector<dim>::project_to_all_faces(quadrature);
+          q_projector =
+            QProjector<dim>::project_to_all_faces(ReferenceCell::Type::Quad,
+                                                  quadrature);
         }
 
       for (const auto &cell : tria.active_cell_iterators())

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.