49 std::vector<
bool>(dim, true)))
76 ref_case < RefinementCase<dim>::isotropic_refinement + 1;
79 const unsigned int nc =
82 for (
unsigned int i = 0; i < nc; ++i)
90 for (
unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
92 FETools::compute_face_embedding_matrices<dim, double>(*
this,
98 unsigned int target_row = 0;
99 for (
unsigned int d = 0; d < GeometryInfo<dim>::max_children_per_face; ++
d)
100 for (
unsigned int i = 0; i < face_embeddings[
d].
m(); ++i)
102 for (
unsigned int j = 0; j < face_embeddings[
d].
n(); ++j)
125 std::ostringstream namebuf;
126 namebuf <<
"FE_RaviartThomasNodal<" << dim <<
">(" << this->degree - 1 <<
")";
128 return namebuf.str();
133 std::unique_ptr<FiniteElement<dim, dim>>
136 return std_cxx14::make_unique<FE_RaviartThomasNodal<dim>>(*this);
150 this->generalized_support_points.resize(this->dofs_per_cell);
151 this->generalized_face_support_points.resize(this->dofs_per_face);
154 unsigned int current = 0;
164 QGauss<dim - 1> face_points(deg + 1);
166 for (
unsigned int k = 0; k < this->dofs_per_face; ++k)
167 this->generalized_face_support_points[k] = face_points.point(k);
170 for (
unsigned int k = 0;
173 this->generalized_support_points[k] =
175 0,
true,
false,
false, this->dofs_per_face));
188 for (
unsigned int d = 0;
d < dim; ++
d)
190 std::unique_ptr<QAnisotropic<dim>> quadrature;
194 quadrature = std_cxx14::make_unique<QAnisotropic<dim>>(high);
197 quadrature = std_cxx14::make_unique<QAnisotropic<dim>>(
198 ((
d == 0) ? low : high), ((
d == 1) ? low : high));
202 std_cxx14::make_unique<QAnisotropic<dim>>(((
d == 0) ? low : high),
203 ((
d == 1) ? low : high),
211 for (
unsigned int k = 0; k < quadrature->size(); ++k)
212 this->generalized_support_points[current++] = quadrature->point(k);
220 std::vector<unsigned int>
225 unsigned int dofs_per_face = 1;
226 for (
unsigned int d = 1;
d < dim; ++
d)
227 dofs_per_face *= deg + 1;
230 const unsigned int interior_dofs = dim * deg * dofs_per_face;
232 std::vector<unsigned int> dpo(dim + 1);
233 dpo[dim - 1] = dofs_per_face;
234 dpo[dim] = interior_dofs;
246 return std::vector<bool>();
255 const unsigned int dofs_per_cell =
257 unsigned int dofs_per_face = deg + 1;
258 for (
unsigned int d = 2;
d < dim; ++
d)
259 dofs_per_face *= deg + 1;
265 std::vector<bool> ret_val(dofs_per_cell,
false);
278 const unsigned int shape_index,
279 const unsigned int face_index)
const
287 const unsigned int support_face = shape_index / this->degree;
308 std::vector<double> & nodal_values)
const
310 Assert(support_point_values.size() == this->generalized_support_points.size(),
312 this->generalized_support_points.size()));
313 Assert(nodal_values.size() == this->dofs_per_cell,
323 unsigned int fbase = 0;
325 for (; f < GeometryInfo<dim>::faces_per_cell;
326 ++f, fbase += this->dofs_per_face)
328 for (
unsigned int i = 0; i < this->dofs_per_face; ++i)
330 nodal_values[fbase + i] = support_point_values[fbase + i](
337 const unsigned int istep = (this->dofs_per_cell - fbase) / dim;
341 while (fbase < this->dofs_per_cell)
343 for (
unsigned int i = 0; i < istep; ++i)
345 nodal_values[fbase + i] = support_point_values[fbase + i](f);
370 std::vector<std::pair<unsigned int, unsigned int>>
380 return std::vector<std::pair<unsigned int, unsigned int>>();
382 return std::vector<std::pair<unsigned int, unsigned int>>();
386 return std::vector<std::pair<unsigned int, unsigned int>>();
393 std::vector<std::pair<unsigned int, unsigned int>>
407 return std::vector<std::pair<unsigned int, unsigned int>>();
428 const unsigned int p = this->degree - 1;
429 const unsigned int q = fe_q_other->degree - 1;
431 std::vector<std::pair<unsigned int, unsigned int>> identities;
434 for (
unsigned int i = 0; i < p + 1; ++i)
435 identities.emplace_back(i, i);
437 else if (p % 2 == 0 && q % 2 == 0)
438 identities.emplace_back(p / 2, q / 2);
446 return std::vector<std::pair<unsigned int, unsigned int>>();
451 return std::vector<std::pair<unsigned int, unsigned int>>();
457 std::vector<std::pair<unsigned int, unsigned int>>
471 return std::vector<std::pair<unsigned int, unsigned int>>();
475 const unsigned int p = this->dofs_per_quad;
476 const unsigned int q = fe_q_other->dofs_per_quad;
478 std::vector<std::pair<unsigned int, unsigned int>> identities;
481 for (
unsigned int i = 0; i < p; ++i)
482 identities.emplace_back(i, i);
484 else if (p % 2 != 0 && q % 2 != 0)
485 identities.emplace_back(p / 2, q / 2);
493 return std::vector<std::pair<unsigned int, unsigned int>>();
498 return std::vector<std::pair<unsigned int, unsigned int>>();
507 const unsigned int codim)
const
517 if (this->degree < fe_rt_nodal_other->degree)
519 else if (this->degree == fe_rt_nodal_other->degree)
527 if (fe_nothing->is_dominating())
575 &x_source_fe) !=
nullptr),
578 Assert(interpolation_matrix.
n() == this->dofs_per_face,
617 double eps = 2
e-13 * this->degree * (dim - 1);
627 const Point<dim> &p = face_projection.point(i);
629 for (
unsigned int j = 0; j < this->dofs_per_face; ++j)
631 double matrix_entry =
632 this->shape_value_component(this->face_to_cell_index(j, 0), p, 0);
645 interpolation_matrix(i, j) = matrix_entry;
658 for (
unsigned int i = 0; i < this->dofs_per_face; ++i)
659 sum += interpolation_matrix(j, i);
671 const unsigned int subface,
679 &x_source_fe) !=
nullptr),
682 Assert(interpolation_matrix.
n() == this->dofs_per_face,
721 double eps = 2
e-13 * this->degree * (dim - 1);
732 const Point<dim> &p = subface_projection.point(i);
734 for (
unsigned int j = 0; j < this->dofs_per_face; ++j)
736 double matrix_entry =
737 this->shape_value_component(this->face_to_cell_index(j, 0), p, 0);
750 interpolation_matrix(i, j) = matrix_entry;
763 for (
unsigned int i = 0; i < this->dofs_per_face; ++i)
764 sum += interpolation_matrix(j, i);
774 #include "fe_raviart_thomas_nodal.inst"