37 std::vector<unsigned int>
38 build_multiplicities(
const std::vector<std::vector<T>> &
functions)
40 std::vector<unsigned int> multiplicities;
41 multiplicities.push_back(1);
42 for (
unsigned int i = 0; i <
functions.size(); i++)
43 multiplicities.push_back(
functions[i].size());
45 return multiplicities;
52 template <
int dim,
int spacedim>
53 std::vector<const FiniteElement<dim, spacedim> *>
58 std::vector<const FiniteElement<dim, spacedim> *> fes;
59 fes.push_back(fe_base);
60 for (
unsigned int i = 0; i < fe_enriched.size(); i++)
61 fes.push_back(fe_enriched[i]);
71 template <
int dim,
int spacedim>
75 const std::vector<unsigned int> & multiplicities,
77 const typename ::Triangulation<dim, spacedim>::cell_iterator
83 "FEs and multiplicities should have the same size"));
91 const unsigned int n_comp_base = fes[0]->n_components();
94 for (
unsigned int fe = 1; fe < fes.size(); fe++)
102 "Only dominating FE_Nothing can be used in FE_Enriched"));
105 fes[fe]->n_components() == n_comp_base,
107 "All elements must have the same number of components"));
117 template <
int dim,
int spacedim>
123 for (
unsigned int fe = 1; fe < fes.size(); fe++)
135 template <
int dim,
int spacedim>
139 FE_Nothing<dim, spacedim>(fe_base.n_components(),
145 template <
int dim,
int spacedim>
154 const typename
Triangulation<dim, spacedim>::cell_iterator &)>>>(
157 const typename
Triangulation<dim, spacedim>::cell_iterator &)>>(
159 [=](const typename
Triangulation<dim, spacedim>::cell_iterator &)
160 -> const
Function<spacedim> * {
return enrichment_function; })))
164 template <
int dim,
int spacedim>
177 template <
int dim,
int spacedim>
180 const std::vector<unsigned int> & multiplicities,
196 Assert(internal::FE_Enriched::consistency_check(fes,
206 for (
unsigned int fe = 1; fe < fes.size(); fe++)
214 for (
unsigned int system_index = 0; system_index < this->
n_dofs_per_cell();
217 const unsigned int base_no =
222 const unsigned int base_m =
227 "Size mismatch for base_no_mult_local_enriched_dofs: " 243 for (
unsigned int base_no = 1;
247 for (
unsigned int m = 0;
251 fes[base_no]->n_dofs_per_cell(),
259 template <
int dim,
int spacedim>
260 const std::vector<std::vector<std::function<const Function<spacedim> *(
268 template <
int dim,
int spacedim>
276 "For enriched finite elements shape_value() can not be defined on the reference element."));
281 template <
int dim,
int spacedim>
282 std::unique_ptr<FiniteElement<dim, spacedim>>
285 std::vector<const FiniteElement<dim, spacedim> *> fes;
286 std::vector<unsigned int> multiplicities;
294 return std::unique_ptr<FE_Enriched<dim, spacedim>>(
299 template <
int dim,
int spacedim>
322 template <
int dim,
int spacedim>
324 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
332 auto update_each_flags = fes_data->update_each;
333 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
334 data_ptr = std::make_unique<InternalData>(std::move(fes_data));
343 const unsigned int n_q_points = quadrature.
size();
351 data.enrichment[base][m].values.resize(n_q_points);
354 data.enrichment[base][m].gradients.resize(n_q_points);
357 data.enrichment[base][m].hessians.resize(n_q_points);
365 template <
int dim,
int spacedim>
366 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
377 fe_system->get_face_data(update_flags, mapping, quadrature, output_data);
386 template <
int dim,
int spacedim>
387 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
397 fe_system->get_subface_data(update_flags, mapping, quadrature, output_data);
406 template <
int dim,
int spacedim>
407 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
415 auto data =
fe_system->get_data(flags, mapping, quadrature, output_data);
424 template <
int dim,
int spacedim>
428 const std::vector<unsigned int> & multiplicities)
430 Assert(fes.size() == multiplicities.size(),
437 for (
unsigned int i = 0; i < fes.size(); i++)
438 if (multiplicities[i] > 0)
454 for (
unsigned int face_no = 0; face_no < this->
n_unique_faces(); ++face_no)
489 fe_system->adjust_line_dof_index_for_line_orientation_table;
494 template <
int dim,
int spacedim>
498 std::ostringstream namebuf;
511 return namebuf.str();
515 template <
int dim,
int spacedim>
523 template <
int dim,
int spacedim>
531 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
538 Assert(dynamic_cast<const InternalData *>(&fe_internal) !=
nullptr,
554 quadrature, fe_data, mapping_data, cell, output_data);
558 template <
int dim,
int spacedim>
562 const unsigned int face_no,
566 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
573 Assert(dynamic_cast<const InternalData *>(&fe_internal) !=
nullptr,
591 quadrature[0], fe_data, mapping_data, cell, output_data);
595 template <
int dim,
int spacedim>
599 const unsigned int face_no,
600 const unsigned int sub_no,
604 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
611 Assert(dynamic_cast<const InternalData *>(&fe_internal) !=
nullptr,
628 quadrature, fe_data, mapping_data, cell, output_data);
632 template <
int dim,
int spacedim>
654 const unsigned int n_q_points = quadrature.
size();
664 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
675 for (
unsigned int system_index = 0;
680 const unsigned int base_index =
691 unsigned int out_index = 0;
692 for (
unsigned int i = 0; i < system_index; ++i)
694 unsigned int in_index = 0;
695 for (
unsigned int i = 0; i < base_index; ++i)
702 for (
unsigned int s = 0;
707 for (
unsigned int q = 0; q < n_q_points; ++q)
709 base_data.shape_values(in_index + s, q);
712 for (
unsigned int q = 0; q < n_q_points; ++q)
714 base_data.shape_gradients[in_index + s][q];
717 for (
unsigned int q = 0; q < n_q_points; ++q)
719 base_data.shape_hessians[in_index + s][q];
728 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
735 for (
unsigned int m = 0;
747 "The pointer to the enrichment function is not set"));
751 "Only scalar-valued enrichment functions are allowed"));
758 fe_data.
enrichment[base_no][m].hessians.size(),
760 for (
unsigned int q = 0; q < n_q_points; q++)
771 fe_data.
enrichment[base_no][m].gradients.size(),
773 for (
unsigned int q = 0; q < n_q_points; q++)
785 for (
unsigned int q = 0; q < n_q_points; q++)
804 for (
unsigned int m = 0;
807 for (
unsigned int i = 0;
811 const unsigned int enriched_dof =
813 for (
unsigned int q = 0; q < n_q_points; ++q)
833 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
835 for (
unsigned int m = 0;
838 for (
unsigned int i = 0;
842 const unsigned int enriched_dof =
844 for (
unsigned int q = 0; q < n_q_points; ++q)
856 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
858 for (
unsigned int m = 0;
861 for (
unsigned int i = 0;
865 const unsigned int enriched_dof =
867 for (
unsigned int q = 0; q < n_q_points; ++q)
877 template <
int dim,
int spacedim>
885 template <
int dim,
int spacedim>
893 template <
int dim,
int spacedim>
898 const unsigned int face_no)
const 903 fe_system->get_face_interpolation_matrix(fe_enr_other->get_fe_system(),
917 template <
int dim,
int spacedim>
921 const unsigned int subface,
923 const unsigned int face_no)
const 928 fe_system->get_subface_interpolation_matrix(fe_enr_other->get_fe_system(),
943 template <
int dim,
int spacedim>
944 std::vector<std::pair<unsigned int, unsigned int>>
951 return fe_system->hp_vertex_dof_identities(fe_enr_other->get_fe_system());
956 return std::vector<std::pair<unsigned int, unsigned int>>();
961 template <
int dim,
int spacedim>
962 std::vector<std::pair<unsigned int, unsigned int>>
969 return fe_system->hp_line_dof_identities(fe_enr_other->get_fe_system());
974 return std::vector<std::pair<unsigned int, unsigned int>>();
979 template <
int dim,
int spacedim>
980 std::vector<std::pair<unsigned int, unsigned int>>
983 const unsigned int face_no)
const 988 return fe_system->hp_quad_dof_identities(fe_enr_other->get_fe_system(),
994 return std::vector<std::pair<unsigned int, unsigned int>>();
999 template <
int dim,
int spacedim>
1003 const unsigned int codim)
const 1022 return fe_system->compare_for_domination(fe_enr_other->get_fe_system(),
1033 template <
int dim,
int spacedim>
1036 const unsigned int child,
1039 return fe_system->get_prolongation_matrix(child, refinement_case);
1043 template <
int dim,
int spacedim>
1046 const unsigned int child,
1049 return fe_system->get_restriction_matrix(child, refinement_case);
1056 template <
int dim,
int spacedim>
1059 : fesystem_data(
std::move(fesystem_data))
1063 template <
int dim,
int spacedim>
1066 const unsigned int base_no)
const 1072 template <
int dim,
int spacedim>
1075 const unsigned int base_no)
const 1085 template <
int dim,
int spacedim>
1093 std::vector<bool> vertices_subdomain_1(
1098 if (predicate_1(cell))
1100 vertices_subdomain_1[cell->vertex_index(v)] =
true;
1104 if (predicate_2(cell))
1106 if (vertices_subdomain_1[cell->vertex_index(v)] ==
true)
1115 template <
int dim,
int spacedim>
1120 std::vector<unsigned int> & predicate_colors)
1122 const unsigned int num_indices = predicates.size();
1127 dsp.
reinit(num_indices, num_indices);
1134 for (
unsigned int i = 0; i < num_indices; ++i)
1135 for (
unsigned int j = i + 1; j < num_indices; ++j)
1147 predicate_colors.resize(num_indices);
1155 template <
int dim,
int spacedim>
1160 const std::vector<unsigned int> & predicate_colors,
1161 std::map<
unsigned int, std::map<unsigned int, unsigned int>>
1162 & cellwise_color_predicate_map,
1163 std::vector<std::set<unsigned int>> &fe_sets)
1167 cellwise_color_predicate_map.clear();
1195 unsigned int map_index = 0;
1199 cell->set_active_fe_index(0);
1208 cell->set_material_id(map_index);
1209 std::set<unsigned int> color_list;
1217 std::map<unsigned int, unsigned int> &cell_map =
1218 cellwise_color_predicate_map[map_index];
1219 for (
unsigned int i = 0; i < predicates.size(); ++i)
1221 if (predicates[i](cell))
1227 auto ret = cell_map.insert(
1228 std::pair<unsigned int, unsigned int>(predicate_colors[i],
1233 "Only one enrichment function per color"));
1235 color_list.insert(predicate_colors[i]);
1252 if (!color_list.empty())
1255 std::find(fe_sets.begin(), fe_sets.end(), color_list);
1257 if (it == fe_sets.end())
1259 fe_sets.push_back(color_list);
1260 cell->set_active_fe_index(fe_sets.size() - 1);
1265 cell->set_active_fe_index(std::distance(fe_sets.begin(), it));
1312 const unsigned int fe_index = cell->active_fe_index();
1313 const std::set<unsigned int> fe_set = fe_sets.at(fe_index);
1320 if (!cell->at_boundary(face) &&
1321 cell->material_id() < cell->neighbor(face)->material_id())
1323 const auto nbr_fe_index =
1324 cell->neighbor(face)->active_fe_index();
1327 const auto nbr_fe_set = fe_sets.at(nbr_fe_index);
1330 std::set<unsigned int> intersection_set;
1331 std::set_intersection(
1336 std::inserter(intersection_set, intersection_set.begin()));
1339 if (!intersection_set.empty())
1341 const auto it = std::find(fe_sets.begin(),
1345 if (it == fe_sets.end())
1347 fe_sets.push_back(intersection_set);
1357 template <
int dim,
int spacedim>
1360 const unsigned int n_colors,
1362 const std::map<
unsigned int, std::map<unsigned int, unsigned int>>
1363 &cellwise_color_predicate_map,
1368 color_enrichments.clear();
1384 color_enrichments.resize(n_colors);
1385 for (
unsigned int i = 0; i < n_colors; ++i)
1387 color_enrichments[i] =
1390 const unsigned int id = cell->material_id();
1397 return enrichments[cellwise_color_predicate_map.at(
id).at(i + 1)]
1405 template <
int dim,
int spacedim>
1408 const unsigned int n_colors,
1409 const std::vector<std::set<unsigned int>> &fe_sets,
1412 & color_enrichments,
1419 const std::function<const Function<spacedim> *(
1425 ExcMessage(
"Called enrichment function for FE_Nothing"));
1433 for (
const auto &fe_set : fe_sets)
1435 std::vector<const FiniteElement<dim, spacedim> *> vec_fe_enriched(
1436 n_colors, &fe_nothing);
1437 std::vector<std::vector<std::function<const Function<spacedim> *(
1441 for (
const unsigned int color_id : fe_set)
1446 const unsigned int ind = color_id - 1;
1456 vec_fe_enriched[ind] = &fe_enriched;
1460 functions[ind][0] = color_enrichments[ind];
1475 template <
int dim,
int spacedim>
1482 , fe_enriched(fe_enriched)
1484 , predicates(predicates)
1491 ExcMessage(
"Number of predicates should be positive"));
1496 template <
int dim,
int spacedim>
1514 internal::make_colorwise_enrichment_functions<dim, spacedim>(
1532 #include "fe_enriched.inst" virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > face_system_to_component_table
bool is_dominating() const
static const unsigned int invalid_unsigned_int
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const override
#define AssertDimension(dim1, dim2)
void set_cellwise_color_set_and_fe_index(DoFHandler< dim, spacedim > &dof_handler, const std::vector< predicate_function< dim, spacedim >> &predicates, const std::vector< unsigned int > &predicate_colors, std::map< unsigned int, std::map< unsigned int, unsigned int >> &cellwise_color_predicate_map, std::vector< std::set< unsigned int >> &fe_sets)
const std::vector< std::vector< std::function< const Function< spacedim > *(const typename Triangulation< dim, spacedim >::cell_iterator &)> > > get_enrichments() const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
const hp::FECollection< dim, spacedim > & build_fe_collection(DoFHandler< dim, spacedim > &dof_handler)
FullMatrix< double > interface_constraints
Contents is actually a matrix.
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override
InternalData(std::unique_ptr< typename FESystem< dim, spacedim >::InternalData > fesystem_data)
unsigned int n_nonzero_components(const unsigned int i) const
constexpr SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
FE_Enriched(const FiniteElement< dim, spacedim > &fe_base, const FiniteElement< dim, spacedim > &fe_enriched, const Function< spacedim > *enrichment_function)
void initialize(const std::vector< const FiniteElement< dim, spacedim > *> &fes, const std::vector< unsigned int > &multiplicities)
const FiniteElement< dim, spacedim > & fe_enriched
void add(const size_type i, const size_type j)
internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > & get_fe_output_object(const unsigned int base_no) const
#define AssertIndexRange(index, range)
bool find_connection_between_subdomains(const DoFHandler< dim, spacedim > &dof_handler, const predicate_function< dim, spacedim > &predicate_1, const predicate_function< dim, spacedim > &predicate_2)
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
std::vector< std::vector< EnrichmentValues > > enrichment
std::function< bool(const typename Triangulation< dim, spacedim >::cell_iterator &)> predicate_function
Transformed quadrature points.
#define AssertThrow(cond, exc)
hp::FECollection< dim, spacedim > fe_collection
static ::ExceptionBase & ExcInterpolationNotImplemented()
void multiply_by_enrichment(const Quadrature< dim_1 > &quadrature, const InternalData &fe_data, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename Triangulation< dim, spacedim >::cell_iterator &cell, internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
std::vector< Point< dim > > unit_support_points
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual std::string get_name() const override
constexpr SymmetricTensor< 4, dim, Number > outer_product(const SymmetricTensor< 2, dim, Number > &t1, const SymmetricTensor< 2, dim, Number > &t2)
const std::vector< std::vector< std::function< const Function< spacedim > *(const typename Triangulation< dim, spacedim >::cell_iterator &)> > > enrichments
const FESystem< dim, spacedim > & get_fe_system() const
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
void push_back(const FiniteElement< dim, spacedim > &new_fe)
std::unique_ptr< typename FESystem< dim, spacedim >::InternalData > fesystem_data
Third derivatives of shape functions.
const std::unique_ptr< const FESystem< dim, spacedim > > fe_system
#define Assert(cond, exc)
unsigned int element_multiplicity(const unsigned int index) const
IteratorRange< active_cell_iterator > active_cell_iterators() const
FiniteElement< dim, spacedim >::InternalDataBase & get_fe_data(const unsigned int base_no) const
void reinit(const size_type m, const size_type n, const IndexSet &rowset=IndexSet())
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
unsigned int size() const
#define DEAL_II_NAMESPACE_CLOSE
void make_fe_collection_from_colored_enrichments(const unsigned int n_colors, const std::vector< std::set< unsigned int >> &fe_sets, const std::vector< std::function< const Function< spacedim > *(const typename Triangulation< dim, spacedim >::cell_iterator &)>> &color_enrichments, const FiniteElement< dim, spacedim > &fe_base, const FiniteElement< dim, spacedim > &fe_enriched, const FE_Nothing< dim, spacedim > &fe_nothing, hp::FECollection< dim, spacedim > &fe_collection)
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
std::vector< std::set< unsigned int > > fe_sets
std::vector< unsigned int > predicate_colors
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
void copy_from(const size_type n_rows, const size_type n_cols, const ForwardIterator begin, const ForwardIterator end)
Second derivatives of shape functions.
unsigned int n_unique_faces() const
void make_colorwise_enrichment_functions(const unsigned int n_colors, const std::vector< std::shared_ptr< Function< spacedim >>> &enrichments, const std::map< unsigned int, std::map< unsigned int, unsigned int >> &cellwise_color_predicate_map, std::vector< std::function< const Function< spacedim > *(const typename Triangulation< dim, spacedim >::cell_iterator &)>> &color_enrichments)
unsigned int color_predicates(const DoFHandler< dim, spacedim > &mesh, const std::vector< predicate_function< dim, spacedim >> &predicates, std::vector< unsigned int > &predicate_colors)
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
std::string dim_string(const int dim, const int spacedim)
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
unsigned int size() const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
const std::vector< predicate_function< dim, spacedim > > predicates
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override
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 typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
std::unique_ptr< To > dynamic_unique_cast(std::unique_ptr< From > &&p)
unsigned int n_components() const
unsigned int n_dofs_per_cell() const
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
const Triangulation< dim, spacedim > & get_triangulation() const
#define DEAL_II_NAMESPACE_OPEN
Shape function gradients.
std::vector< std::vector< std::vector< unsigned int > > > base_no_mult_local_enriched_dofs
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > face_system_to_base_table
void push_back(const size_type size)
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
static ::ExceptionBase & ExcNotImplemented()
std::vector< cell_iterator_function > color_enrichments
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > setup_data(std::unique_ptr< typename FESystem< dim, spacedim >::InternalData > fes_data, const UpdateFlags flags, const Quadrature< dim_1 > &quadrature) const
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
BlockIndices base_to_block_indices
std::map< unsigned int, std::map< unsigned int, unsigned int > > cellwise_color_predicate_map
const std::vector< std::shared_ptr< Function< spacedim > > > enrichments
unsigned int n_base_elements() 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 typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
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 typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
std::vector< int > adjust_line_dof_index_for_line_orientation_table
const FE_Nothing< dim, spacedim > fe_nothing
virtual bool hp_constraints_are_implemented() const override
int(&) functions(const void *v1, const void *v2)
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
static ::ExceptionBase & ExcInternalError()