17 #include <deal.II/fe/fe_enriched.h> 19 #ifdef DEAL_II_WITH_CXX14 21 #include <deal.II/fe/fe_tools.h> 23 DEAL_II_NAMESPACE_OPEN
31 std::vector<unsigned int>
32 build_multiplicities(
const std::vector<std::vector<T > > &functions )
34 std::vector<unsigned int> multiplicities;
35 multiplicities.push_back(1);
36 for (
unsigned int i = 0; i < functions.size(); i++)
37 multiplicities.push_back(functions[i].size());
39 return multiplicities;
46 template <
int dim,
int spacedim>
47 std::vector< const FiniteElement< dim, spacedim > * >
51 std::vector< const FiniteElement< dim, spacedim > * > fes;
52 fes.push_back(fe_base);
53 for (
unsigned int i = 0; i < fe_enriched.size(); i++)
54 fes.push_back(fe_enriched[i]);
64 template <
int dim,
int spacedim>
67 const std::vector< unsigned int > &multiplicities,
73 ExcMessage(
"FEs and multiplicities should have the same size"));
81 const unsigned int n_comp_base = fes[0]->n_components();
84 for (
unsigned int fe=1; fe < fes.size(); fe++)
89 ExcMessage(
"Only dominating FE_Nothing can be used in FE_Enriched"));
92 ExcMessage(
"All elements must have the same number of components"));
101 template <
int dim,
int spacedim>
106 for (
unsigned int fe=1; fe < fes.size(); fe++)
116 template <
int dim,
int spacedim>
126 template <
int dim,
int spacedim>
140 return enrichment_function;
147 template <
int dim,
int spacedim>
152 FE_Enriched<dim,spacedim> (build_fes(fe_base,fe_enriched),
153 build_multiplicities(functions),
158 template <
int dim,
int spacedim>
160 const std::vector< unsigned int > &multiplicities,
163 FiniteElement<dim,spacedim> (
FETools::Compositing::multiply_dof_numbers(fes,multiplicities,false),
164 FETools::Compositing::compute_restriction_is_additive_flags(fes,multiplicities),
165 FETools::Compositing::compute_nonzero_components(fes,multiplicities,false)),
171 Assert(consistency_check(fes,multiplicities,functions),
179 for (
unsigned int fe=1; fe < fes.size(); fe++)
187 for (
unsigned int system_index=0; system_index<this->
dofs_per_cell;
197 ExcMessage(
"Size mismatch for base_no_mult_local_enriched_dofs: " 199 "; base_no = " + std::to_string(base_no) +
200 "; base_m = " + std::to_string(base_m) +
201 "; system_index = " + std::to_string(system_index)));
217 fes[base_no]->dofs_per_cell));
222 template <
int dim,
int spacedim>
230 template <
int dim,
int spacedim>
240 template <
int dim,
int spacedim>
244 std::vector< const FiniteElement< dim, spacedim > * > fes;
245 std::vector< unsigned int > multiplicities;
257 template <
int dim,
int spacedim>
281 template <
int dim,
int spacedim>
304 const unsigned int n_q_points = quadrature.
size();
312 data->
enrichment[base][m].values.resize(n_q_points);
315 data->
enrichment[base][m].gradients.resize(n_q_points);
318 data->
enrichment[base][m].hessians.resize(n_q_points);
326 template <
int dim,
int spacedim>
339 template <
int dim,
int spacedim>
352 template <
int dim,
int spacedim>
365 template <
int dim,
int spacedim>
367 const std::vector<unsigned int> &multiplicities)
369 Assert (fes.size() == multiplicities.size(),
375 for (
unsigned int i=0; i<fes.size(); i++)
376 if (multiplicities[i]>0)
421 template <
int dim,
int spacedim>
425 std::ostringstream namebuf;
427 namebuf <<
"FE_Enriched<" 440 return namebuf.str();
444 template <
int dim,
int spacedim>
452 template <
int dim,
int spacedim>
459 const ::internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data,
464 Assert (dynamic_cast<const InternalData *> (&fe_internal) != NULL,
487 template <
int dim,
int spacedim>
491 const unsigned int face_no,
495 const ::internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data,
500 Assert (dynamic_cast<const InternalData *> (&fe_internal) != NULL,
523 template <
int dim,
int spacedim>
527 const unsigned int face_no,
528 const unsigned int sub_no,
532 const ::internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data,
537 Assert (dynamic_cast<const InternalData *> (&fe_internal) != NULL,
562 template <
int dim,
int spacedim>
582 const unsigned int n_q_points = quadrature.
size();
592 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
603 for (
unsigned int system_index=0; system_index<this->
dofs_per_cell;
617 unsigned int out_index = 0;
618 for (
unsigned int i=0; i<system_index; ++i)
620 unsigned int in_index = 0;
621 for (
unsigned int i=0; i<base_index; ++i)
631 for (
unsigned int q=0; q<n_q_points; ++q)
636 for (
unsigned int q=0; q<n_q_points; ++q)
641 for (
unsigned int q=0; q<n_q_points; ++q)
652 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
660 ExcMessage(
"The pointer to the enrichment function is NULL"));
663 ExcMessage(
"Only scalar-valued enrichment functions are allowed"));
670 for (
unsigned int q=0; q<n_q_points; q++)
679 for (
unsigned int q=0; q<n_q_points; q++)
688 for (
unsigned int q=0; q<n_q_points; q++)
701 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
707 for (
unsigned int q=0; q<n_q_points; ++q)
723 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
729 for (
unsigned int q=0; q<n_q_points; ++q)
739 for (
unsigned int base_no = 1; base_no < this->
n_base_elements(); base_no++)
745 for (
unsigned int q=0; q<n_q_points; ++q)
754 template <
int dim,
int spacedim>
763 template <
int dim,
int spacedim>
772 template <
int dim,
int spacedim>
781 fe_system.get_face_interpolation_matrix(fe_enr_other->get_fe_system(),
787 AssertThrow(
false,
typename FEL::ExcInterpolationNotImplemented());
792 template <
int dim,
int spacedim>
796 const unsigned int subface,
802 fe_system.get_subface_interpolation_matrix(fe_enr_other->get_fe_system(),
809 AssertThrow(
false,
typename FEL::ExcInterpolationNotImplemented());
814 template <
int dim,
int spacedim>
815 std::vector<std::pair<unsigned int, unsigned int> >
822 return fe_system.hp_vertex_dof_identities(fe_enr_other->get_fe_system());
827 return std::vector<std::pair<unsigned int, unsigned int> >();
832 template <
int dim,
int spacedim>
833 std::vector<std::pair<unsigned int, unsigned int> >
840 return fe_system.hp_line_dof_identities(fe_enr_other->get_fe_system());
845 return std::vector<std::pair<unsigned int, unsigned int> >();
850 template <
int dim,
int spacedim>
851 std::vector<std::pair<unsigned int, unsigned int> >
858 return fe_system.hp_quad_dof_identities(fe_enr_other->get_fe_system());
863 return std::vector<std::pair<unsigned int, unsigned int> >();
868 template <
int dim,
int spacedim>
886 return fe_system.compare_for_face_domination(fe_enr_other->get_fe_system());
895 template <
int dim,
int spacedim>
900 return fe_system.get_prolongation_matrix(child, refinement_case);
904 template <
int dim,
int spacedim>
909 return fe_system.get_restriction_matrix(child, refinement_case);
916 template <
int dim,
int spacedim>
919 fesystem_data(
std::move(fesystem_data))
923 template <
int dim,
int spacedim>
932 template <
int dim,
int spacedim>
942 #include "fe_enriched.inst" 944 DEAL_II_NAMESPACE_CLOSE
const std::vector< std::vector< std::function< const Function< spacedim > *(const typename Triangulation< dim, spacedim >::cell_iterator &) > > > get_enrichments() const
FiniteElement< dim, spacedim >::InternalDataBase * setup_data(std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fes_data, const UpdateFlags flags, const Quadrature< dim_1 > &quadrature) const
FullMatrix< double > interface_constraints
InternalData(std::unique_ptr< typename FESystem< dim, spacedim >::InternalData > fesystem_data)
virtual FiniteElement< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
FE_Enriched(const FiniteElement< dim, spacedim > &fe_base, const FiniteElement< dim, spacedim > &fe_enriched, const Function< spacedim > *enrichment_function)
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
bool is_dominating() const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
std::vector< std::vector< EnrichmentValues > > enrichment
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const
Transformed quadrature points.
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
virtual FiniteElement< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
#define AssertThrow(cond, exc)
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > face_system_to_base_table
std::unique_ptr< typename FESystem< dim, spacedim >::InternalData > fesystem_data
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
std::vector< Point< dim > > unit_support_points
unsigned int element_multiplicity(const unsigned int index) const
unsigned int size() const
static::ExceptionBase & ExcMessage(std::string arg1)
void initialize(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities)
Third derivatives of shape functions.
virtual FiniteElement< dim, spacedim > * clone() const
#define Assert(cond, exc)
unsigned int n_nonzero_components(const unsigned int i) const
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
std::vector< Point< dim-1 > > unit_face_support_points
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
FESystem< dim, spacedim > fe_system
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
unsigned int n_base_elements() const
unsigned int n_components() const
Second derivatives of shape functions.
virtual bool hp_constraints_are_implemented() const
std::string dim_string(const int dim, const int spacedim)
internal::FEValues::FiniteElementRelatedData< dim, spacedim > & get_fe_output_object(const unsigned int base_no) const
const unsigned int dofs_per_cell
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
virtual FiniteElement< dim, spacedim >::InternalDataBase * get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
void push_back(const size_type size)
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const
const unsigned int dofs_per_face
static::ExceptionBase & ExcNotImplemented()
std::vector< std::pair< unsigned int, unsigned int > > face_system_to_component_table
virtual std::string get_name() const
const FESystem< dim, spacedim > & get_fe_system() const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
BlockIndices base_to_block_indices
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
std::vector< int > adjust_line_dof_index_for_line_orientation_table
const std::vector< std::vector< std::function< const Function< spacedim > *(const typename Triangulation< dim, spacedim >::cell_iterator &) > > > enrichments
FiniteElement< dim, spacedim >::InternalDataBase & get_fe_data(const unsigned int base_no) const
void multiply_by_enrichment(const Quadrature< dim_1 > &quadrature, const InternalData &fe_data, const internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, const typename Triangulation< dim, spacedim >::cell_iterator &cell, internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
static::ExceptionBase & ExcInternalError()