Make FE_SimplexP_Bubbles inherit from FE_SimplexPoly.

This makes the code simpler and should make this element work automatically with
the current matrix-free simplex implementation.

diff --git a/include/deal.II/fe/fe_simplex_p_bubbles.h b/include/deal.II/fe/fe_simplex_p_bubbles.h
index 074411c44eb983f287fcc98611404e4052039c4d..556fa53e8a379b830d816dbe562e6e38e84ff5ae 100644 (file)
--- a/include/deal.II/fe/fe_simplex_p_bubbles.h
+++ b/include/deal.II/fe/fe_simplex_p_bubbles.h
@@ -20,7 +20,7 @@
 
 #include <deal.II/base/polynomials_barycentric.h>
 
-#include <deal.II/fe/fe_poly.h>
+#include <deal.II/fe/fe_simplex_p.h>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -68,7 +68,7 @@ DEAL_II_NAMESPACE_OPEN
  * are not yet implemented in this class.
  */
 template <int dim, int spacedim = dim>
-class FE_SimplexP_Bubbles : public dealii::FE_Poly<dim, spacedim>
+class FE_SimplexP_Bubbles : public FE_SimplexPoly<dim, spacedim>
 {
 public:
   /**
@@ -77,7 +77,8 @@ public:
    * space for this element: see the general documentation of this class for
    * more information.
    *
-   * @note For <code>degree == 1</code> this element is equivalent to FE_P(1).
+   * @note For <code>degree == 1</code> this element is equivalent to
+   * FE_SimplexP(1).
    */
   FE_SimplexP_Bubbles(const unsigned int degree);
 
@@ -96,14 +97,6 @@ public:
   virtual std::string
   get_name() const override;
 
-  /**
-   * @copydoc dealii::FiniteElement::convert_generalized_support_point_values_to_dof_values()
-   */
-  virtual void
-  convert_generalized_support_point_values_to_dof_values(
-    const std::vector<Vector<double>> &support_point_values,
-    std::vector<double> &              nodal_values) const override;
-
 protected:
   /**
    * Degree of the approximation (i.e., the constructor argument).

diff --git a/source/fe/fe_simplex_p_bubbles.cc b/source/fe/fe_simplex_p_bubbles.cc
index 044ce1643a99d1d8c219cec1aa9758a69c57e588..914aa2919ec5a5e3121c1fd51fd4f63406ced189 100644 (file)
--- a/source/fe/fe_simplex_p_bubbles.cc
+++ b/source/fe/fe_simplex_p_bubbles.cc
@@ -26,91 +26,6 @@
 
 DEAL_II_NAMESPACE_OPEN
 
-namespace
-{
-  /**
-   * Set up a vector that contains the unit (reference) cell support points
-   * for FE_Poly and sufficiently similar elements.
-   */
-  template <int dim>
-  std::vector<Point<dim>>
-  unit_support_points_fe_poly_bubbles(const unsigned int degree)
-  {
-    std::vector<Point<dim>> unit_points;
-
-    // Piecewise constants are a special case: use a support point at the
-    // centroid and only the centroid
-    if (degree == 0)
-      {
-        Point<dim> centroid;
-        std::fill(centroid.begin_raw(),
-                  centroid.end_raw(),
-                  1.0 / double(dim + 1));
-        unit_points.emplace_back(centroid);
-        return unit_points;
-      }
-
-    if (dim == 1)
-      {
-        // We don't really have dim = 1 support for simplex elements yet, but
-        // its convenient for populating the face array
-        Assert(degree <= 2, ExcNotImplemented());
-        if (degree >= 1)
-          {
-            unit_points.emplace_back(0.0);
-            unit_points.emplace_back(1.0);
-
-            if (degree == 2)
-              unit_points.emplace_back(0.5);
-          }
-      }
-    else if (dim == 2)
-      {
-        Assert(degree <= 2, ExcNotImplemented());
-        if (degree >= 1)
-          {
-            unit_points.emplace_back(0.0, 0.0);
-            unit_points.emplace_back(1.0, 0.0);
-            unit_points.emplace_back(0.0, 1.0);
-
-            if (degree == 2)
-              {
-                unit_points.emplace_back(0.5, 0.0);
-                unit_points.emplace_back(0.5, 0.5);
-                unit_points.emplace_back(0.0, 0.5);
-              }
-          }
-      }
-    else if (dim == 3)
-      {
-        Assert(degree <= 2, ExcNotImplemented());
-        if (degree >= 1)
-          {
-            unit_points.emplace_back(0.0, 0.0, 0.0);
-            unit_points.emplace_back(1.0, 0.0, 0.0);
-            unit_points.emplace_back(0.0, 1.0, 0.0);
-            unit_points.emplace_back(0.0, 0.0, 1.0);
-
-            if (degree == 2)
-              {
-                unit_points.emplace_back(0.5, 0.0, 0.0);
-                unit_points.emplace_back(0.5, 0.5, 0.0);
-                unit_points.emplace_back(0.0, 0.5, 0.0);
-                unit_points.emplace_back(0.0, 0.0, 0.5);
-                unit_points.emplace_back(0.5, 0.0, 0.5);
-                unit_points.emplace_back(0.0, 0.5, 0.5);
-              }
-          }
-      }
-    else
-      {
-        Assert(false, ExcNotImplemented());
-      }
-
-    return unit_points;
-  }
-} // namespace
-
 namespace FE_P_BubblesImplementation
 {
   template <int dim>
@@ -148,8 +63,9 @@ namespace FE_P_BubblesImplementation
   unit_support_points(const unsigned int degree)
   {
     Assert(degree < 3, ExcNotImplemented());
-    std::vector<Point<dim>> points =
-      unit_support_points_fe_poly_bubbles<dim>(degree);
+    // Start with the points used by FE_SimplexP, and then add bubbles.
+    FE_SimplexP<dim>        fe_p(degree);
+    std::vector<Point<dim>> points = fe_p.get_unit_support_points();
 
     Point<dim> centroid;
     std::fill(centroid.begin_raw(), centroid.end_raw(), 1.0 / double(dim + 1));
@@ -186,6 +102,17 @@ namespace FE_P_BubblesImplementation
 
 
 
+  template <>
+  std::vector<Point<0>>
+  unit_support_points<0>(const unsigned int /*degree*/)
+  {
+    std::vector<Point<0>> points;
+    points.emplace_back();
+    return points;
+  }
+
+
+
   template <int dim>
   BarycentricPolynomials<dim>
   get_basis(const unsigned int degree)
@@ -314,24 +241,15 @@ namespace FE_P_BubblesImplementation
 template <int dim, int spacedim>
 FE_SimplexP_Bubbles<dim, spacedim>::FE_SimplexP_Bubbles(
   const unsigned int degree)
-  : dealii::FE_Poly<dim, spacedim>(
+  : FE_SimplexPoly<dim, spacedim>(
       FE_P_BubblesImplementation::get_basis<dim>(degree),
       FE_P_BubblesImplementation::get_fe_data<dim>(degree),
-      std::vector<bool>(
-        FE_P_BubblesImplementation::get_fe_data<dim>(degree).dofs_per_cell,
-        true),
-      std::vector<ComponentMask>(
-        FE_P_BubblesImplementation::get_fe_data<dim>(degree).dofs_per_cell,
-        std::vector<bool>(1, true)))
+      FE_P_BubblesImplementation::unit_support_points<dim>(degree),
+      {FE_P_BubblesImplementation::unit_support_points<dim - 1>(degree)},
+      // Interface contraints are not yet implemented
+      FullMatrix<double>())
   , approximation_degree(degree)
-{
-  this->unit_support_points =
-    FE_P_BubblesImplementation::unit_support_points<dim>(degree);
-
-  // TODO
-  // this->unit_face_support_points =
-  //   unit_face_support_points_fe_poly<dim>(degree);
-}
+{}
 
 
 
@@ -345,28 +263,6 @@ FE_SimplexP_Bubbles<dim, spacedim>::get_name() const
 
 
 
-template <int dim, int spacedim>
-void
-FE_SimplexP_Bubbles<dim, spacedim>::
-  convert_generalized_support_point_values_to_dof_values(
-    const std::vector<Vector<double>> &support_point_values,
-    std::vector<double> &              nodal_values) const
-{
-  AssertDimension(support_point_values.size(),
-                  this->get_unit_support_points().size());
-  AssertDimension(support_point_values.size(), nodal_values.size());
-  AssertDimension(this->dofs_per_cell, nodal_values.size());
-
-  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
-    {
-      AssertDimension(support_point_values[i].size(), 1);
-
-      nodal_values[i] = support_point_values[i](0);
-    }
-}
-
-
-
 template <int dim, int spacedim>
 std::unique_ptr<FiniteElement<dim, spacedim>>
 FE_SimplexP_Bubbles<dim, spacedim>::clone() const