39 count_nonzeros(
const std::vector<unsigned int> &vec)
41 return std::count_if(vec.begin(), vec.end(), [](
const unsigned int i) {
54 template <
int dim,
int spacedim = dim>
57 const unsigned int base_no)
65 unsigned int out_index = 0;
66 for (
unsigned int system_index = 0; system_index < fe.
n_dofs_per_cell();
73 const unsigned int base_component =
75 const unsigned int base_index =
80 table[base_component][base_index] = out_index;
91 template <
int dim,
int spacedim = dim>
92 std::vector<typename FESystem<dim, spacedim>::BaseOffsets>
94 const unsigned int base_no)
96 std::vector<typename FESystem<dim, spacedim>::BaseOffsets> table;
99 unsigned int out_index = 0;
100 for (
unsigned int system_index = 0; system_index < fe.
n_dofs_per_cell();
105 const unsigned int base_index =
108 table.emplace_back();
110 table.back().n_nonzero_components =
112 unsigned int in_index = 0;
113 for (
unsigned int i = 0; i < base_index; ++i)
116 table.back().in_index = in_index;
117 table.back().out_index = out_index;
132 template <
int dim,
int spacedim = dim>
136 const unsigned int base_no,
146 const unsigned int n_dofs_per_cell =
149 auto copy_row = [](
const auto row_in,
auto row_out) {
150 std::copy(row_in.begin(), row_in.end(), row_out.begin());
154 for (
unsigned int component = 0; component < n_components; ++component)
155 for (
unsigned int b = 0; b < n_dofs_per_cell; ++b)
158 output_data.
shape_values[base_to_system_table[component][b]]);
161 for (
unsigned int component = 0; component < n_components; ++component)
162 for (
unsigned int b = 0; b < n_dofs_per_cell; ++b)
168 for (
unsigned int component = 0; component < n_components; ++component)
169 for (
unsigned int b = 0; b < n_dofs_per_cell; ++b)
175 for (
unsigned int component = 0; component < n_components; ++component)
176 for (
unsigned int b = 0; b < n_dofs_per_cell; ++b)
187 template <
int dim,
int spacedim = dim>
191 const unsigned int base_no,
192 const unsigned int n_q_points,
204 for (
const auto &offset : offsets)
207 for (
unsigned int s = 0; s < offset.n_nonzero_components; ++s)
208 for (
unsigned int q = 0; q < n_q_points; ++q)
213 for (
unsigned int s = 0; s < offset.n_nonzero_components; ++s)
214 for (
unsigned int q = 0; q < n_q_points; ++q)
219 for (
unsigned int s = 0; s < offset.n_nonzero_components; ++s)
220 for (
unsigned int q = 0; q < n_q_points; ++q)
225 for (
unsigned int s = 0; s < offset.n_nonzero_components; ++s)
226 for (
unsigned int q = 0; q < n_q_points; ++q)
236template <
int dim,
int spacedim>
238 const unsigned int n_base_elements)
239 : base_fe_datas(n_base_elements)
240 , base_fe_output_objects(n_base_elements)
245template <
int dim,
int spacedim>
254template <
int dim,
int spacedim>
257 const unsigned int base_no)
const
260 return *base_fe_datas[base_no];
265template <
int dim,
int spacedim>
268 const unsigned int base_no,
272 base_fe_datas[base_no] = std::move(ptr);
277template <
int dim,
int spacedim>
280 const unsigned int base_no)
const
283 return base_fe_output_objects[base_no];
291template <
int dim,
int spacedim>
295template <
int dim,
int spacedim>
297 const unsigned int n_elements)
305 std::vector<const FiniteElement<dim, spacedim> *> fes;
307 std::vector<unsigned int> multiplicities;
308 multiplicities.push_back(n_elements);
309 initialize(fes, multiplicities);
314template <
int dim,
int spacedim>
316 const unsigned int n1,
318 const unsigned int n2)
326 , base_elements(static_cast<
int>(n1 > 0) + static_cast<
int>(n2 > 0))
328 std::vector<const FiniteElement<dim, spacedim> *> fes;
331 std::vector<unsigned int> multiplicities;
332 multiplicities.push_back(n1);
333 multiplicities.push_back(n2);
334 initialize(fes, multiplicities);
339template <
int dim,
int spacedim>
341 const unsigned int n1,
343 const unsigned int n2,
345 const unsigned int n3)
360 , base_elements(static_cast<
int>(n1 > 0) + static_cast<
int>(n2 > 0) +
361 static_cast<
int>(n3 > 0))
363 std::vector<const FiniteElement<dim, spacedim> *> fes;
367 std::vector<unsigned int> multiplicities;
368 multiplicities.push_back(n1);
369 multiplicities.push_back(n2);
370 multiplicities.push_back(n3);
371 initialize(fes, multiplicities);
376template <
int dim,
int spacedim>
378 const unsigned int n1,
380 const unsigned int n2,
382 const unsigned int n3,
384 const unsigned int n4)
410 , base_elements(static_cast<
int>(n1 > 0) + static_cast<
int>(n2 > 0) +
411 static_cast<
int>(n3 > 0) + static_cast<
int>(n4 > 0))
413 std::vector<const FiniteElement<dim, spacedim> *> fes;
418 std::vector<unsigned int> multiplicities;
419 multiplicities.push_back(n1);
420 multiplicities.push_back(n2);
421 multiplicities.push_back(n3);
422 multiplicities.push_back(n4);
423 initialize(fes, multiplicities);
428template <
int dim,
int spacedim>
430 const unsigned int n1,
432 const unsigned int n2,
434 const unsigned int n3,
436 const unsigned int n4,
438 const unsigned int n5)
441 multiply_dof_numbers(&fe1, n1, &fe2, n2, &fe3, n3, &fe4, n4, &fe5, n5),
462 , base_elements(static_cast<
int>(n1 > 0) + static_cast<
int>(n2 > 0) +
463 static_cast<
int>(n3 > 0) + static_cast<
int>(n4 > 0) +
464 static_cast<
int>(n5 > 0))
466 std::vector<const FiniteElement<dim, spacedim> *> fes;
472 std::vector<unsigned int> multiplicities;
473 multiplicities.push_back(n1);
474 multiplicities.push_back(n2);
475 multiplicities.push_back(n3);
476 multiplicities.push_back(n4);
477 multiplicities.push_back(n5);
478 initialize(fes, multiplicities);
483template <
int dim,
int spacedim>
486 const std::vector<unsigned int> &multiplicities)
493 , base_elements(count_nonzeros(multiplicities))
495 initialize(fes, multiplicities);
500template <
int dim,
int spacedim>
511 std::ostringstream namebuf;
514 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
516 namebuf << base_element(i).get_name();
517 if (this->element_multiplicity(i) != 1)
518 namebuf <<
'^' << this->element_multiplicity(i);
519 if (i != this->n_base_elements() - 1)
524 return namebuf.str();
529template <
int dim,
int spacedim>
530std::unique_ptr<FiniteElement<dim, spacedim>>
533 std::vector<const FiniteElement<dim, spacedim> *> fes;
534 std::vector<unsigned int> multiplicities;
536 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
538 fes.push_back(&base_element(i));
539 multiplicities.push_back(this->element_multiplicity(i));
541 return std::make_unique<FESystem<dim, spacedim>>(fes, multiplicities);
546template <
int dim,
int spacedim>
549 const unsigned int first_component,
550 const unsigned int n_selected_components)
const
552 Assert(first_component + n_selected_components <= this->n_components(),
553 ExcMessage(
"Invalid arguments (not a part of this FiniteElement)."));
555 const unsigned int base_index =
556 this->component_to_base_table[first_component].first.first;
557 const unsigned int component_in_base =
558 this->component_to_base_table[first_component].first.second;
559 const unsigned int base_components =
560 this->base_element(base_index).n_components();
564 if (n_selected_components <= base_components)
565 return this->base_element(base_index)
566 .get_sub_fe(component_in_base, n_selected_components);
568 Assert(n_selected_components == this->n_components(),
569 ExcMessage(
"You can not select a part of a FiniteElement."));
575template <
int dim,
int spacedim>
581 Assert(this->is_primitive(i),
585 return (base_element(this->system_to_base_table[i].
first.first)
586 .shape_value(this->system_to_base_table[i].second, p));
591template <
int dim,
int spacedim>
594 const unsigned int i,
596 const unsigned int component)
const
603 if (this->nonzero_components[i][component] ==
false)
611 const unsigned int base = this->component_to_base_index(component).first;
612 const unsigned int component_in_base =
613 this->component_to_base_index(component).second;
621 return (base_element(base).shape_value_component(
622 this->system_to_base_table[i].
second, p, component_in_base));
627template <
int dim,
int spacedim>
633 Assert(this->is_primitive(i),
637 return (base_element(this->system_to_base_table[i].
first.first)
638 .shape_grad(this->system_to_base_table[i].second, p));
643template <
int dim,
int spacedim>
646 const unsigned int i,
648 const unsigned int component)
const
654 if (this->nonzero_components[i][component] ==
false)
659 const unsigned int base = this->component_to_base_index(component).first;
660 const unsigned int component_in_base =
661 this->component_to_base_index(component).second;
666 return (base_element(base).shape_grad_component(
667 this->system_to_base_table[i].
second, p, component_in_base));
672template <
int dim,
int spacedim>
678 Assert(this->is_primitive(i),
682 return (base_element(this->system_to_base_table[i].
first.first)
683 .shape_grad_grad(this->system_to_base_table[i].second, p));
688template <
int dim,
int spacedim>
691 const unsigned int i,
693 const unsigned int component)
const
699 if (this->nonzero_components[i][component] ==
false)
704 const unsigned int base = this->component_to_base_index(component).first;
705 const unsigned int component_in_base =
706 this->component_to_base_index(component).second;
711 return (base_element(base).shape_grad_grad_component(
712 this->system_to_base_table[i].
second, p, component_in_base));
717template <
int dim,
int spacedim>
723 Assert(this->is_primitive(i),
727 return (base_element(this->system_to_base_table[i].
first.first)
728 .shape_3rd_derivative(this->system_to_base_table[i].second, p));
733template <
int dim,
int spacedim>
736 const unsigned int i,
738 const unsigned int component)
const
744 if (this->nonzero_components[i][component] ==
false)
749 const unsigned int base = this->component_to_base_index(component).first;
750 const unsigned int component_in_base =
751 this->component_to_base_index(component).second;
756 return (base_element(base).shape_3rd_derivative_component(
757 this->system_to_base_table[i].
second, p, component_in_base));
762template <
int dim,
int spacedim>
768 Assert(this->is_primitive(i),
772 return (base_element(this->system_to_base_table[i].
first.first)
773 .shape_4th_derivative(this->system_to_base_table[i].second, p));
778template <
int dim,
int spacedim>
781 const unsigned int i,
783 const unsigned int component)
const
789 if (this->nonzero_components[i][component] ==
false)
794 const unsigned int base = this->component_to_base_index(component).first;
795 const unsigned int component_in_base =
796 this->component_to_base_index(component).second;
801 return (base_element(base).shape_4th_derivative_component(
802 this->system_to_base_table[i].
second, p, component_in_base));
807template <
int dim,
int spacedim>
817 Assert((interpolation_matrix.
m() == this->n_dofs_per_cell()) ||
820 this->n_dofs_per_cell()));
822 (this->n_dofs_per_cell() == 0),
832 (x_source_fe.
get_name().find(
"FESystem<") == 0) ||
846 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
850 spacedim>::ExcInterpolationNotImplemented()));
859 std::vector<FullMatrix<double>> base_matrices(this->n_base_elements());
860 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
862 base_matrices[i].reinit(base_element(i).n_dofs_per_cell(),
864 base_element(i).get_interpolation_matrix(source_fe.
base_element(i),
871 interpolation_matrix = 0;
872 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
874 if (this->system_to_base_table[i].
first ==
876 interpolation_matrix(i, j) =
877 (base_matrices[this->system_to_base_table[i].first.first](
878 this->system_to_base_table[i].second,
884template <
int dim,
int spacedim>
887 const unsigned int child,
894 "Restriction matrices are only available for refined cells!"));
900 if (this->restriction[refinement_case - 1][child].n() == 0)
902 std::lock_guard<std::mutex> lock(restriction_matrix_mutex);
905 if (this->restriction[refinement_case - 1][child].n() ==
906 this->n_dofs_per_cell())
907 return this->restriction[refinement_case - 1][child];
910 std::vector<const FullMatrix<double> *> base_matrices(
911 this->n_base_elements());
913 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
916 &base_element(i).get_restriction_matrix(child, refinement_case);
918 Assert(base_matrices[i]->n() == base_element(i).n_dofs_per_cell(),
923 this->n_dofs_per_cell());
933 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
934 for (
unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
941 if (this->system_to_base_table[i].
first !=
942 this->system_to_base_table[j].
first)
946 const unsigned int base = this->system_to_base_table[i].first.first;
948 const unsigned int base_index_i =
949 this->system_to_base_table[i].second,
951 this->system_to_base_table[j].second;
956 (*base_matrices[base])(base_index_i, base_index_j);
960 this->restriction[refinement_case - 1][child]) = std::move(restriction);
963 return this->restriction[refinement_case - 1][child];
968template <
int dim,
int spacedim>
971 const unsigned int child,
978 "Restriction matrices are only available for refined cells!"));
983 if (this->prolongation[refinement_case - 1][child].n() == 0)
985 std::lock_guard<std::mutex> lock(prolongation_matrix_mutex);
987 if (this->prolongation[refinement_case - 1][child].n() ==
988 this->n_dofs_per_cell())
989 return this->prolongation[refinement_case - 1][child];
991 std::vector<const FullMatrix<double> *> base_matrices(
992 this->n_base_elements());
993 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
996 &base_element(i).get_prolongation_matrix(child, refinement_case);
998 Assert(base_matrices[i]->n() == base_element(i).n_dofs_per_cell(),
1003 this->n_dofs_per_cell());
1005 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
1006 for (
unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
1008 if (this->system_to_base_table[i].
first !=
1009 this->system_to_base_table[j].
first)
1011 const unsigned int base = this->system_to_base_table[i].first.first;
1013 const unsigned int base_index_i =
1014 this->system_to_base_table[i].second,
1016 this->system_to_base_table[j].second;
1018 (*base_matrices[base])(base_index_i, base_index_j);
1022 this->prolongation[refinement_case - 1][child]) = std::move(prolongate);
1025 return this->prolongation[refinement_case - 1][child];
1029template <
int dim,
int spacedim>
1032 const unsigned int face_dof_index,
1033 const unsigned int face,
1034 const unsigned char combined_orientation)
const
1039 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int>
1040 face_base_index = this->face_system_to_base_index(face_dof_index, face);
1042 const unsigned int base_face_to_cell_index =
1043 this->base_element(face_base_index.first.first)
1044 .face_to_cell_index(face_base_index.second, face, combined_orientation);
1051 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int> target =
1052 std::make_pair(face_base_index.first, base_face_to_cell_index);
1053 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
1054 if (this->system_to_base_index(i) == target)
1069template <
int dim,
int spacedim>
1076 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1077 out |= base_element(base_no).requires_update_flags(flags);
1083template <
int dim,
int spacedim>
1084std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1100 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1101 data_ptr = std::make_unique<InternalData>(this->n_base_elements());
1102 auto &data =
dynamic_cast<InternalData &
>(*data_ptr);
1103 data.update_each = requires_update_flags(flags);
1117 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1120 &base_fe_output_object = data.get_fe_output_object(base_no);
1123 base_element(base_no),
1124 flags | base_element(base_no).requires_update_flags(flags));
1132 auto base_fe_data = base_element(base_no).get_data(flags,
1135 base_fe_output_object);
1137 data.set_fe_data(base_no, std::move(base_fe_data));
1146template <
int dim,
int spacedim>
1147std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1163 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1164 data_ptr = std::make_unique<InternalData>(this->n_base_elements());
1165 auto &data =
dynamic_cast<InternalData &
>(*data_ptr);
1166 data.update_each = requires_update_flags(flags);
1180 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1183 &base_fe_output_object = data.get_fe_output_object(base_no);
1186 base_element(base_no),
1187 flags | base_element(base_no).requires_update_flags(flags));
1195 auto base_fe_data = base_element(base_no).get_face_data(
1196 flags, mapping, quadrature, base_fe_output_object);
1198 data.set_fe_data(base_no, std::move(base_fe_data));
1209template <
int dim,
int spacedim>
1210std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1226 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1227 data_ptr = std::make_unique<InternalData>(this->n_base_elements());
1228 auto &data =
dynamic_cast<InternalData &
>(*data_ptr);
1230 data.update_each = requires_update_flags(flags);
1244 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1247 &base_fe_output_object = data.get_fe_output_object(base_no);
1250 base_element(base_no),
1251 flags | base_element(base_no).requires_update_flags(flags));
1259 auto base_fe_data = base_element(base_no).get_subface_data(
1260 flags, mapping, quadrature, base_fe_output_object);
1262 data.set_fe_data(base_no, std::move(base_fe_data));
1270template <
int dim,
int spacedim>
1285 compute_fill(mapping,
1287 invalid_face_number,
1288 invalid_face_number,
1299template <
int dim,
int spacedim>
1303 const unsigned int face_no,
1314 compute_fill(mapping,
1317 invalid_face_number,
1328template <
int dim,
int spacedim>
1332 const unsigned int face_no,
1333 const unsigned int sub_no,
1344 compute_fill(mapping,
1358template <
int dim,
int spacedim>
1359template <
class Q_or_QC>
1364 const unsigned int face_no,
1365 const unsigned int sub_no,
1366 const Q_or_QC &quadrature,
1379 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
1381 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
1398 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1402 fe_data.get_fe_data(base_no);
1405 &base_data = fe_data.get_fe_output_object(base_no);
1411 const Quadrature<dim - 1> *sub_face_quadrature =
nullptr;
1415 if (face_no == invalid_face_number)
1420 n_q_points = cell_quadrature->
size();
1422 else if (sub_no == invalid_face_number)
1431 (*face_quadrature)[face_quadrature->size() == 1 ? 0 : face_no]
1436 sub_face_quadrature =
1437 dynamic_cast<const Quadrature<dim - 1
> *>(&quadrature);
1440 n_q_points = sub_face_quadrature->size();
1452 if (face_no == invalid_face_number)
1461 else if (sub_no == invalid_face_number)
1474 *sub_face_quadrature,
1492 primitive_offset_tables[base_no],
1503 nonprimitive_offset_tables[base_no],
1511template <
int dim,
int spacedim>
1519 for (
unsigned int base = 0; base < this->n_base_elements(); ++base)
1520 if (base_element(base).constraints_are_implemented() ==
false)
1526 const unsigned int face_no = 0;
1528 this->interface_constraints.TableBase<2, double>::reinit(
1529 this->interface_constraints_size());
1539 for (
unsigned int n = 0; n < this->interface_constraints.n(); ++n)
1540 for (
unsigned int m = 0; m < this->interface_constraints.m(); ++m)
1549 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int>
1550 n_index = this->face_system_to_base_table[face_no][n];
1554 std::pair<std::pair<unsigned int, unsigned int>,
unsigned int> m_index;
1571 if (m < this->n_dofs_per_vertex())
1572 m_index = this->system_to_base_table[m];
1576 const unsigned int index_in_line =
1577 (m - this->n_dofs_per_vertex()) % this->n_dofs_per_line();
1578 const unsigned int sub_line =
1579 (m - this->n_dofs_per_vertex()) / this->n_dofs_per_line();
1586 const unsigned int tmp1 =
1587 2 * this->n_dofs_per_vertex() + index_in_line;
1589 this->face_system_to_base_table[face_no][tmp1].first;
1601 this->face_system_to_base_table[face_no][tmp1].
second >=
1603 base_element(m_index.first.first).n_dofs_per_vertex(),
1605 const unsigned int tmp2 =
1606 this->face_system_to_base_table[face_no][tmp1].second -
1607 2 * base_element(m_index.first.first).n_dofs_per_vertex();
1608 Assert(tmp2 < base_element(m_index.first.first)
1612 base_element(m_index.first.first).n_dofs_per_vertex() +
1613 base_element(m_index.first.first).n_dofs_per_line() *
1627 if (m < 5 * this->n_dofs_per_vertex())
1628 m_index = this->system_to_base_table[m];
1631 if (m < 5 * this->n_dofs_per_vertex() +
1632 12 * this->n_dofs_per_line())
1635 const unsigned int index_in_line =
1636 (m - 5 * this->n_dofs_per_vertex()) %
1637 this->n_dofs_per_line();
1638 const unsigned int sub_line =
1639 (m - 5 * this->n_dofs_per_vertex()) /
1640 this->n_dofs_per_line();
1643 const unsigned int tmp1 =
1644 4 * this->n_dofs_per_vertex() + index_in_line;
1646 this->face_system_to_base_table[face_no][tmp1].first;
1649 this->face_system_to_base_table[face_no][tmp1].
second >=
1650 4 * base_element(m_index.first.first)
1651 .n_dofs_per_vertex(),
1653 const unsigned int tmp2 =
1654 this->face_system_to_base_table[face_no][tmp1].second -
1656 base_element(m_index.first.first).n_dofs_per_vertex();
1657 Assert(tmp2 < base_element(m_index.first.first)
1661 5 * base_element(m_index.first.first)
1662 .n_dofs_per_vertex() +
1663 base_element(m_index.first.first).n_dofs_per_line() *
1671 const unsigned int index_in_quad =
1672 (m - 5 * this->n_dofs_per_vertex() -
1673 12 * this->n_dofs_per_line()) %
1674 this->n_dofs_per_quad(face_no);
1675 Assert(index_in_quad < this->n_dofs_per_quad(face_no),
1677 const unsigned int sub_quad =
1678 ((m - 5 * this->n_dofs_per_vertex() -
1679 12 * this->n_dofs_per_line()) /
1680 this->n_dofs_per_quad(face_no));
1683 const unsigned int tmp1 = 4 * this->n_dofs_per_vertex() +
1684 4 * this->n_dofs_per_line() +
1687 this->face_system_to_base_table[face_no].size(),
1690 this->face_system_to_base_table[face_no][tmp1].first;
1693 this->face_system_to_base_table[face_no][tmp1].
second >=
1694 4 * base_element(m_index.first.first)
1695 .n_dofs_per_vertex() +
1696 4 * base_element(m_index.first.first)
1699 const unsigned int tmp2 =
1700 this->face_system_to_base_table[face_no][tmp1].second -
1701 4 * base_element(m_index.first.first)
1702 .n_dofs_per_vertex() -
1703 4 * base_element(m_index.first.first).n_dofs_per_line();
1704 Assert(tmp2 < base_element(m_index.first.first)
1705 .n_dofs_per_quad(face_no),
1708 5 * base_element(m_index.first.first)
1709 .n_dofs_per_vertex() +
1711 base_element(m_index.first.first).n_dofs_per_line() +
1712 base_element(m_index.first.first)
1713 .n_dofs_per_quad(face_no) *
1728 if (n_index.first == m_index.first)
1729 this->interface_constraints(m, n) =
1730 (base_element(n_index.first.first)
1731 .constraints()(m_index.second, n_index.second));
1737template <
int dim,
int spacedim>
1741 const std::vector<unsigned int> &multiplicities)
1743 Assert(fes.size() == multiplicities.size(),
1746 ExcMessage(
"Need to pass at least one finite element."));
1747 Assert(count_nonzeros(multiplicities) > 0,
1748 ExcMessage(
"You only passed FiniteElements with multiplicity 0."));
1753 this->base_to_block_indices.reinit(0, 0);
1755 for (
unsigned int i = 0; i < fes.size(); ++i)
1756 if (multiplicities[i] > 0)
1757 this->base_to_block_indices.push_back(multiplicities[i]);
1762 unsigned int ind = 0;
1763 for (
unsigned int i = 0; i < fes.size(); ++i)
1764 if (multiplicities[i] > 0)
1767 base_elements[ind] = {fes[i]->clone(), multiplicities[i]};
1781 this->system_to_component_table.resize(this->n_dofs_per_cell());
1784 this->system_to_component_table,
1785 this->component_to_base_table,
1788 this->face_system_to_component_table.resize(this->n_unique_faces());
1790 for (
unsigned int face_no = 0; face_no < this->n_unique_faces(); ++face_no)
1792 this->face_system_to_component_table[face_no].resize(
1793 this->n_dofs_per_face(face_no));
1796 this->face_system_to_base_table[face_no],
1797 this->face_system_to_component_table[face_no],
1818 for (
unsigned int base_el = 0; base_el < this->n_base_elements(); ++base_el)
1819 if (!base_element(base_el).has_support_points() &&
1820 base_element(base_el).n_dofs_per_cell() != 0)
1829 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
1831 const unsigned int base = this->system_to_base_table[i].first.first,
1832 base_index = this->system_to_base_table[i].second;
1837 base_element(base).unit_support_points[base_index];
1842 primitive_offset_tables.resize(this->n_base_elements());
1844 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1845 if (base_element(base_no).is_primitive())
1846 primitive_offset_tables[base_no] =
1851 nonprimitive_offset_tables.resize(this->n_base_elements());
1853 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1854 if (!base_element(base_no).is_primitive())
1855 nonprimitive_offset_tables[base_no] =
1862 for (
unsigned int face_no = 0; face_no < this->n_unique_faces();
1875 bool flag_has_no_support_points =
false;
1877 for (
unsigned int base_el = 0; base_el < this->n_base_elements();
1879 if (!base_element(base_el).has_support_points() &&
1880 (base_element(base_el).n_dofs_per_face(face_no) > 0))
1882 this->unit_face_support_points[face_no].resize(0);
1883 flag_has_no_support_points =
true;
1888 if (flag_has_no_support_points)
1893 this->unit_face_support_points[face_no].resize(
1894 this->n_dofs_per_face(face_no));
1896 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
1898 const unsigned int base_i =
1899 this->face_system_to_base_table[face_no][i].first.first;
1900 const unsigned int index_in_base =
1901 this->face_system_to_base_table[face_no][i].second;
1905 base_element(base_i).unit_face_support_points[face_no].size(),
1908 this->unit_face_support_points[face_no][i] =
1909 base_element(base_i)
1910 .unit_face_support_points[face_no][index_in_base];
1923 generalized_support_points_index_table.resize(this->n_base_elements());
1925 for (
unsigned int base = 0; base < this->n_base_elements(); ++base)
1936 if (!base_element(base).has_generalized_support_points())
1939 for (
const auto &point :
1940 base_element(base).get_generalized_support_points())
1944 std::find(std::begin(this->generalized_support_points),
1945 std::end(this->generalized_support_points),
1948 if (p == std::end(this->generalized_support_points))
1951 const auto n = this->generalized_support_points.size();
1952 generalized_support_points_index_table[base].push_back(n);
1953 this->generalized_support_points.push_back(point);
1958 const auto n =
p - std::begin(this->generalized_support_points);
1959 generalized_support_points_index_table[base].push_back(n);
1966 for (
unsigned int i = 0; i < base_elements.size(); ++i)
1968 if (!base_element(i).has_generalized_support_points())
1971 const auto &points =
1972 base_elements[i].first->get_generalized_support_points();
1973 for (
unsigned int j = 0; j < points.size(); ++j)
1975 const auto n = generalized_support_points_index_table[i][j];
1976 Assert(this->generalized_support_points[n] == points[j],
1986 for (
unsigned int face_no = 0; face_no < this->n_unique_faces();
1991 Assert(this->adjust_quad_dof_index_for_face_orientation_table[face_no]
1994 this->n_dofs_per_quad(face_no),
1999 unsigned int index = 0;
2000 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2003 this->base_element(b)
2004 .adjust_quad_dof_index_for_face_orientation_table[face_no];
2005 for (
unsigned int c = 0; c < this->element_multiplicity(b); ++c)
2007 for (
unsigned int i = 0; i < temp.size(0); ++i)
2008 for (
unsigned int j = 0;
2012 this->adjust_quad_dof_index_for_face_orientation_table
2013 [face_no](index + i, j) = temp(i, j);
2014 index += temp.size(0);
2021 Assert(this->adjust_line_dof_index_for_line_orientation_table.size() ==
2022 this->n_dofs_per_line(),
2024 unsigned int index = 0;
2025 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2027 const std::vector<int> &temp2 =
2028 this->base_element(b)
2029 .adjust_line_dof_index_for_line_orientation_table;
2030 for (
unsigned int c = 0; c < this->element_multiplicity(b); ++c)
2035 this->adjust_line_dof_index_for_line_orientation_table.begin() +
2037 index += temp2.size();
2049template <
int dim,
int spacedim>
2053 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2054 if (base_element(b).hp_constraints_are_implemented() ==
false)
2062template <
int dim,
int spacedim>
2067 const unsigned int face_no)
const
2069 Assert(interpolation_matrix.
n() == this->n_dofs_per_face(face_no),
2071 this->n_dofs_per_face(face_no)));
2086 if (
const auto *fe_other_system =
2090 interpolation_matrix = 0;
2094 unsigned int base_index = 0, base_index_other = 0;
2095 unsigned int multiplicity = 0, multiplicity_other = 0;
2103 fe_other_system->base_element(
2110 base_to_base_interpolation.reinit(base_other.n_dofs_per_face(face_no),
2113 base_to_base_interpolation,
2119 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
2120 if (this->face_system_to_base_index(i, face_no).first ==
2121 std::make_pair(base_index, multiplicity))
2122 for (
unsigned int j = 0;
2123 j < fe_other_system->n_dofs_per_face(face_no);
2125 if (fe_other_system->face_system_to_base_index(j, face_no)
2127 std::make_pair(base_index_other, multiplicity_other))
2128 interpolation_matrix(j, i) = base_to_base_interpolation(
2129 fe_other_system->face_system_to_base_index(j, face_no)
2131 this->face_system_to_base_index(i, face_no).second);
2137 if (multiplicity == this->element_multiplicity(base_index))
2142 ++multiplicity_other;
2143 if (multiplicity_other ==
2144 fe_other_system->element_multiplicity(base_index_other))
2146 multiplicity_other = 0;
2152 if (base_index == this->n_base_elements())
2154 Assert(base_index_other == fe_other_system->n_base_elements(),
2161 Assert(base_index_other != fe_other_system->n_base_elements(),
2172 spacedim>::ExcInterpolationNotImplemented()));
2178template <
int dim,
int spacedim>
2182 const unsigned int subface,
2184 const unsigned int face_no)
const
2187 (x_source_fe.
get_name().find(
"FESystem<") == 0) ||
2191 Assert(interpolation_matrix.
n() == this->n_dofs_per_face(face_no),
2193 this->n_dofs_per_face(face_no)));
2210 if (fe_other_system !=
nullptr)
2213 interpolation_matrix = 0;
2217 unsigned int base_index = 0, base_index_other = 0;
2218 unsigned int multiplicity = 0, multiplicity_other = 0;
2233 base_to_base_interpolation.reinit(base_other.n_dofs_per_face(face_no),
2237 base_to_base_interpolation,
2243 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
2244 if (this->face_system_to_base_index(i, face_no).first ==
2245 std::make_pair(base_index, multiplicity))
2246 for (
unsigned int j = 0;
2251 std::make_pair(base_index_other, multiplicity_other))
2252 interpolation_matrix(j, i) = base_to_base_interpolation(
2255 this->face_system_to_base_index(i, face_no).second);
2261 if (multiplicity == this->element_multiplicity(base_index))
2266 ++multiplicity_other;
2267 if (multiplicity_other ==
2270 multiplicity_other = 0;
2276 if (base_index == this->n_base_elements())
2293 fe_other_system !=
nullptr,
2295 spacedim>::ExcInterpolationNotImplemented()));
2301template <
int dim,
int spacedim>
2302template <
int structdim>
2303std::vector<std::pair<unsigned int, unsigned int>>
2306 const unsigned int face_no)
const
2325 unsigned int base_index = 0, base_index_other = 0;
2326 unsigned int multiplicity = 0, multiplicity_other = 0;
2330 unsigned int dof_offset = 0, dof_offset_other = 0;
2332 std::vector<std::pair<unsigned int, unsigned int>> identities;
2346 std::vector<std::pair<unsigned int, unsigned int>> base_identities;
2363 for (
const auto &base_identity : base_identities)
2364 identities.emplace_back(base_identity.
first + dof_offset,
2365 base_identity.
second + dof_offset_other);
2368 dof_offset += base.template n_dofs_per_object<structdim>();
2370 base_other.template n_dofs_per_object<structdim>();
2376 if (multiplicity == this->element_multiplicity(base_index))
2381 ++multiplicity_other;
2382 if (multiplicity_other ==
2385 multiplicity_other = 0;
2391 if (base_index == this->n_base_elements())
2409 return std::vector<std::pair<unsigned int, unsigned int>>();
2415template <
int dim,
int spacedim>
2416std::vector<std::pair<unsigned int, unsigned int>>
2420 return hp_object_dof_identities<0>(fe_other);
2423template <
int dim,
int spacedim>
2424std::vector<std::pair<unsigned int, unsigned int>>
2428 return hp_object_dof_identities<1>(fe_other);
2433template <
int dim,
int spacedim>
2434std::vector<std::pair<unsigned int, unsigned int>>
2437 const unsigned int face_no)
const
2439 return hp_object_dof_identities<2>(fe_other, face_no);
2444template <
int dim,
int spacedim>
2448 const unsigned int codim)
const
2459 Assert(this->n_components() == fe_sys_other->n_components(),
2461 Assert(this->n_base_elements() == fe_sys_other->n_base_elements(),
2468 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2470 Assert(this->base_element(b).n_components() ==
2471 fe_sys_other->base_element(b).n_components(),
2473 Assert(this->element_multiplicity(b) ==
2474 fe_sys_other->element_multiplicity(b),
2480 (this->base_element(b).compare_for_domination(
2481 fe_sys_other->base_element(b), codim));
2482 domination = domination & base_domination;
2494template <
int dim,
int spacedim>
2499 return *base_elements[
index].first;
2504template <
int dim,
int spacedim>
2507 const unsigned int shape_index,
2508 const unsigned int face_index)
const
2510 return (base_element(this->system_to_base_index(shape_index).
first.first)
2511 .has_support_on_face(this->system_to_base_index(shape_index).second,
2517template <
int dim,
int spacedim>
2523 (this->unit_support_points.empty()),
2532 return (base_element(this->system_to_base_index(index).
first.first)
2533 .unit_support_point(this->system_to_base_index(index).second));
2538template <
int dim,
int spacedim>
2541 const unsigned int index,
2542 const unsigned int face_no)
const
2546 (this->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2547 .size() == this->n_dofs_per_face(face_no)) ||
2548 (this->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2550 (typename
FiniteElement<dim, spacedim>::ExcFEHasNoSupportPoints()));
2553 if (this->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2556 ->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2562 base_element(this->face_system_to_base_index(index, face_no).
first.first)
2563 .unit_face_support_point(
2564 this->face_system_to_base_index(index, face_no).second, face_no));
2569template <
int dim,
int spacedim>
2570std::pair<Table<2, bool>, std::vector<unsigned int>>
2576 Table<2, bool> constant_modes(this->n_components(), this->n_dofs_per_cell());
2577 std::vector<unsigned int> components;
2578 for (
unsigned int i = 0; i < base_elements.size(); ++i)
2580 const std::pair<Table<2, bool>, std::vector<unsigned int>> base_table =
2581 base_elements[i].first->get_constant_modes();
2583 const unsigned int element_multiplicity = this->element_multiplicity(i);
2588 const unsigned int comp = components.size();
2589 if (constant_modes.n_rows() <
2590 comp + base_table.first.n_rows() * element_multiplicity)
2592 Table<2, bool> new_constant_modes(comp + base_table.first.n_rows() *
2593 element_multiplicity,
2594 constant_modes.n_cols());
2595 for (
unsigned int r = 0; r < comp; ++r)
2596 for (
unsigned int c = 0; c < this->n_dofs_per_cell(); ++c)
2597 new_constant_modes(r, c) = constant_modes(r, c);
2599 constant_modes = std::move(new_constant_modes);
2604 for (
unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
2606 std::pair<std::pair<unsigned int, unsigned int>,
unsigned int> ind =
2607 this->system_to_base_index(k);
2608 if (ind.first.first == i)
2609 for (
unsigned int c = 0; c < base_table.first.n_rows(); ++c)
2610 constant_modes(comp +
2611 ind.first.second * base_table.first.n_rows() + c,
2612 k) = base_table.first(c, ind.second);
2614 for (
unsigned int r = 0; r < element_multiplicity; ++r)
2615 for (
const unsigned int c : base_table.
second)
2616 components.push_back(
2617 comp + r * this->base_elements[i].
first->n_components() + c);
2620 return std::pair<Table<2, bool>, std::vector<unsigned int>>(constant_modes,
2626template <
int dim,
int spacedim>
2630 std::vector<double> &dof_values)
const
2632 Assert(this->has_generalized_support_points(),
2633 ExcMessage(
"The FESystem does not have generalized support points"));
2636 this->get_generalized_support_points().size());
2639 std::vector<double> base_dof_values;
2640 std::vector<Vector<double>> base_point_values;
2645 unsigned int current_vector_component = 0;
2646 for (
unsigned int base = 0; base < base_elements.size(); ++base)
2651 const auto &base_element = this->base_element(base);
2652 const unsigned int multiplicity = this->element_multiplicity(base);
2653 const unsigned int n_base_dofs = base_element.n_dofs_per_cell();
2654 const unsigned int n_base_components = base_element.n_components();
2659 if (n_base_dofs == 0)
2661 current_vector_component += multiplicity * n_base_components;
2665 if (base_element.has_generalized_support_points())
2667 const unsigned int n_base_points =
2668 base_element.get_generalized_support_points().size();
2670 base_dof_values.resize(n_base_dofs);
2671 base_point_values.resize(n_base_points);
2673 for (
unsigned int m = 0; m < multiplicity;
2674 ++m, current_vector_component += n_base_components)
2678 for (
unsigned int j = 0; j < base_point_values.size(); ++j)
2680 base_point_values[j].reinit(n_base_components,
false);
2683 generalized_support_points_index_table[base][j];
2687 const auto *
const begin =
2688 std::begin(point_values[n]) + current_vector_component;
2689 const auto *
const end =
begin + n_base_components;
2690 std::copy(begin, end, std::begin(base_point_values[j]));
2694 .convert_generalized_support_point_values_to_dof_values(
2695 base_point_values, base_dof_values);
2704 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
2705 if (this->system_to_base_index(i).first ==
2706 std::make_pair(base, m))
2708 base_dof_values[this->system_to_base_index(i).second];
2720 for (
unsigned int m = 0; m < multiplicity; ++m)
2721 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
2722 if (this->system_to_base_index(i).first ==
2723 std::make_pair(base, m))
2724 dof_values[i] = std::numeric_limits<double>::signaling_NaN();
2726 current_vector_component += multiplicity * n_base_components;
2733template <
int dim,
int spacedim>
2741 sizeof(base_elements));
2742 for (
unsigned int i = 0; i < base_elements.size(); ++i)
2750#include "fe_system.inst"
void set_fe_data(const unsigned int base_no, std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase >)
InternalData(const unsigned int n_base_elements)
std::vector< std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > > base_fe_datas
FiniteElement< dim, spacedim >::InternalDataBase & get_fe_data(const unsigned int base_no) const
internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > & get_fe_output_object(const unsigned int base_no) const
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const override
virtual Point< dim > unit_support_point(const unsigned int index) const override
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
std::vector< std::pair< std::unique_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 override
virtual std::string get_name() const override
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override
virtual const FiniteElement< dim, spacedim > & get_sub_fe(const unsigned int first_component, const unsigned int n_selected_components) const override
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
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
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
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 Q_or_QC &quadrature, const CellSimilarity::Similarity cell_similarity, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_data, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const override
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_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 bool hp_constraints_are_implemented() const override
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
void build_interface_constraints()
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const unsigned char combined_orientation=ReferenceCell::default_combined_face_orientation()) const override
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index, const unsigned int face_no=0) const override
void initialize(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities)
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double > > &support_point_values, std::vector< double > &dof_values) const override
std::vector< std::pair< unsigned int, unsigned int > > hp_object_dof_identities(const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) 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 typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual double shape_value(const unsigned int i, const Point< dim > &p) 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 Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) 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::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
virtual std::size_t memory_consumption() const override
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const override
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
unsigned int n_dofs_per_cell() const
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_components() const
virtual std::string get_name() const =0
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
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
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
bool is_primitive() 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 void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > face_system_to_base_index(const unsigned int index, const unsigned int face_no=0) const
unsigned int element_multiplicity(const unsigned int index) const
unsigned int n_nonzero_components(const unsigned int i) const
unsigned int n_base_elements() const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const
virtual std::size_t memory_consumption() const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Abstract base class for mapping classes.
unsigned int size() const
unsigned int max_n_quadrature_points() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
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_3rd_derivatives
Third derivatives of shape functions.
@ update_gradients
Shape function gradients.
@ update_default
No update.
Task< RT > new_task(const std::function< RT()> &function)
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
@ neither_element_dominates
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
std::enable_if_t< std::is_fundamental_v< T >, std::size_t > memory_consumption(const T &t)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
std::string dim_string(const int dim, const int spacedim)
std::vector< Point< dim > > unit_support_points(const std::vector< Point< 1 > > &line_support_points, const std::vector< unsigned int > &renumbering)
std::vector< typename FESystem< dim, spacedim >::BaseOffsets > setup_nonprimitive_offset_table(const FESystem< dim, spacedim > &fe, const unsigned int base_no)
Table< 2, unsigned int > setup_primitive_offset_table(const FESystem< dim, spacedim > &fe, const unsigned int base_no)
void copy_nonprimitive_base_element_values(const FESystem< dim, spacedim > &fe, const unsigned int base_no, const unsigned int n_q_points, const UpdateFlags base_flags, const std::vector< typename FESystem< dim, spacedim >::BaseOffsets > &offsets, const FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &base_data, FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data)
void copy_primitive_base_element_values(const FESystem< dim, spacedim > &fe, const unsigned int base_no, const UpdateFlags base_flags, const Table< 2, unsigned int > &base_to_system_table, const FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &base_data, FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data)
static const unsigned int invalid_unsigned_int
const InputIterator OutputIterator out
const Iterator const std_cxx20::type_identity_t< Iterator > & end