42 template <
int spacedim>
44 get_dof_sign_change_h_div(
45 const typename ::Triangulation<1, spacedim>::cell_iterator &,
47 const std::vector<MappingKind> &,
48 std::vector<double> &)
61 template <
int spacedim>
63 get_dof_sign_change_h_div(
64 const typename ::Triangulation<2, spacedim>::cell_iterator &cell,
66 const std::vector<MappingKind> &mapping_kind,
67 std::vector<double> &face_sign)
69 const unsigned int dim = 2;
72 const CellId this_cell_id = cell->id();
74 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
76 typename ::Triangulation<dim, spacedim>::face_iterator face =
78 if (!face->at_boundary())
80 const unsigned int nn = cell->neighbor_face_no(f);
81 const typename ::Triangulation<dim,
82 spacedim>::cell_iterator
83 neighbor_cell_at_face = cell->neighbor(f);
84 const CellId neighbor_cell_id = neighbor_cell_at_face->id();
88 if (((nn + f) % 2 == 0) && this_cell_id < neighbor_cell_id)
95 Assert(mapping_kind.size() == 1 ||
96 cell_j < mapping_kind.size(),
101 if ((mapping_kind.size() > 1 ?
102 mapping_kind[cell_j] :
112 template <
int spacedim>
114 get_dof_sign_change_h_div(
115 const typename ::Triangulation<3, spacedim>::cell_iterator
118 const std::vector<MappingKind> & ,
119 std::vector<double> & )
125 template <
int spacedim>
127 get_dof_sign_change_nedelec(
128 const typename ::Triangulation<1, spacedim>::cell_iterator
131 const std::vector<MappingKind> & ,
132 std::vector<double> & )
137 template <
int spacedim>
139 get_dof_sign_change_nedelec(
140 const typename ::Triangulation<2, spacedim>::cell_iterator &cell,
142 const std::vector<MappingKind> &mapping_kind,
143 std::vector<double> &line_dof_sign)
168 for (
unsigned int l = 0; l < GeometryInfo<2>::lines_per_cell; ++l)
169 if (cell->line_orientation(l) !=
177 line_dof_sign[l] = -1.0;
196 for (
unsigned int local_line_dof = 0;
197 local_line_dof < (k + 1);
199 if (local_line_dof % 2 == 0)
200 line_dof_sign[local_line_dof + l * (k + 1)] = -1.0;
205 template <
int spacedim>
207 get_dof_sign_change_nedelec(
208 const typename ::Triangulation<3, spacedim>::cell_iterator &cell,
210 const std::vector<MappingKind> &mapping_kind,
211 std::vector<double> &line_dof_sign)
239 for (
unsigned int l = 0; l < GeometryInfo<3>::lines_per_cell; ++l)
240 if (cell->line_orientation(l) !=
248 line_dof_sign[l] = -1.0;
267 for (
unsigned int local_line_dof = 0;
268 local_line_dof < (k + 1);
270 if (local_line_dof % 2 == 0)
271 line_dof_sign[local_line_dof + l * (k + 1)] = -1.0;
281template <
int dim,
int spacedim>
285 const std::vector<bool> &restriction_is_additive_flags,
286 const std::vector<ComponentMask> &nonzero_components)
288 restriction_is_additive_flags,
291 , poly_space(polynomials.clone())
293 cached_point[0] = -1;
301 for (
unsigned int comp = 0; comp < this->n_components(); ++comp)
302 this->component_to_base_table[comp].
first.second = comp;
306 adjust_quad_dof_sign_for_face_orientation_table.resize(
307 this->n_unique_2d_subobjects());
309 for (
unsigned int f = 0; f < this->n_unique_2d_subobjects(); ++f)
311 adjust_quad_dof_sign_for_face_orientation_table[f] =
314 adjust_quad_dof_sign_for_face_orientation_table[f].fill(
false);
321template <
int dim,
int spacedim>
324 , mapping_kind(fe.mapping_kind)
325 , adjust_quad_dof_sign_for_face_orientation_table(
326 fe.adjust_quad_dof_sign_for_face_orientation_table)
327 , poly_space(fe.poly_space->clone())
328 , inverse_node_matrix(fe.inverse_node_matrix)
333template <
int dim,
int spacedim>
337 return mapping_kind.size() == 1;
341template <
int dim,
int spacedim>
344 const unsigned int index,
345 const unsigned int face,
346 const bool face_orientation,
347 const bool face_flip,
348 const bool face_rotation)
const
361 Assert(adjust_quad_dof_sign_for_face_orientation_table
362 [this->n_unique_2d_subobjects() == 1 ? 0 : face]
365 this->n_dofs_per_quad(face),
368 return adjust_quad_dof_sign_for_face_orientation_table
369 [this->n_unique_2d_subobjects() == 1 ? 0 : face](
377template <
int dim,
int spacedim>
381 if (single_mapping_kind())
382 return mapping_kind[0];
385 return mapping_kind[i];
390template <
int dim,
int spacedim>
402template <
int dim,
int spacedim>
405 const unsigned int i,
407 const unsigned int component)
const
412 std::lock_guard<std::mutex> lock(cache_mutex);
414 if (cached_point != p || cached_values.empty())
417 cached_values.resize(poly_space->n());
419 std::vector<Tensor<4, dim>> dummy1;
420 std::vector<Tensor<5, dim>> dummy2;
421 poly_space->evaluate(
422 p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2);
426 if (inverse_node_matrix.n_cols() == 0)
427 return cached_values[i][component];
429 for (
unsigned int j = 0; j < inverse_node_matrix.n_cols(); ++j)
430 s += inverse_node_matrix(j, i) * cached_values[j][component];
436template <
int dim,
int spacedim>
447template <
int dim,
int spacedim>
450 const unsigned int i,
452 const unsigned int component)
const
457 std::lock_guard<std::mutex> lock(cache_mutex);
459 if (cached_point != p || cached_grads.empty())
462 cached_grads.resize(poly_space->n());
464 std::vector<Tensor<4, dim>> dummy1;
465 std::vector<Tensor<5, dim>> dummy2;
466 poly_space->evaluate(
467 p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2);
471 if (inverse_node_matrix.n_cols() == 0)
472 return cached_grads[i][component];
474 for (
unsigned int j = 0; j < inverse_node_matrix.n_cols(); ++j)
475 s += inverse_node_matrix(j, i) * cached_grads[j][component];
482template <
int dim,
int spacedim>
493template <
int dim,
int spacedim>
496 const unsigned int i,
498 const unsigned int component)
const
503 std::lock_guard<std::mutex> lock(cache_mutex);
505 if (cached_point != p || cached_grad_grads.empty())
508 cached_grad_grads.resize(poly_space->n());
510 std::vector<Tensor<4, dim>> dummy1;
511 std::vector<Tensor<5, dim>> dummy2;
512 poly_space->evaluate(
513 p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2);
517 if (inverse_node_matrix.n_cols() == 0)
518 return cached_grad_grads[i][component];
520 for (
unsigned int j = 0; j < inverse_node_matrix.n_cols(); ++j)
521 s += inverse_node_matrix(i, j) * cached_grad_grads[j][component];
531template <
int dim,
int spacedim>
550 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
552 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
554 const unsigned int n_q_points = quadrature.
size();
557 fe_data.shape_values.size()[0] == this->n_dofs_per_cell(),
559 this->n_dofs_per_cell()));
561 fe_data.shape_values.size()[1] == n_q_points,
569 std::fill(fe_data.dof_sign_change.begin(),
570 fe_data.dof_sign_change.end(),
573 internal::FE_PolyTensor::get_dof_sign_change_nedelec(
574 cell, *
this, this->mapping_kind, fe_data.dof_sign_change);
581 internal::FE_PolyTensor::get_dof_sign_change_h_div(cell,
584 fe_data.dof_sign_change);
587 const unsigned int first_quad_index = this->get_first_quad_index();
589 const unsigned int n_dofs_per_quad = this->n_dofs_per_quad();
590 const unsigned int n_quad_dofs =
593 for (
unsigned int dof_index = 0; dof_index < this->n_dofs_per_cell();
600 const bool is_quad_dof =
602 (first_quad_index <= dof_index) &&
603 (dof_index < first_quad_index + n_quad_dofs));
613 double dof_sign = 1.0;
616 dof_sign = fe_data.dof_sign_change[dof_index];
624 const unsigned int face_index_from_dof_index =
625 (dof_index - first_quad_index) / (n_dofs_per_quad);
627 const unsigned int local_quad_dof_index = dof_index % n_dofs_per_quad;
630 if (adjust_quad_dof_sign_for_face_orientation(
631 local_quad_dof_index,
632 face_index_from_dof_index,
633 cell->face_orientation(face_index_from_dof_index),
634 cell->face_flip(face_index_from_dof_index),
635 cell->face_rotation(face_index_from_dof_index)))
639 const MappingKind mapping_kind = get_mapping_kind(dof_index);
641 const unsigned int first =
642 output_data.shape_function_to_row_table
643 [dof_index * this->n_components() +
644 this->get_nonzero_components(dof_index).first_selected_component()];
658 switch (mapping_kind)
662 for (
unsigned int k = 0; k < n_q_points; ++k)
663 for (
unsigned int d = 0;
d < dim; ++
d)
664 output_data.shape_values(
first + d, k) =
665 fe_data.shape_values[dof_index][k][
d];
678 for (
unsigned int k = 0; k < n_q_points; ++k)
679 for (
unsigned int d = 0;
d < dim; ++
d)
680 output_data.shape_values(
first + d, k) =
681 fe_data.transformed_shape_values[k][
d];
694 for (
unsigned int k = 0; k < n_q_points; ++k)
695 for (
unsigned int d = 0;
d < dim; ++
d)
696 output_data.shape_values(
first + d, k) =
697 dof_sign * fe_data.transformed_shape_values[k][
d];
709 for (
unsigned int k = 0; k < n_q_points; ++k)
710 for (
unsigned int d = 0;
d < dim; ++
d)
711 output_data.shape_values(
first + d, k) =
712 dof_sign * fe_data.transformed_shape_values[k][
d];
730 switch (mapping_kind)
739 for (
unsigned int k = 0; k < n_q_points; ++k)
740 for (
unsigned int d = 0;
d < dim; ++
d)
741 output_data.shape_gradients[
first + d][k] =
742 fe_data.transformed_shape_grads[k][d];
753 for (
unsigned int k = 0; k < n_q_points; ++k)
754 for (
unsigned int d = 0;
d < spacedim; ++
d)
755 for (
unsigned int n = 0; n < spacedim; ++n)
756 fe_data.transformed_shape_grads[k][d] -=
757 output_data.shape_values(
first + n, k) *
760 for (
unsigned int k = 0; k < n_q_points; ++k)
761 for (
unsigned int d = 0;
d < dim; ++
d)
762 output_data.shape_gradients[
first + d][k] =
763 fe_data.transformed_shape_grads[k][d];
769 for (
unsigned int k = 0; k < n_q_points; ++k)
770 fe_data.untransformed_shape_grads[k] =
771 fe_data.shape_grads[dof_index][k];
778 for (
unsigned int k = 0; k < n_q_points; ++k)
779 for (
unsigned int d = 0;
d < spacedim; ++
d)
780 for (
unsigned int n = 0; n < spacedim; ++n)
781 fe_data.transformed_shape_grads[k][d] +=
782 output_data.shape_values(
first + n, k) *
786 for (
unsigned int k = 0; k < n_q_points; ++k)
787 for (
unsigned int d = 0;
d < dim; ++
d)
788 output_data.shape_gradients[
first + d][k] =
789 fe_data.transformed_shape_grads[k][d];
796 for (
unsigned int k = 0; k < n_q_points; ++k)
797 fe_data.untransformed_shape_grads[k] =
798 fe_data.shape_grads[dof_index][k];
805 for (
unsigned int k = 0; k < n_q_points; ++k)
806 for (
unsigned int d = 0;
d < spacedim; ++
d)
807 for (
unsigned int n = 0; n < spacedim; ++n)
808 fe_data.transformed_shape_grads[k][d] +=
809 (output_data.shape_values(
first + n, k) *
812 (output_data.shape_values(
first + d, k) *
815 for (
unsigned int k = 0; k < n_q_points; ++k)
816 for (
unsigned int d = 0;
d < dim; ++
d)
817 output_data.shape_gradients[
first + d][k] =
818 dof_sign * fe_data.transformed_shape_grads[k][d];
835 for (
unsigned int k = 0; k < n_q_points; ++k)
836 fe_data.untransformed_shape_grads[k] =
837 fe_data.shape_grads[dof_index][k];
845 for (
unsigned int k = 0; k < n_q_points; ++k)
846 for (
unsigned int d = 0;
d < spacedim; ++
d)
847 for (
unsigned int n = 0; n < spacedim; ++n)
848 fe_data.transformed_shape_grads[k][d] -=
849 output_data.shape_values(
first + n, k) *
852 for (
unsigned int k = 0; k < n_q_points; ++k)
853 for (
unsigned int d = 0;
d < dim; ++
d)
854 output_data.shape_gradients[
first + d][k] =
855 dof_sign * fe_data.transformed_shape_grads[k][d];
873 switch (mapping_kind)
883 for (
unsigned int k = 0; k < n_q_points; ++k)
884 for (
unsigned int d = 0;
d < spacedim; ++
d)
885 for (
unsigned int n = 0; n < spacedim; ++n)
886 fe_data.transformed_shape_hessians[k][d] -=
887 output_data.shape_gradients[
first + d][k][n] *
890 for (
unsigned int k = 0; k < n_q_points; ++k)
891 for (
unsigned int d = 0;
d < dim; ++
d)
892 output_data.shape_hessians[
first + d][k] =
893 fe_data.transformed_shape_hessians[k][d];
899 for (
unsigned int k = 0; k < n_q_points; ++k)
900 fe_data.untransformed_shape_hessian_tensors[k] =
901 fe_data.shape_grad_grads[dof_index][k];
905 fe_data.untransformed_shape_hessian_tensors),
910 for (
unsigned int k = 0; k < n_q_points; ++k)
911 for (
unsigned int d = 0;
d < spacedim; ++
d)
912 for (
unsigned int n = 0; n < spacedim; ++n)
913 for (
unsigned int i = 0; i < spacedim; ++i)
914 for (
unsigned int j = 0; j < spacedim; ++j)
916 fe_data.transformed_shape_hessians[k][
d][i][j] -=
917 (output_data.shape_values(
first + n, k) *
921 (output_data.shape_gradients[
first + d][k][n] *
924 (output_data.shape_gradients[
first + n][k][i] *
927 (output_data.shape_gradients[
first + n][k][j] *
932 for (
unsigned int k = 0; k < n_q_points; ++k)
933 for (
unsigned int d = 0;
d < dim; ++
d)
934 output_data.shape_hessians[
first + d][k] =
935 fe_data.transformed_shape_hessians[k][d];
941 for (
unsigned int k = 0; k < n_q_points; ++k)
942 fe_data.untransformed_shape_hessian_tensors[k] =
943 fe_data.shape_grad_grads[dof_index][k];
947 fe_data.untransformed_shape_hessian_tensors),
952 for (
unsigned int k = 0; k < n_q_points; ++k)
953 for (
unsigned int d = 0;
d < spacedim; ++
d)
954 for (
unsigned int n = 0; n < spacedim; ++n)
955 for (
unsigned int i = 0; i < spacedim; ++i)
956 for (
unsigned int j = 0; j < spacedim; ++j)
958 fe_data.transformed_shape_hessians[k][
d][i][j] +=
959 (output_data.shape_values(
first + n, k) *
963 (output_data.shape_gradients[
first + n][k][i] *
966 (output_data.shape_gradients[
first + n][k][j] *
969 (output_data.shape_gradients[
first + d][k][n] *
972 for (
unsigned int m = 0; m < spacedim; ++m)
974 .transformed_shape_hessians[k][d][i][j] -=
976 .jacobian_pushed_forward_grads[k][d][i]
979 .jacobian_pushed_forward_grads[k][m][n]
981 output_data.shape_values(
first + n, k)) +
983 .jacobian_pushed_forward_grads[k][d][m]
986 .jacobian_pushed_forward_grads[k][m][i]
988 output_data.shape_values(
first + n, k));
991 for (
unsigned int k = 0; k < n_q_points; ++k)
992 for (
unsigned int d = 0;
d < dim; ++
d)
993 output_data.shape_hessians[
first + d][k] =
994 fe_data.transformed_shape_hessians[k][d];
1001 for (
unsigned int k = 0; k < n_q_points; ++k)
1002 fe_data.untransformed_shape_hessian_tensors[k] =
1003 fe_data.shape_grad_grads[dof_index][k];
1007 fe_data.untransformed_shape_hessian_tensors),
1012 for (
unsigned int k = 0; k < n_q_points; ++k)
1013 for (
unsigned int d = 0;
d < spacedim; ++
d)
1014 for (
unsigned int n = 0; n < spacedim; ++n)
1015 for (
unsigned int i = 0; i < spacedim; ++i)
1016 for (
unsigned int j = 0; j < spacedim; ++j)
1018 fe_data.transformed_shape_hessians[k][
d][i][j] +=
1019 (output_data.shape_values(
first + n, k) *
1023 (output_data.shape_gradients[
first + n][k][i] *
1026 (output_data.shape_gradients[
first + n][k][j] *
1029 (output_data.shape_gradients[
first + d][k][n] *
1033 fe_data.transformed_shape_hessians[k][
d][i][j] -=
1034 (output_data.shape_values(
first + d, k) *
1038 (output_data.shape_gradients[
first + d][k][i] *
1041 (output_data.shape_gradients[
first +
d][k][j] *
1045 for (
unsigned int m = 0; m < spacedim; ++m)
1048 .transformed_shape_hessians[k][
d][i][j] -=
1055 output_data.shape_values(
first + n, k)) +
1057 .jacobian_pushed_forward_grads[k][d][m]
1060 .jacobian_pushed_forward_grads[k][m][i]
1062 output_data.shape_values(
first + n, k));
1065 .transformed_shape_hessians[k][
d][i][j] +=
1072 output_data.shape_values(
first + d, k)) +
1074 .jacobian_pushed_forward_grads[k][n][m]
1077 .jacobian_pushed_forward_grads[k][m][i]
1079 output_data.shape_values(
first + d, k));
1083 for (
unsigned int k = 0; k < n_q_points; ++k)
1084 for (
unsigned int d = 0;
d < dim; ++
d)
1085 output_data.shape_hessians[
first + d][k] =
1086 dof_sign * fe_data.transformed_shape_hessians[k][d];
1093 for (
unsigned int k = 0; k < n_q_points; ++k)
1094 fe_data.untransformed_shape_hessian_tensors[k] =
1095 fe_data.shape_grad_grads[dof_index][k];
1099 fe_data.untransformed_shape_hessian_tensors),
1104 for (
unsigned int k = 0; k < n_q_points; ++k)
1105 for (
unsigned int d = 0;
d < spacedim; ++
d)
1106 for (
unsigned int n = 0; n < spacedim; ++n)
1107 for (
unsigned int i = 0; i < spacedim; ++i)
1108 for (
unsigned int j = 0; j < spacedim; ++j)
1110 fe_data.transformed_shape_hessians[k][
d][i][j] -=
1111 (output_data.shape_values(
first + n, k) *
1115 (output_data.shape_gradients[
first + d][k][n] *
1118 (output_data.shape_gradients[
first + n][k][i] *
1121 (output_data.shape_gradients[
first + n][k][j] *
1126 for (
unsigned int k = 0; k < n_q_points; ++k)
1127 for (
unsigned int d = 0;
d < dim; ++
d)
1128 output_data.shape_hessians[
first + d][k] =
1129 dof_sign * fe_data.transformed_shape_hessians[k][d];
1153template <
int dim,
int spacedim>
1157 const unsigned int face_no,
1174 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
1176 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
1178 const unsigned int n_q_points = quadrature[0].
size();
1206 std::fill(fe_data.dof_sign_change.begin(),
1207 fe_data.dof_sign_change.end(),
1210 internal::FE_PolyTensor::get_dof_sign_change_nedelec(
1211 cell, *
this, this->mapping_kind, fe_data.dof_sign_change);
1218 internal::FE_PolyTensor::get_dof_sign_change_h_div(cell,
1221 fe_data.dof_sign_change);
1224 const unsigned int first_quad_index = this->get_first_quad_index();
1226 const unsigned int n_dofs_per_quad = this->n_dofs_per_quad();
1227 const unsigned int n_quad_dofs =
1230 for (
unsigned int dof_index = 0; dof_index < this->n_dofs_per_cell();
1237 const bool is_quad_dof =
1239 (first_quad_index <= dof_index) &&
1240 (dof_index < first_quad_index + n_quad_dofs));
1250 double dof_sign = 1.0;
1253 dof_sign = fe_data.dof_sign_change[dof_index];
1261 unsigned int face_index_from_dof_index =
1262 (dof_index - first_quad_index) / (n_dofs_per_quad);
1264 unsigned int local_quad_dof_index = dof_index % n_dofs_per_quad;
1267 if (adjust_quad_dof_sign_for_face_orientation(
1268 local_quad_dof_index,
1269 face_index_from_dof_index,
1270 cell->face_orientation(face_index_from_dof_index),
1271 cell->face_flip(face_index_from_dof_index),
1272 cell->face_rotation(face_index_from_dof_index)))
1276 const MappingKind mapping_kind = get_mapping_kind(dof_index);
1278 const unsigned int first =
1279 output_data.shape_function_to_row_table
1280 [dof_index * this->n_components() +
1281 this->get_nonzero_components(dof_index).first_selected_component()];
1285 switch (mapping_kind)
1289 for (
unsigned int k = 0; k < n_q_points; ++k)
1290 for (
unsigned int d = 0;
d < dim; ++
d)
1291 output_data.shape_values(
first + d, k) =
1292 fe_data.shape_values[dof_index][k + offset][
d];
1300 transformed_shape_values =
1310 transformed_shape_values);
1312 for (
unsigned int k = 0; k < n_q_points; ++k)
1313 for (
unsigned int d = 0;
d < dim; ++
d)
1314 output_data.shape_values(
first + d, k) =
1315 transformed_shape_values[k][
d];
1323 transformed_shape_values =
1333 transformed_shape_values);
1334 for (
unsigned int k = 0; k < n_q_points; ++k)
1335 for (
unsigned int d = 0;
d < dim; ++
d)
1336 output_data.shape_values(
first + d, k) =
1337 dof_sign * transformed_shape_values[k][
d];
1344 transformed_shape_values =
1354 transformed_shape_values);
1356 for (
unsigned int k = 0; k < n_q_points; ++k)
1357 for (
unsigned int d = 0;
d < dim; ++
d)
1358 output_data.shape_values(
first + d, k) =
1359 dof_sign * transformed_shape_values[k][
d];
1371 switch (mapping_kind)
1385 transformed_shape_grads);
1386 for (
unsigned int k = 0; k < n_q_points; ++k)
1387 for (
unsigned int d = 0;
d < dim; ++
d)
1388 output_data.shape_gradients[
first + d][k] =
1389 transformed_shape_grads[k][d];
1405 transformed_shape_grads);
1407 for (
unsigned int k = 0; k < n_q_points; ++k)
1408 for (
unsigned int d = 0;
d < spacedim; ++
d)
1409 for (
unsigned int n = 0; n < spacedim; ++n)
1410 transformed_shape_grads[k][d] -=
1411 output_data.shape_values(
first + n, k) *
1414 for (
unsigned int k = 0; k < n_q_points; ++k)
1415 for (
unsigned int d = 0;
d < dim; ++
d)
1416 output_data.shape_gradients[
first + d][k] =
1417 transformed_shape_grads[k][d];
1427 for (
unsigned int k = 0; k < n_q_points; ++k)
1428 fe_data.untransformed_shape_grads[k + offset] =
1429 fe_data.shape_grads[dof_index][k + offset];
1436 transformed_shape_grads);
1438 for (
unsigned int k = 0; k < n_q_points; ++k)
1439 for (
unsigned int d = 0;
d < spacedim; ++
d)
1440 for (
unsigned int n = 0; n < spacedim; ++n)
1441 transformed_shape_grads[k][d] +=
1442 output_data.shape_values(
first + n, k) *
1445 for (
unsigned int k = 0; k < n_q_points; ++k)
1446 for (
unsigned int d = 0;
d < dim; ++
d)
1447 output_data.shape_gradients[
first + d][k] =
1448 transformed_shape_grads[k][d];
1460 for (
unsigned int k = 0; k < n_q_points; ++k)
1461 fe_data.untransformed_shape_grads[k + offset] =
1462 fe_data.shape_grads[dof_index][k + offset];
1469 transformed_shape_grads);
1471 for (
unsigned int k = 0; k < n_q_points; ++k)
1472 for (
unsigned int d = 0;
d < spacedim; ++
d)
1473 for (
unsigned int n = 0; n < spacedim; ++n)
1474 transformed_shape_grads[k][d] +=
1475 (output_data.shape_values(
first + n, k) *
1478 (output_data.shape_values(
first + d, k) *
1481 for (
unsigned int k = 0; k < n_q_points; ++k)
1482 for (
unsigned int d = 0;
d < dim; ++
d)
1483 output_data.shape_gradients[
first + d][k] =
1484 dof_sign * transformed_shape_grads[k][d];
1501 for (
unsigned int k = 0; k < n_q_points; ++k)
1502 fe_data.untransformed_shape_grads[k + offset] =
1503 fe_data.shape_grads[dof_index][k + offset];
1515 transformed_shape_grads);
1517 for (
unsigned int k = 0; k < n_q_points; ++k)
1518 for (
unsigned int d = 0;
d < spacedim; ++
d)
1519 for (
unsigned int n = 0; n < spacedim; ++n)
1520 transformed_shape_grads[k][d] -=
1521 output_data.shape_values(
first + n, k) *
1524 for (
unsigned int k = 0; k < n_q_points; ++k)
1525 for (
unsigned int d = 0;
d < dim; ++
d)
1526 output_data.shape_gradients[
first + d][k] =
1527 dof_sign * transformed_shape_grads[k][d];
1539 switch (mapping_kind)
1544 transformed_shape_hessians =
1554 transformed_shape_hessians);
1556 for (
unsigned int k = 0; k < n_q_points; ++k)
1557 for (
unsigned int d = 0;
d < spacedim; ++
d)
1558 for (
unsigned int n = 0; n < spacedim; ++n)
1559 transformed_shape_hessians[k][d] -=
1560 output_data.shape_gradients[
first + d][k][n] *
1563 for (
unsigned int k = 0; k < n_q_points; ++k)
1564 for (
unsigned int d = 0;
d < dim; ++
d)
1565 output_data.shape_hessians[
first + d][k] =
1566 transformed_shape_hessians[k][d];
1572 for (
unsigned int k = 0; k < n_q_points; ++k)
1573 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1574 fe_data.shape_grad_grads[dof_index][k + offset];
1577 transformed_shape_hessians =
1587 transformed_shape_hessians);
1589 for (
unsigned int k = 0; k < n_q_points; ++k)
1590 for (
unsigned int d = 0;
d < spacedim; ++
d)
1591 for (
unsigned int n = 0; n < spacedim; ++n)
1592 for (
unsigned int i = 0; i < spacedim; ++i)
1593 for (
unsigned int j = 0; j < spacedim; ++j)
1595 transformed_shape_hessians[k][
d][i][j] -=
1596 (output_data.shape_values(
first + n, k) *
1600 (output_data.shape_gradients[
first + d][k][n] *
1603 (output_data.shape_gradients[
first + n][k][i] *
1606 (output_data.shape_gradients[
first + n][k][j] *
1611 for (
unsigned int k = 0; k < n_q_points; ++k)
1612 for (
unsigned int d = 0;
d < dim; ++
d)
1613 output_data.shape_hessians[
first + d][k] =
1614 transformed_shape_hessians[k][d];
1621 for (
unsigned int k = 0; k < n_q_points; ++k)
1622 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1623 fe_data.shape_grad_grads[dof_index][k + offset];
1626 transformed_shape_hessians =
1636 transformed_shape_hessians);
1638 for (
unsigned int k = 0; k < n_q_points; ++k)
1639 for (
unsigned int d = 0;
d < spacedim; ++
d)
1640 for (
unsigned int n = 0; n < spacedim; ++n)
1641 for (
unsigned int i = 0; i < spacedim; ++i)
1642 for (
unsigned int j = 0; j < spacedim; ++j)
1644 transformed_shape_hessians[k][
d][i][j] +=
1645 (output_data.shape_values(
first + n, k) *
1649 (output_data.shape_gradients[
first + n][k][i] *
1652 (output_data.shape_gradients[
first + n][k][j] *
1655 (output_data.shape_gradients[
first + d][k][n] *
1658 for (
unsigned int m = 0; m < spacedim; ++m)
1659 transformed_shape_hessians[k][d][i][j] -=
1661 .jacobian_pushed_forward_grads[k][d][i]
1664 .jacobian_pushed_forward_grads[k][m][n]
1666 output_data.shape_values(
first + n, k)) +
1668 .jacobian_pushed_forward_grads[k][d][m]
1671 .jacobian_pushed_forward_grads[k][m][i]
1673 output_data.shape_values(
first + n, k));
1676 for (
unsigned int k = 0; k < n_q_points; ++k)
1677 for (
unsigned int d = 0;
d < dim; ++
d)
1678 output_data.shape_hessians[
first + d][k] =
1679 transformed_shape_hessians[k][d];
1687 for (
unsigned int k = 0; k < n_q_points; ++k)
1688 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1689 fe_data.shape_grad_grads[dof_index][k + offset];
1692 transformed_shape_hessians =
1702 transformed_shape_hessians);
1704 for (
unsigned int k = 0; k < n_q_points; ++k)
1705 for (
unsigned int d = 0;
d < spacedim; ++
d)
1706 for (
unsigned int n = 0; n < spacedim; ++n)
1707 for (
unsigned int i = 0; i < spacedim; ++i)
1708 for (
unsigned int j = 0; j < spacedim; ++j)
1710 transformed_shape_hessians[k][
d][i][j] +=
1711 (output_data.shape_values(
first + n, k) *
1715 (output_data.shape_gradients[
first + n][k][i] *
1718 (output_data.shape_gradients[
first + n][k][j] *
1721 (output_data.shape_gradients[
first + d][k][n] *
1725 transformed_shape_hessians[k][
d][i][j] -=
1726 (output_data.shape_values(
first + d, k) *
1730 (output_data.shape_gradients[
first + d][k][i] *
1733 (output_data.shape_gradients[
first +
d][k][j] *
1737 for (
unsigned int m = 0; m < spacedim; ++m)
1739 transformed_shape_hessians[k][
d][i][j] -=
1746 output_data.shape_values(
first + n, k)) +
1748 .jacobian_pushed_forward_grads[k][d][m]
1751 .jacobian_pushed_forward_grads[k][m][i]
1753 output_data.shape_values(
first + n, k));
1755 transformed_shape_hessians[k][
d][i][j] +=
1762 output_data.shape_values(
first + d, k)) +
1764 .jacobian_pushed_forward_grads[k][n][m]
1767 .jacobian_pushed_forward_grads[k][m][i]
1769 output_data.shape_values(
first + d, k));
1773 for (
unsigned int k = 0; k < n_q_points; ++k)
1774 for (
unsigned int d = 0;
d < dim; ++
d)
1775 output_data.shape_hessians[
first + d][k] =
1776 dof_sign * transformed_shape_hessians[k][d];
1783 for (
unsigned int k = 0; k < n_q_points; ++k)
1784 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1785 fe_data.shape_grad_grads[dof_index][k + offset];
1788 transformed_shape_hessians =
1798 transformed_shape_hessians);
1800 for (
unsigned int k = 0; k < n_q_points; ++k)
1801 for (
unsigned int d = 0;
d < spacedim; ++
d)
1802 for (
unsigned int n = 0; n < spacedim; ++n)
1803 for (
unsigned int i = 0; i < spacedim; ++i)
1804 for (
unsigned int j = 0; j < spacedim; ++j)
1806 transformed_shape_hessians[k][
d][i][j] -=
1807 (output_data.shape_values(
first + n, k) *
1811 (output_data.shape_gradients[
first + d][k][n] *
1814 (output_data.shape_gradients[
first + n][k][i] *
1817 (output_data.shape_gradients[
first + n][k][j] *
1822 for (
unsigned int k = 0; k < n_q_points; ++k)
1823 for (
unsigned int d = 0;
d < dim; ++
d)
1824 output_data.shape_hessians[
first + d][k] =
1825 dof_sign * transformed_shape_hessians[k][d];
1845template <
int dim,
int spacedim>
1849 const unsigned int face_no,
1850 const unsigned int sub_no,
1865 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
1867 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
1869 const unsigned int n_q_points = quadrature.
size();
1878 cell->combined_face_orientation(
1881 cell->subface_case(face_no));
1890 std::fill(fe_data.dof_sign_change.begin(),
1891 fe_data.dof_sign_change.end(),
1894 internal::FE_PolyTensor::get_dof_sign_change_nedelec(
1895 cell, *
this, this->mapping_kind, fe_data.dof_sign_change);
1902 internal::FE_PolyTensor::get_dof_sign_change_h_div(cell,
1905 fe_data.dof_sign_change);
1908 const unsigned int first_quad_index = this->get_first_quad_index();
1910 const unsigned int n_dofs_per_quad = this->n_dofs_per_quad();
1911 const unsigned int n_quad_dofs =
1914 for (
unsigned int dof_index = 0; dof_index < this->n_dofs_per_cell();
1921 const bool is_quad_dof =
1923 (first_quad_index <= dof_index) &&
1924 (dof_index < first_quad_index + n_quad_dofs));
1934 double dof_sign = 1.0;
1937 dof_sign = fe_data.dof_sign_change[dof_index];
1945 unsigned int face_index_from_dof_index =
1946 (dof_index - first_quad_index) / (n_dofs_per_quad);
1948 unsigned int local_quad_dof_index = dof_index % n_dofs_per_quad;
1951 if (adjust_quad_dof_sign_for_face_orientation(
1952 local_quad_dof_index,
1953 face_index_from_dof_index,
1954 cell->face_orientation(face_index_from_dof_index),
1955 cell->face_flip(face_index_from_dof_index),
1956 cell->face_rotation(face_index_from_dof_index)))
1960 const MappingKind mapping_kind = get_mapping_kind(dof_index);
1962 const unsigned int first =
1963 output_data.shape_function_to_row_table
1964 [dof_index * this->n_components() +
1965 this->get_nonzero_components(dof_index).first_selected_component()];
1969 switch (mapping_kind)
1973 for (
unsigned int k = 0; k < n_q_points; ++k)
1974 for (
unsigned int d = 0;
d < dim; ++
d)
1975 output_data.shape_values(
first + d, k) =
1976 fe_data.shape_values[dof_index][k + offset][
d];
1984 transformed_shape_values =
1994 transformed_shape_values);
1996 for (
unsigned int k = 0; k < n_q_points; ++k)
1997 for (
unsigned int d = 0;
d < dim; ++
d)
1998 output_data.shape_values(
first + d, k) =
1999 transformed_shape_values[k][
d];
2008 transformed_shape_values =
2019 transformed_shape_values);
2020 for (
unsigned int k = 0; k < n_q_points; ++k)
2021 for (
unsigned int d = 0;
d < dim; ++
d)
2022 output_data.shape_values(
first + d, k) =
2023 dof_sign * transformed_shape_values[k][
d];
2030 transformed_shape_values =
2041 transformed_shape_values);
2043 for (
unsigned int k = 0; k < n_q_points; ++k)
2044 for (
unsigned int d = 0;
d < dim; ++
d)
2045 output_data.shape_values(
first + d, k) =
2046 dof_sign * transformed_shape_values[k][
d];
2062 switch (mapping_kind)
2072 transformed_shape_grads);
2073 for (
unsigned int k = 0; k < n_q_points; ++k)
2074 for (
unsigned int d = 0;
d < dim; ++
d)
2075 output_data.shape_gradients[
first + d][k] =
2076 transformed_shape_grads[k][d];
2088 transformed_shape_grads);
2090 for (
unsigned int k = 0; k < n_q_points; ++k)
2091 for (
unsigned int d = 0;
d < spacedim; ++
d)
2092 for (
unsigned int n = 0; n < spacedim; ++n)
2093 transformed_shape_grads[k][d] -=
2094 output_data.shape_values(
first + n, k) *
2097 for (
unsigned int k = 0; k < n_q_points; ++k)
2098 for (
unsigned int d = 0;
d < dim; ++
d)
2099 output_data.shape_gradients[
first + d][k] =
2100 transformed_shape_grads[k][d];
2107 for (
unsigned int k = 0; k < n_q_points; ++k)
2108 fe_data.untransformed_shape_grads[k + offset] =
2109 fe_data.shape_grads[dof_index][k + offset];
2117 transformed_shape_grads);
2119 for (
unsigned int k = 0; k < n_q_points; ++k)
2120 for (
unsigned int d = 0;
d < spacedim; ++
d)
2121 for (
unsigned int n = 0; n < spacedim; ++n)
2122 transformed_shape_grads[k][d] +=
2123 output_data.shape_values(
first + n, k) *
2126 for (
unsigned int k = 0; k < n_q_points; ++k)
2127 for (
unsigned int d = 0;
d < dim; ++
d)
2128 output_data.shape_gradients[
first + d][k] =
2129 transformed_shape_grads[k][d];
2137 for (
unsigned int k = 0; k < n_q_points; ++k)
2138 fe_data.untransformed_shape_grads[k + offset] =
2139 fe_data.shape_grads[dof_index][k + offset];
2147 transformed_shape_grads);
2149 for (
unsigned int k = 0; k < n_q_points; ++k)
2150 for (
unsigned int d = 0;
d < spacedim; ++
d)
2151 for (
unsigned int n = 0; n < spacedim; ++n)
2152 transformed_shape_grads[k][d] +=
2153 (output_data.shape_values(
first + n, k) *
2156 (output_data.shape_values(
first + d, k) *
2159 for (
unsigned int k = 0; k < n_q_points; ++k)
2160 for (
unsigned int d = 0;
d < dim; ++
d)
2161 output_data.shape_gradients[
first + d][k] =
2162 dof_sign * transformed_shape_grads[k][d];
2178 for (
unsigned int k = 0; k < n_q_points; ++k)
2179 fe_data.untransformed_shape_grads[k + offset] =
2180 fe_data.shape_grads[dof_index][k + offset];
2188 transformed_shape_grads);
2190 for (
unsigned int k = 0; k < n_q_points; ++k)
2191 for (
unsigned int d = 0;
d < spacedim; ++
d)
2192 for (
unsigned int n = 0; n < spacedim; ++n)
2193 transformed_shape_grads[k][d] -=
2194 output_data.shape_values(
first + n, k) *
2197 for (
unsigned int k = 0; k < n_q_points; ++k)
2198 for (
unsigned int d = 0;
d < dim; ++
d)
2199 output_data.shape_gradients[
first + d][k] =
2200 dof_sign * transformed_shape_grads[k][d];
2212 switch (mapping_kind)
2217 transformed_shape_hessians =
2227 transformed_shape_hessians);
2229 for (
unsigned int k = 0; k < n_q_points; ++k)
2230 for (
unsigned int d = 0;
d < spacedim; ++
d)
2231 for (
unsigned int n = 0; n < spacedim; ++n)
2232 transformed_shape_hessians[k][d] -=
2233 output_data.shape_gradients[
first + d][k][n] *
2236 for (
unsigned int k = 0; k < n_q_points; ++k)
2237 for (
unsigned int d = 0;
d < dim; ++
d)
2238 output_data.shape_hessians[
first + d][k] =
2239 transformed_shape_hessians[k][d];
2245 for (
unsigned int k = 0; k < n_q_points; ++k)
2246 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2247 fe_data.shape_grad_grads[dof_index][k + offset];
2250 transformed_shape_hessians =
2260 transformed_shape_hessians);
2262 for (
unsigned int k = 0; k < n_q_points; ++k)
2263 for (
unsigned int d = 0;
d < spacedim; ++
d)
2264 for (
unsigned int n = 0; n < spacedim; ++n)
2265 for (
unsigned int i = 0; i < spacedim; ++i)
2266 for (
unsigned int j = 0; j < spacedim; ++j)
2268 transformed_shape_hessians[k][
d][i][j] -=
2269 (output_data.shape_values(
first + n, k) *
2273 (output_data.shape_gradients[
first + d][k][n] *
2276 (output_data.shape_gradients[
first + n][k][i] *
2279 (output_data.shape_gradients[
first + n][k][j] *
2284 for (
unsigned int k = 0; k < n_q_points; ++k)
2285 for (
unsigned int d = 0;
d < dim; ++
d)
2286 output_data.shape_hessians[
first + d][k] =
2287 transformed_shape_hessians[k][d];
2294 for (
unsigned int k = 0; k < n_q_points; ++k)
2295 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2296 fe_data.shape_grad_grads[dof_index][k + offset];
2299 transformed_shape_hessians =
2309 transformed_shape_hessians);
2311 for (
unsigned int k = 0; k < n_q_points; ++k)
2312 for (
unsigned int d = 0;
d < spacedim; ++
d)
2313 for (
unsigned int n = 0; n < spacedim; ++n)
2314 for (
unsigned int i = 0; i < spacedim; ++i)
2315 for (
unsigned int j = 0; j < spacedim; ++j)
2317 transformed_shape_hessians[k][
d][i][j] +=
2318 (output_data.shape_values(
first + n, k) *
2322 (output_data.shape_gradients[
first + n][k][i] *
2325 (output_data.shape_gradients[
first + n][k][j] *
2328 (output_data.shape_gradients[
first + d][k][n] *
2331 for (
unsigned int m = 0; m < spacedim; ++m)
2332 transformed_shape_hessians[k][d][i][j] -=
2334 .jacobian_pushed_forward_grads[k][d][i]
2337 .jacobian_pushed_forward_grads[k][m][n]
2339 output_data.shape_values(
first + n, k)) +
2341 .jacobian_pushed_forward_grads[k][d][m]
2344 .jacobian_pushed_forward_grads[k][m][i]
2346 output_data.shape_values(
first + n, k));
2349 for (
unsigned int k = 0; k < n_q_points; ++k)
2350 for (
unsigned int d = 0;
d < dim; ++
d)
2351 output_data.shape_hessians[
first + d][k] =
2352 transformed_shape_hessians[k][d];
2360 for (
unsigned int k = 0; k < n_q_points; ++k)
2361 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2362 fe_data.shape_grad_grads[dof_index][k + offset];
2365 transformed_shape_hessians =
2375 transformed_shape_hessians);
2377 for (
unsigned int k = 0; k < n_q_points; ++k)
2378 for (
unsigned int d = 0;
d < spacedim; ++
d)
2379 for (
unsigned int n = 0; n < spacedim; ++n)
2380 for (
unsigned int i = 0; i < spacedim; ++i)
2381 for (
unsigned int j = 0; j < spacedim; ++j)
2383 transformed_shape_hessians[k][
d][i][j] +=
2384 (output_data.shape_values(
first + n, k) *
2388 (output_data.shape_gradients[
first + n][k][i] *
2391 (output_data.shape_gradients[
first + n][k][j] *
2394 (output_data.shape_gradients[
first + d][k][n] *
2398 transformed_shape_hessians[k][
d][i][j] -=
2399 (output_data.shape_values(
first + d, k) *
2403 (output_data.shape_gradients[
first + d][k][i] *
2406 (output_data.shape_gradients[
first +
d][k][j] *
2409 for (
unsigned int m = 0; m < spacedim; ++m)
2411 transformed_shape_hessians[k][
d][i][j] -=
2418 output_data.shape_values(
first + n, k)) +
2420 .jacobian_pushed_forward_grads[k][d][m]
2423 .jacobian_pushed_forward_grads[k][m][i]
2425 output_data.shape_values(
first + n, k));
2427 transformed_shape_hessians[k][
d][i][j] +=
2434 output_data.shape_values(
first + d, k)) +
2436 .jacobian_pushed_forward_grads[k][n][m]
2439 .jacobian_pushed_forward_grads[k][m][i]
2441 output_data.shape_values(
first + d, k));
2445 for (
unsigned int k = 0; k < n_q_points; ++k)
2446 for (
unsigned int d = 0;
d < dim; ++
d)
2447 output_data.shape_hessians[
first + d][k] =
2448 dof_sign * transformed_shape_hessians[k][d];
2455 for (
unsigned int k = 0; k < n_q_points; ++k)
2456 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2457 fe_data.shape_grad_grads[dof_index][k + offset];
2460 transformed_shape_hessians =
2470 transformed_shape_hessians);
2472 for (
unsigned int k = 0; k < n_q_points; ++k)
2473 for (
unsigned int d = 0;
d < spacedim; ++
d)
2474 for (
unsigned int n = 0; n < spacedim; ++n)
2475 for (
unsigned int i = 0; i < spacedim; ++i)
2476 for (
unsigned int j = 0; j < spacedim; ++j)
2478 transformed_shape_hessians[k][
d][i][j] -=
2479 (output_data.shape_values(
first + n, k) *
2483 (output_data.shape_gradients[
first + d][k][n] *
2486 (output_data.shape_gradients[
first + n][k][i] *
2489 (output_data.shape_gradients[
first + n][k][j] *
2494 for (
unsigned int k = 0; k < n_q_points; ++k)
2495 for (
unsigned int d = 0;
d < dim; ++
d)
2496 output_data.shape_hessians[
first + d][k] =
2497 dof_sign * transformed_shape_hessians[k][d];
2517template <
int dim,
int spacedim>
2524 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
2526 const MappingKind mapping_kind = get_mapping_kind(i);
2528 switch (mapping_kind)
2619#include "fe/fe_poly_tensor.inst"
ArrayView< std::remove_reference_t< typename std::iterator_traits< Iterator >::reference >, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) 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 double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual double shape_value(const unsigned int i, const Point< dim > &p) 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
bool adjust_quad_dof_sign_for_face_orientation(const unsigned int index, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation) 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
MappingKind get_mapping_kind(const unsigned int i) const
FE_PolyTensor(const TensorPolynomialsBase< dim > &polynomials, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
bool single_mapping_kind() const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int tensor_degree() const
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const types::geometric_orientation combined_orientation=numbers::default_geometric_orientation) const
Abstract base class for mapping classes.
virtual void transform(const ArrayView< const Tensor< 1, dim > > &input, const MappingKind kind, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim > > &output) const =0
static DataSetDescriptor subface(const ReferenceCell &reference_cell, const unsigned int face_no, const unsigned int subface_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points, const internal::SubfaceCase< dim > ref_case=internal::SubfaceCase< dim >::case_isotropic)
static DataSetDescriptor face(const ReferenceCell &reference_cell, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)
unsigned int size() const
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NOT_IMPLEMENTED()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_contravariant_transformation
Contravariant transformation.
@ update_jacobian_pushed_forward_grads
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_3rd_derivatives
Third derivatives of shape functions.
@ update_covariant_transformation
Covariant transformation.
@ update_gradients
Shape function gradients.
@ update_default
No update.
@ update_piola
Values needed for Piola transform.
@ mapping_covariant_gradient
@ mapping_contravariant_hessian
@ mapping_covariant_hessian
@ mapping_contravariant_gradient
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
types::geometric_orientation combined_face_orientation(const bool face_orientation, const bool face_rotation, const bool face_flip)
constexpr types::geometric_orientation default_geometric_orientation