20#include <deal.II/fe/fe_poly_face.templates.h>
33 namespace FE_FaceQImplementation
38 get_QGaussLobatto_points(
const unsigned int degree)
43 return std::vector<Point<1>>(1,
Point<1>(0.5));
33 namespace FE_FaceQImplementation {
…}
49template <
int dim,
int spacedim>
54 internal::FE_FaceQImplementation::get_QGaussLobatto_points(degree))),
62 const unsigned int codim = dim - 1;
64 Utilities::fixed_power<codim>(this->degree + 1));
66 if (this->degree == 0)
67 for (
unsigned int d = 0; d < codim; ++d)
71 std::vector<Point<1>> points =
72 internal::FE_FaceQImplementation::get_QGaussLobatto_points(
degree);
75 for (
unsigned int iz = 0; iz <= ((codim > 2) ? this->degree : 0); ++iz)
76 for (
unsigned int iy = 0; iy <= ((codim > 1) ? this->degree : 0); ++iy)
77 for (
unsigned int ix = 0; ix <= this->
degree; ++ix)
97 for (
unsigned int i = 0; i < n_face_dofs; ++i)
98 for (
unsigned int d = 0; d < dim; ++d)
100 for (
unsigned int e = 0, c = 0; e < dim; ++e)
104 unsigned int renumber = i;
105 if (dim == 3 && d == 1)
119template <
int dim,
int spacedim>
120std::unique_ptr<FiniteElement<dim, spacedim>>
123 return std::make_unique<FE_FaceQ<dim, spacedim>>(this->degree);
128template <
int dim,
int spacedim>
135 std::ostringstream namebuf;
137 << this->degree <<
")";
139 return namebuf.str();
144template <
int dim,
int spacedim>
149 const unsigned int face_no)
const
151 get_subface_interpolation_matrix(source_fe,
153 interpolation_matrix,
159template <
int dim,
int spacedim>
163 const unsigned int subface,
165 const unsigned int face_no)
const
169 Assert(interpolation_matrix.
n() == this->n_dofs_per_face(face_no),
171 this->n_dofs_per_face(face_no)));
186 this->n_dofs_per_face(face_no) <= source_fe->n_dofs_per_face(face_no),
188 spacedim>::ExcInterpolationNotImplemented()));
192 source_fe->get_unit_face_support_points(face_no));
197 const double eps = 2e-13 * (this->degree + 1) * (dim - 1);
201 for (
unsigned int i = 0; i < source_fe->n_dofs_per_face(face_no); ++i)
203 const Point<dim - 1> p =
205 face_quadrature.point(i) :
207 face_quadrature.point(i), subface);
209 for (
unsigned int j = 0; j < this->n_dofs_per_face(face_no); ++j)
211 double matrix_entry = this->poly_space.compute_value(j, p);
216 if (std::fabs(matrix_entry - 1.0) < eps)
218 if (std::fabs(matrix_entry) < eps)
221 interpolation_matrix(i, j) = matrix_entry;
229 for (
unsigned int j = 0; j < source_fe->n_dofs_per_face(face_no); ++j)
233 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
234 sum += interpolation_matrix(j, i);
240 else if (
dynamic_cast<const FE_Nothing<dim> *
>(&x_source_fe) !=
nullptr)
248 spacedim>::ExcInterpolationNotImplemented()));
253template <
int dim,
int spacedim>
256 const unsigned int shape_index,
257 const unsigned int face_index)
const
259 return (face_index == (shape_index / this->n_dofs_per_face(face_index)));
264template <
int dim,
int spacedim>
265std::vector<unsigned int>
268 std::vector<unsigned int> dpo(dim + 1, 0U);
269 dpo[dim - 1] = deg + 1;
270 for (
unsigned int i = 1; i < dim - 1; ++i)
271 dpo[dim - 1] *= deg + 1;
277template <
int dim,
int spacedim>
286template <
int dim,
int spacedim>
287std::vector<std::pair<unsigned int, unsigned int>>
292 return std::vector<std::pair<unsigned int, unsigned int>>();
297template <
int dim,
int spacedim>
298std::vector<std::pair<unsigned int, unsigned int>>
306 return std::vector<std::pair<unsigned int, unsigned int>>();
319 const unsigned int p = this->degree;
320 const unsigned int q = fe_q_other->degree;
322 std::vector<std::pair<unsigned int, unsigned int>> identities;
324 const std::vector<unsigned int> &index_map_inverse =
325 this->poly_space.get_numbering_inverse();
326 const std::vector<unsigned int> &index_map_inverse_other =
327 fe_q_other->poly_space.get_numbering_inverse();
329 for (
unsigned int i = 0; i < p + 1; ++i)
330 for (
unsigned int j = 0; j < q + 1; ++j)
332 this->unit_support_points[index_map_inverse[i]][dim - 1] -
333 fe_q_other->unit_support_points[index_map_inverse_other[j]]
335 identities.emplace_back(i, j);
343 return std::vector<std::pair<unsigned int, unsigned int>>();
355 return std::vector<std::pair<unsigned int, unsigned int>>();
360 return std::vector<std::pair<unsigned int, unsigned int>>();
367template <
int dim,
int spacedim>
368std::vector<std::pair<unsigned int, unsigned int>>
371 const unsigned int)
const
377 return std::vector<std::pair<unsigned int, unsigned int>>();
389 const unsigned int p = this->degree;
390 const unsigned int q = fe_q_other->degree;
392 std::vector<std::pair<unsigned int, unsigned int>> identities;
394 const std::vector<unsigned int> &index_map_inverse =
395 this->poly_space.get_numbering_inverse();
396 const std::vector<unsigned int> &index_map_inverse_other =
397 fe_q_other->poly_space.get_numbering_inverse();
399 std::vector<std::pair<unsigned int, unsigned int>> identities_1d;
401 for (
unsigned int i = 0; i < p + 1; ++i)
402 for (
unsigned int j = 0; j < q + 1; ++j)
404 this->unit_support_points[index_map_inverse[i]][dim - 2] -
405 fe_q_other->unit_support_points[index_map_inverse_other[j]]
407 identities_1d.emplace_back(i, j);
409 for (
unsigned int n1 = 0; n1 < identities_1d.size(); ++n1)
410 for (
unsigned int n2 = 0; n2 < identities_1d.size(); ++n2)
411 identities.emplace_back(identities_1d[n1].first * (p + 1) +
412 identities_1d[n2].first,
413 identities_1d[n1].second * (q + 1) +
414 identities_1d[n2].
second);
422 return std::vector<std::pair<unsigned int, unsigned int>>();
434 return std::vector<std::pair<unsigned int, unsigned int>>();
439 return std::vector<std::pair<unsigned int, unsigned int>>();
446template <
int dim,
int spacedim>
450 const unsigned int codim)
const
459 if (this->degree < fe_faceq_other->degree)
461 else if (this->degree == fe_faceq_other->degree)
469 if (fe_nothing->is_dominating())
482template <
int dim,
int spacedim>
483std::pair<Table<2, bool>, std::vector<unsigned int>>
487 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
488 constant_modes(0, i) =
true;
489 return std::pair<Table<2, bool>, std::vector<unsigned int>>(
490 constant_modes, std::vector<unsigned int>(1, 0));
493template <
int dim,
int spacedim>
497 std::vector<double> &nodal_values)
const
500 this->get_unit_support_points().size());
504 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
508 nodal_values[i] = support_point_values[i](0);
514template <
int spacedim>
534template <
int spacedim>
535std::unique_ptr<FiniteElement<1, spacedim>>
538 return std::make_unique<FE_FaceQ<1, spacedim>>(this->degree);
543template <
int spacedim>
550 std::ostringstream namebuf;
552 << this->degree <<
")";
554 return namebuf.str();
559template <
int spacedim>
564 const unsigned int face_no)
const
566 get_subface_interpolation_matrix(source_fe,
568 interpolation_matrix,
574template <
int spacedim>
580 const unsigned int face_no)
const
582 Assert(interpolation_matrix.
n() == this->n_dofs_per_face(face_no),
584 this->n_dofs_per_face(face_no)));
588 interpolation_matrix(0, 0) = 1.;
593template <
int spacedim>
596 const unsigned int face_index)
const
599 return (face_index == shape_index);
604template <
int spacedim>
605std::vector<unsigned int>
608 std::vector<unsigned int> dpo(2, 0U);
615template <
int spacedim>
622template <
int spacedim>
623std::vector<std::pair<unsigned int, unsigned int>>
628 return std::vector<std::pair<unsigned int, unsigned int>>(1,
635template <
int spacedim>
636std::vector<std::pair<unsigned int, unsigned int>>
641 return std::vector<std::pair<unsigned int, unsigned int>>();
646template <
int spacedim>
647std::vector<std::pair<unsigned int, unsigned int>>
650 const unsigned int)
const
653 return std::vector<std::pair<unsigned int, unsigned int>>();
658template <
int spacedim>
659std::pair<Table<2, bool>, std::vector<unsigned int>>
663 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
664 constant_modes(0, i) =
true;
665 return std::pair<Table<2, bool>, std::vector<unsigned int>>(
666 constant_modes, std::vector<unsigned int>(1, 0));
671template <
int spacedim>
687template <
int spacedim>
706template <
int spacedim>
710 const unsigned int face,
720 const unsigned int foffset = face;
723 for (
unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
724 output_data.shape_values(k, 0) = 0.;
725 output_data.shape_values(foffset, 0) = 1;
730template <
int spacedim>
752template <
int dim,
int spacedim>
756 Polynomials::Legendre::generate_complete_basis(degree)),
766template <
int dim,
int spacedim>
767std::unique_ptr<FiniteElement<dim, spacedim>>
770 return std::make_unique<FE_FaceP<dim, spacedim>>(this->degree);
775template <
int dim,
int spacedim>
782 std::ostringstream namebuf;
784 << this->degree <<
")";
786 return namebuf.str();
791template <
int dim,
int spacedim>
794 const unsigned int shape_index,
795 const unsigned int face_index)
const
797 return (face_index == (shape_index / this->n_dofs_per_face(face_index)));
802template <
int dim,
int spacedim>
803std::vector<unsigned int>
806 std::vector<unsigned int> dpo(dim + 1, 0U);
807 dpo[dim - 1] = deg + 1;
808 for (
unsigned int i = 1; i < dim - 1; ++i)
810 dpo[dim - 1] *= deg + 1 + i;
811 dpo[dim - 1] /= i + 1;
818template <
int dim,
int spacedim>
827template <
int dim,
int spacedim>
831 const unsigned int codim)
const
840 if (this->degree < fe_facep_other->degree)
842 else if (this->degree == fe_facep_other->degree)
850 if (fe_nothing->is_dominating())
865template <
int dim,
int spacedim>
870 const unsigned int face_no)
const
872 get_subface_interpolation_matrix(source_fe,
874 interpolation_matrix,
880template <
int dim,
int spacedim>
884 const unsigned int subface,
886 const unsigned int face_no)
const
890 Assert(interpolation_matrix.
n() == this->n_dofs_per_face(face_no),
892 this->n_dofs_per_face(face_no)));
907 this->n_dofs_per_face(face_no) <= source_fe->n_dofs_per_face(face_no),
909 spacedim>::ExcInterpolationNotImplemented()));
914 const QGauss<dim - 1> face_quadrature(source_fe->degree + 1);
919 const double eps = 2e-13 * (this->degree + 1) * (dim - 1);
922 source_fe->n_dofs_per_face(face_no));
924 for (
unsigned int k = 0; k < face_quadrature.size(); ++k)
926 const Point<dim - 1> p =
928 face_quadrature.point(k) :
930 face_quadrature.point(k), subface);
932 for (
unsigned int j = 0; j < source_fe->n_dofs_per_face(face_no); ++j)
933 mass(k, j) = source_fe->poly_space.compute_value(j, p);
943 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
945 for (
unsigned int k = 0; k < face_quadrature.size(); ++k)
947 const Point<dim - 1> p =
949 face_quadrature.point(k) :
951 face_quadrature.point(k), subface);
952 v_in(k) = this->poly_space.compute_value(i, p);
954 [[maybe_unused]]
const double result = H.
least_squares(v_out, v_in);
957 for (
unsigned int j = 0; j < source_fe->n_dofs_per_face(face_no); ++j)
959 double matrix_entry = v_out(j);
964 if (std::fabs(matrix_entry - 1.0) < eps)
966 if (std::fabs(matrix_entry) < eps)
969 interpolation_matrix(j, i) = matrix_entry;
973 else if (
dynamic_cast<const FE_Nothing<dim> *
>(&x_source_fe) !=
nullptr)
981 spacedim>::ExcInterpolationNotImplemented()));
986template <
int dim,
int spacedim>
987std::pair<Table<2, bool>, std::vector<unsigned int>>
992 constant_modes(0, face * this->n_dofs_per_face(face)) =
true;
993 return std::pair<Table<2, bool>, std::vector<unsigned int>>(
994 constant_modes, std::vector<unsigned int>(1, 0));
999template <
int spacedim>
1006template <
int spacedim>
1013 std::ostringstream namebuf;
1015 << this->degree <<
")";
1017 return namebuf.str();
1023#include "fe/fe_face.inst"
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
virtual bool hp_constraints_are_implemented() const override
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
virtual std::string get_name() const override
static std::vector< unsigned int > get_dpo_vector(const unsigned int deg)
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
virtual std::string get_name() const override
static std::vector< unsigned int > get_dpo_vector(const unsigned int deg)
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
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
FE_FaceQ(const unsigned int p)
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) const override
virtual bool hp_constraints_are_implemented() const override
const unsigned int degree
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_unique_faces() const
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const =0
std::vector< Point< dim > > unit_support_points
double least_squares(Vector< number2 > &dst, const Vector< number2 > &src) const
Abstract base class for mapping classes.
#define DEAL_II_NAMESPACE_OPEN
constexpr bool running_in_debug_mode()
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NOT_IMPLEMENTED()
#define Assert(cond, exc)
static ::ExceptionBase & ExcLeastSquaresError(double arg1)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ update_gradients
Shape function gradients.
@ either_element_can_dominate
@ other_element_dominates
@ neither_element_dominates
std::string dim_string(const int dim, const int spacedim)
constexpr unsigned int invalid_unsigned_int
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices()
static Point< dim > child_to_cell_coordinates(const Point< dim > &p, const unsigned int child_index, const RefinementCase< dim > refine_case=RefinementCase< dim >::isotropic_refinement)