]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Introduce FE_DGQLegendre and FE_DGQHermite
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 11 Jan 2017 17:49:02 +0000 (18:49 +0100)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 11 Jan 2017 21:52:55 +0000 (22:52 +0100)
doc/news/changes/minor/20170111MartinKronbichler-c [new file with mode: 0644]
include/deal.II/fe/fe_dgq.h
source/fe/fe_dgq.cc
source/fe/fe_dgq.inst.in
source/fe/fe_tools.cc
tests/fe/fe_data_test.cc
tests/fe/fe_data_test.output
tests/fe/get_fe_by_name_01.cc
tests/fe/get_fe_by_name_01.output

diff --git a/doc/news/changes/minor/20170111MartinKronbichler-c b/doc/news/changes/minor/20170111MartinKronbichler-c
new file mode 100644 (file)
index 0000000..15a5621
--- /dev/null
@@ -0,0 +1,5 @@
+New: There are two new finite elements FE_DGQLegendre and FE_DGQHermite which
+are similar to FE_DGQ but use the class Polynomials::Legendre and
+Polynomials::HermiteInterpolation as shape functions.
+<br>
+(Martin Kronbichler, 2017/01/11)
index eea42928c5fee0eb36c46b8f1a8b9bbdf26d429f..ce6c865994d23aa72a2a17366c7d2750612d3bd0 100644 (file)
@@ -312,15 +312,14 @@ public:
 
 protected:
   /**
-   * Constructor for tensor product polynomials based on Polynomials::Lagrange
-   * interpolation of the support points in the quadrature rule
-   * <tt>points</tt>. The degree of these polynomials is
-   * <tt>points.size()-1</tt>.
+   * Constructor for tensor product polynomials based on an arbitrary vector
+   * of polynomials. This constructor is used in derived classes to construct
+   * e.g. elements with arbitrary nodes or elements based on Legendre
+   * polynomials.
    *
-   * Note: The FE_DGQ::clone function does not work properly for FE with
-   * arbitrary nodes!
+   * The degree of these polynomials is <tt>polynomials.size()-1</tt>.
    */
-  FE_DGQ (const Quadrature<1> &points);
+  FE_DGQ (const std::vector<Polynomials::Polynomial<double> > &polynomials);
 
   /**
    * @p clone function instead of a copy constructor.
@@ -419,6 +418,104 @@ protected:
 };
 
 
+
+/**
+ * Implementation of scalar, discontinuous tensor product elements based on
+ * Legendre polynomials, described by the tensor product of the polynomial
+ * space Polynomials::Legendre. As opposed to the basic FE_DGQ element, these
+ * elements are not interpolatory and no support points are defined.
+ *
+ * See the base class documentation in FE_DGQ for details.
+ *
+ * @author Martin Kronbichler, 2017
+ */
+template <int dim,int spacedim=dim>
+class FE_DGQLegendre : public FE_DGQ<dim,spacedim>
+{
+public:
+  /**
+   * Constructor for tensor product polynomials based on Polynomials::Legendre
+   * interpolation.
+   */
+  FE_DGQLegendre (const unsigned int degree);
+
+  /**
+   * Return a list of constant modes of the element. For the Legendre basis,
+   * it returns one row where the first element (corresponding to the constant
+   * mode) is set to true and all other elements are set to false.
+   */
+  virtual std::pair<Table<2,bool>, std::vector<unsigned int> >
+  get_constant_modes () const;
+
+  /**
+   * Return a string that uniquely identifies a finite element. This class
+   * returns <tt>FE_DGQLegendre<dim>(degree)</tt> with <tt>dim</tt> and
+   * <tt>degree</tt> replaced by the values given by the template parameter
+   * and the argument passed to the constructor, respectively.
+   */
+  virtual std::string get_name () const;
+
+protected:
+  /**
+   * @p clone function instead of a copy constructor.
+   *
+   * This function is needed by the constructors of @p FESystem.
+   */
+  virtual FiniteElement<dim,spacedim> *clone() const;
+};
+
+
+
+/**
+ * Implementation of scalar, discontinuous tensor product elements based on
+ * Hermite polynomials, described by the polynomial space
+ * Polynomials::HermiteInterpolation. As opposed to the basic FE_DGQ element,
+ * these elements are not interpolatory and no support points are defined.
+ *
+ * This element is only a Hermite polynomials for degrees larger or equal to
+ * three. For degrees zero to two, a usual Lagrange basis is selected.
+ *
+ * See the base class documentation in FE_DGQ for details.
+ *
+ * @author Martin Kronbichler, 2017
+ */
+template <int dim,int spacedim=dim>
+class FE_DGQHermite : public FE_DGQ<dim,spacedim>
+{
+public:
+  /**
+   * Constructor for tensor product polynomials based on
+   * Polynomials::HermiteInterpolation.
+   */
+  FE_DGQHermite (const unsigned int degree);
+
+  /**
+   * Return a list of constant modes of the element. For the Hermite basis of
+   * degree three and larger, it returns one row where the first two elements
+   * (corresponding to the value left and the value at the right) are set to
+   * true and all other elements are set to false.
+   */
+  virtual std::pair<Table<2,bool>, std::vector<unsigned int> >
+  get_constant_modes () const;
+
+  /**
+   * Return a string that uniquely identifies a finite element. This class
+   * returns <tt>FE_DGQHermite<dim>(degree)</tt>, with <tt>dim</tt> and
+   * <tt>degree</tt> replaced by the values given by the template parameter
+   * and the argument passed to the constructor, respectively.
+   */
+  virtual std::string get_name () const;
+
+protected:
+  /**
+   * @p clone function instead of a copy constructor.
+   *
+   * This function is needed by the constructors of @p FESystem.
+   */
+  virtual FiniteElement<dim,spacedim> *clone() const;
+};
+
+
 /*@}*/
 
 DEAL_II_NAMESPACE_CLOSE
index 1eb6709a9d2fad0125cc60d1abee35ce375684f7..16e91bc78bdb9da502466c87b2baf8897eff12f9 100644 (file)
@@ -67,21 +67,16 @@ FE_DGQ<dim, spacedim>::FE_DGQ (const unsigned int degree)
 
 
 template <int dim, int spacedim>
-FE_DGQ<dim, spacedim>::FE_DGQ (const Quadrature<1> &points)
+FE_DGQ<dim, spacedim>::FE_DGQ (const std::vector<Polynomials::Polynomial<double> > &polynomials)
   :
   FE_Poly<TensorProductPolynomials<dim>, dim, spacedim> (
-    TensorProductPolynomials<dim>(Polynomials::generate_complete_Lagrange_basis(points.get_points())),
-    FiniteElementData<dim>(get_dpo_vector(points.size()-1), 1, points.size()-1, FiniteElementData<dim>::L2),
-    std::vector<bool>(FiniteElementData<dim>(get_dpo_vector(points.size()-1),1, points.size()-1).dofs_per_cell, true),
-    std::vector<ComponentMask>(FiniteElementData<dim>(get_dpo_vector(points.size()-1),1, points.size()-1).dofs_per_cell, std::vector<bool>(1,true)))
+    TensorProductPolynomials<dim>(polynomials),
+    FiniteElementData<dim>(get_dpo_vector(polynomials.size()-1), 1, polynomials.size()-1, FiniteElementData<dim>::L2),
+    std::vector<bool>(FiniteElementData<dim>(get_dpo_vector(polynomials.size()-1),1, polynomials.size()-1).dofs_per_cell, true),
+    std::vector<ComponentMask>(FiniteElementData<dim>(get_dpo_vector(polynomials.size()-1),1, polynomials.size()-1).dofs_per_cell, std::vector<bool>(1,true)))
 {
-  // Compute support points, which are the tensor product of the Lagrange
-  // interpolation points in the constructor.
-  Assert (points.size() > 0,
-          (typename FiniteElement<dim, spacedim>::ExcFEHasNoSupportPoints ()));
-  Quadrature<dim> support_quadrature(points);
-  this->unit_support_points = support_quadrature.get_points();
-
+  // No support points can be defined in general. Derived classes might define
+  // support points like the class FE_DGQArbitraryNodes
 
   // do not initialize embedding and restriction here. these matrices are
   // initialized on demand in get_restriction_matrix and
@@ -558,6 +553,12 @@ FE_DGQ<dim, spacedim>::has_support_on_face (const unsigned int shape_index,
 
   unsigned int n = this->degree+1;
 
+  // For DGQ elements that do not define support points, we cannot define
+  // whether they have support at the boundary easily, so return true in any
+  // case
+  if (this->unit_support_points.empty())
+    return true;
+
   // for DGQ(0) elements or arbitrary node DGQ with support points not located
   // at the element boundary, the single shape functions is constant and
   // therefore lives on the boundary
@@ -656,10 +657,17 @@ FE_DGQ<dim, spacedim>::memory_consumption () const
 
 
 
+// ------------------------------ FE_DGQArbitraryNodes -----------------------
+
 template <int dim, int spacedim>
 FE_DGQArbitraryNodes<dim,spacedim>::FE_DGQArbitraryNodes (const Quadrature<1> &points)
-  : FE_DGQ<dim,spacedim>(points)
-{}
+  : FE_DGQ<dim,spacedim>(Polynomials::generate_complete_Lagrange_basis(points.get_points()))
+{
+  Assert (points.size() > 0,
+          (typename FiniteElement<dim, spacedim>::ExcFEHasNoSupportPoints ()));
+  Quadrature<dim> support_quadrature(points);
+  this->unit_support_points = support_quadrature.get_points();
+}
 
 
 
@@ -767,6 +775,103 @@ FE_DGQArbitraryNodes<dim,spacedim>::clone() const
 
 
 
+// ---------------------------------- FE_DGQLegendre -------------------------
+
+template <int dim, int spacedim>
+FE_DGQLegendre<dim,spacedim>::FE_DGQLegendre (const unsigned int degree)
+  : FE_DGQ<dim,spacedim>(Polynomials::Legendre::generate_complete_basis(degree))
+{}
+
+
+
+template <int dim, int spacedim>
+std::pair<Table<2,bool>, std::vector<unsigned int> >
+FE_DGQLegendre<dim,spacedim>::get_constant_modes () const
+{
+  // Legendre represents a constant function by one in the first basis
+  // function and zero in all others
+  Table<2,bool> constant_modes(1, this->dofs_per_cell);
+  constant_modes(0,0) = true;
+  return std::pair<Table<2,bool>, std::vector<unsigned int> >
+         (constant_modes, std::vector<unsigned int>(1, 0));
+}
+
+
+
+template <int dim, int spacedim>
+std::string
+FE_DGQLegendre<dim,spacedim>::get_name () const
+{
+  return "FE_DGQLegendre<" + Utilities::dim_string(dim,spacedim) + ">("
+         + Utilities::int_to_string(this->degree) + ")";
+}
+
+
+
+template <int dim, int spacedim>
+FiniteElement<dim,spacedim> *
+FE_DGQLegendre<dim,spacedim>::clone() const
+{
+  return new FE_DGQLegendre<dim,spacedim>(this->degree);
+}
+
+
+
+// ---------------------------------- FE_DGQHermite --------------------------
+
+template <int dim, int spacedim>
+FE_DGQHermite<dim,spacedim>::FE_DGQHermite (const unsigned int degree)
+  : FE_DGQ<dim,spacedim>(degree < 3 ?
+                         Polynomials::generate_complete_Lagrange_basis(get_QGaussLobatto_points(degree))
+                         :
+                         Polynomials::HermiteInterpolation::generate_complete_basis(degree))
+{}
+
+
+
+template <int dim, int spacedim>
+std::pair<Table<2,bool>, std::vector<unsigned int> >
+FE_DGQHermite<dim,spacedim>::get_constant_modes () const
+{
+  if (this->degree < 3)
+    return this->FE_DGQ<dim,spacedim>::get_constant_modes();
+  else
+    {
+      // The first two basis functions in the Hermite polynomials represent
+      // the value 1 in the left and right end point of the element. Expand
+      // them into the tensor product.
+      AssertThrow(dim<=3, ExcNotImplemented());
+      Table<2,bool> constant_modes(1, this->dofs_per_cell);
+      for (unsigned int i=0; i<(dim>2?2:1); ++i)
+        for (unsigned int j=0; j<(dim>1?2:1); ++j)
+          for (unsigned int k=0; k<2; ++k)
+            constant_modes(0,i*(this->degree+1)*(this->degree+1)+j*(this->degree+1)+k) = true;
+      return std::pair<Table<2,bool>, std::vector<unsigned int> >
+             (constant_modes, std::vector<unsigned int>(1, 0));
+    }
+}
+
+
+
+template <int dim, int spacedim>
+std::string
+FE_DGQHermite<dim,spacedim>::get_name () const
+{
+  return "FE_DGQHermite<" + Utilities::dim_string(dim,spacedim) + ">("
+         + Utilities::int_to_string(this->degree) + ")";
+}
+
+
+
+template <int dim, int spacedim>
+FiniteElement<dim,spacedim> *
+FE_DGQHermite<dim,spacedim>::clone() const
+{
+  return new FE_DGQHermite<dim,spacedim>(this->degree);
+}
+
+
+
 // explicit instantiations
 #include "fe_dgq.inst"
 
index 27bee60a59f8522595b5664a34de4c3e29baf58d..ab6ef74dffd67524ec7436ce33db05c846ac9331 100644 (file)
@@ -20,5 +20,7 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
 #if deal_II_dimension <= deal_II_space_dimension
     template class FE_DGQ<deal_II_dimension, deal_II_space_dimension>;
     template class FE_DGQArbitraryNodes<deal_II_dimension, deal_II_space_dimension>;
+    template class FE_DGQLegendre<deal_II_dimension, deal_II_space_dimension>;
+    template class FE_DGQHermite<deal_II_dimension, deal_II_space_dimension>;
 #endif
 }
index 0d0b73553564f418b7953b8bec6de4db6808beb0..ddeb7db4c13eadd8d41021970bb436d6ac9ca9a9 100644 (file)
@@ -1095,6 +1095,10 @@ namespace
       = FEFactoryPointer(new FETools::FEFactory<FE_DGQ<dim> >);
     result["FE_DGQArbitraryNodes"]
       = FEFactoryPointer(new FETools::FEFactory<FE_DGQ<dim> >);
+    result["FE_DGQLegendre"]
+      = FEFactoryPointer(new FETools::FEFactory<FE_DGQLegendre<dim> >);
+    result["FE_DGQHermite"]
+      = FEFactoryPointer(new FETools::FEFactory<FE_DGQHermite<dim> >);
     result["FE_FaceQ"]
       = FEFactoryPointer(new FETools::FEFactory<FE_FaceQ<dim> >);
     result["FE_FaceP"]
@@ -1134,6 +1138,10 @@ namespace
       = FEFactoryPointer(new FETools::FEFactory<FE_Nothing<dim,spacedim> >);
     result["FE_DGQArbitraryNodes"]
       = FEFactoryPointer(new FETools::FEFactory<FE_DGQ<dim,spacedim> >);
+    result["FE_DGQLegendre"]
+      = FEFactoryPointer(new FETools::FEFactory<FE_DGQLegendre<dim,spacedim> >);
+    result["FE_DGQHermite"]
+      = FEFactoryPointer(new FETools::FEFactory<FE_DGQHermite<dim,spacedim> >);
     result["FE_Q_Bubbles"]
       = FEFactoryPointer(new FETools::FEFactory<FE_Q_Bubbles<dim,spacedim> >);
     result["FE_Q_DG0"]
index 9b72b18ff30b7f920ea6ee5714981d0e5ee1468b..28ef7b52535ade28e6cc914527cc1b056c90cd4a 100644 (file)
@@ -121,6 +121,12 @@ void test_fe_datas()
   deallog << (*fe_datas.rbegin())->get_name() << std::endl;
   fe_datas.push_back(new FE_DGQArbitraryNodes<dim> (QGauss<1>(3)));
   deallog << (*fe_datas.rbegin())->get_name() << std::endl;
+  fe_datas.push_back(new FE_DGQLegendre<dim> (1));
+  deallog << (*fe_datas.rbegin())->get_name() << std::endl;
+  fe_datas.push_back(new FE_DGQLegendre<dim> (2));
+  deallog << (*fe_datas.rbegin())->get_name() << std::endl;
+  fe_datas.push_back(new FE_DGQHermite<dim> (3));
+  deallog << (*fe_datas.rbegin())->get_name() << std::endl;
   fe_datas.push_back(new FE_DGP<dim> (1));
   deallog << (*fe_datas.rbegin())->get_name() << std::endl;
   fe_datas.push_back(new FE_DGP<dim> (2));
index a6f7857d171d20d653a0c5ab345ff6268050b716..d51493c3122af9bbf23f1f1cf71cf4734d4dee42 100644 (file)
@@ -10,6 +10,9 @@ DEAL::FE_DGQ<1>(2)
 DEAL::FE_DGQArbitraryNodes<1>(QIterated(QTrapez(),4))
 DEAL::FE_DGQ<1>(4)
 DEAL::FE_DGQArbitraryNodes<1>(QGauss(3))
+DEAL::FE_DGQLegendre<1>(1)
+DEAL::FE_DGQLegendre<1>(2)
+DEAL::FE_DGQHermite<1>(3)
 DEAL::FE_DGP<1>(1)
 DEAL::FE_DGP<1>(2)
 DEAL::FESystem<1>[FE_Q<1>(2)^2]
@@ -326,7 +329,7 @@ DEAL::face_to_cell_index:
 DEAL::support on face 0:       0       1       2
 DEAL::support on face 1:       0       1       2
 DEAL::
-DEAL::fe_data[11]:FE_DGP<1>(1)
+DEAL::fe_data[11]:FE_DGQLegendre<1>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=0
@@ -353,7 +356,7 @@ DEAL::face_to_cell_index:
 DEAL::support on face 0:       0       1
 DEAL::support on face 1:       0       1
 DEAL::
-DEAL::fe_data[12]:FE_DGP<1>(2)
+DEAL::fe_data[12]:FE_DGQLegendre<1>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=3
 DEAL::dofs_per_quad=0
@@ -380,7 +383,88 @@ DEAL::face_to_cell_index:
 DEAL::support on face 0:       0       1       2
 DEAL::support on face 1:       0       1       2
 DEAL::
-DEAL::fe_data[13]:FESystem<1>[FE_Q<1>(2)^2]
+DEAL::fe_data[13]:FE_DGQHermite<1>(3)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=4
+DEAL::dofs_per_quad=0
+DEAL::dofs_per_hex=0
+DEAL::first_line_index=0
+DEAL::first_quad_index=4
+DEAL::first_hex_index=4
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=4
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[4]->4
+DEAL::degree=3
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2       3
+DEAL::support on face 1:       0       1       2       3
+DEAL::
+DEAL::fe_data[14]:FE_DGP<1>(1)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=2
+DEAL::dofs_per_quad=0
+DEAL::dofs_per_hex=0
+DEAL::first_line_index=0
+DEAL::first_quad_index=2
+DEAL::first_hex_index=2
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=2
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[2]->2
+DEAL::degree=1
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1
+DEAL::support on face 1:       0       1
+DEAL::
+DEAL::fe_data[15]:FE_DGP<1>(2)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=3
+DEAL::dofs_per_quad=0
+DEAL::dofs_per_hex=0
+DEAL::first_line_index=0
+DEAL::first_quad_index=3
+DEAL::first_hex_index=3
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=3
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[3]->3
+DEAL::degree=2
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2
+DEAL::support on face 1:       0       1       2
+DEAL::
+DEAL::fe_data[16]:FESystem<1>[FE_Q<1>(2)^2]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=0
@@ -407,7 +491,7 @@ DEAL::face_to_cell_index: 2 3
 DEAL::support on face 0:       0       1
 DEAL::support on face 1:       2       3
 DEAL::
-DEAL::fe_data[14]:FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)]
+DEAL::fe_data[17]:FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)]
 DEAL::dofs_per_vertex=3
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=0
@@ -434,7 +518,7 @@ DEAL::face_to_cell_index: 3 4 5
 DEAL::support on face 0:       0       1       2
 DEAL::support on face 1:       3       4       5
 DEAL::
-DEAL::fe_data[15]:FE_FaceQ<1>(0)
+DEAL::fe_data[18]:FE_FaceQ<1>(0)
 DEAL::dofs_per_vertex=1
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -461,7 +545,7 @@ DEAL::face_to_cell_index: 1
 DEAL::support on face 0:       0
 DEAL::support on face 1:       1
 DEAL::
-DEAL::fe_data[16]:FE_FaceQ<1>(1)
+DEAL::fe_data[19]:FE_FaceQ<1>(1)
 DEAL::dofs_per_vertex=1
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -488,7 +572,7 @@ DEAL::face_to_cell_index: 1
 DEAL::support on face 0:       0
 DEAL::support on face 1:       1
 DEAL::
-DEAL::fe_data[17]:FE_FaceQ<1>(3)
+DEAL::fe_data[20]:FE_FaceQ<1>(3)
 DEAL::dofs_per_vertex=1
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -515,7 +599,7 @@ DEAL::face_to_cell_index: 1
 DEAL::support on face 0:       0
 DEAL::support on face 1:       1
 DEAL::
-DEAL::fe_data[18]:FE_FaceP<1>(0)
+DEAL::fe_data[21]:FE_FaceP<1>(0)
 DEAL::dofs_per_vertex=1
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -542,7 +626,7 @@ DEAL::face_to_cell_index: 1
 DEAL::support on face 0:       0
 DEAL::support on face 1:       1
 DEAL::
-DEAL::fe_data[19]:FE_FaceP<1>(1)
+DEAL::fe_data[22]:FE_FaceP<1>(1)
 DEAL::dofs_per_vertex=1
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -569,7 +653,7 @@ DEAL::face_to_cell_index: 1
 DEAL::support on face 0:       0
 DEAL::support on face 1:       1
 DEAL::
-DEAL::fe_data[20]:FE_FaceP<1>(3)
+DEAL::fe_data[23]:FE_FaceP<1>(3)
 DEAL::dofs_per_vertex=1
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -596,7 +680,7 @@ DEAL::face_to_cell_index: 1
 DEAL::support on face 0:       0
 DEAL::support on face 1:       1
 DEAL::
-DEAL::fe_data[21]:FESystem<1>[FE_Q<1>(3)^2]
+DEAL::fe_data[24]:FESystem<1>[FE_Q<1>(3)^2]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=4
 DEAL::dofs_per_quad=0
@@ -623,7 +707,7 @@ DEAL::face_to_cell_index: 2 3
 DEAL::support on face 0:       0       1
 DEAL::support on face 1:       2       3
 DEAL::
-DEAL::fe_data[22]:FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(3)]
+DEAL::fe_data[25]:FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(3)]
 DEAL::dofs_per_vertex=3
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=0
@@ -650,7 +734,7 @@ DEAL::face_to_cell_index: 3 4 5
 DEAL::support on face 0:       0       1       2
 DEAL::support on face 1:       3       4       5
 DEAL::
-DEAL::fe_data[23]:FESystem<1>[FE_Q<1>(4)^2]
+DEAL::fe_data[26]:FESystem<1>[FE_Q<1>(4)^2]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=6
 DEAL::dofs_per_quad=0
@@ -677,7 +761,7 @@ DEAL::face_to_cell_index: 2 3
 DEAL::support on face 0:       0       1
 DEAL::support on face 1:       2       3
 DEAL::
-DEAL::fe_data[24]:FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]^2]
+DEAL::fe_data[27]:FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]^2]
 DEAL::dofs_per_vertex=4
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -704,7 +788,7 @@ DEAL::face_to_cell_index: 4 5 6 7
 DEAL::support on face 0:       0       1       2       3
 DEAL::support on face 1:       4       5       6       7
 DEAL::
-DEAL::fe_data[25]:FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]-FESystem<1>[FE_DGQ<1>(1)^2]]
+DEAL::fe_data[28]:FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]-FESystem<1>[FE_DGQ<1>(1)^2]]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=4
 DEAL::dofs_per_quad=0
@@ -731,7 +815,7 @@ DEAL::face_to_cell_index: 2 3
 DEAL::support on face 0:       0       1       4       6
 DEAL::support on face 1:       2       3       5       7
 DEAL::
-DEAL::fe_data[26]:FESystem<1>[FESystem<1>[FE_Q<1>(1)-FE_Q<1>(2)]-FESystem<1>[FE_Q<1>(2)^2]-FESystem<1>[FE_DGQ<1>(2)^2]]
+DEAL::fe_data[29]:FESystem<1>[FESystem<1>[FE_Q<1>(1)-FE_Q<1>(2)]-FESystem<1>[FE_Q<1>(2)^2]-FESystem<1>[FE_DGQ<1>(2)^2]]
 DEAL::dofs_per_vertex=4
 DEAL::dofs_per_line=9
 DEAL::dofs_per_quad=0
@@ -758,7 +842,7 @@ DEAL::face_to_cell_index: 4 5 6 7
 DEAL::support on face 0:       0       1       2       3       11      14
 DEAL::support on face 1:       4       5       6       7       13      16
 DEAL::
-DEAL::fe_data[27]:FESystem<1>[FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]^2]^2-FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]-FESystem<1>[FE_DGQ<1>(1)^2]]-FESystem<1>[FESystem<1>[FE_Q<1>(1)-FE_Q<1>(2)]-FESystem<1>[FE_Q<1>(2)^2]-FESystem<1>[FE_DGQ<1>(2)^2]]^2]
+DEAL::fe_data[30]:FESystem<1>[FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]^2]^2-FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]-FESystem<1>[FE_DGQ<1>(1)^2]]-FESystem<1>[FESystem<1>[FE_Q<1>(1)-FE_Q<1>(2)]-FESystem<1>[FE_Q<1>(2)^2]-FESystem<1>[FE_DGQ<1>(2)^2]]^2]
 DEAL::dofs_per_vertex=18
 DEAL::dofs_per_line=22
 DEAL::dofs_per_quad=0
@@ -796,6 +880,9 @@ DEAL::FE_DGQ<2>(2)
 DEAL::FE_DGQArbitraryNodes<2>(QIterated(QTrapez(),4))
 DEAL::FE_DGQ<2>(4)
 DEAL::FE_DGQArbitraryNodes<2>(QGauss(3))
+DEAL::FE_DGQLegendre<2>(1)
+DEAL::FE_DGQLegendre<2>(2)
+DEAL::FE_DGQHermite<2>(3)
 DEAL::FE_DGP<2>(1)
 DEAL::FE_DGP<2>(2)
 DEAL::FESystem<2>[FE_Q<2>(2)^2]
@@ -1175,7 +1262,100 @@ DEAL::support on face 1:        0       1       2       3       4       5       6       7       8
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8
 DEAL::
-DEAL::fe_data[11]:FE_DGP<2>(1)
+DEAL::fe_data[11]:FE_DGQLegendre<2>(1)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=0
+DEAL::dofs_per_quad=4
+DEAL::dofs_per_hex=0
+DEAL::first_line_index=0
+DEAL::first_quad_index=0
+DEAL::first_hex_index=4
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=4
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[4]->4
+DEAL::degree=1
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2       3
+DEAL::support on face 1:       0       1       2       3
+DEAL::support on face 2:       0       1       2       3
+DEAL::support on face 3:       0       1       2       3
+DEAL::
+DEAL::fe_data[12]:FE_DGQLegendre<2>(2)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=0
+DEAL::dofs_per_quad=9
+DEAL::dofs_per_hex=0
+DEAL::first_line_index=0
+DEAL::first_quad_index=0
+DEAL::first_hex_index=9
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=9
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[9]->9
+DEAL::degree=2
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2       3       4       5       6       7       8
+DEAL::support on face 1:       0       1       2       3       4       5       6       7       8
+DEAL::support on face 2:       0       1       2       3       4       5       6       7       8
+DEAL::support on face 3:       0       1       2       3       4       5       6       7       8
+DEAL::
+DEAL::fe_data[13]:FE_DGQHermite<2>(3)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=0
+DEAL::dofs_per_quad=16
+DEAL::dofs_per_hex=0
+DEAL::first_line_index=0
+DEAL::first_quad_index=0
+DEAL::first_hex_index=16
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=16
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[16]->16
+DEAL::degree=3
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15
+DEAL::support on face 1:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15
+DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15
+DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15
+DEAL::
+DEAL::fe_data[14]:FE_DGP<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=3
@@ -1206,7 +1386,7 @@ DEAL::support on face 1:  0       1       2
 DEAL::support on face 2:       0       1       2
 DEAL::support on face 3:       0       1       2
 DEAL::
-DEAL::fe_data[12]:FE_DGP<2>(2)
+DEAL::fe_data[15]:FE_DGP<2>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=6
@@ -1237,7 +1417,7 @@ DEAL::support on face 1:  0       1       2       3       4       5
 DEAL::support on face 2:       0       1       2       3       4       5
 DEAL::support on face 3:       0       1       2       3       4       5
 DEAL::
-DEAL::fe_data[13]:FESystem<2>[FE_Q<2>(2)^2]
+DEAL::fe_data[16]:FESystem<2>[FE_Q<2>(2)^2]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=2
@@ -1268,7 +1448,7 @@ DEAL::support on face 1:  2       3       6       7       10      11
 DEAL::support on face 2:       0       1       2       3       12      13
 DEAL::support on face 3:       4       5       6       7       14      15
 DEAL::
-DEAL::fe_data[14]:FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)]
+DEAL::fe_data[17]:FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)]
 DEAL::dofs_per_vertex=3
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=1
@@ -1299,7 +1479,7 @@ DEAL::support on face 1:  3       4       5       9       10      11      13
 DEAL::support on face 2:       0       1       2       3       4       5       14
 DEAL::support on face 3:       6       7       8       9       10      11      15
 DEAL::
-DEAL::fe_data[15]:FE_FaceQ<2>(0)
+DEAL::fe_data[18]:FE_FaceQ<2>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=0
@@ -1330,7 +1510,7 @@ DEAL::support on face 1:  1
 DEAL::support on face 2:       2
 DEAL::support on face 3:       3
 DEAL::
-DEAL::fe_data[16]:FE_FaceQ<2>(1)
+DEAL::fe_data[19]:FE_FaceQ<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=0
@@ -1361,7 +1541,7 @@ DEAL::support on face 1:  2       3
 DEAL::support on face 2:       4       5
 DEAL::support on face 3:       6       7
 DEAL::
-DEAL::fe_data[17]:FE_FaceQ<2>(3)
+DEAL::fe_data[20]:FE_FaceQ<2>(3)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=4
 DEAL::dofs_per_quad=0
@@ -1392,7 +1572,7 @@ DEAL::support on face 1:  4       5       6       7
 DEAL::support on face 2:       8       9       10      11
 DEAL::support on face 3:       12      13      14      15
 DEAL::
-DEAL::fe_data[18]:FE_FaceP<2>(0)
+DEAL::fe_data[21]:FE_FaceP<2>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=0
@@ -1423,7 +1603,7 @@ DEAL::support on face 1:  1
 DEAL::support on face 2:       2
 DEAL::support on face 3:       3
 DEAL::
-DEAL::fe_data[19]:FE_FaceP<2>(1)
+DEAL::fe_data[22]:FE_FaceP<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=0
@@ -1454,7 +1634,7 @@ DEAL::support on face 1:  2       3
 DEAL::support on face 2:       4       5
 DEAL::support on face 3:       6       7
 DEAL::
-DEAL::fe_data[20]:FE_FaceP<2>(3)
+DEAL::fe_data[23]:FE_FaceP<2>(3)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=4
 DEAL::dofs_per_quad=0
@@ -1485,7 +1665,7 @@ DEAL::support on face 1:  4       5       6       7
 DEAL::support on face 2:       8       9       10      11
 DEAL::support on face 3:       12      13      14      15
 DEAL::
-DEAL::fe_data[21]:FE_DGRaviartThomas<2>(0)
+DEAL::fe_data[24]:FE_DGRaviartThomas<2>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=4
@@ -1516,7 +1696,7 @@ DEAL::support on face 1:  0       1       2       3
 DEAL::support on face 2:       0       1       2       3
 DEAL::support on face 3:       0       1       2       3
 DEAL::
-DEAL::fe_data[22]:FE_DGRaviartThomas<2>(1)
+DEAL::fe_data[25]:FE_DGRaviartThomas<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=12
@@ -1547,7 +1727,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::
-DEAL::fe_data[23]:FE_DGBDM<2>(1)
+DEAL::fe_data[26]:FE_DGBDM<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=8
@@ -1578,7 +1758,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7
 DEAL::support on face 2:       0       1       2       3       4       5       6       7
 DEAL::support on face 3:       0       1       2       3       4       5       6       7
 DEAL::
-DEAL::fe_data[24]:FE_DGBDM<2>(2)
+DEAL::fe_data[27]:FE_DGBDM<2>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=14
@@ -1609,7 +1789,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::
-DEAL::fe_data[25]:FE_DGNedelec<2>(0)
+DEAL::fe_data[28]:FE_DGNedelec<2>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=4
@@ -1640,7 +1820,7 @@ DEAL::support on face 1:  0       1       2       3
 DEAL::support on face 2:       0       1       2       3
 DEAL::support on face 3:       0       1       2       3
 DEAL::
-DEAL::fe_data[26]:FE_DGNedelec<2>(1)
+DEAL::fe_data[29]:FE_DGNedelec<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=12
@@ -1671,7 +1851,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::
-DEAL::fe_data[27]:FE_RaviartThomas<2>(0)
+DEAL::fe_data[30]:FE_RaviartThomas<2>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=0
@@ -1702,7 +1882,7 @@ DEAL::support on face 1:  1       2       3
 DEAL::support on face 2:       0       1       2
 DEAL::support on face 3:       0       1       3
 DEAL::
-DEAL::fe_data[28]:FE_RaviartThomas<2>(1)
+DEAL::fe_data[31]:FE_RaviartThomas<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=4
@@ -1733,7 +1913,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::
-DEAL::fe_data[29]:FE_RaviartThomas<2>(2)
+DEAL::fe_data[32]:FE_RaviartThomas<2>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=3
 DEAL::dofs_per_quad=12
@@ -1764,7 +1944,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23
 DEAL::
-DEAL::fe_data[30]:FESystem<2>[FE_RaviartThomas<2>(1)-FE_DGQ<2>(1)]
+DEAL::fe_data[33]:FESystem<2>[FE_RaviartThomas<2>(1)-FE_DGQ<2>(1)]
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=8
@@ -1795,7 +1975,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11      13      15
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      14      15
 DEAL::
-DEAL::fe_data[31]:FE_BDM<2>(1)
+DEAL::fe_data[34]:FE_BDM<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=0
@@ -1826,7 +2006,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7
 DEAL::support on face 2:       0       1       2       3       4       5       6       7
 DEAL::support on face 3:       0       1       2       3       4       5       6       7
 DEAL::
-DEAL::fe_data[32]:FE_BDM<2>(2)
+DEAL::fe_data[35]:FE_BDM<2>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=3
 DEAL::dofs_per_quad=2
@@ -1857,7 +2037,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::
-DEAL::fe_data[33]:FE_Nedelec<2>(0)
+DEAL::fe_data[36]:FE_Nedelec<2>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=0
@@ -1888,7 +2068,7 @@ DEAL::support on face 1:  1       2       3
 DEAL::support on face 2:       0       1       2
 DEAL::support on face 3:       0       1       3
 DEAL::
-DEAL::fe_data[34]:FE_Nedelec<2>(1)
+DEAL::fe_data[37]:FE_Nedelec<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=4
@@ -1919,7 +2099,7 @@ DEAL::support on face 1:  2       3       4       5       6       7
 DEAL::support on face 2:       0       1       2       3       4       5
 DEAL::support on face 3:       0       1       2       3       6       7       8       9       10      11
 DEAL::
-DEAL::fe_data[35]:FE_DGBDM<2>(1)
+DEAL::fe_data[38]:FE_DGBDM<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=8
@@ -1950,7 +2130,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7
 DEAL::support on face 2:       0       1       2       3       4       5       6       7
 DEAL::support on face 3:       0       1       2       3       4       5       6       7
 DEAL::
-DEAL::fe_data[36]:FE_DGBDM<2>(2)
+DEAL::fe_data[39]:FE_DGBDM<2>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=14
@@ -1981,7 +2161,7 @@ DEAL::support on face 1:  0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      12      13
 DEAL::
-DEAL::fe_data[37]:FE_RaviartThomasNodal<2>(0)
+DEAL::fe_data[40]:FE_RaviartThomasNodal<2>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=0
@@ -2012,7 +2192,7 @@ DEAL::support on face 1:  1       2       3
 DEAL::support on face 2:       0       1       2
 DEAL::support on face 3:       0       1       3
 DEAL::
-DEAL::fe_data[38]:FE_RaviartThomasNodal<2>(1)
+DEAL::fe_data[41]:FE_RaviartThomasNodal<2>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=4
@@ -2043,7 +2223,7 @@ DEAL::support on face 1:  2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 2:       0       1       2       3       4       5       8       9       10      11
 DEAL::support on face 3:       0       1       2       3       6       7       8       9       10      11
 DEAL::
-DEAL::fe_data[39]:FESystem<2>[FE_RaviartThomasNodal<2>(1)-FE_DGQ<2>(1)]
+DEAL::fe_data[42]:FESystem<2>[FE_RaviartThomasNodal<2>(1)-FE_DGQ<2>(1)]
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=8
@@ -2074,7 +2254,7 @@ DEAL::support on face 1:  2       3       4       5       6       7       8       9       10      11      13      15
 DEAL::support on face 2:       0       1       2       3       4       5       8       9       10      11      12      13
 DEAL::support on face 3:       0       1       2       3       6       7       8       9       10      11      14      15
 DEAL::
-DEAL::fe_data[40]:FESystem<2>[FE_Q<2>(3)^2]
+DEAL::fe_data[43]:FESystem<2>[FE_Q<2>(3)^2]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=4
 DEAL::dofs_per_quad=8
@@ -2105,7 +2285,7 @@ DEAL::support on face 1:  2       3       6       7       12      13      14      15
 DEAL::support on face 2:       0       1       2       3       16      17      18      19
 DEAL::support on face 3:       4       5       6       7       20      21      22      23
 DEAL::
-DEAL::fe_data[41]:FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(3)]
+DEAL::fe_data[44]:FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(3)]
 DEAL::dofs_per_vertex=3
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=4
@@ -2136,7 +2316,7 @@ DEAL::support on face 1:  3       4       5       9       10      11      14      15
 DEAL::support on face 2:       0       1       2       3       4       5       16      17
 DEAL::support on face 3:       6       7       8       9       10      11      18      19
 DEAL::
-DEAL::fe_data[42]:FESystem<2>[FE_Q<2>(4)^2]
+DEAL::fe_data[45]:FESystem<2>[FE_Q<2>(4)^2]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=6
 DEAL::dofs_per_quad=18
@@ -2167,7 +2347,7 @@ DEAL::support on face 1:  2       3       6       7       14      15      16      17      18      19
 DEAL::support on face 2:       0       1       2       3       20      21      22      23      24      25
 DEAL::support on face 3:       4       5       6       7       26      27      28      29      30      31
 DEAL::
-DEAL::fe_data[43]:FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]^2]
+DEAL::fe_data[46]:FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]^2]
 DEAL::dofs_per_vertex=4
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -2198,7 +2378,7 @@ DEAL::support on face 1:  4       5       6       7       12      13      14      15
 DEAL::support on face 2:       0       1       2       3       4       5       6       7
 DEAL::support on face 3:       8       9       10      11      12      13      14      15
 DEAL::
-DEAL::fe_data[44]:FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]-FESystem<2>[FE_DGQ<2>(1)^2]]
+DEAL::fe_data[47]:FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]-FESystem<2>[FE_DGQ<2>(1)^2]]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=8
@@ -2229,7 +2409,7 @@ DEAL::support on face 1:  2       3       6       7       9       11      13      15
 DEAL::support on face 2:       0       1       2       3       8       9       12      13
 DEAL::support on face 3:       4       5       6       7       10      11      14      15
 DEAL::
-DEAL::fe_data[45]:FESystem<2>[FESystem<2>[FE_Q<2>(1)-FE_Q<2>(2)]-FESystem<2>[FE_Q<2>(2)^2]-FESystem<2>[FE_DGQ<2>(2)^2]]
+DEAL::fe_data[48]:FESystem<2>[FESystem<2>[FE_Q<2>(1)-FE_Q<2>(2)]-FESystem<2>[FE_Q<2>(2)^2]-FESystem<2>[FE_DGQ<2>(2)^2]]
 DEAL::dofs_per_vertex=4
 DEAL::dofs_per_line=3
 DEAL::dofs_per_quad=21
@@ -2260,7 +2440,7 @@ DEAL::support on face 1:  4       5       6       7       12      13      14      15      19      20      21      33      36      39      42      45      48
 DEAL::support on face 2:       0       1       2       3       4       5       6       7       22      23      24      31      32      33      40      41      42
 DEAL::support on face 3:       8       9       10      11      12      13      14      15      25      26      27      37      38      39      46      47      48
 DEAL::
-DEAL::fe_data[46]:FESystem<2>[FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]^2]^2-FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]-FESystem<2>[FE_DGQ<2>(1)^2]]-FESystem<2>[FESystem<2>[FE_Q<2>(1)-FE_Q<2>(2)]-FESystem<2>[FE_Q<2>(2)^2]-FESystem<2>[FE_DGQ<2>(2)^2]]^2]
+DEAL::fe_data[49]:FESystem<2>[FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]^2]^2-FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]-FESystem<2>[FE_DGQ<2>(1)^2]]-FESystem<2>[FESystem<2>[FE_Q<2>(1)-FE_Q<2>(2)]-FESystem<2>[FE_Q<2>(2)^2]-FESystem<2>[FE_DGQ<2>(2)^2]]^2]
 DEAL::dofs_per_vertex=18
 DEAL::dofs_per_line=6
 DEAL::dofs_per_quad=50
@@ -2302,6 +2482,9 @@ DEAL::FE_DGQ<3>(2)
 DEAL::FE_DGQArbitraryNodes<3>(QIterated(QTrapez(),4))
 DEAL::FE_DGQ<3>(4)
 DEAL::FE_DGQArbitraryNodes<3>(QGauss(3))
+DEAL::FE_DGQLegendre<3>(1)
+DEAL::FE_DGQLegendre<3>(2)
+DEAL::FE_DGQHermite<3>(3)
 DEAL::FE_DGP<3>(1)
 DEAL::FE_DGP<3>(2)
 DEAL::FESystem<3>[FE_Q<3>(2)^2]
@@ -2720,7 +2903,112 @@ DEAL::support on face 3:        0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
 DEAL::
-DEAL::fe_data[11]:FE_DGP<3>(1)
+DEAL::fe_data[11]:FE_DGQLegendre<3>(1)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=0
+DEAL::dofs_per_quad=0
+DEAL::dofs_per_hex=8
+DEAL::first_line_index=0
+DEAL::first_quad_index=0
+DEAL::first_hex_index=0
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=8
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[8]->8
+DEAL::degree=1
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2       3       4       5       6       7
+DEAL::support on face 1:       0       1       2       3       4       5       6       7
+DEAL::support on face 2:       0       1       2       3       4       5       6       7
+DEAL::support on face 3:       0       1       2       3       4       5       6       7
+DEAL::support on face 4:       0       1       2       3       4       5       6       7
+DEAL::support on face 5:       0       1       2       3       4       5       6       7
+DEAL::
+DEAL::fe_data[12]:FE_DGQLegendre<3>(2)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=0
+DEAL::dofs_per_quad=0
+DEAL::dofs_per_hex=27
+DEAL::first_line_index=0
+DEAL::first_quad_index=0
+DEAL::first_hex_index=0
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=27
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[27]->27
+DEAL::degree=2
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
+DEAL::support on face 1:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
+DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
+DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
+DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
+DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26
+DEAL::
+DEAL::fe_data[13]:FE_DGQHermite<3>(3)
+DEAL::dofs_per_vertex=0
+DEAL::dofs_per_line=0
+DEAL::dofs_per_quad=0
+DEAL::dofs_per_hex=64
+DEAL::first_line_index=0
+DEAL::first_quad_index=0
+DEAL::first_hex_index=0
+DEAL::first_face_line_index=0
+DEAL::first_face_quad_index=0
+DEAL::dofs_per_face=0
+DEAL::dofs_per_cell=64
+DEAL::primitive=yes
+DEAL::components=1
+DEAL::blocks=1:[64]->64
+DEAL::degree=3
+DEAL::conformity= L2
+DEAL::unit_support_points=0
+DEAL::unit_face_support_points=0
+DEAL::generalized_support_points=0
+DEAL::generalized_face_support_points=0
+DEAL::face_to_equivalent_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::face_to_cell_index:
+DEAL::support on face 0:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63
+DEAL::support on face 1:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63
+DEAL::support on face 2:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63
+DEAL::support on face 3:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63
+DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63
+DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63
+DEAL::
+DEAL::fe_data[14]:FE_DGP<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -2755,7 +3043,7 @@ DEAL::support on face 3:  0       1       2       3
 DEAL::support on face 4:       0       1       2       3
 DEAL::support on face 5:       0       1       2       3
 DEAL::
-DEAL::fe_data[12]:FE_DGP<3>(2)
+DEAL::fe_data[15]:FE_DGP<3>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -2790,7 +3078,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9
 DEAL::
-DEAL::fe_data[13]:FESystem<3>[FE_Q<3>(2)^2]
+DEAL::fe_data[16]:FESystem<3>[FE_Q<3>(2)^2]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=2
@@ -2825,7 +3113,7 @@ DEAL::support on face 3:  4       5       6       7       12      13      14      15      22      23      30      31      36      37      38      39      46      47
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       16      17      18      19      20      21      22      23      48      49
 DEAL::support on face 5:       8       9       10      11      12      13      14      15      24      25      26      27      28      29      30      31      50      51
 DEAL::
-DEAL::fe_data[14]:FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)]
+DEAL::fe_data[17]:FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)]
 DEAL::dofs_per_vertex=3
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=1
@@ -2860,7 +3148,7 @@ DEAL::support on face 3:  6       7       8       9       10      11      18      19      20      21      22      23      27      31      34      35      39
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      24      25      26      27      40
 DEAL::support on face 5:       12      13      14      15      16      17      18      19      20      21      22      23      28      29      30      31      41
 DEAL::
-DEAL::fe_data[15]:FE_FaceQ<3>(0)
+DEAL::fe_data[18]:FE_FaceQ<3>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=1
@@ -2895,7 +3183,7 @@ DEAL::support on face 3:  3
 DEAL::support on face 4:       4
 DEAL::support on face 5:       5
 DEAL::
-DEAL::fe_data[16]:FE_FaceQ<3>(1)
+DEAL::fe_data[19]:FE_FaceQ<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=4
@@ -2930,7 +3218,7 @@ DEAL::support on face 3:  12      13      14      15
 DEAL::support on face 4:       16      17      18      19
 DEAL::support on face 5:       20      21      22      23
 DEAL::
-DEAL::fe_data[17]:FE_FaceQ<3>(3)
+DEAL::fe_data[20]:FE_FaceQ<3>(3)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=16
@@ -2965,7 +3253,7 @@ DEAL::support on face 3:  48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63
 DEAL::support on face 4:       64      65      66      67      68      69      70      71      72      73      74      75      76      77      78      79
 DEAL::support on face 5:       80      81      82      83      84      85      86      87      88      89      90      91      92      93      94      95
 DEAL::
-DEAL::fe_data[18]:FE_FaceP<3>(0)
+DEAL::fe_data[21]:FE_FaceP<3>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=1
@@ -3000,7 +3288,7 @@ DEAL::support on face 3:  3
 DEAL::support on face 4:       4
 DEAL::support on face 5:       5
 DEAL::
-DEAL::fe_data[19]:FE_FaceP<3>(1)
+DEAL::fe_data[22]:FE_FaceP<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=3
@@ -3035,7 +3323,7 @@ DEAL::support on face 3:  9       10      11
 DEAL::support on face 4:       12      13      14
 DEAL::support on face 5:       15      16      17
 DEAL::
-DEAL::fe_data[20]:FE_FaceP<3>(3)
+DEAL::fe_data[23]:FE_FaceP<3>(3)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=10
@@ -3070,7 +3358,7 @@ DEAL::support on face 3:  30      31      32      33      34      35      36      37      38      39
 DEAL::support on face 4:       40      41      42      43      44      45      46      47      48      49
 DEAL::support on face 5:       50      51      52      53      54      55      56      57      58      59
 DEAL::
-DEAL::fe_data[21]:FE_DGRaviartThomas<3>(0)
+DEAL::fe_data[24]:FE_DGRaviartThomas<3>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3105,7 +3393,7 @@ DEAL::support on face 3:  0       1       2       3       4       5
 DEAL::support on face 4:       0       1       2       3       4       5
 DEAL::support on face 5:       0       1       2       3       4       5
 DEAL::
-DEAL::fe_data[22]:FE_DGRaviartThomas<3>(1)
+DEAL::fe_data[25]:FE_DGRaviartThomas<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3140,7 +3428,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35
 DEAL::
-DEAL::fe_data[23]:FE_DGBDM<3>(1)
+DEAL::fe_data[26]:FE_DGBDM<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3175,7 +3463,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17
 DEAL::
-DEAL::fe_data[24]:FE_DGBDM<3>(2)
+DEAL::fe_data[27]:FE_DGBDM<3>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3210,7 +3498,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38
 DEAL::
-DEAL::fe_data[25]:FE_DGNedelec<3>(0)
+DEAL::fe_data[28]:FE_DGNedelec<3>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3245,7 +3533,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11
 DEAL::
-DEAL::fe_data[26]:FE_DGNedelec<3>(1)
+DEAL::fe_data[29]:FE_DGNedelec<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3280,7 +3568,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53
 DEAL::
-DEAL::fe_data[27]:FE_RaviartThomas<3>(0)
+DEAL::fe_data[30]:FE_RaviartThomas<3>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=1
@@ -3315,7 +3603,7 @@ DEAL::support on face 3:  0       1       2       3       4       5
 DEAL::support on face 4:       0       1       2       3       4       5
 DEAL::support on face 5:       0       1       2       3       4       5
 DEAL::
-DEAL::fe_data[28]:FE_RaviartThomas<3>(1)
+DEAL::fe_data[31]:FE_RaviartThomas<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=4
@@ -3350,7 +3638,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35
 DEAL::
-DEAL::fe_data[29]:FE_RaviartThomas<3>(2)
+DEAL::fe_data[32]:FE_RaviartThomas<3>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=9
@@ -3385,7 +3673,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63      64      65      66      67      68      69      70      71      72      73      74      75      76      77      78      79      80      81      82      83      84      85      86      87      88      89      90      91      92      93      94      95      96      97      98      99      100     101     102     103     104     105     106     107
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39      40      41      42      43      44      45      46      47      48      49      50      51      52      53      54      55      56      57      58      59      60      61      62      63      64      65      66      67      68      69      70      71      72      73      74      75      76      77      78      79      80      81      82      83      84      85      86      87      88      89      90      91      92      93      94      95      96      97      98      99      100     101     102     103     104     105     106     107
 DEAL::
-DEAL::fe_data[30]:FESystem<3>[FE_RaviartThomas<3>(1)-FE_DGQ<3>(1)]
+DEAL::fe_data[33]:FESystem<3>[FE_RaviartThomas<3>(1)-FE_DGQ<3>(1)]
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=4
@@ -3420,7 +3708,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      40      41      42      43
 DEAL::
-DEAL::fe_data[31]:FE_BDM<3>(1)
+DEAL::fe_data[34]:FE_BDM<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=3
@@ -3455,7 +3743,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17
 DEAL::
-DEAL::fe_data[32]:FE_BDM<3>(2)
+DEAL::fe_data[35]:FE_BDM<3>(2)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=6
@@ -3490,7 +3778,7 @@ DEAL::support on face 3:  0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38
 DEAL::
-DEAL::fe_data[33]:FE_Nedelec<3>(0)
+DEAL::fe_data[36]:FE_Nedelec<3>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=1
 DEAL::dofs_per_quad=0
@@ -3525,7 +3813,7 @@ DEAL::support on face 3:  0       1       3       4       5       7       10      11
 DEAL::support on face 4:       0       1       2       3       8       9       10      11
 DEAL::support on face 5:       4       5       6       7       8       9       10      11
 DEAL::
-DEAL::fe_data[34]:FE_Nedelec<3>(1)
+DEAL::fe_data[37]:FE_Nedelec<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=2
 DEAL::dofs_per_quad=4
@@ -3560,7 +3848,7 @@ DEAL::support on face 3:  0       1       2       3       6       7       8       9       10      11      14      15      20      21      22      23      26      27      30      31      36
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       16      17      18      19      20      21      22      23      24      25      28      29      34      35      38      39      40      41      42      43
 DEAL::support on face 5:       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      28      29      34      35      38      39      44      45      46      47
 DEAL::
-DEAL::fe_data[35]:FE_RaviartThomasNodal<3>(0)
+DEAL::fe_data[38]:FE_RaviartThomasNodal<3>(0)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=1
@@ -3595,7 +3883,7 @@ DEAL::support on face 3:  0       1       3       4       5
 DEAL::support on face 4:       0       1       2       3       4
 DEAL::support on face 5:       0       1       2       3       5
 DEAL::
-DEAL::fe_data[36]:FE_RaviartThomasNodal<3>(1)
+DEAL::fe_data[39]:FE_RaviartThomasNodal<3>(1)
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=4
@@ -3630,7 +3918,7 @@ DEAL::support on face 3:  0       1       2       3       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35
 DEAL::
-DEAL::fe_data[37]:FESystem<3>[FE_RaviartThomasNodal<3>(1)-FE_DGQ<3>(1)]
+DEAL::fe_data[40]:FESystem<3>[FE_RaviartThomasNodal<3>(1)-FE_DGQ<3>(1)]
 DEAL::dofs_per_vertex=0
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=4
@@ -3665,7 +3953,7 @@ DEAL::support on face 3:  0       1       2       3       6       7       8       9       10      11      12      13      14      15      16      17      18      19      20      21      22
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      36      37      38      39
 DEAL::support on face 5:       0       1       2       3       4       5       6       7       10      11      12      13      14      15      16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      32      33      34      35      40      41      42      43
 DEAL::
-DEAL::fe_data[38]:FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]^2]
+DEAL::fe_data[41]:FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]^2]
 DEAL::dofs_per_vertex=4
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3700,7 +3988,7 @@ DEAL::support on face 3:  8       9       10      11      12      13      14      15      24      25      26      27      28      29      30      31
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15
 DEAL::support on face 5:       16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31
 DEAL::
-DEAL::fe_data[39]:FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]-FESystem<3>[FE_DGQ<3>(1)^2]]
+DEAL::fe_data[42]:FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]-FESystem<3>[FE_DGQ<3>(1)^2]]
 DEAL::dofs_per_vertex=2
 DEAL::dofs_per_line=0
 DEAL::dofs_per_quad=0
@@ -3735,7 +4023,7 @@ DEAL::support on face 3:  4       5       6       7       12      13      14      15      18      19      22      23      26      27      30      31
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       16      17      18      19      24      25      26      27
 DEAL::support on face 5:       8       9       10      11      12      13      14      15      20      21      22      23      28      29      30      31
 DEAL::
-DEAL::fe_data[40]:FESystem<3>[FESystem<3>[FE_Q<3>(1)-FE_Q<3>(2)]-FESystem<3>[FE_Q<3>(2)^2]-FESystem<3>[FE_DGQ<3>(2)^2]]
+DEAL::fe_data[43]:FESystem<3>[FESystem<3>[FE_Q<3>(1)-FE_Q<3>(2)]-FESystem<3>[FE_Q<3>(2)^2]-FESystem<3>[FE_DGQ<3>(2)^2]]
 DEAL::dofs_per_vertex=4
 DEAL::dofs_per_line=3
 DEAL::dofs_per_quad=3
@@ -3770,7 +4058,7 @@ DEAL::support on face 3:  8       9       10      11      12      13      14      15      24      25      26      27      28      29      30      31      41      42      43
 DEAL::support on face 4:       0       1       2       3       4       5       6       7       8       9       10      11      12      13      14      15      32      33      34      35      36      37      38      39      40      41      42      43      80      81      82      89      90      91      92      93      94      95      96      97      116     117     118     119     120     121     122     123     124
 DEAL::support on face 5:       16      17      18      19      20      21      22      23      24      25      26      27      28      29      30      31      44      45      46      47      48      49      50      51      52      53      54      55      83      84      85      107     108     109     110     111     112     113     114     115     134     135     136     137     138     139     140     141     142
 DEAL::
-DEAL::fe_data[41]:FESystem<3>[FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]^2]^2-FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]-FESystem<3>[FE_DGQ<3>(1)^2]]-FESystem<3>[FESystem<3>[FE_Q<3>(1)-FE_Q<3>(2)]-FESystem<3>[FE_Q<3>(2)^2]-FESystem<3>[FE_DGQ<3>(2)^2]]^2]
+DEAL::fe_data[44]:FESystem<3>[FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]^2]^2-FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]-FESystem<3>[FE_DGQ<3>(1)^2]]-FESystem<3>[FESystem<3>[FE_Q<3>(1)-FE_Q<3>(2)]-FESystem<3>[FE_Q<3>(2)^2]-FESystem<3>[FE_DGQ<3>(2)^2]]^2]
 DEAL::dofs_per_vertex=18
 DEAL::dofs_per_line=6
 DEAL::dofs_per_quad=6
index f01b5c2c3823d67252801882a788bc2e12795a04..4e5a0111526c888eef41e1f508e21088b17e6bf7 100644 (file)
@@ -119,6 +119,7 @@ int main ()
   gen.generate_all_codim("FE_Bernstein(1)");
   gen.generate_all_codim("FE_DGP(1)");
   gen.generate_all_codim("FE_DGQ(1)");
+  gen.generate_all_codim("FE_DGQLegendre(1)");
   gen.generate_all_codim("FE_Nothing()");
   gen.generate_all_codim("FE_DGQArbitraryNodes(QGauss(2))");
   gen.generate_all_codim("FE_Q_Bubbles(1)");
@@ -134,6 +135,7 @@ int main ()
   gen.generate_all_codim("FE_Q_iso_Q1(2)");
   gen.generate_all_codim("FE_Q(2)");
   gen.generate_all_codim("FE_Bernstein(2)");
+  gen.generate_all_codim("FE_DGQHermite(3)");
 
   //systems
   gen.generate_all_dim("FESystem[FE_Q_Hierarchical(1)^2-FE_Q_Hierarchical(1)]");
index 55458b2cd171c7d37e72a1201632339959883f50..8bb30049ca42de7187daea72ea1d1ed7cf9fdb15 100644 (file)
@@ -236,6 +236,24 @@ DEAL::FE_DGQ<2,3>(1)
 DEAL::Read FE_DGQ(1)
 DEAL::Generated :
 DEAL::FE_DGQ<3>(1)
+DEAL::Read FE_DGQLegendre(1)
+DEAL::Generated :
+DEAL::FE_DGQLegendre<1>(1)
+DEAL::Read FE_DGQLegendre(1)
+DEAL::Generated :
+DEAL::FE_DGQLegendre<1,2>(1)
+DEAL::Read FE_DGQLegendre(1)
+DEAL::Generated :
+DEAL::FE_DGQLegendre<2>(1)
+DEAL::Read FE_DGQLegendre(1)
+DEAL::Generated :
+DEAL::FE_DGQLegendre<1,3>(1)
+DEAL::Read FE_DGQLegendre(1)
+DEAL::Generated :
+DEAL::FE_DGQLegendre<2,3>(1)
+DEAL::Read FE_DGQLegendre(1)
+DEAL::Generated :
+DEAL::FE_DGQLegendre<3>(1)
 DEAL::Read FE_Nothing()
 DEAL::Generated :
 DEAL::FE_Nothing<1>()
@@ -506,6 +524,24 @@ DEAL::FE_Bernstein<2>(2)
 DEAL::Read FE_Bernstein(2)
 DEAL::Generated :
 DEAL::FE_Bernstein<3>(2)
+DEAL::Read FE_DGQHermite(3)
+DEAL::Generated :
+DEAL::FE_DGQHermite<1>(3)
+DEAL::Read FE_DGQHermite(3)
+DEAL::Generated :
+DEAL::FE_DGQHermite<1,2>(3)
+DEAL::Read FE_DGQHermite(3)
+DEAL::Generated :
+DEAL::FE_DGQHermite<2>(3)
+DEAL::Read FE_DGQHermite(3)
+DEAL::Generated :
+DEAL::FE_DGQHermite<1,3>(3)
+DEAL::Read FE_DGQHermite(3)
+DEAL::Generated :
+DEAL::FE_DGQHermite<2,3>(3)
+DEAL::Read FE_DGQHermite(3)
+DEAL::Generated :
+DEAL::FE_DGQHermite<3>(3)
 DEAL::Read FESystem[FE_Q_Hierarchical(1)^2-FE_Q_Hierarchical(1)]
 DEAL::Generated :
 DEAL::FESystem<1>[FE_Q_Hierarchical<1>(1)^2-FE_Q_Hierarchical<1>(1)]

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.