--- /dev/null
+Fixed: Previously, FiniteElement::has_generalized_support_points() and
+FiniteElement::has_face_support_points() returned false for the FE_Nothing element. Now, FE_Nothing::has_generalized_support_points() and
+FE_Nothing::has_face_support_points() correctly return true, as the empty
+arrays returned by FE_Nothing::get_generalized_support_points() and FE_Nothing::get_unit_face_support_points() accurately describe the support
+points of the element(i.e., they don't have any, as there are no
+degrees of freedom).
+<br>
+(Oreste Marquis, 2024/08/21)
/**
* Return whether a finite element has defined support points on faces. If
* the result is true, then a call to the get_unit_face_support_points()
- * yields a non-empty vector.
+ * yields an array with `dofs_per_face` entries. Note that the function's name
+ * is poorly chosen: The function does not return *whether* an element has
+ * support points, but whether its implementation is able to *provide*
+ * a list of support points via get_unit_face_support_points().
*
* For more information, see the documentation for the has_support_points()
* function.
/**
* Return whether a finite element has defined generalized support
* points. If the result is true, then a call to the
- * get_generalized_support_points() yields a non-empty vector.
+ * get_generalized_support_points() yields an
+ * array with minimal set of *unique* support points. Note that the function's
+ * name is poorly chosen: The function does not return *whether* an element
+ * has support points, but whether its implementation is able to *provide*
+ * a list of support points via get_unit_support_points().
*
* See the
* @ref GlossGeneralizedSupport "glossary entry on generalized support points"
- * for more information.
+ * or the documentation for the has_support_points() function for more
+ * information.
*/
bool
has_generalized_support_points() const;
bool
FiniteElement<dim, spacedim>::has_generalized_support_points() const
{
- return (get_generalized_support_points().size() != 0);
+ if (this->dofs_per_cell > 0)
+ return (get_generalized_support_points().size() != 0);
+ else
+ {
+ // If the FE has no DoFs, the array size should be zero:
+ AssertDimension(get_generalized_support_points().size(), 0);
+
+ // A finite element without DoFs *has* generalized support points
+ // (which is then an empty array)
+ return true;
+ }
}
FiniteElement<dim, spacedim>::has_face_support_points(
const unsigned int face_no) const
{
- return (unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
- .size() != 0);
+ const unsigned int face_index = this->n_unique_faces() == 1 ? 0 : face_no;
+ if (this->n_dofs_per_face(face_index) > 0)
+ return (unit_face_support_points[face_index].size() != 0);
+ else
+ {
+ // If the FE has no DoFs on face, the array size should be zero
+ AssertDimension(unit_face_support_points[face_index].size(), 0);
+
+ // A finite element without DoFs *has* face support points
+ // (which is then an empty array)
+ return true;
+ }
}
Assert(source_fe.n_components() == this->n_components(),
ExcDimensionMismatch(source_fe.n_components(), this->n_components()));
- if (source_fe.has_face_support_points(face_no))
+ if ((source_fe.has_face_support_points(face_no)) &&
+ (source_fe.n_dofs_per_face(face_no) > 0))
{
// have this test in here since a table of size 2x0 reports its size as
// 0x0
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2001 - 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+// Test FiniteElement::has_support_points() for various elements,
+// including FE_Nothing.
+
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/fe/fe_bdm.h>
+#include <deal.II/fe/fe_dg_vector.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_face.h>
+#include <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_q_hierarchical.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <deal.II/grid/tria_iterator.h>
+
+#include <iostream>
+
+#include "../tests.h"
+
+template <int dim>
+void
+test_2d_3d(std::vector<FiniteElement<dim> *> &finite_elements)
+{
+ // Vector DG elements
+ finite_elements.push_back(new FE_DGRaviartThomas<dim>(0));
+ finite_elements.push_back(new FE_DGRaviartThomas<dim>(1));
+ finite_elements.push_back(new FE_DGBDM<dim>(1));
+ finite_elements.push_back(new FE_DGBDM<dim>(2));
+ finite_elements.push_back(new FE_DGNedelec<dim>(0));
+ finite_elements.push_back(new FE_DGNedelec<dim>(1));
+
+ // Hdiv elements
+ FE_RaviartThomas<dim> *rt0 = new FE_RaviartThomas<dim>(0);
+ finite_elements.push_back(rt0);
+
+ FE_RaviartThomas<dim> *rt1 = new FE_RaviartThomas<dim>(1);
+ finite_elements.push_back(rt1);
+
+ finite_elements.push_back(new FE_RaviartThomas<dim>(2));
+ finite_elements.push_back(new FESystem<dim>(*rt1, 1, FE_DGQ<dim>(1), 1));
+
+ finite_elements.push_back(new FE_BDM<dim>(1));
+ finite_elements.push_back(new FE_BDM<dim>(2));
+
+ // Hcurl elements
+ FE_Nedelec<dim> *ned0 = new FE_Nedelec<dim>(0);
+ finite_elements.push_back(ned0);
+ FE_Nedelec<dim> *ned1 = new FE_Nedelec<dim>(1);
+ finite_elements.push_back(ned1);
+}
+
+void test_2d_3d(std::vector<FiniteElement<1> *> /*&
+ finite_elements*/)
+{}
+
+template <int dim>
+void
+test_finite_elements()
+{
+ std::vector<FiniteElement<dim> *> finite_elements;
+ finite_elements.push_back(new FE_Nothing<dim>());
+ finite_elements.push_back(new FE_Q<dim>(1));
+ finite_elements.push_back(new FE_Q<dim>(2));
+ finite_elements.push_back(new FE_Q<dim>(4));
+ finite_elements.push_back(new FE_Q_Hierarchical<dim>(1));
+ finite_elements.push_back(new FE_Q_Hierarchical<dim>(2));
+ finite_elements.push_back(new FE_Q_Hierarchical<dim>(4));
+ finite_elements.push_back(new FE_DGQ<dim>(1));
+ finite_elements.push_back(new FE_DGQ<dim>(2));
+ finite_elements.push_back(
+ new FE_DGQArbitraryNodes<dim>(QIterated<1>(QTrapezoid<1>(), 4)));
+ finite_elements.push_back(new FE_DGQ<dim>(4));
+ finite_elements.push_back(new FE_DGQArbitraryNodes<dim>(QGauss<1>(3)));
+ finite_elements.push_back(new FE_DGQLegendre<dim>(1));
+ finite_elements.push_back(new FE_DGQLegendre<dim>(2));
+ finite_elements.push_back(new FE_DGQHermite<dim>(3));
+ finite_elements.push_back(new FE_DGP<dim>(1));
+ finite_elements.push_back(new FE_DGP<dim>(2));
+ finite_elements.push_back(new FESystem<dim>(FE_Q<dim>(2), 2));
+ finite_elements.push_back(
+ new FESystem<dim>(FE_Q<dim>(1), 2, FE_Q<dim>(2), 1));
+ finite_elements.push_back(
+ new FESystem<dim>(FE_Q<dim>(1), 2, FE_Q<dim>(2), 1, FE_Nothing<dim>(), 2));
+
+ // Face Q elements
+ finite_elements.push_back(new FE_FaceQ<dim>(0));
+ finite_elements.push_back(new FE_FaceQ<dim>(1));
+ finite_elements.push_back(new FE_FaceQ<dim>(3));
+ // Face P elements
+ finite_elements.push_back(new FE_FaceP<dim>(0));
+ finite_elements.push_back(new FE_FaceP<dim>(1));
+ finite_elements.push_back(new FE_FaceP<dim>(3));
+
+ // Check vector elements in 2d and higher only
+ test_2d_3d(finite_elements);
+
+ if (dim == 2)
+ {
+ finite_elements.push_back(new FE_DGBDM<dim>(1));
+ finite_elements.push_back(new FE_DGBDM<dim>(2));
+ }
+ if (dim > 1)
+ {
+ FE_RaviartThomasNodal<dim> *rt0 = new FE_RaviartThomasNodal<dim>(0);
+ FE_RaviartThomasNodal<dim> *rt1 = new FE_RaviartThomasNodal<dim>(1);
+ finite_elements.push_back(rt0);
+ finite_elements.push_back(rt1);
+ finite_elements.push_back(new FESystem<dim>(*rt1, 1, FE_DGQ<dim>(1), 1));
+ finite_elements.push_back(
+ new FESystem<dim>(*rt1, 1, FE_DGQ<dim>(1), 1, FE_Nothing<dim>(), 1));
+ }
+
+ finite_elements.push_back(new FESystem<dim>(FE_Q<dim>(3), 2));
+ finite_elements.push_back(
+ new FESystem<dim>(FE_Q<dim>(1), 2, FE_Q<dim>(3), 1));
+ finite_elements.push_back(new FESystem<dim>(FE_Q<dim>(4), 2));
+
+ // have systems of systems, and
+ // construct hierarchies of
+ // subsequently weirder elements by
+ // taking each of them in turn as
+ // basis of others
+ finite_elements.push_back(
+ new FESystem<dim>(FESystem<dim>(FE_Q<dim>(1), 2), 2));
+ finite_elements.push_back(new FESystem<dim>(
+ FESystem<dim>(FE_Q<dim>(1), 2), 1, FESystem<dim>(FE_DGQ<dim>(1), 2), 1));
+ finite_elements.push_back(
+ new FESystem<dim>(FESystem<dim>(FE_Q<dim>(1), 1, FE_Q<dim>(2), 1),
+ 1,
+ FESystem<dim>(FE_Q<dim>(2), 2),
+ 1,
+ FESystem<dim>(FE_DGQ<dim>(2), 2),
+ 1));
+ finite_elements.push_back(
+ new FESystem<dim>(*finite_elements[finite_elements.size() - 3],
+ 2,
+ *finite_elements[finite_elements.size() - 2],
+ 1,
+ *finite_elements[finite_elements.size() - 1],
+ 2));
+
+ deallog << std::endl << "dim=" << dim << std::endl;
+ for (unsigned int n = 0; n < finite_elements.size(); ++n)
+ {
+ FiniteElement<dim> *fe_data = finite_elements[n];
+ deallog << " " << fe_data->get_name() << std::endl;
+ deallog << " has_face_support_points="
+ << (fe_data->has_face_support_points() ? "true" : "false")
+ << std::endl;
+ }
+
+ // delete all FiniteElement objects
+ for (unsigned int i = 0; i < finite_elements.size(); ++i)
+ delete finite_elements[i];
+}
+
+int
+main()
+{
+ initlog();
+
+ test_finite_elements<1>();
+ test_finite_elements<2>();
+ test_finite_elements<3>();
+}
--- /dev/null
+
+DEAL::
+DEAL::dim=1
+DEAL:: FE_Nothing<1>()
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<1>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<1>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<1>(4)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q_Hierarchical<1>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Q_Hierarchical<1>(2)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Q_Hierarchical<1>(4)
+DEAL:: has_face_support_points=false
+DEAL:: FE_DGQ<1>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQ<1>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQArbitraryNodes<1>(QIterated(QTrapezoid(),4))
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQ<1>(4)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQArbitraryNodes<1>(QGauss(3))
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQLegendre<1>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQLegendre<1>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQHermite<1>(3)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGP<1>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGP<1>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<1>[FE_Q<1>(2)^2]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_Nothing<1>()^2]
+DEAL:: has_face_support_points=false
+DEAL:: FE_FaceQ<1>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<1>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<1>(3)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceP<1>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceP<1>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceP<1>(3)
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<1>[FE_Q<1>(3)^2]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(3)]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<1>[FE_Q<1>(4)^2]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]^2]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]-FESystem<1>[FE_DGQ<1>(1)^2]]
+DEAL:: has_face_support_points=false
+DEAL:: 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:: has_face_support_points=false
+DEAL:: 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:: has_face_support_points=false
+DEAL::
+DEAL::dim=2
+DEAL:: FE_Nothing<2>()
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<2>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<2>(4)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q_Hierarchical<2>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Q_Hierarchical<2>(2)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Q_Hierarchical<2>(4)
+DEAL:: has_face_support_points=false
+DEAL:: FE_DGQ<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQ<2>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQArbitraryNodes<2>(QIterated(QTrapezoid(),4))
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQ<2>(4)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQArbitraryNodes<2>(QGauss(3))
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQLegendre<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQLegendre<2>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQHermite<2>(3)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGP<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGP<2>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(2)^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_Nothing<2>()^2]
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<2>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<2>(3)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceP<2>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_FaceP<2>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_FaceP<2>(3)
+DEAL:: has_face_support_points=false
+DEAL:: FE_DGRaviartThomas<2>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGRaviartThomas<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGBDM<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGBDM<2>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGNedelec<2>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGNedelec<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_RaviartThomas<2>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_RaviartThomas<2>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_RaviartThomas<2>(2)
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<2>[FE_RaviartThomas<2>(1)-FE_DGQ<2>(1)]
+DEAL:: has_face_support_points=false
+DEAL:: FE_BDM<2>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_BDM<2>(2)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Nedelec<2>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Nedelec<2>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_DGBDM<2>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGBDM<2>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_RaviartThomasNodal<2>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_RaviartThomasNodal<2>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<2>[FE_RaviartThomasNodal<2>(1)-FE_DGQ<2>(1)]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<2>[FE_RaviartThomasNodal<2>(1)-FE_DGQ<2>(1)-FE_Nothing<2>()]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<2>[FE_Q<2>(3)^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(3)]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(4)^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]-FESystem<2>[FE_DGQ<2>(1)^2]]
+DEAL:: has_face_support_points=true
+DEAL:: 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:: has_face_support_points=true
+DEAL:: 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:: has_face_support_points=true
+DEAL::
+DEAL::dim=3
+DEAL:: FE_Nothing<3>()
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<3>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q<3>(4)
+DEAL:: has_face_support_points=true
+DEAL:: FE_Q_Hierarchical<3>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Q_Hierarchical<3>(2)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Q_Hierarchical<3>(4)
+DEAL:: has_face_support_points=false
+DEAL:: FE_DGQ<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQ<3>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQArbitraryNodes<3>(QIterated(QTrapezoid(),4))
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQ<3>(4)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQArbitraryNodes<3>(QGauss(3))
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQLegendre<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQLegendre<3>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGQHermite<3>(3)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGP<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGP<3>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(2)^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_Nothing<3>()^2]
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<3>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceQ<3>(3)
+DEAL:: has_face_support_points=true
+DEAL:: FE_FaceP<3>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_FaceP<3>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_FaceP<3>(3)
+DEAL:: has_face_support_points=false
+DEAL:: FE_DGRaviartThomas<3>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGRaviartThomas<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGBDM<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGBDM<3>(2)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGNedelec<3>(0)
+DEAL:: has_face_support_points=true
+DEAL:: FE_DGNedelec<3>(1)
+DEAL:: has_face_support_points=true
+DEAL:: FE_RaviartThomas<3>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_RaviartThomas<3>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_RaviartThomas<3>(2)
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<3>[FE_RaviartThomas<3>(1)-FE_DGQ<3>(1)]
+DEAL:: has_face_support_points=false
+DEAL:: FE_BDM<3>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_BDM<3>(2)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Nedelec<3>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_Nedelec<3>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FE_RaviartThomasNodal<3>(0)
+DEAL:: has_face_support_points=false
+DEAL:: FE_RaviartThomasNodal<3>(1)
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<3>[FE_RaviartThomasNodal<3>(1)-FE_DGQ<3>(1)]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<3>[FE_RaviartThomasNodal<3>(1)-FE_DGQ<3>(1)-FE_Nothing<3>()]
+DEAL:: has_face_support_points=false
+DEAL:: FESystem<3>[FE_Q<3>(3)^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(3)]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(4)^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]^2]
+DEAL:: has_face_support_points=true
+DEAL:: FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]-FESystem<3>[FE_DGQ<3>(1)^2]]
+DEAL:: has_face_support_points=true
+DEAL:: 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:: has_face_support_points=true
+DEAL:: 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:: has_face_support_points=true
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2001 - 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+// Test FiniteElement::has_support_points() for various elements,
+// including FE_Nothing.
+
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/fe/fe_bdm.h>
+#include <deal.II/fe/fe_dg_vector.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_face.h>
+#include <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_q_hierarchical.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <deal.II/grid/tria_iterator.h>
+
+#include <iostream>
+
+#include "../tests.h"
+
+template <int dim>
+void
+test_2d_3d(std::vector<FiniteElement<dim> *> &finite_elements)
+{
+ // Vector DG elements
+ finite_elements.push_back(new FE_DGRaviartThomas<dim>(0));
+ finite_elements.push_back(new FE_DGRaviartThomas<dim>(1));
+ finite_elements.push_back(new FE_DGBDM<dim>(1));
+ finite_elements.push_back(new FE_DGBDM<dim>(2));
+ finite_elements.push_back(new FE_DGNedelec<dim>(0));
+ finite_elements.push_back(new FE_DGNedelec<dim>(1));
+
+ // Hdiv elements
+ FE_RaviartThomas<dim> *rt0 = new FE_RaviartThomas<dim>(0);
+ finite_elements.push_back(rt0);
+
+ FE_RaviartThomas<dim> *rt1 = new FE_RaviartThomas<dim>(1);
+ finite_elements.push_back(rt1);
+
+ finite_elements.push_back(new FE_RaviartThomas<dim>(2));
+ finite_elements.push_back(new FESystem<dim>(*rt1, 1, FE_DGQ<dim>(1), 1));
+
+ finite_elements.push_back(new FE_BDM<dim>(1));
+ finite_elements.push_back(new FE_BDM<dim>(2));
+
+ // Hcurl elements
+ FE_Nedelec<dim> *ned0 = new FE_Nedelec<dim>(0);
+ finite_elements.push_back(ned0);
+ FE_Nedelec<dim> *ned1 = new FE_Nedelec<dim>(1);
+ finite_elements.push_back(ned1);
+}
+
+void
+test_2d_3d(std::vector<FiniteElement<1> *> /*&finite_elements*/)
+{}
+
+template <int dim>
+void
+test_finite_elements()
+{
+ std::vector<FiniteElement<dim> *> finite_elements;
+ finite_elements.push_back(new FE_Nothing<dim>());
+ finite_elements.push_back(new FE_Q<dim>(1));
+ finite_elements.push_back(new FE_Q<dim>(2));
+ finite_elements.push_back(new FE_Q<dim>(4));
+ finite_elements.push_back(new FE_Q_Hierarchical<dim>(1));
+ finite_elements.push_back(new FE_Q_Hierarchical<dim>(2));
+ finite_elements.push_back(new FE_Q_Hierarchical<dim>(4));
+ finite_elements.push_back(new FE_DGQ<dim>(1));
+ finite_elements.push_back(new FE_DGQ<dim>(2));
+ finite_elements.push_back(
+ new FE_DGQArbitraryNodes<dim>(QIterated<1>(QTrapezoid<1>(), 4)));
+ finite_elements.push_back(new FE_DGQ<dim>(4));
+ finite_elements.push_back(new FE_DGQArbitraryNodes<dim>(QGauss<1>(3)));
+ finite_elements.push_back(new FE_DGQLegendre<dim>(1));
+ finite_elements.push_back(new FE_DGQLegendre<dim>(2));
+ finite_elements.push_back(new FE_DGQHermite<dim>(3));
+ finite_elements.push_back(new FE_DGP<dim>(1));
+ finite_elements.push_back(new FE_DGP<dim>(2));
+ finite_elements.push_back(new FESystem<dim>(FE_Q<dim>(2), 2));
+ finite_elements.push_back(
+ new FESystem<dim>(FE_Q<dim>(1), 2, FE_Q<dim>(2), 1));
+ finite_elements.push_back(
+ new FESystem<dim>(FE_Q<dim>(1), 2, FE_Q<dim>(2), 1, FE_Nothing<dim>(), 2));
+
+ // Face Q elements
+ finite_elements.push_back(new FE_FaceQ<dim>(0));
+ finite_elements.push_back(new FE_FaceQ<dim>(1));
+ finite_elements.push_back(new FE_FaceQ<dim>(3));
+ // Face P elements
+ finite_elements.push_back(new FE_FaceP<dim>(0));
+ finite_elements.push_back(new FE_FaceP<dim>(1));
+ finite_elements.push_back(new FE_FaceP<dim>(3));
+
+ // Check vector elements in 2d and higher only
+ test_2d_3d(finite_elements);
+
+ if (dim == 2)
+ {
+ finite_elements.push_back(new FE_DGBDM<dim>(1));
+ finite_elements.push_back(new FE_DGBDM<dim>(2));
+ }
+ if (dim > 1)
+ {
+ FE_RaviartThomasNodal<dim> *rt0 = new FE_RaviartThomasNodal<dim>(0);
+ FE_RaviartThomasNodal<dim> *rt1 = new FE_RaviartThomasNodal<dim>(1);
+ finite_elements.push_back(rt0);
+ finite_elements.push_back(rt1);
+ finite_elements.push_back(new FESystem<dim>(*rt1, 1, FE_DGQ<dim>(1), 1));
+ finite_elements.push_back(
+ new FESystem<dim>(*rt1, 1, FE_DGQ<dim>(1), 1, FE_Nothing<dim>(), 1));
+ }
+
+ finite_elements.push_back(new FESystem<dim>(FE_Q<dim>(3), 2));
+ finite_elements.push_back(
+ new FESystem<dim>(FE_Q<dim>(1), 2, FE_Q<dim>(3), 1));
+ finite_elements.push_back(new FESystem<dim>(FE_Q<dim>(4), 2));
+
+ // have systems of systems, and
+ // construct hierarchies of
+ // subsequently weirder elements by
+ // taking each of them in turn as
+ // basis of others
+ finite_elements.push_back(
+ new FESystem<dim>(FESystem<dim>(FE_Q<dim>(1), 2), 2));
+ finite_elements.push_back(new FESystem<dim>(
+ FESystem<dim>(FE_Q<dim>(1), 2), 1, FESystem<dim>(FE_DGQ<dim>(1), 2), 1));
+ finite_elements.push_back(
+ new FESystem<dim>(FESystem<dim>(FE_Q<dim>(1), 1, FE_Q<dim>(2), 1),
+ 1,
+ FESystem<dim>(FE_Q<dim>(2), 2),
+ 1,
+ FESystem<dim>(FE_DGQ<dim>(2), 2),
+ 1));
+ finite_elements.push_back(
+ new FESystem<dim>(*finite_elements[finite_elements.size() - 3],
+ 2,
+ *finite_elements[finite_elements.size() - 2],
+ 1,
+ *finite_elements[finite_elements.size() - 1],
+ 2));
+
+ deallog << std::endl << "dim=" << dim << std::endl;
+ for (unsigned int n = 0; n < finite_elements.size(); ++n)
+ {
+ FiniteElement<dim> *fe_data = finite_elements[n];
+ deallog << " " << fe_data->get_name() << std::endl;
+ deallog << " has_generalized_support_points="
+ << (fe_data->has_generalized_support_points() ? "true" : "false")
+ << std::endl;
+ }
+
+ // delete all FiniteElement objects
+ for (unsigned int i = 0; i < finite_elements.size(); ++i)
+ delete finite_elements[i];
+}
+
+int
+main()
+{
+ initlog();
+
+ test_finite_elements<1>();
+ test_finite_elements<2>();
+ test_finite_elements<3>();
+}
--- /dev/null
+
+DEAL::
+DEAL::dim=1
+DEAL:: FE_Nothing<1>()
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<1>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<1>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<1>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<1>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<1>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<1>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<1>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<1>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQArbitraryNodes<1>(QIterated(QTrapezoid(),4))
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<1>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQArbitraryNodes<1>(QGauss(3))
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQLegendre<1>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGQLegendre<1>(2)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGQHermite<1>(3)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGP<1>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGP<1>(2)
+DEAL:: has_generalized_support_points=false
+DEAL:: FESystem<1>[FE_Q<1>(2)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_Nothing<1>()^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<1>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<1>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<1>(3)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceP<1>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceP<1>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceP<1>(3)
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<1>[FE_Q<1>(3)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(3)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<1>[FE_Q<1>(4)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<1>[FESystem<1>[FE_Q<1>(1)^2]-FESystem<1>[FE_DGQ<1>(1)^2]]
+DEAL:: has_generalized_support_points=true
+DEAL:: 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:: has_generalized_support_points=true
+DEAL:: 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:: has_generalized_support_points=true
+DEAL::
+DEAL::dim=2
+DEAL:: FE_Nothing<2>()
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<2>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<2>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<2>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<2>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<2>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQArbitraryNodes<2>(QIterated(QTrapezoid(),4))
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<2>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQArbitraryNodes<2>(QGauss(3))
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQLegendre<2>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGQLegendre<2>(2)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGQHermite<2>(3)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGP<2>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGP<2>(2)
+DEAL:: has_generalized_support_points=false
+DEAL:: FESystem<2>[FE_Q<2>(2)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_Nothing<2>()^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<2>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<2>(3)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceP<2>(0)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_FaceP<2>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_FaceP<2>(3)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGRaviartThomas<2>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGRaviartThomas<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGBDM<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGBDM<2>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGNedelec<2>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGNedelec<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomas<2>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomas<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomas<2>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_RaviartThomas<2>(1)-FE_DGQ<2>(1)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_BDM<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_BDM<2>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Nedelec<2>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Nedelec<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGBDM<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGBDM<2>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomasNodal<2>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomasNodal<2>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_RaviartThomasNodal<2>(1)-FE_DGQ<2>(1)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_RaviartThomasNodal<2>(1)-FE_DGQ<2>(1)-FE_Nothing<2>()]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(3)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(3)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FE_Q<2>(4)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<2>[FESystem<2>[FE_Q<2>(1)^2]-FESystem<2>[FE_DGQ<2>(1)^2]]
+DEAL:: has_generalized_support_points=true
+DEAL:: 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:: has_generalized_support_points=true
+DEAL:: 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:: has_generalized_support_points=true
+DEAL::
+DEAL::dim=3
+DEAL:: FE_Nothing<3>()
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<3>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q<3>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<3>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Q_Hierarchical<3>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<3>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQArbitraryNodes<3>(QIterated(QTrapezoid(),4))
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQ<3>(4)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQArbitraryNodes<3>(QGauss(3))
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGQLegendre<3>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGQLegendre<3>(2)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGQHermite<3>(3)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGP<3>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGP<3>(2)
+DEAL:: has_generalized_support_points=false
+DEAL:: FESystem<3>[FE_Q<3>(2)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_Nothing<3>()^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<3>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceQ<3>(3)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_FaceP<3>(0)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_FaceP<3>(1)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_FaceP<3>(3)
+DEAL:: has_generalized_support_points=false
+DEAL:: FE_DGRaviartThomas<3>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGRaviartThomas<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGBDM<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGBDM<3>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGNedelec<3>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_DGNedelec<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomas<3>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomas<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomas<3>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_RaviartThomas<3>(1)-FE_DGQ<3>(1)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_BDM<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_BDM<3>(2)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Nedelec<3>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_Nedelec<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomasNodal<3>(0)
+DEAL:: has_generalized_support_points=true
+DEAL:: FE_RaviartThomasNodal<3>(1)
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_RaviartThomasNodal<3>(1)-FE_DGQ<3>(1)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_RaviartThomasNodal<3>(1)-FE_DGQ<3>(1)-FE_Nothing<3>()]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(3)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(3)]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FE_Q<3>(4)^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]^2]
+DEAL:: has_generalized_support_points=true
+DEAL:: FESystem<3>[FESystem<3>[FE_Q<3>(1)^2]-FESystem<3>[FE_DGQ<3>(1)^2]]
+DEAL:: has_generalized_support_points=true
+DEAL:: 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:: has_generalized_support_points=true
+DEAL:: 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:: has_generalized_support_points=true