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!"));
898 if (this->restriction[refinement_case - 1][child].n() == 0)
900 std::lock_guard<std::mutex> lock(restriction_matrix_mutex);
903 if (this->restriction[refinement_case - 1][child].n() ==
904 this->n_dofs_per_cell())
905 return this->restriction[refinement_case - 1][child];
908 std::vector<const FullMatrix<double> *> base_matrices(
909 this->n_base_elements());
911 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
914 &base_element(i).get_restriction_matrix(child, refinement_case);
916 Assert(base_matrices[i]->n() == base_element(i).n_dofs_per_cell(),
921 this->n_dofs_per_cell());
931 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
932 for (
unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
939 if (this->system_to_base_table[i].
first !=
940 this->system_to_base_table[j].
first)
944 const unsigned int base = this->system_to_base_table[i].first.first;
946 const unsigned int base_index_i =
947 this->system_to_base_table[i].second,
949 this->system_to_base_table[j].second;
954 (*base_matrices[base])(base_index_i, base_index_j);
958 this->restriction[refinement_case - 1][child]) = std::move(restriction);
961 return this->restriction[refinement_case - 1][child];
966template <
int dim,
int spacedim>
969 const unsigned int child,
976 "Restriction matrices are only available for refined cells!"));
981 if (this->prolongation[refinement_case - 1][child].n() == 0)
983 std::lock_guard<std::mutex> lock(prolongation_matrix_mutex);
985 if (this->prolongation[refinement_case - 1][child].n() ==
986 this->n_dofs_per_cell())
987 return this->prolongation[refinement_case - 1][child];
989 std::vector<const FullMatrix<double> *> base_matrices(
990 this->n_base_elements());
991 for (
unsigned int i = 0; i < this->n_base_elements(); ++i)
994 &base_element(i).get_prolongation_matrix(child, refinement_case);
996 Assert(base_matrices[i]->n() == base_element(i).n_dofs_per_cell(),
1001 this->n_dofs_per_cell());
1003 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
1004 for (
unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
1006 if (this->system_to_base_table[i].
first !=
1007 this->system_to_base_table[j].
first)
1009 const unsigned int base = this->system_to_base_table[i].first.first;
1011 const unsigned int base_index_i =
1012 this->system_to_base_table[i].second,
1014 this->system_to_base_table[j].second;
1016 (*base_matrices[base])(base_index_i, base_index_j);
1020 this->prolongation[refinement_case - 1][child]) = std::move(prolongate);
1023 return this->prolongation[refinement_case - 1][child];
1027template <
int dim,
int spacedim>
1030 const unsigned int face_dof_index,
1031 const unsigned int face,
1032 const unsigned char combined_orientation)
const
1037 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int>
1038 face_base_index = this->face_system_to_base_index(face_dof_index, face);
1040 const unsigned int base_face_to_cell_index =
1041 this->base_element(face_base_index.first.first)
1042 .face_to_cell_index(face_base_index.second, face, combined_orientation);
1049 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int> target =
1050 std::make_pair(face_base_index.first, base_face_to_cell_index);
1051 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
1052 if (this->system_to_base_index(i) == target)
1067template <
int dim,
int spacedim>
1074 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1075 out |= base_element(base_no).requires_update_flags(flags);
1081template <
int dim,
int spacedim>
1082std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1098 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1099 data_ptr = std::make_unique<InternalData>(this->n_base_elements());
1100 auto &data =
dynamic_cast<InternalData &
>(*data_ptr);
1101 data.update_each = requires_update_flags(flags);
1115 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1118 &base_fe_output_object = data.get_fe_output_object(base_no);
1121 base_element(base_no),
1122 flags | base_element(base_no).requires_update_flags(flags));
1130 auto base_fe_data = base_element(base_no).get_data(flags,
1133 base_fe_output_object);
1135 data.set_fe_data(base_no, std::move(base_fe_data));
1144template <
int dim,
int spacedim>
1145std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1161 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1162 data_ptr = std::make_unique<InternalData>(this->n_base_elements());
1163 auto &data =
dynamic_cast<InternalData &
>(*data_ptr);
1164 data.update_each = requires_update_flags(flags);
1178 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1181 &base_fe_output_object = data.get_fe_output_object(base_no);
1184 base_element(base_no),
1185 flags | base_element(base_no).requires_update_flags(flags));
1193 auto base_fe_data = base_element(base_no).get_face_data(
1194 flags, mapping, quadrature, base_fe_output_object);
1196 data.set_fe_data(base_no, std::move(base_fe_data));
1207template <
int dim,
int spacedim>
1208std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1224 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1225 data_ptr = std::make_unique<InternalData>(this->n_base_elements());
1226 auto &data =
dynamic_cast<InternalData &
>(*data_ptr);
1228 data.update_each = requires_update_flags(flags);
1242 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1245 &base_fe_output_object = data.get_fe_output_object(base_no);
1248 base_element(base_no),
1249 flags | base_element(base_no).requires_update_flags(flags));
1257 auto base_fe_data = base_element(base_no).get_subface_data(
1258 flags, mapping, quadrature, base_fe_output_object);
1260 data.set_fe_data(base_no, std::move(base_fe_data));
1268template <
int dim,
int spacedim>
1283 compute_fill(mapping,
1285 invalid_face_number,
1286 invalid_face_number,
1297template <
int dim,
int spacedim>
1301 const unsigned int face_no,
1312 compute_fill(mapping,
1315 invalid_face_number,
1326template <
int dim,
int spacedim>
1330 const unsigned int face_no,
1331 const unsigned int sub_no,
1342 compute_fill(mapping,
1356template <
int dim,
int spacedim>
1357template <
class Q_or_QC>
1362 const unsigned int face_no,
1363 const unsigned int sub_no,
1364 const Q_or_QC &quadrature,
1377 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
1379 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
1396 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1400 fe_data.get_fe_data(base_no);
1403 &base_data = fe_data.get_fe_output_object(base_no);
1409 const Quadrature<dim - 1> *sub_face_quadrature =
nullptr;
1413 if (face_no == invalid_face_number)
1418 n_q_points = cell_quadrature->
size();
1420 else if (sub_no == invalid_face_number)
1429 (*face_quadrature)[face_quadrature->size() == 1 ? 0 : face_no]
1434 sub_face_quadrature =
1435 dynamic_cast<const Quadrature<dim - 1
> *>(&quadrature);
1438 n_q_points = sub_face_quadrature->size();
1450 if (face_no == invalid_face_number)
1459 else if (sub_no == invalid_face_number)
1472 *sub_face_quadrature,
1490 primitive_offset_tables[base_no],
1501 nonprimitive_offset_tables[base_no],
1509template <
int dim,
int spacedim>
1517 for (
unsigned int base = 0; base < this->n_base_elements(); ++base)
1518 if (base_element(base).constraints_are_implemented() ==
false)
1524 const unsigned int face_no = 0;
1526 this->interface_constraints.TableBase<2, double>::reinit(
1527 this->interface_constraints_size());
1537 for (
unsigned int n = 0; n < this->interface_constraints.n(); ++n)
1538 for (
unsigned int m = 0; m < this->interface_constraints.m(); ++m)
1547 const std::pair<std::pair<unsigned int, unsigned int>,
unsigned int>
1548 n_index = this->face_system_to_base_table[face_no][n];
1552 std::pair<std::pair<unsigned int, unsigned int>,
unsigned int> m_index;
1569 if (m < this->n_dofs_per_vertex())
1570 m_index = this->system_to_base_table[m];
1574 const unsigned int index_in_line =
1575 (m - this->n_dofs_per_vertex()) % this->n_dofs_per_line();
1576 const unsigned int sub_line =
1577 (m - this->n_dofs_per_vertex()) / this->n_dofs_per_line();
1584 const unsigned int tmp1 =
1585 2 * this->n_dofs_per_vertex() + index_in_line;
1587 this->face_system_to_base_table[face_no][tmp1].first;
1599 this->face_system_to_base_table[face_no][tmp1].
second >=
1601 base_element(m_index.first.first).n_dofs_per_vertex(),
1603 const unsigned int tmp2 =
1604 this->face_system_to_base_table[face_no][tmp1].second -
1605 2 * base_element(m_index.first.first).n_dofs_per_vertex();
1606 Assert(tmp2 < base_element(m_index.first.first)
1610 base_element(m_index.first.first).n_dofs_per_vertex() +
1611 base_element(m_index.first.first).n_dofs_per_line() *
1625 if (m < 5 * this->n_dofs_per_vertex())
1626 m_index = this->system_to_base_table[m];
1629 if (m < 5 * this->n_dofs_per_vertex() +
1630 12 * this->n_dofs_per_line())
1633 const unsigned int index_in_line =
1634 (m - 5 * this->n_dofs_per_vertex()) %
1635 this->n_dofs_per_line();
1636 const unsigned int sub_line =
1637 (m - 5 * this->n_dofs_per_vertex()) /
1638 this->n_dofs_per_line();
1641 const unsigned int tmp1 =
1642 4 * this->n_dofs_per_vertex() + index_in_line;
1644 this->face_system_to_base_table[face_no][tmp1].first;
1647 this->face_system_to_base_table[face_no][tmp1].
second >=
1648 4 * base_element(m_index.first.first)
1649 .n_dofs_per_vertex(),
1651 const unsigned int tmp2 =
1652 this->face_system_to_base_table[face_no][tmp1].second -
1654 base_element(m_index.first.first).n_dofs_per_vertex();
1655 Assert(tmp2 < base_element(m_index.first.first)
1659 5 * base_element(m_index.first.first)
1660 .n_dofs_per_vertex() +
1661 base_element(m_index.first.first).n_dofs_per_line() *
1669 const unsigned int index_in_quad =
1670 (m - 5 * this->n_dofs_per_vertex() -
1671 12 * this->n_dofs_per_line()) %
1672 this->n_dofs_per_quad(face_no);
1673 Assert(index_in_quad < this->n_dofs_per_quad(face_no),
1675 const unsigned int sub_quad =
1676 ((m - 5 * this->n_dofs_per_vertex() -
1677 12 * this->n_dofs_per_line()) /
1678 this->n_dofs_per_quad(face_no));
1681 const unsigned int tmp1 = 4 * this->n_dofs_per_vertex() +
1682 4 * this->n_dofs_per_line() +
1685 this->face_system_to_base_table[face_no].size(),
1688 this->face_system_to_base_table[face_no][tmp1].first;
1691 this->face_system_to_base_table[face_no][tmp1].
second >=
1692 4 * base_element(m_index.first.first)
1693 .n_dofs_per_vertex() +
1694 4 * base_element(m_index.first.first)
1697 const unsigned int tmp2 =
1698 this->face_system_to_base_table[face_no][tmp1].second -
1699 4 * base_element(m_index.first.first)
1700 .n_dofs_per_vertex() -
1701 4 * base_element(m_index.first.first).n_dofs_per_line();
1702 Assert(tmp2 < base_element(m_index.first.first)
1703 .n_dofs_per_quad(face_no),
1706 5 * base_element(m_index.first.first)
1707 .n_dofs_per_vertex() +
1709 base_element(m_index.first.first).n_dofs_per_line() +
1710 base_element(m_index.first.first)
1711 .n_dofs_per_quad(face_no) *
1726 if (n_index.first == m_index.first)
1727 this->interface_constraints(m, n) =
1728 (base_element(n_index.first.first)
1729 .constraints()(m_index.second, n_index.second));
1735template <
int dim,
int spacedim>
1739 const std::vector<unsigned int> &multiplicities)
1741 Assert(fes.size() == multiplicities.size(),
1744 ExcMessage(
"Need to pass at least one finite element."));
1745 Assert(count_nonzeros(multiplicities) > 0,
1746 ExcMessage(
"You only passed FiniteElements with multiplicity 0."));
1751 this->base_to_block_indices.reinit(0, 0);
1753 for (
unsigned int i = 0; i < fes.size(); ++i)
1754 if (multiplicities[i] > 0)
1755 this->base_to_block_indices.push_back(multiplicities[i]);
1760 unsigned int ind = 0;
1761 for (
unsigned int i = 0; i < fes.size(); ++i)
1762 if (multiplicities[i] > 0)
1765 base_elements[ind] = {fes[i]->clone(), multiplicities[i]};
1779 this->system_to_component_table.resize(this->n_dofs_per_cell());
1782 this->system_to_component_table,
1783 this->component_to_base_table,
1786 this->face_system_to_component_table.resize(this->n_unique_faces());
1788 for (
unsigned int face_no = 0; face_no < this->n_unique_faces(); ++face_no)
1790 this->face_system_to_component_table[face_no].resize(
1791 this->n_dofs_per_face(face_no));
1794 this->face_system_to_base_table[face_no],
1795 this->face_system_to_component_table[face_no],
1816 for (
unsigned int base_el = 0; base_el < this->n_base_elements(); ++base_el)
1817 if (!base_element(base_el).has_support_points() &&
1818 base_element(base_el).n_dofs_per_cell() != 0)
1827 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
1829 const unsigned int base = this->system_to_base_table[i].first.first,
1830 base_index = this->system_to_base_table[i].second;
1835 base_element(base).unit_support_points[base_index];
1840 primitive_offset_tables.resize(this->n_base_elements());
1842 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1843 if (base_element(base_no).is_primitive())
1844 primitive_offset_tables[base_no] =
1849 nonprimitive_offset_tables.resize(this->n_base_elements());
1851 for (
unsigned int base_no = 0; base_no < this->n_base_elements(); ++base_no)
1852 if (!base_element(base_no).is_primitive())
1853 nonprimitive_offset_tables[base_no] =
1860 for (
unsigned int face_no = 0; face_no < this->n_unique_faces();
1873 bool flag_has_no_support_points =
false;
1875 for (
unsigned int base_el = 0; base_el < this->n_base_elements();
1877 if (!base_element(base_el).has_support_points() &&
1878 (base_element(base_el).n_dofs_per_face(face_no) > 0))
1880 this->unit_face_support_points[face_no].resize(0);
1881 flag_has_no_support_points =
true;
1886 if (flag_has_no_support_points)
1891 this->unit_face_support_points[face_no].resize(
1892 this->n_dofs_per_face(face_no));
1894 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
1896 const unsigned int base_i =
1897 this->face_system_to_base_table[face_no][i].first.first;
1898 const unsigned int index_in_base =
1899 this->face_system_to_base_table[face_no][i].second;
1903 base_element(base_i).unit_face_support_points[face_no].size(),
1906 this->unit_face_support_points[face_no][i] =
1907 base_element(base_i)
1908 .unit_face_support_points[face_no][index_in_base];
1921 generalized_support_points_index_table.resize(this->n_base_elements());
1923 for (
unsigned int base = 0; base < this->n_base_elements(); ++base)
1934 if (!base_element(base).has_generalized_support_points())
1937 for (
const auto &point :
1938 base_element(base).get_generalized_support_points())
1942 std::find(std::begin(this->generalized_support_points),
1943 std::end(this->generalized_support_points),
1946 if (p == std::end(this->generalized_support_points))
1949 const auto n = this->generalized_support_points.size();
1950 generalized_support_points_index_table[base].push_back(n);
1951 this->generalized_support_points.push_back(point);
1956 const auto n =
p - std::begin(this->generalized_support_points);
1957 generalized_support_points_index_table[base].push_back(n);
1964 for (
unsigned int i = 0; i < base_elements.size(); ++i)
1966 if (!base_element(i).has_generalized_support_points())
1969 const auto &points =
1970 base_elements[i].first->get_generalized_support_points();
1971 for (
unsigned int j = 0; j < points.size(); ++j)
1973 const auto n = generalized_support_points_index_table[i][j];
1974 Assert(this->generalized_support_points[n] == points[j],
1984 for (
unsigned int face_no = 0; face_no < this->n_unique_faces();
1989 Assert(this->adjust_quad_dof_index_for_face_orientation_table[face_no]
1992 this->n_dofs_per_quad(face_no),
1997 unsigned int index = 0;
1998 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2001 this->base_element(b)
2002 .adjust_quad_dof_index_for_face_orientation_table[face_no];
2003 for (
unsigned int c = 0; c < this->element_multiplicity(b); ++c)
2005 for (
unsigned int i = 0; i < temp.size(0); ++i)
2006 for (
unsigned int j = 0;
2010 this->adjust_quad_dof_index_for_face_orientation_table
2011 [face_no](index + i, j) = temp(i, j);
2012 index += temp.size(0);
2019 Assert(this->adjust_line_dof_index_for_line_orientation_table.size() ==
2020 this->n_dofs_per_line(),
2022 unsigned int index = 0;
2023 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2025 const std::vector<int> &temp2 =
2026 this->base_element(b)
2027 .adjust_line_dof_index_for_line_orientation_table;
2028 for (
unsigned int c = 0; c < this->element_multiplicity(b); ++c)
2033 this->adjust_line_dof_index_for_line_orientation_table.begin() +
2035 index += temp2.size();
2047template <
int dim,
int spacedim>
2051 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2052 if (base_element(b).hp_constraints_are_implemented() ==
false)
2060template <
int dim,
int spacedim>
2065 const unsigned int face_no)
const
2067 Assert(interpolation_matrix.
n() == this->n_dofs_per_face(face_no),
2069 this->n_dofs_per_face(face_no)));
2084 if (
const auto *fe_other_system =
2088 interpolation_matrix = 0;
2092 unsigned int base_index = 0, base_index_other = 0;
2093 unsigned int multiplicity = 0, multiplicity_other = 0;
2101 fe_other_system->base_element(
2108 base_to_base_interpolation.reinit(base_other.n_dofs_per_face(face_no),
2111 base_to_base_interpolation,
2117 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
2118 if (this->face_system_to_base_index(i, face_no).first ==
2119 std::make_pair(base_index, multiplicity))
2120 for (
unsigned int j = 0;
2121 j < fe_other_system->n_dofs_per_face(face_no);
2123 if (fe_other_system->face_system_to_base_index(j, face_no)
2125 std::make_pair(base_index_other, multiplicity_other))
2126 interpolation_matrix(j, i) = base_to_base_interpolation(
2127 fe_other_system->face_system_to_base_index(j, face_no)
2129 this->face_system_to_base_index(i, face_no).second);
2135 if (multiplicity == this->element_multiplicity(base_index))
2140 ++multiplicity_other;
2141 if (multiplicity_other ==
2142 fe_other_system->element_multiplicity(base_index_other))
2144 multiplicity_other = 0;
2150 if (base_index == this->n_base_elements())
2152 Assert(base_index_other == fe_other_system->n_base_elements(),
2159 Assert(base_index_other != fe_other_system->n_base_elements(),
2170 spacedim>::ExcInterpolationNotImplemented()));
2176template <
int dim,
int spacedim>
2180 const unsigned int subface,
2182 const unsigned int face_no)
const
2185 (x_source_fe.
get_name().find(
"FESystem<") == 0) ||
2189 Assert(interpolation_matrix.
n() == this->n_dofs_per_face(face_no),
2191 this->n_dofs_per_face(face_no)));
2208 if (fe_other_system !=
nullptr)
2211 interpolation_matrix = 0;
2215 unsigned int base_index = 0, base_index_other = 0;
2216 unsigned int multiplicity = 0, multiplicity_other = 0;
2231 base_to_base_interpolation.reinit(base_other.n_dofs_per_face(face_no),
2235 base_to_base_interpolation,
2241 for (
unsigned int i = 0; i < this->n_dofs_per_face(face_no); ++i)
2242 if (this->face_system_to_base_index(i, face_no).first ==
2243 std::make_pair(base_index, multiplicity))
2244 for (
unsigned int j = 0;
2249 std::make_pair(base_index_other, multiplicity_other))
2250 interpolation_matrix(j, i) = base_to_base_interpolation(
2253 this->face_system_to_base_index(i, face_no).second);
2259 if (multiplicity == this->element_multiplicity(base_index))
2264 ++multiplicity_other;
2265 if (multiplicity_other ==
2268 multiplicity_other = 0;
2274 if (base_index == this->n_base_elements())
2291 fe_other_system !=
nullptr,
2293 spacedim>::ExcInterpolationNotImplemented()));
2299template <
int dim,
int spacedim>
2300template <
int structdim>
2301std::vector<std::pair<unsigned int, unsigned int>>
2304 const unsigned int face_no)
const
2323 unsigned int base_index = 0, base_index_other = 0;
2324 unsigned int multiplicity = 0, multiplicity_other = 0;
2328 unsigned int dof_offset = 0, dof_offset_other = 0;
2330 std::vector<std::pair<unsigned int, unsigned int>> identities;
2344 std::vector<std::pair<unsigned int, unsigned int>> base_identities;
2361 for (
const auto &base_identity : base_identities)
2362 identities.emplace_back(base_identity.
first + dof_offset,
2363 base_identity.
second + dof_offset_other);
2366 dof_offset += base.template n_dofs_per_object<structdim>();
2368 base_other.template n_dofs_per_object<structdim>();
2374 if (multiplicity == this->element_multiplicity(base_index))
2379 ++multiplicity_other;
2380 if (multiplicity_other ==
2383 multiplicity_other = 0;
2389 if (base_index == this->n_base_elements())
2407 return std::vector<std::pair<unsigned int, unsigned int>>();
2413template <
int dim,
int spacedim>
2414std::vector<std::pair<unsigned int, unsigned int>>
2418 return hp_object_dof_identities<0>(fe_other);
2421template <
int dim,
int spacedim>
2422std::vector<std::pair<unsigned int, unsigned int>>
2426 return hp_object_dof_identities<1>(fe_other);
2431template <
int dim,
int spacedim>
2432std::vector<std::pair<unsigned int, unsigned int>>
2435 const unsigned int face_no)
const
2437 return hp_object_dof_identities<2>(fe_other, face_no);
2442template <
int dim,
int spacedim>
2446 const unsigned int codim)
const
2457 Assert(this->n_components() == fe_sys_other->n_components(),
2459 Assert(this->n_base_elements() == fe_sys_other->n_base_elements(),
2466 for (
unsigned int b = 0;
b < this->n_base_elements(); ++
b)
2468 Assert(this->base_element(b).n_components() ==
2469 fe_sys_other->base_element(b).n_components(),
2471 Assert(this->element_multiplicity(b) ==
2472 fe_sys_other->element_multiplicity(b),
2478 (this->base_element(b).compare_for_domination(
2479 fe_sys_other->base_element(b), codim));
2480 domination = domination & base_domination;
2492template <
int dim,
int spacedim>
2497 return *base_elements[
index].first;
2502template <
int dim,
int spacedim>
2505 const unsigned int shape_index,
2506 const unsigned int face_index)
const
2508 return (base_element(this->system_to_base_index(shape_index).
first.first)
2509 .has_support_on_face(this->system_to_base_index(shape_index).second,
2515template <
int dim,
int spacedim>
2521 (this->unit_support_points.empty()),
2530 return (base_element(this->system_to_base_index(index).
first.first)
2531 .unit_support_point(this->system_to_base_index(index).second));
2536template <
int dim,
int spacedim>
2539 const unsigned int index,
2540 const unsigned int face_no)
const
2544 (this->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2545 .size() == this->n_dofs_per_face(face_no)) ||
2546 (this->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2548 (typename
FiniteElement<dim, spacedim>::ExcFEHasNoSupportPoints()));
2551 if (this->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2554 ->unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
2560 base_element(this->face_system_to_base_index(index, face_no).
first.first)
2561 .unit_face_support_point(
2562 this->face_system_to_base_index(index, face_no).second, face_no));
2567template <
int dim,
int spacedim>
2568std::pair<Table<2, bool>, std::vector<unsigned int>>
2574 Table<2, bool> constant_modes(this->n_components(), this->n_dofs_per_cell());
2575 std::vector<unsigned int> components;
2576 for (
unsigned int i = 0; i < base_elements.size(); ++i)
2578 const std::pair<Table<2, bool>, std::vector<unsigned int>> base_table =
2579 base_elements[i].first->get_constant_modes();
2581 const unsigned int element_multiplicity = this->element_multiplicity(i);
2586 const unsigned int comp = components.size();
2587 if (constant_modes.n_rows() <
2588 comp + base_table.first.n_rows() * element_multiplicity)
2590 Table<2, bool> new_constant_modes(comp + base_table.first.n_rows() *
2591 element_multiplicity,
2592 constant_modes.n_cols());
2593 for (
unsigned int r = 0; r < comp; ++r)
2594 for (
unsigned int c = 0; c < this->n_dofs_per_cell(); ++c)
2595 new_constant_modes(r, c) = constant_modes(r, c);
2597 constant_modes = std::move(new_constant_modes);
2602 for (
unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
2604 std::pair<std::pair<unsigned int, unsigned int>,
unsigned int> ind =
2605 this->system_to_base_index(k);
2606 if (ind.first.first == i)
2607 for (
unsigned int c = 0; c < base_table.first.n_rows(); ++c)
2608 constant_modes(comp +
2609 ind.first.second * base_table.first.n_rows() + c,
2610 k) = base_table.first(c, ind.second);
2612 for (
unsigned int r = 0; r < element_multiplicity; ++r)
2613 for (
const unsigned int c : base_table.
second)
2614 components.push_back(
2615 comp + r * this->base_elements[i].
first->n_components() + c);
2618 return std::pair<Table<2, bool>, std::vector<unsigned int>>(constant_modes,
2624template <
int dim,
int spacedim>
2628 std::vector<double> &dof_values)
const
2630 Assert(this->has_generalized_support_points(),
2631 ExcMessage(
"The FESystem does not have generalized support points"));
2634 this->get_generalized_support_points().size());
2637 std::vector<double> base_dof_values;
2638 std::vector<Vector<double>> base_point_values;
2643 unsigned int current_vector_component = 0;
2644 for (
unsigned int base = 0; base < base_elements.size(); ++base)
2649 const auto &base_element = this->base_element(base);
2650 const unsigned int multiplicity = this->element_multiplicity(base);
2651 const unsigned int n_base_dofs = base_element.n_dofs_per_cell();
2652 const unsigned int n_base_components = base_element.n_components();
2657 if (n_base_dofs == 0)
2659 current_vector_component += multiplicity * n_base_components;
2663 if (base_element.has_generalized_support_points())
2665 const unsigned int n_base_points =
2666 base_element.get_generalized_support_points().size();
2668 base_dof_values.resize(n_base_dofs);
2669 base_point_values.resize(n_base_points);
2671 for (
unsigned int m = 0; m < multiplicity;
2672 ++m, current_vector_component += n_base_components)
2676 for (
unsigned int j = 0; j < base_point_values.size(); ++j)
2678 base_point_values[j].reinit(n_base_components,
false);
2681 generalized_support_points_index_table[base][j];
2685 const auto *
const begin =
2686 std::begin(point_values[n]) + current_vector_component;
2687 const auto *
const end =
begin + n_base_components;
2688 std::copy(begin, end, std::begin(base_point_values[j]));
2692 .convert_generalized_support_point_values_to_dof_values(
2693 base_point_values, base_dof_values);
2702 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
2703 if (this->system_to_base_index(i).first ==
2704 std::make_pair(base, m))
2706 base_dof_values[this->system_to_base_index(i).second];
2718 for (
unsigned int m = 0; m < multiplicity; ++m)
2719 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
2720 if (this->system_to_base_index(i).first ==
2721 std::make_pair(base, m))
2722 dof_values[i] = std::numeric_limits<double>::signaling_NaN();
2724 current_vector_component += multiplicity * n_base_components;
2731template <
int dim,
int spacedim>
2739 sizeof(base_elements));
2740 for (
unsigned int i = 0; i < base_elements.size(); ++i)
2748#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