author = {S. Bartels and C. Carstensen and G. Dolzmann},
title = {Inhomogeneous Dirichlet conditions in a priori and a posteriori finite element error analysis},
journal = {Numerische Mathematik}
-}
\ No newline at end of file
+}
+
+@article{fried1975finite,
+ title={Finite element mass matrix lumping by numerical integration with no convergence rate loss},
+ author={Fried, Isaac and Malkus, David S},
+ journal={International Journal of Solids and Structures},
+ volume={11},
+ number={4},
+ pages={461--466},
+ year={1975},
+ publisher={Elsevier}
+}
+
+@article{Geevers_2018,
+ title={New Higher-Order Mass-Lumped Tetrahedral Elements for Wave Propagation Modelling},
+ volume={40},
+ ISSN={1095-7197},
+ url={http://dx.doi.org/10.1137/18M1175549},
+ DOI={10.1137/18m1175549},
+ number={5},
+ journal={SIAM Journal on Scientific Computing},
+ publisher={Society for Industrial & Applied Mathematics (SIAM)},
+ author={Geevers, S. and Mulder, W. A. and van der Vegt, J. J. W.},
+ year={2018},
+ month={Jan},
+ pages={A2830–A2857}
+}
get_name() const override;
};
+ /**
+ * @brief Enriched version of FE_P that can be used with nodal quadrature.
+ *
+ * Many explicit time integration schemes require solving a mass matrix at
+ * each time step. There are various ways around this requirement - for
+ * example, step-48 replaces the mass matrix with a diagonal approximation,
+ * which makes the solution step trivial. In step-48, and also commonly for
+ * tensor-product elements, this is done by computing the mass matrix with a
+ * lower-order quadrature point based on the nodes of the finite element
+ * (i.e., the nodal quadrature rule one obtains by using the shape functions
+ * as an interpolatory basis).
+ *
+ * A major drawback of standard simplex-based finite elements is that they
+ * cannot be used with nodal quadrature since some of the quadrature weights
+ * end up being either zero or negative, resulting in either an unsolvable or
+ * unstable approximation to the mass matrix. For example: the shape functions
+ * of FE_P<2>(2) with support points at vertices have mean values of zero so
+ * that element cannot be used with mass lumping.
+
+ * This element avoids this issue by replacing the shape functions of FE_P
+ * with an augmented space amendable to the construction of nodal quadrature
+ * rules. For example, on the triangle a single basis function is added
+ * corresponding to interpolation at the centroid (and all other basis
+ * functions are updated to preserve the partition of unity property). This
+ * results in shape functions with positive means (i.e., a valid nodal
+ * quadrature formula). Similarly, in 3D, the polynomial space of FE_P<3>(2)
+ * is enriched with five additional degrees of freedom (where four have
+ * support points at face centroids and one has a support point at the
+ * centroid) to enable construction of valid nodal quadrature rule.
+ *
+ * Since this FE space includes bubbles (i.e., extra functions which are
+ * nonzero only on element interiors), the polynomial degrees of the component
+ * basis functions are higher than the actual approximation degree of the
+ * element. For example, with a constructor argument <code>degree = 2</code>
+ * in 3D, the polynomials are in fact cubic (degree 3) but the order of the
+ * approximation is the same as if we were using quadratic (degree 2) finite
+ * elements.
+ *
+ * The 2D quadratic element was first described in @cite fried1975finite. The
+ * 3D quadratic element implemented here was first described in
+ * @cite Geevers_2018. Higher degree elements amendable to lumping exist but
+ * are not yet implemented in this class.
+ */
+ template <int dim, int spacedim = dim>
+ class FE_P_Bubbles : public dealii::FE_Poly<dim, spacedim>
+ {
+ public:
+ /**
+ * Constructor, taking the approximation degree as an argument. The
+ * polynomial space is typically one degree higher than the approximation
+ * 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).
+ */
+ FE_P_Bubbles(const unsigned int degree);
+
+ /**
+ * @copydoc dealii::FiniteElement::clone()
+ */
+ virtual std::unique_ptr<FiniteElement<dim, spacedim>>
+ clone() const override;
+ /**
+ * Return a string that uniquely identifies a finite element. This class
+ * returns <tt>Simplex::FE_F_Bubbles<dim,spacedim>(degree)</tt>, with
+ * @p dim, @p spacedim, and @p degree replaced by appropriate values. As
+ * usual, @p spacedim is omitted in the codimension zero case.
+ */
+ 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).
+ */
+ unsigned int approximation_degree;
+ };
} // namespace Simplex
DEAL_II_NAMESPACE_CLOSE
unit_support_points_fe_poly(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
std::vector<std::vector<Point<dim - 1>>>
unit_face_support_points_fe_poly(const unsigned int degree)
{
- Assert(dim == 2 || dim == 3, ExcNotImplemented());
+ // this concept doesn't exist in 1D so just return an empty vector
+ if (dim == 1)
+ return {};
+
const auto &info = ReferenceCell::internal::Info::get_cell(
dim == 2 ? ReferenceCell::Type::Tri : ReferenceCell::Type::Tet);
std::vector<std::vector<Point<dim - 1>>> unit_face_points;
}
+ namespace FE_P_BubblesImplementation
+ {
+ template <int dim>
+ std::vector<unsigned int>
+ get_dpo_vector(const unsigned int degree)
+ {
+ std::vector<unsigned int> dpo(dim + 1);
+ if (degree == 0)
+ {
+ dpo[dim] = 1; // single interior dof
+ }
+ else
+ {
+ Assert(degree == 1 || degree == 2, ExcNotImplemented());
+ dpo[0] = 1; // vertex dofs
+
+ if (degree == 2)
+ {
+ dpo[1] = 1; // line dofs
+
+ if (dim > 1)
+ dpo[dim] = 1; // the internal bubble function
+ if (dim == 3)
+ dpo[dim - 1] = 1; // face bubble functions
+ }
+ }
+
+ return dpo;
+ }
+
+
+
+ template <int dim>
+ std::vector<Point<dim>>
+ unit_support_points(const unsigned int degree)
+ {
+ Assert(degree < 3, ExcNotImplemented());
+ std::vector<Point<dim>> points = unit_support_points_fe_poly<dim>(degree);
+
+ Point<dim> centroid;
+ std::fill(centroid.begin_raw(),
+ centroid.end_raw(),
+ 1.0 / double(dim + 1));
+
+ switch (dim)
+ {
+ case 1:
+ // nothing more to do
+ return points;
+ case 2:
+ {
+ if (degree == 2)
+ points.push_back(centroid);
+ return points;
+ }
+ case 3:
+ {
+ if (degree == 2)
+ {
+ const double q13 = 1.0 / 3.0;
+ points.emplace_back(q13, q13, 0.0);
+ points.emplace_back(q13, 0.0, q13);
+ points.emplace_back(0.0, q13, q13);
+ points.emplace_back(q13, q13, q13);
+ points.push_back(centroid);
+ }
+ return points;
+ }
+ default:
+ Assert(false, ExcNotImplemented());
+ }
+ return points;
+ }
+
+
+
+ template <int dim>
+ BarycentricPolynomials<dim>
+ get_basis(const unsigned int degree)
+ {
+ Point<dim> centroid;
+ std::fill(centroid.begin_raw(),
+ centroid.end_raw(),
+ 1.0 / double(dim + 1));
+
+ auto M = [](const unsigned int d) {
+ return BarycentricPolynomial<dim, double>::monomial(d);
+ };
+
+ switch (degree)
+ {
+ // we don't need to add bubbles to P0 or P1
+ case 0:
+ case 1:
+ return BarycentricPolynomials<dim>::get_fe_p_basis(degree);
+ case 2:
+ {
+ const auto fe_p =
+ BarycentricPolynomials<dim>::get_fe_p_basis(degree);
+ // no further work is needed in 1D
+ if (dim == 1)
+ return fe_p;
+
+ // in 2D and 3D we add a centroid bubble function
+ auto c_bubble = BarycentricPolynomial<dim>() + 1;
+ for (unsigned int d = 0; d < dim + 1; ++d)
+ c_bubble = c_bubble * M(d);
+ c_bubble = c_bubble / c_bubble.value(centroid);
+
+ std::vector<BarycentricPolynomial<dim>> bubble_functions;
+ if (dim == 2)
+ {
+ bubble_functions.push_back(c_bubble);
+ }
+ else if (dim == 3)
+ {
+ // need 'face bubble' functions in addition to the centroid.
+ // Furthermore we need to subtract them off from the other
+ // functions so that we end up with an interpolatory basis
+ auto b0 = 27 * M(0) * M(1) * M(2);
+ bubble_functions.push_back(b0 -
+ b0.value(centroid) * c_bubble);
+ auto b1 = 27 * M(0) * M(1) * M(3);
+ bubble_functions.push_back(b1 -
+ b1.value(centroid) * c_bubble);
+ auto b2 = 27 * M(0) * M(2) * M(3);
+ bubble_functions.push_back(b2 -
+ b2.value(centroid) * c_bubble);
+ auto b3 = 27 * M(1) * M(2) * M(3);
+ bubble_functions.push_back(b3 -
+ b3.value(centroid) * c_bubble);
+
+ bubble_functions.push_back(c_bubble);
+ }
+
+ // Extract out the support points for the extra bubble (both
+ // volume and face) functions:
+ const std::vector<Point<dim>> support_points =
+ unit_support_points<dim>(degree);
+ const std::vector<Point<dim>> bubble_support_points(
+ support_points.begin() + fe_p.n(), support_points.end());
+ Assert(bubble_support_points.size() == bubble_functions.size(),
+ ExcInternalError());
+ const unsigned int n_bubbles = bubble_support_points.size();
+
+ // Assemble the final basis:
+ std::vector<BarycentricPolynomial<dim>> lump_polys;
+ for (unsigned int i = 0; i < fe_p.n(); ++i)
+ {
+ BarycentricPolynomial<dim> p = fe_p[i];
+
+ for (unsigned int j = 0; j < n_bubbles; ++j)
+ {
+ p = p - p.value(bubble_support_points[j]) *
+ bubble_functions[j];
+ }
+
+ lump_polys.push_back(p);
+ }
+
+ for (auto &p : bubble_functions)
+ lump_polys.push_back(std::move(p));
+
+ // Sanity check:
+#ifdef DEBUG
+ BarycentricPolynomial<dim> unity;
+ for (const auto &p : lump_polys)
+ unity = unity + p;
+
+ Point<dim> test;
+ for (unsigned int d = 0; d < dim; ++d)
+ test[d] = 2.0;
+ Assert(std::abs(unity.value(test) - 1.0) < 1e-10,
+ ExcInternalError());
+#endif
+
+ return BarycentricPolynomials<dim>(lump_polys);
+ }
+ default:
+ Assert(degree < 3, ExcNotImplemented());
+ }
+
+ Assert(degree < 3, ExcNotImplemented());
+ // bogus return to placate compilers
+ return BarycentricPolynomials<dim>::get_fe_p_basis(degree);
+ }
+
+
+
+ template <int dim>
+ FiniteElementData<dim>
+ get_fe_data(const unsigned int degree)
+ {
+ // It's not efficient, but delegate computation of the degree of the
+ // finite element (which is different from the input argument) to the
+ // basis.
+ const auto polys = get_basis<dim>(degree);
+ return FiniteElementData<dim>(get_dpo_vector<dim>(degree),
+ ReferenceCell::Type::get_simplex<dim>(),
+ 1, // n_components
+ polys.degree(),
+ FiniteElementData<dim>::H1);
+ }
+ } // namespace FE_P_BubblesImplementation
+
+
+
+ template <int dim, int spacedim>
+ FE_P_Bubbles<dim, spacedim>::FE_P_Bubbles(const unsigned int degree)
+ : dealii::FE_Poly<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)))
+ , 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);
+ }
+
+
+
+ template <int dim, int spacedim>
+ std::string
+ FE_P_Bubbles<dim, spacedim>::get_name() const
+ {
+ return "Simplex::FE_P_Bubbles<" + Utilities::dim_string(dim, spacedim) +
+ ">" + "(" + std::to_string(approximation_degree) + ")";
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ FE_P_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_P_Bubbles<dim, spacedim>::clone() const
+ {
+ return std::make_unique<FE_P_Bubbles<dim, spacedim>>(*this);
+ }
} // namespace Simplex
// explicit instantiations
deal_II_space_dimension>;
template class Simplex::FE_PyramidDGP<deal_II_dimension,
deal_II_space_dimension>;
+ template class Simplex::FE_P_Bubbles<deal_II_dimension,
+ deal_II_space_dimension>;
#endif
}
} // namespace
+
template <int dim>
double
ScalarWedgePolynomial<dim>::compute_value(const unsigned int i,
template class ScalarPyramidPolynomial<1>;
template class ScalarPyramidPolynomial<2>;
template class ScalarPyramidPolynomial<3>;
-
} // namespace Simplex
DEAL_II_NAMESPACE_CLOSE
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Evaluate Simplex::FE_P_Bubbles.
+
+
+#include <deal.II/simplex/fe_lib.h>
+#include <deal.II/simplex/quadrature_lib.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim, int spacedim>
+void
+test(const FiniteElement<dim, spacedim> &fe, const Quadrature<dim> &quad)
+{
+ deallog.push(fe.get_name());
+
+ for (const auto &point : quad.get_points())
+ {
+ deallog << point << " : " << std::endl;
+ for (unsigned int i = 0; i < fe.n_dofs_per_cell(); i++)
+ deallog << fe.shape_value(i, point) << " ";
+ for (unsigned int i = 0; i < fe.n_dofs_per_cell(); i++)
+ deallog << fe.shape_grad(i, point) << " ";
+ for (unsigned int i = 0; i < fe.n_dofs_per_cell(); i++)
+ deallog << fe.shape_grad_grad(i, point) << " ";
+ deallog << std::endl;
+ }
+ deallog << std::endl;
+
+ deallog.pop();
+}
+
+
+
+template <int dim, int spacedim = dim>
+void
+test_unit_support_points()
+{
+ deallog << "Test support points for dim = " << dim
+ << " and spacedim = " << spacedim << std::endl;
+ for (unsigned int degree = 0; degree < 3; ++degree)
+ {
+ deallog << "approximation degree = " << degree << std::endl;
+ Simplex::FE_P_Bubbles<dim, spacedim> fe(degree);
+ deallog << "element tensor degree = " << fe.tensor_degree() << std::endl;
+ Quadrature<dim> quad(fe.get_unit_support_points());
+ test(fe, quad);
+ }
+}
+
+
+
+int
+main()
+{
+ initlog();
+
+ test_unit_support_points<1>();
+ test_unit_support_points<1, 2>();
+ test_unit_support_points<1, 3>();
+
+ test_unit_support_points<2>();
+ test_unit_support_points<2, 3>();
+
+ test_unit_support_points<3>();
+}
--- /dev/null
+
+DEAL::Test support points for dim = 1 and spacedim = 1
+DEAL::approximation degree = 0
+DEAL::element tensor degree = 0
+DEAL:Simplex::FE_P_Bubbles<1>(0)::0.500000 :
+DEAL:Simplex::FE_P_Bubbles<1>(0)::1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1>(0)::
+DEAL::approximation degree = 1
+DEAL::element tensor degree = 1
+DEAL:Simplex::FE_P_Bubbles<1>(1)::0.00000 :
+DEAL:Simplex::FE_P_Bubbles<1>(1)::1.00000 0.00000 -1.00000 1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1>(1)::1.00000 :
+DEAL:Simplex::FE_P_Bubbles<1>(1)::0.00000 1.00000 -1.00000 1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1>(1)::
+DEAL::approximation degree = 2
+DEAL::element tensor degree = 2
+DEAL:Simplex::FE_P_Bubbles<1>(2)::0.00000 :
+DEAL:Simplex::FE_P_Bubbles<1>(2)::1.00000 0.00000 0.00000 -3.00000 -1.00000 4.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1>(2)::1.00000 :
+DEAL:Simplex::FE_P_Bubbles<1>(2)::0.00000 1.00000 0.00000 1.00000 3.00000 -4.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1>(2)::0.500000 :
+DEAL:Simplex::FE_P_Bubbles<1>(2)::0.00000 0.00000 1.00000 -1.00000 1.00000 0.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1>(2)::
+DEAL::Test support points for dim = 1 and spacedim = 2
+DEAL::approximation degree = 0
+DEAL::element tensor degree = 0
+DEAL:Simplex::FE_P_Bubbles<1,2>(0)::0.500000 :
+DEAL:Simplex::FE_P_Bubbles<1,2>(0)::1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1,2>(0)::
+DEAL::approximation degree = 1
+DEAL::element tensor degree = 1
+DEAL:Simplex::FE_P_Bubbles<1,2>(1)::0.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,2>(1)::1.00000 0.00000 -1.00000 1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1,2>(1)::1.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,2>(1)::0.00000 1.00000 -1.00000 1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1,2>(1)::
+DEAL::approximation degree = 2
+DEAL::element tensor degree = 2
+DEAL:Simplex::FE_P_Bubbles<1,2>(2)::0.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,2>(2)::1.00000 0.00000 0.00000 -3.00000 -1.00000 4.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1,2>(2)::1.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,2>(2)::0.00000 1.00000 0.00000 1.00000 3.00000 -4.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1,2>(2)::0.500000 :
+DEAL:Simplex::FE_P_Bubbles<1,2>(2)::0.00000 0.00000 1.00000 -1.00000 1.00000 0.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1,2>(2)::
+DEAL::Test support points for dim = 1 and spacedim = 3
+DEAL::approximation degree = 0
+DEAL::element tensor degree = 0
+DEAL:Simplex::FE_P_Bubbles<1,3>(0)::0.500000 :
+DEAL:Simplex::FE_P_Bubbles<1,3>(0)::1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1,3>(0)::
+DEAL::approximation degree = 1
+DEAL::element tensor degree = 1
+DEAL:Simplex::FE_P_Bubbles<1,3>(1)::0.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,3>(1)::1.00000 0.00000 -1.00000 1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1,3>(1)::1.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,3>(1)::0.00000 1.00000 -1.00000 1.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<1,3>(1)::
+DEAL::approximation degree = 2
+DEAL::element tensor degree = 2
+DEAL:Simplex::FE_P_Bubbles<1,3>(2)::0.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,3>(2)::1.00000 0.00000 0.00000 -3.00000 -1.00000 4.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1,3>(2)::1.00000 :
+DEAL:Simplex::FE_P_Bubbles<1,3>(2)::0.00000 1.00000 0.00000 1.00000 3.00000 -4.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1,3>(2)::0.500000 :
+DEAL:Simplex::FE_P_Bubbles<1,3>(2)::0.00000 0.00000 1.00000 -1.00000 1.00000 0.00000 4.00000 4.00000 -8.00000
+DEAL:Simplex::FE_P_Bubbles<1,3>(2)::
+DEAL::Test support points for dim = 2 and spacedim = 2
+DEAL::approximation degree = 0
+DEAL::element tensor degree = 0
+DEAL:Simplex::FE_P_Bubbles<2>(0)::0.333333 0.333333 :
+DEAL:Simplex::FE_P_Bubbles<2>(0)::1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2>(0)::
+DEAL::approximation degree = 1
+DEAL::element tensor degree = 1
+DEAL:Simplex::FE_P_Bubbles<2>(1)::0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2>(1)::1.00000 0.00000 0.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2>(1)::1.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2>(1)::0.00000 1.00000 0.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2>(1)::0.00000 1.00000 :
+DEAL:Simplex::FE_P_Bubbles<2>(1)::0.00000 0.00000 1.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2>(1)::
+DEAL::approximation degree = 2
+DEAL::element tensor degree = 3
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2>(2)::1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -3.00000 -3.00000 -1.00000 0.00000 0.00000 -1.00000 4.00000 0.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 4.00000 7.00000 7.00000 4.00000 4.00000 3.00000 3.00000 0.00000 0.00000 3.00000 3.00000 4.00000 -8.00000 -16.0000 -16.0000 0.00000 0.00000 -8.00000 -8.00000 0.00000 0.00000 -16.0000 -16.0000 -8.00000 0.00000 27.0000 27.0000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2>(2)::1.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 3.00000 0.00000 0.00000 -1.00000 -4.00000 -4.00000 0.00000 4.00000 0.00000 0.00000 0.00000 0.00000 4.00000 1.00000 1.00000 -2.00000 4.00000 -3.00000 -3.00000 -6.00000 0.00000 -3.00000 -3.00000 -2.00000 -8.00000 8.00000 8.00000 24.0000 0.00000 16.0000 16.0000 24.0000 0.00000 8.00000 8.00000 16.0000 0.00000 -27.0000 -27.0000 -54.0000
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 1.00000 :
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 -1.00000 0.00000 0.00000 3.00000 0.00000 0.00000 4.00000 0.00000 -4.00000 -4.00000 0.00000 0.00000 -2.00000 1.00000 1.00000 4.00000 -2.00000 -3.00000 -3.00000 0.00000 -6.00000 -3.00000 -3.00000 4.00000 16.0000 8.00000 8.00000 0.00000 24.0000 16.0000 16.0000 0.00000 24.0000 8.00000 8.00000 -8.00000 -54.0000 -27.0000 -27.0000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.500000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 -1.00000 -0.250000 1.00000 0.750000 0.00000 -0.250000 0.00000 -5.00000 0.00000 -1.00000 0.00000 -1.00000 0.00000 6.75000 4.00000 4.00000 4.00000 1.00000 4.00000 0.00000 0.00000 -3.00000 0.00000 0.00000 0.00000 1.00000 -8.00000 -4.00000 -4.00000 12.0000 0.00000 4.00000 4.00000 12.0000 0.00000 -4.00000 -4.00000 4.00000 0.00000 0.00000 0.00000 -27.0000
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.500000 0.500000 :
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.250000 0.250000 0.250000 -0.750000 -0.750000 0.250000 1.00000 1.00000 5.00000 5.00000 1.00000 1.00000 -6.75000 -6.75000 1.00000 1.00000 1.00000 1.00000 1.00000 -3.00000 -3.00000 -3.00000 -3.00000 -3.00000 -3.00000 1.00000 4.00000 8.00000 8.00000 12.0000 12.0000 16.0000 16.0000 12.0000 12.0000 8.00000 8.00000 4.00000 -27.0000 -27.0000 -27.0000 -27.0000
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 0.500000 :
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 -0.250000 -1.00000 -0.250000 0.00000 0.750000 1.00000 -1.00000 0.00000 -1.00000 0.00000 -5.00000 0.00000 6.75000 0.00000 1.00000 4.00000 4.00000 4.00000 1.00000 0.00000 0.00000 0.00000 -3.00000 0.00000 0.00000 4.00000 4.00000 -4.00000 -4.00000 0.00000 12.0000 4.00000 4.00000 0.00000 12.0000 -4.00000 -4.00000 -8.00000 -27.0000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.333333 0.333333 :
+DEAL:Simplex::FE_P_Bubbles<2>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 -0.333333 -0.333333 0.333333 1.11022e-16 1.11022e-16 0.333333 2.22045e-16 -1.33333 1.33333 1.33333 -1.33333 2.22045e-16 4.44089e-16 4.44089e-16 2.00000 3.00000 3.00000 2.00000 2.00000 -1.00000 -1.00000 -2.00000 -2.00000 -1.00000 -1.00000 2.00000 0.00000 -8.88178e-16 -8.88178e-16 8.00000 8.00000 8.00000 8.00000 8.00000 8.00000 -8.88178e-16 -8.88178e-16 0.00000 -18.0000 -9.00000 -9.00000 -18.0000
+DEAL:Simplex::FE_P_Bubbles<2>(2)::
+DEAL::Test support points for dim = 2 and spacedim = 3
+DEAL::approximation degree = 0
+DEAL::element tensor degree = 0
+DEAL:Simplex::FE_P_Bubbles<2,3>(0)::0.333333 0.333333 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(0)::1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2,3>(0)::
+DEAL::approximation degree = 1
+DEAL::element tensor degree = 1
+DEAL:Simplex::FE_P_Bubbles<2,3>(1)::0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(1)::1.00000 0.00000 0.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2,3>(1)::1.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(1)::0.00000 1.00000 0.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2,3>(1)::0.00000 1.00000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(1)::0.00000 0.00000 1.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2,3>(1)::
+DEAL::approximation degree = 2
+DEAL::element tensor degree = 3
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -3.00000 -3.00000 -1.00000 0.00000 0.00000 -1.00000 4.00000 0.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 4.00000 7.00000 7.00000 4.00000 4.00000 3.00000 3.00000 0.00000 0.00000 3.00000 3.00000 4.00000 -8.00000 -16.0000 -16.0000 0.00000 0.00000 -8.00000 -8.00000 0.00000 0.00000 -16.0000 -16.0000 -8.00000 0.00000 27.0000 27.0000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::1.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 3.00000 0.00000 0.00000 -1.00000 -4.00000 -4.00000 0.00000 4.00000 0.00000 0.00000 0.00000 0.00000 4.00000 1.00000 1.00000 -2.00000 4.00000 -3.00000 -3.00000 -6.00000 0.00000 -3.00000 -3.00000 -2.00000 -8.00000 8.00000 8.00000 24.0000 0.00000 16.0000 16.0000 24.0000 0.00000 8.00000 8.00000 16.0000 0.00000 -27.0000 -27.0000 -54.0000
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 1.00000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 -1.00000 0.00000 0.00000 3.00000 0.00000 0.00000 4.00000 0.00000 -4.00000 -4.00000 0.00000 0.00000 -2.00000 1.00000 1.00000 4.00000 -2.00000 -3.00000 -3.00000 0.00000 -6.00000 -3.00000 -3.00000 4.00000 16.0000 8.00000 8.00000 0.00000 24.0000 16.0000 16.0000 0.00000 24.0000 8.00000 8.00000 -8.00000 -54.0000 -27.0000 -27.0000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.500000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 -1.00000 -0.250000 1.00000 0.750000 0.00000 -0.250000 0.00000 -5.00000 0.00000 -1.00000 0.00000 -1.00000 0.00000 6.75000 4.00000 4.00000 4.00000 1.00000 4.00000 0.00000 0.00000 -3.00000 0.00000 0.00000 0.00000 1.00000 -8.00000 -4.00000 -4.00000 12.0000 0.00000 4.00000 4.00000 12.0000 0.00000 -4.00000 -4.00000 4.00000 0.00000 0.00000 0.00000 -27.0000
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.500000 0.500000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.250000 0.250000 0.250000 -0.750000 -0.750000 0.250000 1.00000 1.00000 5.00000 5.00000 1.00000 1.00000 -6.75000 -6.75000 1.00000 1.00000 1.00000 1.00000 1.00000 -3.00000 -3.00000 -3.00000 -3.00000 -3.00000 -3.00000 1.00000 4.00000 8.00000 8.00000 12.0000 12.0000 16.0000 16.0000 12.0000 12.0000 8.00000 8.00000 4.00000 -27.0000 -27.0000 -27.0000 -27.0000
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 0.500000 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 -0.250000 -1.00000 -0.250000 0.00000 0.750000 1.00000 -1.00000 0.00000 -1.00000 0.00000 -5.00000 0.00000 6.75000 0.00000 1.00000 4.00000 4.00000 4.00000 1.00000 0.00000 0.00000 0.00000 -3.00000 0.00000 0.00000 4.00000 4.00000 -4.00000 -4.00000 0.00000 12.0000 4.00000 4.00000 0.00000 12.0000 -4.00000 -4.00000 -8.00000 -27.0000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.333333 0.333333 :
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 -0.333333 -0.333333 0.333333 1.11022e-16 1.11022e-16 0.333333 2.22045e-16 -1.33333 1.33333 1.33333 -1.33333 2.22045e-16 4.44089e-16 4.44089e-16 2.00000 3.00000 3.00000 2.00000 2.00000 -1.00000 -1.00000 -2.00000 -2.00000 -1.00000 -1.00000 2.00000 0.00000 -8.88178e-16 -8.88178e-16 8.00000 8.00000 8.00000 8.00000 8.00000 8.00000 -8.88178e-16 -8.88178e-16 0.00000 -18.0000 -9.00000 -9.00000 -18.0000
+DEAL:Simplex::FE_P_Bubbles<2,3>(2)::
+DEAL::Test support points for dim = 3 and spacedim = 3
+DEAL::approximation degree = 0
+DEAL::element tensor degree = 0
+DEAL:Simplex::FE_P_Bubbles<3>(0)::0.250000 0.250000 0.250000 :
+DEAL:Simplex::FE_P_Bubbles<3>(0)::1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(0)::
+DEAL::approximation degree = 1
+DEAL::element tensor degree = 1
+DEAL:Simplex::FE_P_Bubbles<3>(1)::0.00000 0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(1)::1.00000 0.00000 0.00000 0.00000 -1.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(1)::1.00000 0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(1)::0.00000 1.00000 0.00000 0.00000 -1.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(1)::0.00000 1.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(1)::0.00000 0.00000 1.00000 0.00000 -1.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(1)::0.00000 0.00000 1.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(1)::0.00000 0.00000 0.00000 1.00000 -1.00000 -1.00000 -1.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(1)::
+DEAL::approximation degree = 2
+DEAL::element tensor degree = 3
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -3.00000 -3.00000 -3.00000 -1.00000 0.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 -1.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 4.00000 7.00000 7.00000 7.00000 4.00000 7.00000 7.00000 7.00000 4.00000 4.00000 3.00000 3.00000 3.00000 0.00000 0.00000 3.00000 0.00000 0.00000 0.00000 3.00000 0.00000 3.00000 4.00000 3.00000 0.00000 3.00000 0.00000 0.00000 0.00000 3.00000 0.00000 0.00000 3.00000 3.00000 3.00000 4.00000 -8.00000 -16.0000 -16.0000 -16.0000 0.00000 0.00000 -16.0000 0.00000 0.00000 0.00000 -8.00000 0.00000 -8.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -16.0000 0.00000 -16.0000 -8.00000 -16.0000 0.00000 -16.0000 0.00000 0.00000 0.00000 -16.0000 0.00000 0.00000 -16.0000 -16.0000 -16.0000 -8.00000 0.00000 0.00000 -8.00000 0.00000 0.00000 0.00000 -8.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -8.00000 0.00000 -8.00000 0.00000 0.00000 27.0000 0.00000 27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 27.0000 0.00000 0.00000 0.00000 27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 27.0000 0.00000 27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::1.00000 0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 3.00000 0.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 -1.00000 -4.00000 -4.00000 -4.00000 0.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 4.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 1.00000 -2.00000 -2.00000 4.00000 -3.00000 -3.00000 -3.00000 -6.00000 -3.00000 -3.00000 -3.00000 -6.00000 0.00000 -3.00000 0.00000 -3.00000 -2.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -3.00000 0.00000 0.00000 0.00000 -3.00000 0.00000 -2.00000 -8.00000 8.00000 8.00000 8.00000 24.0000 24.0000 8.00000 24.0000 24.0000 0.00000 16.0000 0.00000 16.0000 24.0000 0.00000 0.00000 0.00000 0.00000 0.00000 8.00000 0.00000 8.00000 16.0000 8.00000 0.00000 8.00000 0.00000 0.00000 0.00000 8.00000 0.00000 0.00000 8.00000 8.00000 8.00000 16.0000 0.00000 0.00000 16.0000 0.00000 0.00000 0.00000 16.0000 0.00000 24.0000 0.00000 0.00000 0.00000 0.00000 0.00000 -8.00000 0.00000 -8.00000 0.00000 0.00000 -27.0000 0.00000 -27.0000 -54.0000 -27.0000 0.00000 -27.0000 0.00000 0.00000 0.00000 -27.0000 0.00000 0.00000 -27.0000 -27.0000 -27.0000 -54.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 27.0000 0.00000 27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 1.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 -1.00000 0.00000 0.00000 0.00000 3.00000 0.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 -4.00000 -4.00000 -4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -2.00000 1.00000 -2.00000 1.00000 4.00000 1.00000 -2.00000 1.00000 -2.00000 -2.00000 -3.00000 0.00000 -3.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -6.00000 -3.00000 -3.00000 -3.00000 4.00000 -3.00000 -3.00000 -3.00000 -6.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -3.00000 0.00000 -3.00000 -2.00000 16.0000 8.00000 8.00000 8.00000 0.00000 0.00000 8.00000 0.00000 0.00000 24.0000 16.0000 0.00000 16.0000 0.00000 0.00000 0.00000 0.00000 0.00000 24.0000 8.00000 24.0000 8.00000 -8.00000 8.00000 24.0000 8.00000 24.0000 0.00000 0.00000 8.00000 0.00000 0.00000 8.00000 8.00000 8.00000 16.0000 0.00000 0.00000 -8.00000 0.00000 0.00000 0.00000 -8.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 16.0000 0.00000 16.0000 24.0000 -54.0000 -27.0000 -27.0000 -27.0000 0.00000 0.00000 -27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -27.0000 0.00000 0.00000 -27.0000 -27.0000 -27.0000 -54.0000 0.00000 0.00000 27.0000 0.00000 0.00000 0.00000 27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 1.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 -1.00000 0.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 3.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -4.00000 -4.00000 -4.00000 4.00000 0.00000 0.00000 0.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -2.00000 -2.00000 1.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 4.00000 -2.00000 0.00000 -3.00000 0.00000 0.00000 0.00000 -3.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -2.00000 -3.00000 0.00000 -3.00000 0.00000 -6.00000 -3.00000 -3.00000 -3.00000 -6.00000 -3.00000 -3.00000 -3.00000 4.00000 16.0000 8.00000 8.00000 8.00000 0.00000 0.00000 8.00000 0.00000 0.00000 0.00000 -8.00000 0.00000 -8.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 8.00000 0.00000 8.00000 16.0000 8.00000 0.00000 8.00000 0.00000 24.0000 24.0000 8.00000 24.0000 24.0000 8.00000 8.00000 8.00000 -8.00000 24.0000 0.00000 16.0000 0.00000 0.00000 0.00000 16.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 24.0000 16.0000 0.00000 16.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -54.0000 -27.0000 -27.0000 -27.0000 0.00000 0.00000 -27.0000 0.00000 0.00000 0.00000 -27.0000 0.00000 -27.0000 -54.0000 -27.0000 0.00000 -27.0000 0.00000 0.00000 27.0000 0.00000 27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.500000 0.00000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -1.00000 -0.250000 -0.250000 1.00000 0.750000 0.750000 0.00000 -0.250000 0.00000 0.00000 0.00000 -0.250000 0.00000 -5.00000 -5.00000 0.00000 -1.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 0.00000 6.75000 0.00000 0.00000 0.00000 6.75000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 4.00000 4.00000 4.00000 4.00000 1.00000 1.50000 4.00000 1.50000 1.00000 4.00000 0.00000 0.00000 0.00000 -3.00000 -2.50000 0.00000 -2.50000 -3.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.500000 0.00000 0.500000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.500000 0.00000 0.500000 1.00000 -8.00000 -4.00000 -4.00000 -4.00000 12.0000 20.0000 -4.00000 20.0000 12.0000 0.00000 4.00000 0.00000 4.00000 12.0000 8.00000 0.00000 8.00000 0.00000 0.00000 -4.00000 0.00000 -4.00000 4.00000 4.00000 0.00000 4.00000 0.00000 0.00000 0.00000 -4.00000 0.00000 0.00000 4.00000 -4.00000 4.00000 4.00000 0.00000 0.00000 4.00000 0.00000 0.00000 8.00000 4.00000 8.00000 12.0000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -27.0000 -40.5000 0.00000 -40.5000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -40.5000 0.00000 -40.5000 -27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 -13.5000 0.00000 -13.5000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -13.5000 0.00000 -13.5000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 64.0000 0.00000 64.0000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.500000 0.500000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.250000 0.250000 0.250000 0.250000 -0.750000 0.00000 -0.750000 0.250000 0.00000 0.00000 0.00000 -0.250000 1.00000 1.00000 1.00000 5.00000 5.00000 0.00000 1.00000 1.00000 1.00000 0.00000 0.00000 -1.66533e-15 0.00000 0.00000 -1.00000 0.00000 0.00000 -1.00000 -6.75000 -6.75000 -6.75000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 6.75000 0.00000 0.00000 0.00000 1.00000 1.00000 0.500000 1.00000 1.00000 0.500000 0.500000 0.500000 8.88178e-15 1.00000 -3.00000 -0.500000 -3.00000 -3.00000 -0.500000 -0.500000 -0.500000 -1.00000 -3.00000 -3.00000 -0.500000 -3.00000 1.00000 -0.500000 -0.500000 -0.500000 -1.00000 0.00000 0.00000 -0.500000 0.00000 0.00000 -0.500000 -0.500000 -0.500000 0.00000 4.00000 8.00000 -1.24345e-14 8.00000 12.0000 4.00000 -1.24345e-14 4.00000 -4.00000 12.0000 16.0000 -8.00000 16.0000 12.0000 -8.00000 -8.00000 -8.00000 -16.0000 12.0000 8.00000 4.00000 8.00000 4.00000 -1.24345e-14 4.00000 -1.24345e-14 -4.00000 0.00000 0.00000 -1.24345e-14 0.00000 0.00000 -1.24345e-14 -1.24345e-14 -1.24345e-14 -1.77636e-14 0.00000 0.00000 -4.00000 0.00000 0.00000 -8.00000 -4.00000 -8.00000 -4.00000 0.00000 0.00000 -8.00000 0.00000 0.00000 -4.00000 -8.00000 -4.00000 -4.00000 -27.0000 -27.0000 13.5000 -27.0000 -27.0000 13.5000 13.5000 13.5000 54.0000 0.00000 0.00000 13.5000 0.00000 0.00000 13.5000 13.5000 13.5000 27.0000 0.00000 0.00000 13.5000 0.00000 0.00000 13.5000 13.5000 13.5000 27.0000 0.00000 0.00000 40.5000 0.00000 0.00000 40.5000 40.5000 40.5000 54.0000 0.00000 0.00000 -64.0000 0.00000 0.00000 -64.0000 -64.0000 -64.0000 -128.000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.500000 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -0.250000 -1.00000 -0.250000 -0.250000 0.00000 0.00000 0.750000 1.00000 0.750000 0.00000 0.00000 -0.250000 -1.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 -5.00000 0.00000 -5.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -1.00000 6.75000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 6.75000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 4.00000 1.50000 4.00000 4.00000 4.00000 1.50000 4.00000 1.00000 1.00000 0.00000 0.500000 0.00000 0.00000 0.00000 0.500000 0.00000 0.00000 -3.00000 0.00000 -2.50000 0.00000 4.00000 0.00000 -2.50000 0.00000 -3.00000 0.00000 0.00000 0.500000 0.00000 0.00000 0.00000 0.500000 0.00000 1.00000 4.00000 -4.00000 4.00000 -4.00000 0.00000 0.00000 4.00000 0.00000 0.00000 12.0000 4.00000 8.00000 4.00000 0.00000 0.00000 8.00000 0.00000 0.00000 12.0000 -4.00000 20.0000 -4.00000 -8.00000 -4.00000 20.0000 -4.00000 12.0000 0.00000 0.00000 4.00000 0.00000 0.00000 -4.00000 4.00000 -4.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 8.00000 0.00000 0.00000 4.00000 8.00000 4.00000 12.0000 -27.0000 0.00000 -40.5000 0.00000 0.00000 0.00000 -40.5000 0.00000 0.00000 0.00000 0.00000 -13.5000 0.00000 0.00000 0.00000 -13.5000 0.00000 0.00000 0.00000 0.00000 -40.5000 0.00000 0.00000 0.00000 -40.5000 0.00000 -27.0000 0.00000 0.00000 -13.5000 0.00000 0.00000 0.00000 -13.5000 0.00000 0.00000 0.00000 0.00000 64.0000 0.00000 0.00000 0.00000 64.0000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.500000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -0.250000 -0.250000 -1.00000 -0.250000 0.00000 0.00000 0.00000 -0.250000 0.00000 0.750000 0.750000 1.00000 -1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -1.00000 0.00000 -5.00000 -5.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 0.00000 6.75000 0.00000 0.00000 0.00000 6.75000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.50000 4.00000 1.50000 1.00000 4.00000 4.00000 4.00000 4.00000 1.00000 0.500000 0.00000 0.500000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.500000 0.00000 0.500000 1.00000 0.00000 0.00000 0.00000 0.00000 -3.00000 -2.50000 0.00000 -2.50000 -3.00000 0.00000 0.00000 0.00000 4.00000 4.00000 4.00000 -4.00000 4.00000 0.00000 0.00000 -4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 4.00000 0.00000 4.00000 4.00000 -4.00000 0.00000 -4.00000 0.00000 12.0000 20.0000 -4.00000 20.0000 12.0000 -4.00000 -4.00000 -4.00000 -8.00000 12.0000 8.00000 4.00000 8.00000 0.00000 0.00000 4.00000 0.00000 0.00000 0.00000 8.00000 0.00000 8.00000 12.0000 4.00000 0.00000 4.00000 0.00000 0.00000 -13.5000 0.00000 -13.5000 0.00000 0.00000 0.00000 0.00000 0.00000 -27.0000 -40.5000 0.00000 -40.5000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -40.5000 0.00000 -40.5000 -27.0000 0.00000 0.00000 0.00000 0.00000 0.00000 -13.5000 0.00000 -13.5000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 64.0000 0.00000 64.0000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.500000 0.00000 0.500000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.250000 0.250000 0.250000 0.250000 0.00000 -0.750000 0.00000 -0.250000 0.00000 -0.750000 0.00000 0.250000 1.00000 1.00000 1.00000 0.00000 -1.00000 0.00000 0.00000 -1.66533e-15 0.00000 1.00000 1.00000 1.00000 5.00000 0.00000 5.00000 0.00000 -1.00000 0.00000 0.00000 0.00000 0.00000 -6.75000 -6.75000 -6.75000 0.00000 0.00000 0.00000 0.00000 6.75000 0.00000 0.00000 0.00000 0.00000 1.00000 0.500000 1.00000 0.500000 8.88178e-15 0.500000 1.00000 0.500000 1.00000 1.00000 -0.500000 -3.00000 -0.500000 -1.00000 -0.500000 -3.00000 -0.500000 -3.00000 0.00000 -0.500000 0.00000 -0.500000 0.00000 -0.500000 0.00000 -0.500000 0.00000 -3.00000 -0.500000 -3.00000 -0.500000 -1.00000 -0.500000 -3.00000 -0.500000 1.00000 4.00000 -1.24345e-14 8.00000 -1.24345e-14 -4.00000 4.00000 8.00000 4.00000 12.0000 0.00000 -4.00000 0.00000 -4.00000 -4.00000 -8.00000 0.00000 -8.00000 0.00000 0.00000 -1.24345e-14 0.00000 -1.24345e-14 -1.77636e-14 -1.24345e-14 0.00000 -1.24345e-14 0.00000 12.0000 4.00000 8.00000 4.00000 -4.00000 -1.24345e-14 8.00000 -1.24345e-14 4.00000 12.0000 -8.00000 16.0000 -8.00000 -16.0000 -8.00000 16.0000 -8.00000 12.0000 0.00000 -8.00000 0.00000 -8.00000 -4.00000 -4.00000 0.00000 -4.00000 0.00000 0.00000 13.5000 0.00000 13.5000 27.0000 13.5000 0.00000 13.5000 0.00000 -27.0000 13.5000 -27.0000 13.5000 54.0000 13.5000 -27.0000 13.5000 -27.0000 0.00000 13.5000 0.00000 13.5000 27.0000 13.5000 0.00000 13.5000 0.00000 0.00000 40.5000 0.00000 40.5000 54.0000 40.5000 0.00000 40.5000 0.00000 0.00000 -64.0000 0.00000 -64.0000 -128.000 -64.0000 0.00000 -64.0000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.500000 0.500000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.250000 0.250000 0.250000 -0.250000 0.00000 0.00000 0.00000 0.250000 -0.750000 0.00000 -0.750000 0.250000 -1.66533e-15 0.00000 0.00000 -1.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 -1.00000 0.00000 0.00000 0.00000 5.00000 5.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 -6.75000 -6.75000 -6.75000 6.75000 0.00000 0.00000 0.00000 0.00000 0.00000 8.88178e-15 0.500000 0.500000 0.500000 1.00000 1.00000 0.500000 1.00000 1.00000 0.00000 -0.500000 -0.500000 -0.500000 0.00000 0.00000 -0.500000 0.00000 0.00000 -1.00000 -0.500000 -0.500000 -0.500000 1.00000 -3.00000 -0.500000 -3.00000 -3.00000 -1.00000 -0.500000 -0.500000 -0.500000 -3.00000 -3.00000 -0.500000 -3.00000 1.00000 -1.77636e-14 -1.24345e-14 -1.24345e-14 -1.24345e-14 0.00000 0.00000 -1.24345e-14 0.00000 0.00000 -4.00000 -4.00000 -8.00000 -4.00000 0.00000 0.00000 -8.00000 0.00000 0.00000 -4.00000 -1.24345e-14 4.00000 -1.24345e-14 4.00000 8.00000 4.00000 8.00000 12.0000 -4.00000 4.00000 -1.24345e-14 4.00000 12.0000 8.00000 -1.24345e-14 8.00000 4.00000 -4.00000 -8.00000 -4.00000 -8.00000 0.00000 0.00000 -4.00000 0.00000 0.00000 -16.0000 -8.00000 -8.00000 -8.00000 12.0000 16.0000 -8.00000 16.0000 12.0000 27.0000 13.5000 13.5000 13.5000 0.00000 0.00000 13.5000 0.00000 0.00000 27.0000 13.5000 13.5000 13.5000 0.00000 0.00000 13.5000 0.00000 0.00000 54.0000 13.5000 13.5000 13.5000 -27.0000 -27.0000 13.5000 -27.0000 -27.0000 54.0000 40.5000 40.5000 40.5000 0.00000 0.00000 40.5000 0.00000 0.00000 -128.000 -64.0000 -64.0000 -64.0000 0.00000 0.00000 -64.0000 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.333333 0.333333 0.00000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::2.77556e-17 2.77556e-17 2.77556e-17 0.00000 -5.55112e-17 -1.11022e-16 -5.55112e-17 0.00000 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 0.00000 -0.333333 -0.333333 -0.148148 0.333333 5.55112e-17 0.185185 5.55112e-17 0.333333 0.185185 0.00000 0.00000 -0.148148 2.22045e-16 -1.33333 -0.148148 1.33333 1.33333 1.18519 -1.33333 2.22045e-16 -0.148148 0.00000 0.00000 -0.148148 0.00000 0.00000 -0.148148 0.00000 0.00000 -0.148148 4.44089e-16 4.44089e-16 -7.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 9.48148 2.00000 3.00000 2.00000 3.00000 2.00000 2.00000 2.00000 2.00000 0.888889 2.00000 -1.00000 1.11022e-16 -1.00000 -2.00000 -1.00000 1.11022e-16 -1.00000 -1.11111 -2.00000 -1.00000 -1.00000 -1.00000 2.00000 1.11022e-16 -1.00000 1.11022e-16 -1.11111 0.00000 0.00000 1.11022e-16 0.00000 0.00000 1.11022e-16 1.11022e-16 1.11022e-16 0.888889 0.00000 -1.77636e-15 -2.66454e-15 -1.77636e-15 8.00000 8.00000 -2.66454e-15 8.00000 0.888889 8.00000 8.00000 8.88178e-16 8.00000 8.00000 8.88178e-16 8.88178e-16 8.88178e-16 -7.11111 8.00000 -1.77636e-15 8.00000 -1.77636e-15 0.00000 -2.66454e-15 8.00000 -2.66454e-15 0.888889 0.00000 0.00000 -2.66454e-15 0.00000 0.00000 -2.66454e-15 -2.66454e-15 -2.66454e-15 0.888889 0.00000 0.00000 0.00000 0.00000 0.00000 8.88178e-16 0.00000 8.88178e-16 0.888889 0.00000 0.00000 8.88178e-16 0.00000 0.00000 0.00000 8.88178e-16 0.00000 0.888889 -18.0000 -9.00000 -9.00000 -9.00000 -18.0000 -9.00000 -9.00000 -9.00000 24.0000 0.00000 0.00000 0.00000 0.00000 0.00000 -9.00000 0.00000 -9.00000 6.00000 0.00000 0.00000 -9.00000 0.00000 0.00000 0.00000 -9.00000 0.00000 6.00000 0.00000 0.00000 9.00000 0.00000 0.00000 9.00000 9.00000 9.00000 24.0000 0.00000 0.00000 7.10543e-15 0.00000 0.00000 7.10543e-15 7.10543e-15 7.10543e-15 -56.8889
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.333333 0.00000 0.333333 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::2.77556e-17 2.77556e-17 0.00000 2.77556e-17 -5.55112e-17 0.00000 0.00000 -5.55112e-17 -1.11022e-16 0.00000 0.00000 1.00000 0.00000 0.00000 0.00000 -0.333333 -0.148148 -0.333333 0.333333 0.185185 5.55112e-17 0.00000 -0.148148 0.00000 5.55112e-17 0.185185 0.333333 2.22045e-16 -0.148148 -1.33333 0.00000 -0.148148 0.00000 0.00000 -0.148148 0.00000 -1.33333 -0.148148 2.22045e-16 1.33333 1.18519 1.33333 0.00000 -0.148148 0.00000 0.00000 -1.00000 0.00000 4.44089e-16 -7.00000 4.44089e-16 0.00000 -1.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 9.48148 0.00000 2.00000 2.00000 3.00000 2.00000 0.888889 2.00000 3.00000 2.00000 2.00000 2.00000 1.11022e-16 -1.00000 1.11022e-16 -1.11111 -1.00000 -1.00000 -1.00000 -2.00000 0.00000 1.11022e-16 0.00000 1.11022e-16 0.888889 1.11022e-16 0.00000 1.11022e-16 0.00000 -2.00000 -1.00000 -1.00000 -1.00000 -1.11111 1.11022e-16 -1.00000 1.11022e-16 2.00000 0.00000 -2.66454e-15 -1.77636e-15 -2.66454e-15 0.888889 8.00000 -1.77636e-15 8.00000 8.00000 0.00000 0.00000 0.00000 0.00000 0.888889 8.88178e-16 0.00000 8.88178e-16 0.00000 0.00000 -2.66454e-15 0.00000 -2.66454e-15 0.888889 -2.66454e-15 0.00000 -2.66454e-15 0.00000 8.00000 8.00000 -1.77636e-15 8.00000 0.888889 -2.66454e-15 -1.77636e-15 -2.66454e-15 0.00000 8.00000 8.88178e-16 8.00000 8.88178e-16 -7.11111 8.88178e-16 8.00000 8.88178e-16 8.00000 0.00000 8.88178e-16 0.00000 8.88178e-16 0.888889 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 6.00000 -9.00000 0.00000 -9.00000 0.00000 -18.0000 -9.00000 -9.00000 -9.00000 24.0000 -9.00000 -9.00000 -9.00000 -18.0000 0.00000 -9.00000 0.00000 -9.00000 6.00000 0.00000 0.00000 0.00000 0.00000 0.00000 9.00000 0.00000 9.00000 24.0000 9.00000 0.00000 9.00000 0.00000 0.00000 7.10543e-15 0.00000 7.10543e-15 -56.8889 7.10543e-15 0.00000 7.10543e-15 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.333333 0.333333 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::2.77556e-17 0.00000 2.77556e-17 2.77556e-17 0.00000 0.00000 -5.55112e-17 -5.55112e-17 0.00000 -1.11022e-16 0.00000 0.00000 1.00000 0.00000 0.00000 -0.148148 -0.333333 -0.333333 -0.148148 0.00000 0.00000 0.185185 0.333333 5.55112e-17 0.185185 5.55112e-17 0.333333 -0.148148 0.00000 0.00000 -0.148148 0.00000 0.00000 -0.148148 2.22045e-16 -1.33333 -0.148148 -1.33333 2.22045e-16 -0.148148 0.00000 0.00000 1.18519 1.33333 1.33333 -1.00000 0.00000 0.00000 -1.00000 0.00000 0.00000 -7.00000 4.44089e-16 4.44089e-16 -1.00000 0.00000 0.00000 9.48148 0.00000 0.00000 0.888889 2.00000 2.00000 2.00000 2.00000 3.00000 2.00000 3.00000 2.00000 0.888889 1.11022e-16 1.11022e-16 1.11022e-16 0.00000 0.00000 1.11022e-16 0.00000 0.00000 -1.11111 1.11022e-16 -1.00000 1.11022e-16 2.00000 -1.00000 -1.00000 -1.00000 -2.00000 -1.11111 -1.00000 1.11022e-16 -1.00000 -2.00000 -1.00000 1.11022e-16 -1.00000 2.00000 0.888889 -2.66454e-15 -2.66454e-15 -2.66454e-15 0.00000 0.00000 -2.66454e-15 0.00000 0.00000 0.888889 0.00000 8.88178e-16 0.00000 0.00000 0.00000 8.88178e-16 0.00000 0.00000 0.888889 -2.66454e-15 8.00000 -2.66454e-15 0.00000 -1.77636e-15 8.00000 -1.77636e-15 8.00000 0.888889 8.00000 -2.66454e-15 8.00000 8.00000 -1.77636e-15 -2.66454e-15 -1.77636e-15 0.00000 0.888889 8.88178e-16 0.00000 8.88178e-16 0.00000 0.00000 0.00000 0.00000 0.00000 -7.11111 8.88178e-16 8.88178e-16 8.88178e-16 8.00000 8.00000 8.88178e-16 8.00000 8.00000 6.00000 0.00000 -9.00000 0.00000 0.00000 0.00000 -9.00000 0.00000 0.00000 6.00000 -9.00000 0.00000 -9.00000 0.00000 0.00000 0.00000 0.00000 0.00000 24.0000 -9.00000 -9.00000 -9.00000 -18.0000 -9.00000 -9.00000 -9.00000 -18.0000 24.0000 9.00000 9.00000 9.00000 0.00000 0.00000 9.00000 0.00000 0.00000 -56.8889 7.10543e-15 7.10543e-15 7.10543e-15 0.00000 0.00000 7.10543e-15 0.00000 0.00000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.333333 0.333333 0.333333 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::1.31582e-16 5.75671e-17 5.75671e-17 5.75671e-17 -2.63164e-16 -1.64477e-17 -2.63164e-16 -2.63164e-16 -1.64477e-17 -1.64477e-17 -1.11022e-16 -1.11022e-16 -1.11022e-16 1.00000 1.05266e-15 0.148148 0.148148 0.148148 0.148148 -0.185185 -0.185185 -0.185185 0.148148 -0.185185 -0.185185 -0.185185 0.148148 0.148148 0.148148 0.148148 0.148148 0.148148 -1.18519 0.148148 0.148148 0.148148 0.148148 0.148148 0.148148 0.148148 -1.18519 0.148148 -1.18519 0.148148 0.148148 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 7.00000 7.00000 7.00000 -9.48148 -9.48148 -9.48148 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 -1.11111 -1.11111 -1.11111 -1.11111 -0.111111 -1.11111 -0.111111 -1.11111 -1.11111 -1.11111 -0.111111 -1.11111 0.888889 -1.11111 -0.111111 -1.11111 -1.11111 -1.11111 -0.111111 -1.11111 -0.111111 -1.11111 -1.11111 -1.11111 -1.11111 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 -7.11111 0.888889 0.888889 -7.11111 -7.11111 -7.11111 -7.11111 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 0.888889 -7.11111 0.888889 -7.11111 -7.11111 -7.11111 0.888889 -7.11111 0.888889 -7.11111 -7.11111 -7.11111 -7.11111 0.888889 0.888889 -7.11111 0.888889 0.888889 6.00000 6.00000 15.0000 6.00000 6.00000 15.0000 15.0000 15.0000 24.0000 6.00000 15.0000 6.00000 15.0000 24.0000 15.0000 6.00000 15.0000 6.00000 24.0000 15.0000 15.0000 15.0000 6.00000 6.00000 15.0000 6.00000 6.00000 24.0000 33.0000 33.0000 33.0000 24.0000 33.0000 33.0000 33.0000 24.0000 -56.8889 -56.8889 -56.8889 -56.8889 -56.8889 -56.8889 -56.8889 -56.8889 -56.8889
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.250000 0.250000 0.250000 :
+DEAL:Simplex::FE_P_Bubbles<3>(2)::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 -0.187500 -0.187500 -0.187500 0.187500 0.00000 0.00000 0.00000 0.187500 0.00000 0.00000 0.00000 0.187500 -4.44089e-16 -0.250000 -0.250000 0.250000 0.250000 0.00000 -0.250000 -4.44089e-16 -0.250000 -0.250000 -0.250000 -4.44089e-16 0.250000 0.00000 0.250000 0.00000 0.250000 0.250000 0.00000 0.00000 -1.68750 0.00000 -1.68750 0.00000 -1.68750 0.00000 0.00000 1.68750 1.68750 1.68750 0.00000 0.00000 0.00000 1.50000 2.00000 2.00000 2.00000 1.50000 2.00000 2.00000 2.00000 1.50000 1.50000 -0.500000 -0.500000 -0.500000 -1.00000 -0.500000 -0.500000 -0.500000 -1.00000 -1.00000 -0.500000 -0.500000 -0.500000 1.50000 -0.500000 -0.500000 -0.500000 -1.00000 -1.00000 -0.500000 -0.500000 -0.500000 -1.00000 -0.500000 -0.500000 -0.500000 1.50000 -3.55271e-15 -3.55271e-15 -3.55271e-15 -3.55271e-15 2.00000 4.00000 -3.55271e-15 4.00000 2.00000 2.00000 2.00000 -2.00000 2.00000 2.00000 -2.00000 -2.00000 -2.00000 -4.00000 2.00000 -3.55271e-15 4.00000 -3.55271e-15 -3.55271e-15 -3.10862e-15 4.00000 -3.10862e-15 2.00000 2.00000 4.00000 -3.55271e-15 4.00000 2.00000 -3.10862e-15 -3.55271e-15 -3.10862e-15 -3.55271e-15 2.00000 -2.00000 2.00000 -2.00000 -4.00000 -2.00000 2.00000 -2.00000 2.00000 -4.00000 -2.00000 -2.00000 -2.00000 2.00000 2.00000 -2.00000 2.00000 2.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 13.5000 0.00000 0.00000 0.00000 0.00000 13.5000 0.00000 0.00000 0.00000 0.00000 13.5000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 13.5000 13.5000 13.5000 13.5000 13.5000 13.5000 13.5000 13.5000 13.5000 -32.0000 -16.0000 -16.0000 -16.0000 -32.0000 -16.0000 -16.0000 -16.0000 -32.0000
+DEAL:Simplex::FE_P_Bubbles<3>(2)::
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Verify that FE_P_Bubbles can be used with a lumped mass matrix by computing a
+// convergence rate.
+
+#include <deal.II/base/function_lib.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_values.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/reference_cell.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools_integrate_difference.h>
+#include <deal.II/numerics/vector_tools_interpolate.h>
+
+#include <deal.II/simplex/fe_lib.h>
+#include <deal.II/simplex/quadrature_lib.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim, int spacedim = dim>
+Quadrature<dim>
+compute_nodal_quadrature(const FiniteElement<dim, spacedim> &fe)
+{
+ Assert(fe.n_blocks() == 1, ExcNotImplemented());
+ Assert(fe.n_components() == 1, ExcNotImplemented());
+
+ const ReferenceCell::Type type = fe.reference_cell_type();
+
+ const Quadrature<dim> q_gauss =
+ ReferenceCell::get_gauss_type_quadrature<dim>(type, fe.tensor_degree() + 1);
+ Triangulation<dim, spacedim> tria;
+ ReferenceCell::make_triangulation(type, tria);
+ const Mapping<dim, spacedim> &mapping =
+ ReferenceCell::get_default_linear_mapping<dim, spacedim>(type);
+
+ auto cell = tria.begin_active();
+ FEValues<dim, spacedim> fe_values(mapping,
+ fe,
+ q_gauss,
+ update_values | update_JxW_values);
+ fe_values.reinit(cell);
+
+ std::vector<Point<dim>> nodal_quad_points = fe.get_unit_support_points();
+ std::vector<double> nodal_quad_weights(nodal_quad_points.size());
+ Assert(nodal_quad_points.size() > 0, ExcNotImplemented());
+ for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
+ {
+ double integral = 0.0;
+ for (unsigned int q = 0; q < q_gauss.size(); ++q)
+ integral += fe_values.shape_value(i, q) * fe_values.JxW(q);
+ nodal_quad_weights[i] = integral;
+ }
+ return {nodal_quad_points, nodal_quad_weights};
+}
+
+
+
+template <int dim, int spacedim = dim>
+void
+test_interpolate()
+{
+ deallog << std::endl << "test interpolation" << std::endl;
+ Functions::CosineFunction<dim> func(1);
+
+ const unsigned int n_refinements = 7 - dim;
+ std::vector<Triangulation<dim, spacedim>> trias;
+ for (unsigned int refinement_n = 0; refinement_n < n_refinements;
+ ++refinement_n)
+ {
+ Triangulation<dim, spacedim> tria;
+ Triangulation<dim, spacedim> tria_hypercube;
+ GridGenerator::subdivided_hyper_cube(tria_hypercube,
+ std::pow(2, refinement_n));
+ if (dim == spacedim)
+ GridTools::distort_random(0.2, tria_hypercube);
+ GridGenerator::convert_hypercube_to_simplex_mesh(tria_hypercube, tria);
+ trias.emplace_back(std::move(tria));
+ }
+
+ deallog << "dim = " << dim << std::endl;
+ for (unsigned int degree = 0; degree < 3; ++degree)
+ {
+ deallog << "degree = " << degree << std::endl;
+ double old_error = -1.0;
+ for (unsigned int refinement_n = 0; refinement_n < n_refinements - degree;
+ ++refinement_n)
+ {
+ const Triangulation<dim, spacedim> &tria = trias[refinement_n];
+ deallog << "number of cells = " << tria.n_active_cells() << std::endl;
+
+ Simplex::FE_P_Bubbles<dim, spacedim> fe(degree);
+
+ const ReferenceCell::Type type = fe.reference_cell_type();
+ DoFHandler<dim, spacedim> dh(tria);
+ dh.distribute_dofs(fe);
+ deallog << "number of dofs = " << dh.n_dofs() << std::endl;
+
+ const Mapping<dim, spacedim> &map =
+ ReferenceCell::get_default_linear_mapping<dim, spacedim>(type);
+
+ Vector<double> solution(dh.n_dofs());
+ VectorTools::interpolate(map, dh, func, solution);
+
+ Simplex::QGauss<dim> error_quad(4);
+ Vector<double> out_l2(tria.n_active_cells());
+ VectorTools::integrate_difference(
+ map, dh, solution, func, out_l2, error_quad, VectorTools::L2_norm);
+ const double new_error =
+ VectorTools::compute_global_error(tria,
+ out_l2,
+ VectorTools::L2_norm);
+ deallog << "error = " << new_error << std::endl;
+ if (old_error != -1.0)
+ deallog << "ratio = " << old_error / new_error << std::endl;
+ old_error = new_error;
+ }
+ }
+}
+
+
+
+template <int dim, int spacedim = dim>
+void
+test_lumped_project()
+{
+ deallog << std::endl << "test lumped project" << std::endl;
+ const auto func = Functions::CosineFunction<dim>(1);
+
+ const unsigned int n_refinements = 7 - dim;
+ std::vector<Triangulation<dim, spacedim>> trias;
+ for (unsigned int refinement_n = 0; refinement_n < n_refinements;
+ ++refinement_n)
+ {
+ Triangulation<dim, spacedim> tria;
+ Triangulation<dim, spacedim> tria_hypercube;
+ GridGenerator::subdivided_hyper_cube(tria_hypercube,
+ std::pow(2, refinement_n));
+ if (dim == spacedim)
+ GridTools::distort_random(0.2, tria_hypercube);
+ GridGenerator::convert_hypercube_to_simplex_mesh(tria_hypercube, tria);
+ trias.emplace_back(std::move(tria));
+ }
+
+
+ deallog << "dim = " << dim << std::endl;
+ for (unsigned int degree = 0; degree < 3; ++degree)
+ {
+ deallog << "degree = " << degree << std::endl;
+ double old_error = -1.0;
+ for (unsigned int refinement_n = 0; refinement_n < n_refinements;
+ ++refinement_n)
+ {
+ const Triangulation<dim, spacedim> &tria = trias[refinement_n];
+ deallog << "number of cells = " << tria.n_active_cells() << std::endl;
+
+ Simplex::FE_P_Bubbles<dim, spacedim> fe(degree);
+
+ const ReferenceCell::Type type = fe.reference_cell_type();
+ DoFHandler<dim, spacedim> dh(tria);
+ dh.distribute_dofs(fe);
+ deallog << "number of dofs = " << dh.n_dofs() << std::endl;
+ const Quadrature<dim> nodal_quad = compute_nodal_quadrature(fe);
+ const Quadrature<dim> cell_quad = Simplex::QGauss<dim>(
+ std::max<unsigned int>(fe.tensor_degree() + 1, 2));
+
+ Vector<double> lumped_mass(dh.n_dofs());
+ Vector<double> consistent_rhs(dh.n_dofs());
+
+ const Mapping<dim, spacedim> &map =
+ ReferenceCell::get_default_linear_mapping<dim, spacedim>(type);
+
+ FEValues<dim> lumped_fev(map,
+ fe,
+ nodal_quad,
+ update_values | update_JxW_values);
+ FEValues<dim> consistent_fev(map,
+ fe,
+ cell_quad,
+ update_quadrature_points |
+ update_values | update_JxW_values);
+
+ std::vector<types::global_dof_index> dofs(fe.dofs_per_cell);
+ for (const auto &cell : dh.active_cell_iterators())
+ {
+ lumped_fev.reinit(cell);
+ consistent_fev.reinit(cell);
+ cell->get_dof_indices(dofs);
+
+ for (unsigned int q = 0; q < nodal_quad.size(); ++q)
+ {
+ for (unsigned int i = 0; i < fe.dofs_per_cell; ++i)
+ {
+ const double v = lumped_fev.shape_value(i, q);
+ Assert(std::abs(v - double(i == q)) < 1e-14,
+ ExcInternalError());
+ lumped_mass[dofs[i]] +=
+ lumped_fev.shape_value(i, q) * lumped_fev.JxW(q);
+ }
+ }
+
+ for (unsigned int q = 0; q < cell_quad.size(); ++q)
+ for (unsigned int i = 0; i < fe.dofs_per_cell; ++i)
+ consistent_rhs[dofs[i]] +=
+ consistent_fev.shape_value(i, q) *
+ func.value(consistent_fev.quadrature_point(q)) *
+ consistent_fev.JxW(q);
+ }
+
+ Vector<double> solution(dh.n_dofs());
+ for (std::size_t i = 0; i < solution.size(); ++i)
+ solution[i] = consistent_rhs[i] / lumped_mass[i];
+
+ Simplex::QGauss<dim> error_quad(4);
+ Vector<double> out_l2(tria.n_active_cells());
+ VectorTools::integrate_difference(
+ map, dh, solution, func, out_l2, error_quad, VectorTools::L2_norm);
+
+ const auto new_error =
+ VectorTools::compute_global_error(tria,
+ out_l2,
+ VectorTools::L2_norm);
+ deallog << "error = " << new_error << std::endl;
+ if (old_error != -1.0)
+ deallog << "ratio = " << old_error / new_error << std::endl;
+ old_error = new_error;
+
+#if 0
+ static unsigned int counter = 0;
+ for (unsigned int n_subdivisions = 1; n_subdivisions <= 2;
+ ++n_subdivisions)
+ {
+ DataOut<dim> data_out;
+
+ data_out.attach_dof_handler(dh);
+ data_out.add_data_vector(solution, "solution",
+ DataOut<dim>::type_dof_data);
+ data_out.build_patches(map, n_subdivisions);
+ std::ofstream output("test." + std::to_string(dim) + "." +
+ std::to_string(counter++) + ".vtk");
+ data_out.write_vtk(output);
+ }
+#endif
+ }
+ }
+}
+
+
+
+int
+main()
+{
+ initlog();
+
+ test_interpolate<2, 2>();
+ test_interpolate<3, 3>();
+
+ test_lumped_project<2, 2>();
+ test_lumped_project<3, 3>();
+}
--- /dev/null
+
+DEAL::
+DEAL::test interpolation
+DEAL::dim = 2
+DEAL::degree = 0
+DEAL::number of cells = 8
+DEAL::number of dofs = 8
+DEAL::error = 0.129785
+DEAL::number of cells = 32
+DEAL::number of dofs = 32
+DEAL::error = 0.0657918
+DEAL::ratio = 1.97266
+DEAL::number of cells = 128
+DEAL::number of dofs = 128
+DEAL::error = 0.0340434
+DEAL::ratio = 1.93258
+DEAL::number of cells = 512
+DEAL::number of dofs = 512
+DEAL::error = 0.0169138
+DEAL::ratio = 2.01276
+DEAL::number of cells = 2048
+DEAL::number of dofs = 2048
+DEAL::error = 0.00849451
+DEAL::ratio = 1.99115
+DEAL::degree = 1
+DEAL::number of cells = 8
+DEAL::number of dofs = 9
+DEAL::error = 0.0600352
+DEAL::number of cells = 32
+DEAL::number of dofs = 25
+DEAL::error = 0.0162206
+DEAL::ratio = 3.70117
+DEAL::number of cells = 128
+DEAL::number of dofs = 81
+DEAL::error = 0.00410053
+DEAL::ratio = 3.95573
+DEAL::number of cells = 512
+DEAL::number of dofs = 289
+DEAL::error = 0.00108554
+DEAL::ratio = 3.77741
+DEAL::degree = 2
+DEAL::number of cells = 8
+DEAL::number of dofs = 33
+DEAL::error = 0.00424328
+DEAL::number of cells = 32
+DEAL::number of dofs = 113
+DEAL::error = 0.000564602
+DEAL::ratio = 7.51554
+DEAL::number of cells = 128
+DEAL::number of dofs = 417
+DEAL::error = 7.95939e-05
+DEAL::ratio = 7.09353
+DEAL::
+DEAL::test interpolation
+DEAL::dim = 3
+DEAL::degree = 0
+DEAL::number of cells = 24
+DEAL::number of dofs = 24
+DEAL::error = 0.112063
+DEAL::number of cells = 192
+DEAL::number of dofs = 192
+DEAL::error = 0.0562048
+DEAL::ratio = 1.99384
+DEAL::number of cells = 1536
+DEAL::number of dofs = 1536
+DEAL::error = 0.0288430
+DEAL::ratio = 1.94865
+DEAL::number of cells = 12288
+DEAL::number of dofs = 12288
+DEAL::error = 0.0145439
+DEAL::ratio = 1.98317
+DEAL::degree = 1
+DEAL::number of cells = 24
+DEAL::number of dofs = 14
+DEAL::error = 0.0754385
+DEAL::number of cells = 192
+DEAL::number of dofs = 63
+DEAL::error = 0.0202341
+DEAL::ratio = 3.72829
+DEAL::number of cells = 1536
+DEAL::number of dofs = 365
+DEAL::error = 0.00523449
+DEAL::ratio = 3.86553
+DEAL::degree = 2
+DEAL::number of cells = 24
+DEAL::number of dofs = 147
+DEAL::error = 0.00517967
+DEAL::number of cells = 192
+DEAL::number of dofs = 989
+DEAL::error = 0.000673567
+DEAL::ratio = 7.68990
+DEAL::
+DEAL::test lumped project
+DEAL::dim = 2
+DEAL::degree = 0
+DEAL::number of cells = 8
+DEAL::number of dofs = 8
+DEAL::error = 0.128417
+DEAL::number of cells = 32
+DEAL::number of dofs = 32
+DEAL::error = 0.0656044
+DEAL::ratio = 1.95745
+DEAL::number of cells = 128
+DEAL::number of dofs = 128
+DEAL::error = 0.0340204
+DEAL::ratio = 1.92838
+DEAL::number of cells = 512
+DEAL::number of dofs = 512
+DEAL::error = 0.0169106
+DEAL::ratio = 2.01178
+DEAL::number of cells = 2048
+DEAL::number of dofs = 2048
+DEAL::error = 0.00849411
+DEAL::ratio = 1.99086
+DEAL::degree = 1
+DEAL::number of cells = 8
+DEAL::number of dofs = 9
+DEAL::error = 0.112188
+DEAL::number of cells = 32
+DEAL::number of dofs = 25
+DEAL::error = 0.0427299
+DEAL::ratio = 2.62551
+DEAL::number of cells = 128
+DEAL::number of dofs = 81
+DEAL::error = 0.0170503
+DEAL::ratio = 2.50610
+DEAL::number of cells = 512
+DEAL::number of dofs = 289
+DEAL::error = 0.00659734
+DEAL::ratio = 2.58443
+DEAL::number of cells = 2048
+DEAL::number of dofs = 1089
+DEAL::error = 0.00268244
+DEAL::ratio = 2.45946
+DEAL::degree = 2
+DEAL::number of cells = 8
+DEAL::number of dofs = 33
+DEAL::error = 0.0112526
+DEAL::number of cells = 32
+DEAL::number of dofs = 113
+DEAL::error = 0.00294719
+DEAL::ratio = 3.81808
+DEAL::number of cells = 128
+DEAL::number of dofs = 417
+DEAL::error = 0.000739778
+DEAL::ratio = 3.98389
+DEAL::number of cells = 512
+DEAL::number of dofs = 1601
+DEAL::error = 0.000195335
+DEAL::ratio = 3.78722
+DEAL::number of cells = 2048
+DEAL::number of dofs = 6273
+DEAL::error = 4.89802e-05
+DEAL::ratio = 3.98804
+DEAL::
+DEAL::test lumped project
+DEAL::dim = 3
+DEAL::degree = 0
+DEAL::number of cells = 24
+DEAL::number of dofs = 24
+DEAL::error = 0.110613
+DEAL::number of cells = 192
+DEAL::number of dofs = 192
+DEAL::error = 0.0560060
+DEAL::ratio = 1.97503
+DEAL::number of cells = 1536
+DEAL::number of dofs = 1536
+DEAL::error = 0.0288177
+DEAL::ratio = 1.94346
+DEAL::number of cells = 12288
+DEAL::number of dofs = 12288
+DEAL::error = 0.0145407
+DEAL::ratio = 1.98187
+DEAL::degree = 1
+DEAL::number of cells = 24
+DEAL::number of dofs = 14
+DEAL::error = 0.117148
+DEAL::number of cells = 192
+DEAL::number of dofs = 63
+DEAL::error = 0.0467977
+DEAL::ratio = 2.50329
+DEAL::number of cells = 1536
+DEAL::number of dofs = 365
+DEAL::error = 0.0170371
+DEAL::ratio = 2.74682
+DEAL::number of cells = 12288
+DEAL::number of dofs = 2457
+DEAL::error = 0.00632181
+DEAL::ratio = 2.69496
+DEAL::degree = 2
+DEAL::number of cells = 24
+DEAL::number of dofs = 147
+DEAL::error = 0.0115617
+DEAL::number of cells = 192
+DEAL::number of dofs = 989
+DEAL::error = 0.00300755
+DEAL::ratio = 3.84423
+DEAL::number of cells = 1536
+DEAL::number of dofs = 7257
+DEAL::error = 0.000784223
+DEAL::ratio = 3.83507
+DEAL::number of cells = 12288
+DEAL::number of dofs = 55601
+DEAL::error = 0.000200317
+DEAL::ratio = 3.91491