16 #include <deal.II/base/memory_consumption.h> 17 #include <deal.II/base/quadrature.h> 18 #include <deal.II/base/qprojector.h> 19 #include <deal.II/grid/tria.h> 20 #include <deal.II/grid/tria_iterator.h> 21 #include <deal.II/dofs/dof_accessor.h> 22 #include <deal.II/fe/mapping.h> 23 #include <deal.II/fe/fe_system.h> 24 #include <deal.II/fe/fe_values.h> 25 #include <deal.II/fe/fe_tools.h> 30 DEAL_II_NAMESPACE_OPEN
34 bool IsNonZero (
unsigned int i)
39 unsigned int count_nonzeros(
const std::vector<unsigned int> &vec)
41 return std::count_if(vec.begin(), vec.end(), IsNonZero);
48 template <
int dim,
int spacedim>
51 base_fe_datas(n_base_elements),
52 base_fe_output_objects(n_base_elements)
57 template <
int dim,
int spacedim>
71 template <
int dim,
int spacedim>
83 template <
int dim,
int spacedim>
96 template <
int dim,
int spacedim>
111 template <
int dim,
int spacedim>
115 template <
int dim,
int spacedim>
117 const unsigned int n_elements) :
119 FETools::Compositing::compute_restriction_is_additive_flags (&fe, n_elements),
120 FETools::Compositing::compute_nonzero_components(&fe, n_elements)),
123 std::vector<const FiniteElement<dim,spacedim>*> fes;
125 std::vector<unsigned int> multiplicities;
126 multiplicities.push_back(n_elements);
132 template <
int dim,
int spacedim>
134 const unsigned int n1,
136 const unsigned int n2) :
138 FETools::Compositing::compute_restriction_is_additive_flags (&fe1, n1,
140 FETools::Compositing::compute_nonzero_components(&fe1, n1,
144 std::vector<const FiniteElement<dim,spacedim>*> fes;
147 std::vector<unsigned int> multiplicities;
148 multiplicities.push_back(n1);
149 multiplicities.push_back(n2);
155 template <
int dim,
int spacedim>
157 const unsigned int n1,
159 const unsigned int n2,
161 const unsigned int n3) :
165 FETools::Compositing::compute_restriction_is_additive_flags (&fe1, n1,
168 FETools::Compositing::compute_nonzero_components(&fe1, n1,
173 std::vector<const FiniteElement<dim,spacedim>*> fes;
177 std::vector<unsigned int> multiplicities;
178 multiplicities.push_back(n1);
179 multiplicities.push_back(n2);
180 multiplicities.push_back(n3);
186 template <
int dim,
int spacedim>
188 const unsigned int n1,
190 const unsigned int n2,
192 const unsigned int n3,
194 const unsigned int n4) :
199 FETools::Compositing::compute_restriction_is_additive_flags (&fe1, n1,
203 FETools::Compositing::compute_nonzero_components(&fe1, n1,
209 std::vector<const FiniteElement<dim,spacedim>*> fes;
214 std::vector<unsigned int> multiplicities;
215 multiplicities.push_back(n1);
216 multiplicities.push_back(n2);
217 multiplicities.push_back(n3);
218 multiplicities.push_back(n4);
224 template <
int dim,
int spacedim>
226 const unsigned int n1,
228 const unsigned int n2,
230 const unsigned int n3,
232 const unsigned int n4,
234 const unsigned int n5) :
240 FETools::Compositing::compute_restriction_is_additive_flags (&fe1, n1,
245 FETools::Compositing::compute_nonzero_components(&fe1, n1,
252 std::vector<const FiniteElement<dim,spacedim>*> fes;
258 std::vector<unsigned int> multiplicities;
259 multiplicities.push_back(n1);
260 multiplicities.push_back(n2);
261 multiplicities.push_back(n3);
262 multiplicities.push_back(n4);
263 multiplicities.push_back(n5);
269 template <
int dim,
int spacedim>
272 const std::vector<unsigned int> &multiplicities)
275 FETools::Compositing::compute_restriction_is_additive_flags (fes, multiplicities),
276 FETools::Compositing::compute_nonzero_components(fes, multiplicities)),
284 template <
int dim,
int spacedim>
290 template <
int dim,
int spacedim>
301 std::ostringstream namebuf;
303 namebuf <<
"FESystem<" 316 return namebuf.str();
321 template <
int dim,
int spacedim>
325 std::vector< const FiniteElement<dim,spacedim>* > fes;
326 std::vector<unsigned int> multiplicities;
338 template <
int dim,
int spacedim>
348 .shape_value(this->system_to_base_table[i].second, p));
353 template <
int dim,
int spacedim>
357 const unsigned int component)
const 390 template <
int dim,
int spacedim>
400 .shape_grad(this->system_to_base_table[i].second, p));
405 template <
int dim,
int spacedim>
409 const unsigned int component)
const 435 template <
int dim,
int spacedim>
445 .shape_grad_grad(this->system_to_base_table[i].second, p));
450 template <
int dim,
int spacedim>
454 const unsigned int component)
const 480 template <
int dim,
int spacedim>
490 .shape_3rd_derivative(this->system_to_base_table[i].second, p));
495 template <
int dim,
int spacedim>
499 const unsigned int component)
const 525 template <
int dim,
int spacedim>
535 .shape_4th_derivative(this->system_to_base_table[i].second, p));
540 template <
int dim,
int spacedim>
544 const unsigned int component)
const 570 template <
int dim,
int spacedim>
610 ExcInterpolationNotImplemented());
615 source_fe.element_multiplicity(i),
617 ExcInterpolationNotImplemented());
626 std::vector<FullMatrix<double> > base_matrices (this->
n_base_elements());
630 source_fe.base_element(i).dofs_per_cell);
631 base_element(i).get_interpolation_matrix (source_fe.base_element(i),
638 interpolation_matrix = 0;
640 for (
unsigned int j=0; j<source_fe.dofs_per_cell; ++j)
642 source_fe.system_to_base_table[j].first)
643 interpolation_matrix(i,j)
646 source_fe.system_to_base_table[j].second));
651 template <
int dim,
int spacedim>
660 ExcMessage(
"Restriction matrices are only available for refined cells!"));
665 if (this->
restriction[refinement_case-1][child].n() == 0)
670 if (this->
restriction[refinement_case-1][child].n() ==
672 return this->
restriction[refinement_case-1][child];
675 bool do_restriction =
true;
678 std::vector<const FullMatrix<double> *>
684 &
base_element(i).get_restriction_matrix(child, refinement_case);
686 do_restriction =
false;
727 restriction(i,j) = (*base_matrices[base])(base_index_i,base_index_j);
731 (this->restriction[refinement_case-1][child]));
735 return this->
restriction[refinement_case-1][child];
740 template <
int dim,
int spacedim>
749 ExcMessage(
"Restriction matrices are only available for refined cells!"));
755 if (this->
prolongation[refinement_case-1][child].n() == 0)
763 bool do_prolongation =
true;
764 std::vector<const FullMatrix<double> *>
769 &
base_element(i).get_prolongation_matrix(child, refinement_case);
771 do_prolongation =
false;
793 prolongate(i,j) = (*base_matrices[base])(base_index_i,base_index_j);
804 template <
int dim,
int spacedim>
808 const unsigned int face,
809 const bool face_orientation,
810 const bool face_flip,
811 const bool face_rotation)
const 816 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int>
820 base_face_to_cell_index
821 = this->
base_element(face_base_index.first.first).face_to_cell_index (face_base_index.second,
832 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int>
833 target = std::make_pair (face_base_index.first, base_face_to_cell_index);
850 template <
int dim,
int spacedim>
857 for (
unsigned int base_no=0; base_no<this->
n_base_elements(); ++base_no)
858 out |=
base_element(base_no).requires_update_flags (flags);
864 template <
int dim,
int spacedim>
894 for (
unsigned int base_no=0; base_no<this->
n_base_elements(); ++base_no)
899 flags |
base_element(base_no).requires_update_flags(flags));
908 base_element(base_no).get_data (flags, mapping, quadrature,
909 base_fe_output_object);
920 template <
int dim,
int spacedim>
950 for (
unsigned int base_no=0; base_no<this->
n_base_elements(); ++base_no)
955 flags |
base_element(base_no).requires_update_flags(flags));
964 base_element(base_no).get_face_data (flags, mapping, quadrature,
965 base_fe_output_object);
978 template <
int dim,
int spacedim>
1008 for (
unsigned int base_no=0; base_no<this->
n_base_elements(); ++base_no)
1013 flags |
base_element(base_no).requires_update_flags(flags));
1022 base_element(base_no).get_subface_data (flags, mapping, quadrature,
1023 base_fe_output_object);
1033 template <
int dim,
int spacedim>
1041 const ::internal::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
1046 quadrature, cell_similarity, mapping_internal, fe_internal,
1047 mapping_data, output_data);
1052 template <
int dim,
int spacedim>
1056 const unsigned int face_no,
1060 const ::internal::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
1066 mapping_data, output_data);
1072 template <
int dim,
int spacedim>
1076 const unsigned int face_no,
1077 const unsigned int sub_no,
1081 const ::internal::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
1085 compute_fill (mapping, cell, face_no, sub_no, quadrature,
1087 mapping_data, output_data);
1092 template <
int dim,
int spacedim>
1093 template <
int dim_1>
1098 const unsigned int face_no,
1099 const unsigned int sub_no,
1133 for (
unsigned int base_no=0; base_no<this->
n_base_elements(); ++base_no)
1145 const Quadrature<dim-1> *face_quadrature = 0;
1146 const unsigned int n_q_points = quadrature.
size();
1150 const Subscriptor *quadrature_base_pointer = &quadrature;
1164 Assert (
dynamic_cast<const Quadrature<dim-1
> *>(quadrature_base_pointer) != 0,
1168 =
static_cast<const Quadrature<dim-1
> *>(quadrature_base_pointer);
1182 mapping, mapping_internal, mapping_data,
1183 base_fe_data, base_data);
1187 mapping, mapping_internal, mapping_data,
1188 base_fe_data, base_data);
1192 mapping, mapping_internal, mapping_data,
1193 base_fe_data, base_data);
1214 for (
unsigned int system_index=0; system_index<this->
dofs_per_cell;
1228 unsigned int out_index = 0;
1229 for (
unsigned int i=0; i<system_index; ++i)
1231 unsigned int in_index = 0;
1232 for (
unsigned int i=0; i<base_index; ++i)
1242 for (
unsigned int q=0; q<n_q_points; ++q)
1248 for (
unsigned int q=0; q<n_q_points; ++q)
1254 for (
unsigned int q=0; q<n_q_points; ++q)
1260 for (
unsigned int q=0; q<n_q_points; ++q)
1269 template <
int dim,
int spacedim>
1277 if (
base_element(base).constraints_are_implemented() ==
false)
1301 const std::pair<std::pair<unsigned int,unsigned int>,
unsigned int> n_index
1306 std::pair<std::pair<unsigned int,unsigned int>,
unsigned int> m_index;
1328 const unsigned int index_in_line
1330 const unsigned int sub_line
1338 const unsigned int tmp1 = 2*this->dofs_per_vertex+index_in_line;
1356 m_index.second =
base_element(m_index.first.first).dofs_per_vertex +
1357 base_element(m_index.first.first).dofs_per_line*sub_line +
1377 const unsigned int index_in_line
1379 const unsigned int sub_line
1383 const unsigned int tmp1 = 4*this->dofs_per_vertex+index_in_line;
1393 m_index.second = 5*
base_element(m_index.first.first).dofs_per_vertex +
1394 base_element(m_index.first.first).dofs_per_line*sub_line +
1401 const unsigned int index_in_quad
1406 const unsigned int sub_quad
1408 this->dofs_per_quad);
1411 const unsigned int tmp1 = 4*this->dofs_per_vertex +
1427 m_index.second = 5*
base_element(m_index.first.first).dofs_per_vertex +
1429 base_element(m_index.first.first).dofs_per_quad*sub_quad +
1443 if (n_index.first == m_index.first)
1445 = (
base_element(n_index.first.first).constraints()(m_index.second,
1452 template <
int dim,
int spacedim>
1454 const std::vector<unsigned int> &multiplicities)
1456 Assert (fes.size() == multiplicities.size(),
1459 ExcMessage (
"Need to pass at least one finite element."));
1460 Assert (count_nonzeros(multiplicities) > 0,
1461 ExcMessage(
"You only passed FiniteElements with multiplicity 0."));
1467 for (
unsigned int i=0; i<fes.size(); i++)
1468 if (multiplicities[i]>0)
1471 std::vector<Threads::Task<FiniteElement<dim,spacedim>*> > clone_base_elements;
1473 for (
unsigned int i=0; i<fes.size(); i++)
1474 if (multiplicities[i]>0)
1479 for (
unsigned int i=0; i<fes.size(); i++)
1481 if (multiplicities[i]>0)
1485 (clone_base_elements[ind].return_value()),
1528 template <
int dim,
int spacedim>
1536 for (
unsigned int base_el=0; base_el<this->
n_base_elements(); ++base_el)
1560 template <
int dim,
int spacedim>
1578 for (
unsigned int base_el=0; base_el<this->
n_base_elements(); ++base_el)
1601 =
base_element(base_i).unit_face_support_points[index_in_base];
1607 template <
int dim,
int spacedim>
1622 unsigned int index = 0;
1626 = this->
base_element(b).adjust_quad_dof_index_for_face_orientation_table;
1629 for (
unsigned int i=0; i<temp.
size(0); ++i)
1630 for (
unsigned int j=0; j<8; ++j)
1633 index += temp.
size(0);
1645 const std::vector<int> &temp2
1646 = this->
base_element(b).adjust_line_dof_index_for_line_orientation_table;
1649 std::copy(temp2.begin(), temp2.end(),
1652 index += temp2.size();
1661 template <
int dim,
int spacedim>
1667 if (
base_element(b).hp_constraints_are_implemented() ==
false)
1675 template <
int dim,
int spacedim>
1709 interpolation_matrix = 0;
1713 unsigned int base_index = 0,
1714 base_index_other = 0;
1715 unsigned int multiplicity = 0,
1716 multiplicity_other = 0;
1724 &base_other = fe_other_system->
base_element(base_index_other);
1730 base_to_base_interpolation.
reinit (base_other.dofs_per_face,
1733 base_to_base_interpolation);
1741 std::make_pair (base_index, multiplicity))
1742 for (
unsigned int j=0; j<fe_other_system->
dofs_per_face; ++j)
1745 std::make_pair (base_index_other, multiplicity_other))
1746 interpolation_matrix(j, i)
1759 ++multiplicity_other;
1760 if (multiplicity_other ==
1763 multiplicity_other = 0;
1785 template <
int dim,
int spacedim>
1789 const unsigned int subface,
1820 interpolation_matrix = 0;
1824 unsigned int base_index = 0,
1825 base_index_other = 0;
1826 unsigned int multiplicity = 0,
1827 multiplicity_other = 0;
1835 &base_other = fe_other_system->
base_element(base_index_other);
1841 base_to_base_interpolation.
reinit (base_other.dofs_per_face,
1845 base_to_base_interpolation);
1853 std::make_pair (base_index, multiplicity))
1854 for (
unsigned int j=0; j<fe_other_system->
dofs_per_face; ++j)
1857 std::make_pair (base_index_other, multiplicity_other))
1858 interpolation_matrix(j, i)
1871 ++multiplicity_other;
1872 if (multiplicity_other ==
1875 multiplicity_other = 0;
1897 template <
int dim,
int spacedim>
1898 template <
int structdim>
1899 std::vector<std::pair<unsigned int, unsigned int> >
1919 unsigned int base_index = 0,
1920 base_index_other = 0;
1921 unsigned int multiplicity = 0,
1922 multiplicity_other = 0;
1926 unsigned int dof_offset = 0,
1927 dof_offset_other = 0;
1929 std::vector<std::pair<unsigned int, unsigned int> > identities;
1935 &base_other = fe_other_system->base_element(base_index_other);
1942 std::vector<std::pair<unsigned int, unsigned int> > base_identities;
1958 for (
unsigned int i=0; i<base_identities.size(); ++i)
1959 identities.push_back
1960 (std::make_pair (base_identities[i].first + dof_offset,
1961 base_identities[i].second + dof_offset_other));
1964 dof_offset += base.template n_dofs_per_object<structdim>();
1965 dof_offset_other += base_other.template n_dofs_per_object<structdim>();
1976 ++multiplicity_other;
1977 if (multiplicity_other ==
1978 fe_other_system->element_multiplicity(base_index_other))
1980 multiplicity_other = 0;
1988 Assert (base_index_other == fe_other_system->n_base_elements(),
1995 Assert (base_index_other != fe_other_system->n_base_elements(),
2004 return std::vector<std::pair<unsigned int, unsigned int> >();
2010 template <
int dim,
int spacedim>
2011 std::vector<std::pair<unsigned int, unsigned int> >
2014 return hp_object_dof_identities<0> (fe_other);
2017 template <
int dim,
int spacedim>
2018 std::vector<std::pair<unsigned int, unsigned int> >
2021 return hp_object_dof_identities<1> (fe_other);
2026 template <
int dim,
int spacedim>
2027 std::vector<std::pair<unsigned int, unsigned int> >
2030 return hp_object_dof_identities<2> (fe_other);
2035 template <
int dim,
int spacedim>
2057 fe_sys_other->base_element(b).n_components(),
2060 fe_sys_other->element_multiplicity(b),
2067 .compare_for_face_domination (fe_sys_other->base_element(b)));
2068 domination = domination & base_domination;
2080 template <
int dim,
int spacedim>
2091 template <
int dim,
int spacedim>
2094 const unsigned int face_index)
const 2097 .has_support_on_face(this->system_to_base_index(shape_index).second,
2103 template <
int dim,
int spacedim>
2112 typename FEL::ExcFEHasNoSupportPoints ());
2126 template <
int dim,
int spacedim>
2135 typename FEL::ExcFEHasNoSupportPoints ());
2149 template <
int dim,
int spacedim>
2150 std::pair<Table<2,bool>, std::vector<unsigned int> >
2160 std::pair<Table<2,bool>, std::vector<unsigned int> >
2168 const unsigned int comp = components.size();
2171 Table<2,bool> new_constant_modes(comp+base_table.first.n_rows()*
2173 constant_modes.n_cols());
2174 for (
unsigned int r=0; r<comp; ++r)
2176 new_constant_modes(r,c) = constant_modes(r,c);
2177 constant_modes.
swap(new_constant_modes);
2184 std::pair<std::pair<unsigned int,unsigned int>,
unsigned int> ind
2186 if (ind.first.first == i)
2187 for (
unsigned int c=0; c<base_table.first.n_rows(); ++c)
2188 constant_modes(comp+ind.first.second*base_table.first.n_rows()+c,k)
2189 = base_table.first(c,ind.second);
2192 for (
unsigned int c=0; c<base_table.second.size(); ++c)
2193 components.push_back(comp+r*this->base_elements[i].first->n_components()
2194 +base_table.second[c]);
2197 return std::pair<Table<2,bool>, std::vector<unsigned int> >(constant_modes,
2203 template <
int dim,
int spacedim>
2221 #include "fe_system.inst" 2223 DEAL_II_NAMESPACE_CLOSE
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
void initialize_unit_support_points()
virtual bool hp_constraints_are_implemented() 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::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual std::size_t memory_consumption() const
static const unsigned int invalid_unsigned_int
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
FiniteElement< dim, spacedim >::InternalDataBase & get_fe_data(const unsigned int base_no) const
std::vector< std::vector< FullMatrix< double > > > restriction
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
void build_interface_constraints()
const unsigned int components
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
#define AssertDimension(dim1, dim2)
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
void swap(TableBase< N, T > &v)
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
FullMatrix< double > interface_constraints
std::vector< std::pair< std_cxx11::shared_ptr< const FiniteElement< dim, spacedim > >, unsigned int > > base_elements
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
virtual FiniteElement< dim, spacedim > * clone() const
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const
const unsigned int dofs_per_quad
Table< 2, int > adjust_quad_dof_index_for_face_orientation_table
void initialize_quad_dof_index_permutation()
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
std::vector< internal::FEValues::FiniteElementRelatedData< dim, spacedim > > base_fe_output_objects
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
void set_fe_data(const unsigned int base_no, typename FiniteElement< dim, spacedim >::InternalDataBase *)
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual std::string get_name() const
TableIndices< 2 > interface_constraints_size() const
#define AssertThrow(cond, exc)
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > face_system_to_base_table
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false) const
bool is_primitive() const
static::ExceptionBase & ExcInterpolationNotImplemented()
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
const unsigned int dofs_per_line
InternalData(const unsigned int n_base_elements)
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
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
std::vector< Point< dim > > unit_support_points
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
unsigned int element_multiplicity(const unsigned int index) const
std::vector< std::pair< unsigned int, unsigned int > > hp_object_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const
unsigned int size() const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
static::ExceptionBase & ExcMessage(std::string arg1)
std::vector< std::vector< FullMatrix< double > > > prolongation
Third derivatives of shape functions.
void reinit(const TableIndices< N > &new_size, const bool omit_default_initialization=false)
#define Assert(cond, exc)
unsigned int n_nonzero_components(const unsigned int i) const
static const unsigned int invalid_face_number
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)
virtual Point< dim > unit_support_point(const unsigned int index) const
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
void initialize_unit_face_support_points()
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual Point< dim-1 > unit_face_support_point(const unsigned int index) const
virtual std::string get_name() const =0
unsigned int n_base_elements() const
unsigned int n_components() const
virtual std::size_t memory_consumption() const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
Second derivatives of shape functions.
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
std::string dim_string(const int dim, const int spacedim)
const unsigned int dofs_per_cell
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
void initialize(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities)
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
std_cxx11::enable_if< std_cxx11::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
size_type n_elements() const
virtual FiniteElement< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
void compute_fill(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim_1 > &quadrature, const CellSimilarity::Similarity cell_similarity, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_data, const internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
std::vector< typename FiniteElement< dim, spacedim >::InternalDataBase * > base_fe_datas
internal::FEValues::FiniteElementRelatedData< dim, spacedim > & get_fe_output_object(const unsigned int base_no) const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > face_system_to_base_index(const unsigned int index) 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::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
Shape function gradients.
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
FESystem(const FiniteElement< dim, spacedim > &fe, const unsigned int n_elements)
void push_back(const size_type size)
const unsigned int dofs_per_face
const std::vector< ComponentMask > nonzero_components
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
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
static::ExceptionBase & ExcNotImplemented()
const unsigned int dofs_per_vertex
std::vector< std::pair< unsigned int, unsigned int > > face_system_to_component_table
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
BlockIndices base_to_block_indices
unsigned int size(const unsigned int i) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Task< RT > new_task(const std_cxx11::function< RT()> &function)
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
static::ExceptionBase & ExcInternalError()