// ---------------------------------------------------------------------
+#include <deal.II/base/polynomial.h>
#include <deal.II/base/qprojector.h>
#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/tensor_product_polynomials.h>
#include <deal.II/dofs/dof_accessor.h>
#include <deal.II/fe/fe_nothing.h>
#include <deal.II/fe/fe_raviart_thomas.h>
#include <deal.II/fe/fe_tools.h>
-#include <deal.II/fe/fe_values.h>
#include <deal.II/fe/mapping.h>
-#include <deal.II/grid/tria.h>
-#include <deal.II/grid/tria_iterator.h>
-
#include <memory>
#include <sstream>
DEAL_II_NAMESPACE_OPEN
-// TODO: implement the adjust_quad_dof_index_for_face_orientation_table and
-// adjust_line_dof_index_for_line_orientation_table fields, and write tests
-// similar to bits/face_orientation_and_fe_q_*
+
+// ---------------- polynomial class for FE_RaviartThomasNodal ---------------
+
+namespace
+{
+ template <int dim>
+ class PolynomialsRaviartThomasNodal : public TensorPolynomialsBase<dim>
+ {
+ public:
+ PolynomialsRaviartThomasNodal(const unsigned int degree);
+
+ /**
+ * Compute the value and derivatives of each Raviart-Thomas polynomial at
+ * @p unit_point.
+ *
+ * The size of the vectors must either be zero or equal <tt>n()</tt>. In
+ * the first case, the function will not compute these values.
+ */
+ void
+ evaluate(const Point<dim> & unit_point,
+ std::vector<Tensor<1, dim>> &values,
+ std::vector<Tensor<2, dim>> &grads,
+ std::vector<Tensor<3, dim>> &grad_grads,
+ std::vector<Tensor<4, dim>> &third_derivatives,
+ std::vector<Tensor<5, dim>> &fourth_derivatives) const override;
+
+ /**
+ * Return the name of the space, which is <tt>PolynomialsRaviartThomas</tt>.
+ */
+ std::string
+ name() const override;
+
+ /**
+ * Return the number of polynomials in the space without requiring to
+ * build an object of PolynomialsRaviartThomas. This is required by the
+ * FiniteElement classes.
+ */
+ static unsigned int
+ n_polynomials(const unsigned int degree);
+
+ const std::vector<unsigned int> &
+ get_renumbering() const;
+
+ /**
+ * @copydoc TensorPolynomialsBase<dim>::clone()
+ */
+ virtual std::unique_ptr<TensorPolynomialsBase<dim>>
+ clone() const override;
+
+ /**
+ * Compute the generalized support points of the associated element in the
+ * ordering of the element. Note that they are not support points in the
+ * classical sense as polynomials of the different components have
+ * different points, which need to be combined in terms of Piola
+ * transforms.
+ */
+ std::vector<Point<dim>>
+ get_polynomial_support_points() const;
+
+ private:
+ /**
+ * The degree variable passed to the constructor.
+ */
+ const unsigned int degree;
+
+ /**
+ * An object representing the polynomial space for a single component. We
+ * can re-use it by rotating the coordinates of the evaluation point.
+ */
+ const AnisotropicPolynomials<dim> polynomial_space;
+
+ /**
+ * Renumbering from lexicographic to hierarchic order.
+ */
+ std::vector<unsigned int> lexicographic_to_hierarchic;
+
+ /**
+ * Renumbering from hierarchic to lexicographic order. Inverse of
+ * lexicographic_to_hierarchic.
+ */
+ std::vector<unsigned int> hierarchic_to_lexicographic;
+
+ /**
+ * Renumbering from shifted polynomial spaces to lexicographic one
+ */
+ std::array<std::vector<unsigned int>, dim> renumber_aniso;
+ };
+
+
+
+ // Create nodal Raviart-Thomas polynomials as the tensor product of Lagrange
+ // polynomials on Gauss-Lobatto points of degree + 2 points in the
+ // continuous direction and degree + 1 points in the discontinuous
+ // directions (we could also choose Lagrange polynomials on Gauss points but
+ // those are slightly more expensive to handle in classes).
+ std::vector<std::vector<Polynomials::Polynomial<double>>>
+ create_rt_polynomials(const unsigned int dim, const unsigned int degree)
+ {
+ std::vector<std::vector<Polynomials::Polynomial<double>>> pols(dim);
+ pols[0] = Polynomials::generate_complete_Lagrange_basis(
+ QGaussLobatto<1>(degree + 2).get_points());
+ if (degree > 0)
+ for (unsigned int d = 1; d < dim; ++d)
+ pols[d] = Polynomials::generate_complete_Lagrange_basis(
+ QGaussLobatto<1>(degree + 1).get_points());
+ else
+ for (unsigned int d = 1; d < dim; ++d)
+ pols[d] = Polynomials::generate_complete_Lagrange_basis(
+ QMidpoint<1>().get_points());
+
+ return pols;
+ }
+
+
+
+ // set up the numbering of the rt polynomials
+ std::vector<unsigned int>
+ compute_rt_hierarchic_to_lexicographic(const unsigned int dim,
+ const unsigned int degree)
+ {
+ const unsigned int n_pols =
+ (degree + 2) * Utilities::pow(degree + 1, dim - 1);
+
+ std::vector<unsigned int> hierarchic_to_lexicographic;
+
+ // dofs on faces
+ for (unsigned int face_no = 0; face_no < 2 * dim; ++face_no)
+ {
+ const unsigned int stride_x = face_no < 2 ? degree + 2 : 1;
+ const unsigned int stride_y =
+ face_no < 4 ? (degree + 2) * (degree + 1) : degree + 1;
+ const unsigned int offset =
+ (face_no % 2) * Utilities::pow(degree + 1, 1 + face_no / 2);
+ for (unsigned int j = 0; j < (dim > 2 ? degree + 1 : 1); ++j)
+ for (unsigned int i = 0; i < degree + 1; ++i)
+ hierarchic_to_lexicographic.push_back(
+ (face_no / 2) * n_pols + offset + i * stride_x + j * stride_y);
+ }
+ // dofs on cells, starting with x component...
+ for (unsigned int k = 0; k < (dim > 2 ? degree + 1 : 1); ++k)
+ for (unsigned int j = 0; j < (dim > 1 ? degree + 1 : 1); ++j)
+ for (unsigned int i = 1; i < degree + 1; ++i)
+ hierarchic_to_lexicographic.push_back(
+ k * (degree + 1) * (degree + 2) + j * (degree + 2) + i);
+ // ... then y component ...
+ if (dim > 1)
+ for (unsigned int k = 0; k < (dim > 2 ? degree + 1 : 1); ++k)
+ for (unsigned int j = 1; j < degree + 1; ++j)
+ for (unsigned int i = 0; i < degree + 1; ++i)
+ hierarchic_to_lexicographic.push_back(
+ n_pols + k * (degree + 1) * (degree + 2) + j * (degree + 1) + i);
+ // ... and finally z component
+ if (dim > 2)
+ for (unsigned int k = 1; k < degree + 1; ++k)
+ for (unsigned int j = 0; j < degree + 1; ++j)
+ for (unsigned int i = 0; i < degree + 1; ++i)
+ hierarchic_to_lexicographic.push_back(
+ 2 * n_pols + k * (degree + 1) * (degree + 1) + j * (degree + 1) +
+ i);
+
+ AssertDimension(hierarchic_to_lexicographic.size(), n_pols * dim);
+
+#ifdef DEBUG
+ // assert that we have a valid permutation
+ std::vector<unsigned int> copy(hierarchic_to_lexicographic);
+ std::sort(copy.begin(), copy.end());
+ for (unsigned int i = 0; i < copy.size(); ++i)
+ AssertDimension(i, copy[i]);
+#endif
+
+ return hierarchic_to_lexicographic;
+ }
+
+
+
+ template <int dim>
+ PolynomialsRaviartThomasNodal<dim>::PolynomialsRaviartThomasNodal(
+ const unsigned int degree)
+ : TensorPolynomialsBase<dim>(degree, n_polynomials(degree))
+ , degree(degree)
+ , polynomial_space(create_rt_polynomials(dim, degree))
+ {
+ // create renumbering of the unknowns from the lexicographic order to the
+ // actual order required by the finite element class with unknowns on
+ // faces placed first
+ const unsigned int n_pols = polynomial_space.n();
+ hierarchic_to_lexicographic =
+ compute_rt_hierarchic_to_lexicographic(dim, degree);
+
+ lexicographic_to_hierarchic =
+ Utilities::invert_permutation(hierarchic_to_lexicographic);
+
+ // since we only store an anisotropic polynomial for the first component,
+ // we set up a second numbering to switch out the actual coordinate
+ // direction
+ renumber_aniso[0].resize(n_pols);
+ for (unsigned int i = 0; i < n_pols; ++i)
+ renumber_aniso[0][i] = i;
+ if (dim > 1)
+ {
+ // switch x and y component (i and j loops)
+ renumber_aniso[1].resize(n_pols);
+ for (unsigned int k = 0; k < (dim > 2 ? degree + 1 : 1); ++k)
+ for (unsigned int j = 0; j < degree + 2; ++j)
+ for (unsigned int i = 0; i < degree + 1; ++i)
+ renumber_aniso[1][(k * (degree + 2) + j) * (degree + 1) + i] =
+ j + i * (degree + 2) + k * (degree + 2) * (degree + 1);
+ }
+ if (dim > 2)
+ {
+ // switch x and z component (i and k loops)
+ renumber_aniso[2].resize(n_pols);
+ for (unsigned int k = 0; k < degree + 2; ++k)
+ for (unsigned int j = 0; j < degree + 1; ++j)
+ for (unsigned int i = 0; i < degree + 1; ++i)
+ renumber_aniso[2][(k * (degree + 1) + j) * (degree + 1) + i] =
+ k + i * (degree + 2) + j * (degree + 2) * (degree + 1);
+ }
+ }
+
+
+
+ template <int dim>
+ void
+ PolynomialsRaviartThomasNodal<dim>::evaluate(
+ const Point<dim> & unit_point,
+ std::vector<Tensor<1, dim>> &values,
+ std::vector<Tensor<2, dim>> &grads,
+ std::vector<Tensor<3, dim>> &grad_grads,
+ std::vector<Tensor<4, dim>> &third_derivatives,
+ std::vector<Tensor<5, dim>> &fourth_derivatives) const
+ {
+ Assert(values.size() == this->n() || values.size() == 0,
+ ExcDimensionMismatch(values.size(), this->n()));
+ Assert(grads.size() == this->n() || grads.size() == 0,
+ ExcDimensionMismatch(grads.size(), this->n()));
+ Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+ ExcDimensionMismatch(grad_grads.size(), this->n()));
+ Assert(third_derivatives.size() == this->n() ||
+ third_derivatives.size() == 0,
+ ExcDimensionMismatch(third_derivatives.size(), this->n()));
+ Assert(fourth_derivatives.size() == this->n() ||
+ fourth_derivatives.size() == 0,
+ ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
+
+ std::vector<double> p_values;
+ std::vector<Tensor<1, dim>> p_grads;
+ std::vector<Tensor<2, dim>> p_grad_grads;
+ std::vector<Tensor<3, dim>> p_third_derivatives;
+ std::vector<Tensor<4, dim>> p_fourth_derivatives;
+
+ const unsigned int n_sub = polynomial_space.n();
+ p_values.resize((values.size() == 0) ? 0 : n_sub);
+ p_grads.resize((grads.size() == 0) ? 0 : n_sub);
+ p_grad_grads.resize((grad_grads.size() == 0) ? 0 : n_sub);
+ p_third_derivatives.resize((third_derivatives.size() == 0) ? 0 : n_sub);
+ p_fourth_derivatives.resize((fourth_derivatives.size() == 0) ? 0 : n_sub);
+
+ for (unsigned int d = 0; d < dim; ++d)
+ {
+ // First we copy the point. The polynomial space for component d
+ // consists of polynomials of degree k in x_d and degree k+1 in the
+ // other variables. in order to simplify this, we use the same
+ // AnisotropicPolynomial space and simply rotate the coordinates
+ // through all directions.
+ Point<dim> p;
+ for (unsigned int c = 0; c < dim; ++c)
+ p(c) = unit_point((c + d) % dim);
+
+ polynomial_space.evaluate(p,
+ p_values,
+ p_grads,
+ p_grad_grads,
+ p_third_derivatives,
+ p_fourth_derivatives);
+
+ for (unsigned int i = 0; i < p_values.size(); ++i)
+ values[lexicographic_to_hierarchic[i + d * n_sub]][d] =
+ p_values[renumber_aniso[d][i]];
+
+ for (unsigned int i = 0; i < p_grads.size(); ++i)
+ for (unsigned int d1 = 0; d1 < dim; ++d1)
+ grads[lexicographic_to_hierarchic[i + d * n_sub]][d]
+ [(d1 + d) % dim] = p_grads[renumber_aniso[d][i]][d1];
+
+ for (unsigned int i = 0; i < p_grad_grads.size(); ++i)
+ for (unsigned int d1 = 0; d1 < dim; ++d1)
+ for (unsigned int d2 = 0; d2 < dim; ++d2)
+ grad_grads[lexicographic_to_hierarchic[i + d * n_sub]][d]
+ [(d1 + d) % dim][(d2 + d) % dim] =
+ p_grad_grads[renumber_aniso[d][i]][d1][d2];
+
+ for (unsigned int i = 0; i < p_third_derivatives.size(); ++i)
+ for (unsigned int d1 = 0; d1 < dim; ++d1)
+ for (unsigned int d2 = 0; d2 < dim; ++d2)
+ for (unsigned int d3 = 0; d3 < dim; ++d3)
+ third_derivatives[lexicographic_to_hierarchic[i + d * n_sub]][d]
+ [(d1 + d) % dim][(d2 + d) % dim]
+ [(d3 + d) % dim] =
+ p_third_derivatives[renumber_aniso[d][i]][d1]
+ [d2][d3];
+
+ for (unsigned int i = 0; i < p_fourth_derivatives.size(); ++i)
+ for (unsigned int d1 = 0; d1 < dim; ++d1)
+ for (unsigned int d2 = 0; d2 < dim; ++d2)
+ for (unsigned int d3 = 0; d3 < dim; ++d3)
+ for (unsigned int d4 = 0; d4 < dim; ++d4)
+ fourth_derivatives[lexicographic_to_hierarchic[i + d * n_sub]]
+ [d][(d1 + d) % dim][(d2 + d) % dim]
+ [(d3 + d) % dim][(d4 + d) % dim] =
+ p_fourth_derivatives[renumber_aniso[d][i]]
+ [d1][d2][d3][d4];
+ }
+ }
+
+
+
+ template <int dim>
+ std::string
+ PolynomialsRaviartThomasNodal<dim>::name() const
+ {
+ return "PolynomialsRaviartThomasNodal";
+ }
+
+
+
+ template <int dim>
+ unsigned int
+ PolynomialsRaviartThomasNodal<dim>::n_polynomials(unsigned int degree)
+ {
+ return dim * (degree + 2) * Utilities::pow(degree + 1, dim - 1);
+ }
+
+
+
+ template <int dim>
+ std::unique_ptr<TensorPolynomialsBase<dim>>
+ PolynomialsRaviartThomasNodal<dim>::clone() const
+ {
+ return std::make_unique<PolynomialsRaviartThomasNodal<dim>>(*this);
+ }
+
+
+
+ template <int dim>
+ std::vector<Point<dim>>
+ PolynomialsRaviartThomasNodal<dim>::get_polynomial_support_points() const
+ {
+ Assert(dim > 0 && dim <= 3, ExcImpossibleInDim(dim));
+ const Quadrature<1> low(
+ degree == 0 ? static_cast<Quadrature<1>>(QMidpoint<1>()) :
+ static_cast<Quadrature<1>>(QGaussLobatto<1>(degree + 1)));
+ const QGaussLobatto<1> high(degree + 2);
+ const QAnisotropic<dim> quad =
+ (dim == 1 ? QAnisotropic<dim>(high) :
+ (dim == 2 ? QAnisotropic<dim>(high, low) :
+ QAnisotropic<dim>(high, low, low)));
+
+ const unsigned int n_sub = polynomial_space.n();
+ std::vector<Point<dim>> points(dim * n_sub);
+ points.resize(n_polynomials(degree));
+ for (unsigned int d = 0; d < dim; ++d)
+ for (unsigned int i = 0; i < n_sub; ++i)
+ points[lexicographic_to_hierarchic[i + d * n_sub]] =
+ quad.point(renumber_aniso[d][i]);
+ return points;
+ }
+
+
+
+ // Return a vector of "dofs per object" where the components of the returned
+ // vector refer to:
+ // 0 = vertex
+ // 1 = edge
+ // 2 = face (which is a cell in 2D)
+ // 3 = cell
+ std::vector<unsigned int>
+ get_rt_dpo_vector(const unsigned int dim, const unsigned int degree)
+ {
+ std::vector<unsigned int> dpo(dim + 1);
+ dpo[0] = 0;
+ dpo[1] = 0;
+ unsigned int dofs_per_face = 1;
+ for (unsigned int d = 1; d < dim; ++d)
+ dofs_per_face *= (degree + 1);
+
+ dpo[dim - 1] = dofs_per_face;
+ dpo[dim] = dim * degree * dofs_per_face;
+
+ return dpo;
+ }
+} // namespace
+
+
+
+// --------------------- actual implementation of element --------------------
template <int dim>
-FE_RaviartThomasNodal<dim>::FE_RaviartThomasNodal(const unsigned int deg)
- : FE_PolyTensor<dim>(PolynomialsRaviartThomas<dim>(deg),
- FiniteElementData<dim>(get_dpo_vector(deg),
+FE_RaviartThomasNodal<dim>::FE_RaviartThomasNodal(const unsigned int degree)
+ : FE_PolyTensor<dim>(PolynomialsRaviartThomasNodal<dim>(degree),
+ FiniteElementData<dim>(get_rt_dpo_vector(dim, degree),
dim,
- deg + 1,
+ degree + 1,
FiniteElementData<dim>::Hdiv),
- get_ria_vector(deg),
+ std::vector<bool>(1, false),
std::vector<ComponentMask>(
- PolynomialsRaviartThomas<dim>::n_polynomials(deg),
+ PolynomialsRaviartThomasNodal<dim>::n_polynomials(
+ degree),
std::vector<bool>(dim, true)))
{
Assert(dim >= 2, ExcImpossibleInDim(dim));
- const unsigned int n_dofs = this->n_dofs_per_cell();
this->mapping_kind = {mapping_raviart_thomas};
- // First, initialize the
- // generalized support points and
- // quadrature weights, since they
- // are required for interpolation.
- initialize_support_points(deg);
-
- // Now compute the inverse node matrix, generating the correct
- // basis functions from the raw ones. For a discussion of what
- // exactly happens here, see FETools::compute_node_matrix.
- const FullMatrix<double> M = FETools::compute_node_matrix(*this);
- this->inverse_node_matrix.reinit(n_dofs, n_dofs);
- this->inverse_node_matrix.invert(M);
- // From now on, the shape functions provided by FiniteElement::shape_value
- // and similar functions will be the correct ones, not
- // the raw shape functions from the polynomial space anymore.
-
- // Reinit the vectors of
- // prolongation matrices to the
- // right sizes. There are no
- // restriction matrices implemented
- for (unsigned int ref_case = RefinementCase<dim>::cut_x;
- ref_case < RefinementCase<dim>::isotropic_refinement + 1;
- ++ref_case)
- {
- const unsigned int nc =
- GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
- for (unsigned int i = 0; i < nc; ++i)
- this->prolongation[ref_case - 1][i].reinit(n_dofs, n_dofs);
- }
+ // First, initialize the generalized support points and quadrature weights,
+ // since they are required for interpolation.
+ this->generalized_support_points =
+ PolynomialsRaviartThomasNodal<dim>(degree).get_polynomial_support_points();
+ AssertDimension(this->generalized_support_points.size(),
+ this->n_dofs_per_cell());
- // TODO: the implementation makes the assumption that all faces have the
- // same number of dofs
- AssertDimension(this->n_unique_faces(), 1);
const unsigned int face_no = 0;
+ if (dim > 1)
+ this->generalized_face_support_points[face_no] =
+ degree == 0 ? QGauss<dim - 1>(1).get_points() :
+ QGaussLobatto<dim - 1>(degree + 1).get_points();
- // Fill prolongation matrices with embedding operators
- FETools::compute_embedding_matrices(*this, this->prolongation);
- // TODO[TL]: for anisotropic refinement we will probably need a table of
- // submatrices with an array for each refine case
FullMatrix<double> face_embeddings[GeometryInfo<dim>::max_children_per_face];
for (unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
face_embeddings[i].reinit(this->n_dofs_per_face(face_no),
face_embeddings,
0,
0);
- this->interface_constraints.reinit((1 << (dim - 1)) *
+ this->interface_constraints.reinit(GeometryInfo<dim>::max_children_per_face *
this->n_dofs_per_face(face_no),
this->n_dofs_per_face(face_no));
unsigned int target_row = 0;
std::string
FE_RaviartThomasNodal<dim>::get_name() const
{
- // note that the
- // FETools::get_fe_by_name
- // function depends on the
- // particular format of the string
- // this function returns, so they
- // have to be kept in synch
-
- // note that this->degree is the maximal
- // polynomial degree and is thus one higher
- // than the argument given to the
- // constructor
- std::ostringstream namebuf;
- namebuf << "FE_RaviartThomasNodal<" << dim << ">(" << this->degree - 1 << ")";
-
- return namebuf.str();
+ // note that the FETools::get_fe_by_name function depends on the particular
+ // format of the string this function returns, so they have to be kept in
+ // synch
+
+ // note that this->degree is the maximal polynomial degree and is thus one
+ // higher than the argument given to the constructor
+ return "FE_RaviartThomasNodal<" + std::to_string(dim) + ">(" +
+ std::to_string(this->degree - 1) + ")";
}
-template <int dim>
-void
-FE_RaviartThomasNodal<dim>::initialize_support_points(const unsigned int deg)
-{
- // TODO: the implementation makes the assumption that all faces have the
- // same number of dofs
- AssertDimension(this->n_unique_faces(), 1);
- const unsigned int face_no = 0;
-
- this->generalized_support_points.resize(this->n_dofs_per_cell());
- this->generalized_face_support_points[face_no].resize(
- this->n_dofs_per_face(face_no));
-
- // Number of the point being entered
- unsigned int current = 0;
-
- // On the faces, we choose as many
- // Gauss points as necessary to
- // determine the normal component
- // uniquely. This is the deg of
- // the Raviart-Thomas element plus
- // one.
- if (dim > 1)
- {
- QGauss<dim - 1> face_points(deg + 1);
- Assert(face_points.size() == this->n_dofs_per_face(face_no),
- ExcInternalError());
- for (unsigned int k = 0; k < this->n_dofs_per_face(face_no); ++k)
- this->generalized_face_support_points[face_no][k] =
- face_points.point(k);
- Quadrature<dim> faces =
- QProjector<dim>::project_to_all_faces(this->reference_cell(),
- face_points);
- for (unsigned int k = 0; k < this->n_dofs_per_face(face_no) *
- GeometryInfo<dim>::faces_per_cell;
- ++k)
- this->generalized_support_points[k] = faces.point(
- k + QProjector<dim>::DataSetDescriptor::face(this->reference_cell(),
- 0,
- true,
- false,
- false,
- this->n_dofs_per_face(
- face_no)));
-
- current =
- this->n_dofs_per_face(face_no) * GeometryInfo<dim>::faces_per_cell;
- }
-
- if (deg == 0)
- return;
- // In the interior, we need
- // anisotropic Gauss quadratures,
- // different for each direction.
- QGauss<1> high(deg + 1);
- QGauss<1> low(deg);
-
- for (unsigned int d = 0; d < dim; ++d)
- {
- std::unique_ptr<QAnisotropic<dim>> quadrature;
- switch (dim)
- {
- case 1:
- quadrature = std::make_unique<QAnisotropic<dim>>(high);
- break;
- case 2:
- quadrature =
- std::make_unique<QAnisotropic<dim>>(((d == 0) ? low : high),
- ((d == 1) ? low : high));
- break;
- case 3:
- quadrature =
- std::make_unique<QAnisotropic<dim>>(((d == 0) ? low : high),
- ((d == 1) ? low : high),
- ((d == 2) ? low : high));
- break;
- default:
- Assert(false, ExcNotImplemented());
- }
-
- for (unsigned int k = 0; k < quadrature->size(); ++k)
- this->generalized_support_points[current++] = quadrature->point(k);
- }
- Assert(current == this->n_dofs_per_cell(), ExcInternalError());
-}
-
-
-
template <int dim>
void
FE_RaviartThomasNodal<
if (dim < 3)
return;
- // TODO: Implement this for this class
- return;
-}
-
-
-
-template <int dim>
-std::vector<unsigned int>
-FE_RaviartThomasNodal<dim>::get_dpo_vector(const unsigned int deg)
-{
- // the element is face-based and we have
- // (deg+1)^(dim-1) DoFs per face
- unsigned int dofs_per_face = 1;
- for (unsigned int d = 1; d < dim; ++d)
- dofs_per_face *= deg + 1;
-
- // and then there are interior dofs
- const unsigned int interior_dofs = dim * deg * dofs_per_face;
-
- std::vector<unsigned int> dpo(dim + 1);
- dpo[dim - 1] = dofs_per_face;
- dpo[dim] = interior_dofs;
-
- return dpo;
-}
-
-
-
-template <>
-std::vector<bool>
-FE_RaviartThomasNodal<1>::get_ria_vector(const unsigned int)
-{
- Assert(false, ExcImpossibleInDim(1));
- return std::vector<bool>();
+ const unsigned int n = this->degree;
+ const unsigned int face_no = 0;
+ Assert(n * n == this->n_dofs_per_quad(face_no), ExcInternalError());
+ for (unsigned int local = 0; local < this->n_dofs_per_quad(face_no); ++local)
+ // face support points are in lexicographic ordering with x running
+ // fastest. invert that (y running fastest)
+ {
+ unsigned int i = local % n, j = local / n;
+
+ // face_orientation=false, face_flip=false, face_rotation=false
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 0) =
+ j + i * n - local;
+ // face_orientation=false, face_flip=false, face_rotation=true
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 1) =
+ i + (n - 1 - j) * n - local;
+ // face_orientation=false, face_flip=true, face_rotation=false
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 2) =
+ (n - 1 - j) + (n - 1 - i) * n - local;
+ // face_orientation=false, face_flip=true, face_rotation=true
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 3) =
+ (n - 1 - i) + j * n - local;
+ // face_orientation=true, face_flip=false, face_rotation=false
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 4) = 0;
+ // face_orientation=true, face_flip=false, face_rotation=true
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 5) =
+ j + (n - 1 - i) * n - local;
+ // face_orientation=true, face_flip=true, face_rotation=false
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 6) =
+ (n - 1 - i) + (n - 1 - j) * n - local;
+ // face_orientation=true, face_flip=true, face_rotation=true
+ this->adjust_quad_dof_index_for_face_orientation_table[face_no](local,
+ 7) =
+ (n - 1 - j) + i * n - local;
+
+ // for face_orientation == false, we need to switch the sign
+ for (unsigned int i = 0; i < 4; ++i)
+ this->adjust_quad_dof_sign_for_face_orientation_table[face_no](local,
+ i) = 1;
+ }
}
-template <int dim>
-std::vector<bool>
-FE_RaviartThomasNodal<dim>::get_ria_vector(const unsigned int deg)
-{
- const unsigned int dofs_per_cell =
- PolynomialsRaviartThomas<dim>::n_polynomials(deg);
- unsigned int dofs_per_face = deg + 1;
- for (unsigned int d = 2; d < dim; ++d)
- dofs_per_face *= deg + 1;
- // all face dofs need to be
- // non-additive, since they have
- // continuity requirements.
- // however, the interior dofs are
- // made additive
- std::vector<bool> ret_val(dofs_per_cell, false);
- for (unsigned int i = GeometryInfo<dim>::faces_per_cell * dofs_per_face;
- i < dofs_per_cell;
- ++i)
- ret_val[i] = true;
-
- return ret_val;
-}
-
-
template <int dim>
bool
FE_RaviartThomasNodal<dim>::has_support_on_face(
AssertIndexRange(shape_index, this->n_dofs_per_cell());
AssertIndexRange(face_index, GeometryInfo<dim>::faces_per_cell);
- // The first degrees of freedom are
- // on the faces and each face has
- // degree degrees.
- const unsigned int support_face = shape_index / this->degree;
+ // The first degrees of freedom are on the faces and each face has degree
+ // degrees.
+ const unsigned int support_face = shape_index / this->n_dofs_per_face();
- // The only thing we know for sure
- // is that shape functions with
- // support on one face are zero on
- // the opposite face.
+ // The only thing we know for sure is that shape functions with support on
+ // one face are zero on the opposite face.
if (support_face < GeometryInfo<dim>::faces_per_cell)
return (face_index != GeometryInfo<dim>::opposite_face[support_face]);
- // In all other cases, return true,
- // which is safe
+ // In all other cases, return true, which is safe
return true;
}
ExcDimensionMismatch(support_point_values[0].size(),
this->n_components()));
- // First do interpolation on
- // faces. There, the component
- // evaluated depends on the face
- // direction and orientation.
+ // First do interpolation on faces. There, the component evaluated depends
+ // on the face direction and orientation.
unsigned int fbase = 0;
unsigned int f = 0;
for (; f < GeometryInfo<dim>::faces_per_cell;
}
}
- // The remaining points form dim
- // chunks, one for each component.
+ // The remaining points form dim chunks, one for each component
const unsigned int istep = (this->n_dofs_per_cell() - fbase) / dim;
Assert((this->n_dofs_per_cell() - fbase) % dim == 0, ExcInternalError());
FE_RaviartThomasNodal<dim>::hp_vertex_dof_identities(
const FiniteElement<dim> &fe_other) const
{
- // we can presently only compute these
- // identities if both FEs are
- // FE_RaviartThomasNodals or the other is FE_Nothing.
- // In either case, no dofs are assigned on the vertex,
- // so we shouldn't be getting here at all.
+ // we can presently only compute these identities if both FEs are
+ // FE_RaviartThomasNodals or the other is FE_Nothing. In either case, no
+ // dofs are assigned on the vertex, so we shouldn't be getting here at all.
if (dynamic_cast<const FE_RaviartThomasNodal<dim> *>(&fe_other) != nullptr)
return std::vector<std::pair<unsigned int, unsigned int>>();
else if (dynamic_cast<const FE_Nothing<dim> *>(&fe_other) != nullptr)
FE_RaviartThomasNodal<dim>::hp_line_dof_identities(
const FiniteElement<dim> &fe_other) const
{
- // we can presently only compute
- // these identities if both FEs are
- // FE_RaviartThomasNodals or if the other
- // one is FE_Nothing
+ // we can presently only compute these identities if both FEs are
+ // FE_RaviartThomasNodals or if the other one is FE_Nothing
if (const FE_RaviartThomasNodal<dim> *fe_q_other =
dynamic_cast<const FE_RaviartThomasNodal<dim> *>(&fe_other))
{
- // dofs are located on faces; these are
- // only lines in 2d
+ // dofs are located on faces; these are only lines in 2d
if (dim != 2)
return std::vector<std::pair<unsigned int, unsigned int>>();
- // dofs are located along lines, so two
- // dofs are identical only if in the
- // following two cases (remember that
- // the face support points are Gauss
- // points):
+ // dofs are located along lines, so two dofs are identical only if in
+ // the following two cases (remember that the face support points are
+ // Gauss points):
// 1. this->degree = fe_q_other->degree,
- // in the case, all the dofs on
- // the line are identical
+ // in the case, all the dofs on the line are identical
// 2. this->degree-1 and fe_q_other->degree-1
- // are both even, i.e. this->dof_per_line
- // and fe_q_other->dof_per_line are both odd,
- // there exists only one point (the middle one)
- // such that dofs are identical on this point
+ // are both even, i.e. this->dof_per_line and fe_q_other->dof_per_line
+ // are both odd, there exists only one point (the middle one) such
+ // that dofs are identical on this point
//
- // to understand this, note that
- // this->degree is the *maximal*
- // polynomial degree, and is thus one
- // higher than the argument given to
+ // to understand this, note that this->degree is the *maximal*
+ // polynomial degree, and is thus one higher than the argument given to
// the constructor
const unsigned int p = this->degree - 1;
const unsigned int q = fe_q_other->degree - 1;
const FiniteElement<dim> &fe_other,
const unsigned int face_no) const
{
- // we can presently only compute
- // these identities if both FEs are
- // FE_RaviartThomasNodals or if the other
- // one is FE_Nothing
+ // we can presently only compute these identities if both FEs are
+ // FE_RaviartThomasNodals or if the other one is FE_Nothing
if (const FE_RaviartThomasNodal<dim> *fe_q_other =
dynamic_cast<const FE_RaviartThomasNodal<dim> *>(&fe_other))
{
- // dofs are located on faces; these are
- // only quads in 3d
+ // dofs are located on faces; these are only quads in 3d
if (dim != 3)
return std::vector<std::pair<unsigned int, unsigned int>>();
- // this works exactly like the line
- // case above
+ // this works exactly like the line case above
const unsigned int p = this->n_dofs_per_quad(face_no);
AssertDimension(fe_q_other->n_unique_faces(), 1);
ExcDimensionMismatch(interpolation_matrix.m(),
x_source_fe.n_dofs_per_face(face_no)));
- // ok, source is a RaviartThomasNodal element, so
- // we will be able to do the work
+ // ok, source is a RaviartThomasNodal element, so we will be able to do the
+ // work
const FE_RaviartThomasNodal<dim> &source_fe =
dynamic_cast<const FE_RaviartThomasNodal<dim> &>(x_source_fe);
- // Make sure, that the element,
- // for which the DoFs should be
- // constrained is the one with
- // the higher polynomial degree.
- // Actually the procedure will work
- // also if this assertion is not
- // satisfied. But the matrices
- // produced in that case might
- // lead to problems in the
- // hp-procedures, which use this
+ // Make sure that the element for which the DoFs should be constrained is
+ // the one with the higher polynomial degree. Actually the procedure will
+ // work also if this assertion is not satisfied. But the matrices produced
+ // in that case might lead to problems in the hp-procedures, which use this
// method.
Assert(this->n_dofs_per_face(face_no) <= source_fe.n_dofs_per_face(face_no),
typename FiniteElement<dim>::ExcInterpolationNotImplemented());
- // generate a quadrature
- // with the generalized support points.
- // This is later based as a
- // basis for the QProjector,
- // which returns the support
- // points on the face.
+ // generate a quadrature with the generalized support points. This is later
+ // based as a basis for the QProjector, which returns the support points on
+ // the face.
Quadrature<dim - 1> quad_face_support(
source_fe.generalized_face_support_points[face_no]);
- // Rule of thumb for FP accuracy,
- // that can be expected for a
- // given polynomial degree.
- // This value is used to cut
- // off values close to zero.
+ // Rule of thumb for FP accuracy, that can be expected for a given
+ // polynomial degree. This value is used to cut off values close to zero.
double eps = 2e-13 * this->degree * (dim - 1);
- // compute the interpolation
- // matrix by simply taking the
- // value at the support points.
+ // compute the interpolation matrix by simply taking the value at the
+ // support points.
const Quadrature<dim> face_projection =
QProjector<dim>::project_to_face(this->reference_cell(),
quad_face_support,
double matrix_entry =
this->shape_value_component(this->face_to_cell_index(j, 0), p, 0);
- // Correct the interpolated
- // value. I.e. if it is close
- // to 1 or 0, make it exactly
- // 1 or 0. Unfortunately, this
- // is required to avoid problems
- // with higher order elements.
+ // Correct the interpolated value. I.e. if it is close to 1 or 0,
+ // make it exactly 1 or 0. Unfortunately, this is required to avoid
+ // problems with higher order elements.
if (std::fabs(matrix_entry - 1.0) < eps)
matrix_entry = 1.0;
if (std::fabs(matrix_entry) < eps)
}
}
- // make sure that the row sum of
- // each of the matrices is 1 at
- // this point. this must be so
- // since the shape functions sum up
- // to 1
+ // make sure that the row sum of each of the matrices is 1 at this
+ // point. this must be so since the shape functions sum up to 1
for (unsigned int j = 0; j < source_fe.n_dofs_per_face(face_no); ++j)
{
double sum = 0.;
FullMatrix<double> & interpolation_matrix,
const unsigned int face_no) const
{
- // this is only implemented, if the
- // source FE is also a
- // RaviartThomasNodal element
+ // this is only implemented, if the source FE is also a RaviartThomasNodal
+ // element
AssertThrow((x_source_fe.get_name().find("FE_RaviartThomasNodal<") == 0) ||
(dynamic_cast<const FE_RaviartThomasNodal<dim> *>(
&x_source_fe) != nullptr),
ExcDimensionMismatch(interpolation_matrix.m(),
x_source_fe.n_dofs_per_face(face_no)));
- // ok, source is a RaviartThomasNodal element, so
- // we will be able to do the work
+ // ok, source is a RaviartThomasNodal element, so we will be able to do the
+ // work
const FE_RaviartThomasNodal<dim> &source_fe =
dynamic_cast<const FE_RaviartThomasNodal<dim> &>(x_source_fe);
- // Make sure, that the element,
- // for which the DoFs should be
- // constrained is the one with
- // the higher polynomial degree.
- // Actually the procedure will work
- // also if this assertion is not
- // satisfied. But the matrices
- // produced in that case might
- // lead to problems in the
- // hp-procedures, which use this
+ // Make sure that the element for which the DoFs should be constrained is
+ // the one with the higher polynomial degree. Actually the procedure will
+ // work also if this assertion is not satisfied. But the matrices produced
+ // in that case might lead to problems in the hp-procedures, which use this
// method.
Assert(this->n_dofs_per_face(face_no) <= source_fe.n_dofs_per_face(face_no),
typename FiniteElement<dim>::ExcInterpolationNotImplemented());
- // generate a quadrature
- // with the generalized support points.
- // This is later based as a
- // basis for the QProjector,
- // which returns the support
- // points on the face.
+ // generate a quadrature with the generalized support points. This is later
+ // based as a basis for the QProjector, which returns the support points on
+ // the face.
Quadrature<dim - 1> quad_face_support(
source_fe.generalized_face_support_points[face_no]);
- // Rule of thumb for FP accuracy,
- // that can be expected for a
- // given polynomial degree.
- // This value is used to cut
- // off values close to zero.
+ // Rule of thumb for FP accuracy, that can be expected for a given
+ // polynomial degree. This value is used to cut off values close to zero.
double eps = 2e-13 * this->degree * (dim - 1);
- // compute the interpolation
- // matrix by simply taking the
- // value at the support points.
-
+ // compute the interpolation matrix by simply taking the value at the
+ // support points.
const Quadrature<dim> subface_projection =
QProjector<dim>::project_to_subface(this->reference_cell(),
quad_face_support,
double matrix_entry =
this->shape_value_component(this->face_to_cell_index(j, 0), p, 0);
- // Correct the interpolated
- // value. I.e. if it is close
- // to 1 or 0, make it exactly
- // 1 or 0. Unfortunately, this
- // is required to avoid problems
- // with higher order elements.
+ // Correct the interpolated value. I.e. if it is close to 1 or 0,
+ // make it exactly 1 or 0. Unfortunately, this is required to avoid
+ // problems with higher order elements.
if (std::fabs(matrix_entry - 1.0) < eps)
matrix_entry = 1.0;
if (std::fabs(matrix_entry) < eps)
}
}
- // make sure that the row sum of
- // each of the matrices is 1 at
- // this point. this must be so
- // since the shape functions sum up
- // to 1
+ // make sure that the row sum of each of the matrices is 1 at this
+ // point. this must be so since the shape functions sum up to 1
for (unsigned int j = 0; j < source_fe.n_dofs_per_face(face_no); ++j)
{
double sum = 0.;
+template <int dim>
+const FullMatrix<double> &
+FE_RaviartThomasNodal<dim>::get_prolongation_matrix(
+ const unsigned int child,
+ const RefinementCase<dim> &refinement_case) const
+{
+ AssertIndexRange(refinement_case,
+ RefinementCase<dim>::isotropic_refinement + 1);
+ Assert(refinement_case != RefinementCase<dim>::no_refinement,
+ ExcMessage(
+ "Prolongation matrices are only available for refined cells!"));
+ AssertIndexRange(child, GeometryInfo<dim>::n_children(refinement_case));
+
+ // initialization upon first request
+ if (this->prolongation[refinement_case - 1][child].n() == 0)
+ {
+ std::lock_guard<std::mutex> lock(this->mutex);
+
+ // if matrix got updated while waiting for the lock
+ if (this->prolongation[refinement_case - 1][child].n() ==
+ this->n_dofs_per_cell())
+ return this->prolongation[refinement_case - 1][child];
+
+ // now do the work. need to get a non-const version of data in order to
+ // be able to modify them inside a const function
+ FE_RaviartThomasNodal<dim> &this_nonconst =
+ const_cast<FE_RaviartThomasNodal<dim> &>(*this);
+ if (refinement_case == RefinementCase<dim>::isotropic_refinement)
+ {
+ std::vector<std::vector<FullMatrix<double>>> isotropic_matrices(
+ RefinementCase<dim>::isotropic_refinement);
+ isotropic_matrices.back().resize(
+ GeometryInfo<dim>::n_children(RefinementCase<dim>(refinement_case)),
+ FullMatrix<double>(this->n_dofs_per_cell(),
+ this->n_dofs_per_cell()));
+ FETools::compute_embedding_matrices(*this, isotropic_matrices, true);
+ this_nonconst.prolongation[refinement_case - 1].swap(
+ isotropic_matrices.back());
+ }
+ else
+ {
+ // must compute both restriction and prolongation matrices because
+ // we only check for their size and the reinit call initializes them
+ // all
+ this_nonconst.reinit_restriction_and_prolongation_matrices();
+ FETools::compute_embedding_matrices(*this,
+ this_nonconst.prolongation);
+ FETools::compute_projection_matrices(*this,
+ this_nonconst.restriction);
+ }
+ }
+
+ // finally return the matrix
+ return this->prolongation[refinement_case - 1][child];
+}
+
+
+
+template <int dim>
+const FullMatrix<double> &
+FE_RaviartThomasNodal<dim>::get_restriction_matrix(
+ const unsigned int child,
+ const RefinementCase<dim> &refinement_case) const
+{
+ AssertIndexRange(refinement_case,
+ RefinementCase<dim>::isotropic_refinement + 1);
+ Assert(refinement_case != RefinementCase<dim>::no_refinement,
+ ExcMessage(
+ "Restriction matrices are only available for refined cells!"));
+ AssertIndexRange(child, GeometryInfo<dim>::n_children(refinement_case));
+
+ // initialization upon first request
+ if (this->restriction[refinement_case - 1][child].n() == 0)
+ {
+ std::lock_guard<std::mutex> lock(this->mutex);
+
+ // if matrix got updated while waiting for the lock...
+ if (this->restriction[refinement_case - 1][child].n() ==
+ this->n_dofs_per_cell())
+ return this->restriction[refinement_case - 1][child];
+
+ // now do the work. need to get a non-const version of data in order to
+ // be able to modify them inside a const function
+ FE_RaviartThomasNodal<dim> &this_nonconst =
+ const_cast<FE_RaviartThomasNodal<dim> &>(*this);
+ if (refinement_case == RefinementCase<dim>::isotropic_refinement)
+ {
+ std::vector<std::vector<FullMatrix<double>>> isotropic_matrices(
+ RefinementCase<dim>::isotropic_refinement);
+ isotropic_matrices.back().resize(
+ GeometryInfo<dim>::n_children(RefinementCase<dim>(refinement_case)),
+ FullMatrix<double>(this->n_dofs_per_cell(),
+ this->n_dofs_per_cell()));
+ FETools::compute_projection_matrices(*this, isotropic_matrices, true);
+ this_nonconst.restriction[refinement_case - 1].swap(
+ isotropic_matrices.back());
+ }
+ else
+ {
+ // must compute both restriction and prolongation matrices because
+ // we only check for their size and the reinit call initializes them
+ // all
+ this_nonconst.reinit_restriction_and_prolongation_matrices();
+ FETools::compute_embedding_matrices(*this,
+ this_nonconst.prolongation);
+ FETools::compute_projection_matrices(*this,
+ this_nonconst.restriction);
+ }
+ }
+
+ // finally return the matrix
+ return this->restriction[refinement_case - 1][child];
+}
+
+
+
// explicit instantiations
#include "fe_raviart_thomas_nodal.inst"
--- /dev/null
+
+DEAL::
+DEAL::******* degree 0 ******* orientation case 0 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 0 ******* orientation case 1 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 0 ******* orientation case 2 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 0 ******* orientation case 3 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 0 ******* orientation case 4 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 0 ******* orientation case 5 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 0 ******* orientation case 6 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 0 ******* orientation case 7 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(0)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 6
+DEAL:: n_dofs_per_face : 1
+DEAL:: n_dofs_per_quad : 1
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 0 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 1 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 2 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 3 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 4 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 5 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 6 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 1 ******* orientation case 7 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(1)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 36
+DEAL:: n_dofs_per_face : 4
+DEAL:: n_dofs_per_quad : 4
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 0 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 1 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 2 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 3 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 4 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 5 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 6 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::
+DEAL::******* degree 2 ******* orientation case 7 *******
+DEAL::Element Info:
+DEAL:: name : FE_RaviartThomasNodal<3>(2)
+DEAL:: is_primitive : 0
+DEAL:: n_dofs_per_cell : 108
+DEAL:: n_dofs_per_face : 9
+DEAL:: n_dofs_per_quad : 9
+DEAL:: n_dofs_per_line : 0
+DEAL:: n_dofs_per_vertex : 0
+DEAL::
+DEAL:: first_line_index : 0
+DEAL:: first_quad_index : 0
+DEAL:: first_face_line_index: 0
+DEAL:: first_face_quad_index: 0
+DEAL::
+DEAL:: n_components : 3
+DEAL:: n_blocks : 1
+DEAL:: n_base_elements : 1
+DEAL::
+DEAL::Normal jumps (at quad points) in cell 0_0: at face 1
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+DEAL::Normal jumps (at quad points) in cell 1_0: at face 0
+DEAL:: interface jumps: 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
projects / dealii.git / commitdiff