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 neigbor_cell_id = neighbor_cell_at_face->id();
88 if (((nn + f) % 2 == 0) && this_cell_id < neigbor_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> & )
139 template <
int spacedim>
141 get_dof_sign_change_nedelec(
142 const typename ::Triangulation<2, spacedim>::cell_iterator &cell,
144 const std::vector<MappingKind> &mapping_kind,
145 std::vector<double> &line_dof_sign)
147 const unsigned int dim = 2;
149 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++l)
150 if (!(cell->line_orientation(l)) &&
152 line_dof_sign[l] = -1.0;
156 template <
int spacedim>
158 get_dof_sign_change_nedelec(
159 const typename ::Triangulation<3, spacedim>::cell_iterator &cell,
161 const std::vector<MappingKind> &mapping_kind,
162 std::vector<double> &line_dof_sign)
164 const unsigned int dim = 3;
167 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++l)
168 if (!(cell->line_orientation(l)) &&
170 line_dof_sign[l] = -1.0;
178template <
int dim,
int spacedim>
182 const std::vector<bool> &restriction_is_additive_flags,
183 const std::vector<ComponentMask> &nonzero_components)
185 restriction_is_additive_flags,
188 , poly_space(polynomials.clone())
190 cached_point[0] = -1;
198 for (
unsigned int comp = 0; comp < this->n_components(); ++comp)
199 this->component_to_base_table[comp].
first.second = comp;
203 adjust_quad_dof_sign_for_face_orientation_table.resize(
204 this->n_unique_2d_subobjects());
206 for (
unsigned int f = 0; f < this->n_unique_2d_subobjects(); ++f)
208 adjust_quad_dof_sign_for_face_orientation_table[f] =
211 adjust_quad_dof_sign_for_face_orientation_table[f].fill(
false);
218template <
int dim,
int spacedim>
221 , mapping_kind(fe.mapping_kind)
222 , adjust_quad_dof_sign_for_face_orientation_table(
223 fe.adjust_quad_dof_sign_for_face_orientation_table)
224 , poly_space(fe.poly_space->clone())
225 , inverse_node_matrix(fe.inverse_node_matrix)
230template <
int dim,
int spacedim>
234 return mapping_kind.size() == 1;
238template <
int dim,
int spacedim>
241 const unsigned int index,
242 const unsigned int face,
243 const bool face_orientation,
244 const bool face_flip,
245 const bool face_rotation)
const
258 Assert(adjust_quad_dof_sign_for_face_orientation_table
259 [this->n_unique_2d_subobjects() == 1 ? 0 : face]
262 this->n_dofs_per_quad(face),
265 return adjust_quad_dof_sign_for_face_orientation_table
266 [this->n_unique_2d_subobjects() == 1 ? 0 : face](
274template <
int dim,
int spacedim>
278 if (single_mapping_kind())
279 return mapping_kind[0];
282 return mapping_kind[i];
287template <
int dim,
int spacedim>
299template <
int dim,
int spacedim>
302 const unsigned int i,
304 const unsigned int component)
const
309 std::lock_guard<std::mutex> lock(cache_mutex);
311 if (cached_point != p || cached_values.empty())
314 cached_values.resize(poly_space->n());
316 std::vector<Tensor<4, dim>> dummy1;
317 std::vector<Tensor<5, dim>> dummy2;
318 poly_space->evaluate(
319 p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2);
323 if (inverse_node_matrix.n_cols() == 0)
324 return cached_values[i][component];
326 for (
unsigned int j = 0; j < inverse_node_matrix.n_cols(); ++j)
327 s += inverse_node_matrix(j, i) * cached_values[j][component];
333template <
int dim,
int spacedim>
344template <
int dim,
int spacedim>
347 const unsigned int i,
349 const unsigned int component)
const
354 std::lock_guard<std::mutex> lock(cache_mutex);
356 if (cached_point != p || cached_grads.empty())
359 cached_grads.resize(poly_space->n());
361 std::vector<Tensor<4, dim>> dummy1;
362 std::vector<Tensor<5, dim>> dummy2;
363 poly_space->evaluate(
364 p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2);
368 if (inverse_node_matrix.n_cols() == 0)
369 return cached_grads[i][component];
371 for (
unsigned int j = 0; j < inverse_node_matrix.n_cols(); ++j)
372 s += inverse_node_matrix(j, i) * cached_grads[j][component];
379template <
int dim,
int spacedim>
390template <
int dim,
int spacedim>
393 const unsigned int i,
395 const unsigned int component)
const
400 std::lock_guard<std::mutex> lock(cache_mutex);
402 if (cached_point != p || cached_grad_grads.empty())
405 cached_grad_grads.resize(poly_space->n());
407 std::vector<Tensor<4, dim>> dummy1;
408 std::vector<Tensor<5, dim>> dummy2;
409 poly_space->evaluate(
410 p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2);
414 if (inverse_node_matrix.n_cols() == 0)
415 return cached_grad_grads[i][component];
417 for (
unsigned int j = 0; j < inverse_node_matrix.n_cols(); ++j)
418 s += inverse_node_matrix(i, j) * cached_grad_grads[j][component];
428template <
int dim,
int spacedim>
447 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
449 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
451 const unsigned int n_q_points = quadrature.
size();
454 fe_data.shape_values.size()[0] == this->n_dofs_per_cell(),
456 this->n_dofs_per_cell()));
458 fe_data.shape_values.size()[1] == n_q_points,
466 std::fill(fe_data.dof_sign_change.begin(),
467 fe_data.dof_sign_change.end(),
470 internal::FE_PolyTensor::get_dof_sign_change_nedelec(
471 cell, *
this, this->mapping_kind, fe_data.dof_sign_change);
478 internal::FE_PolyTensor::get_dof_sign_change_h_div(cell,
481 fe_data.dof_sign_change);
484 const unsigned int first_quad_index = this->get_first_quad_index();
486 const unsigned int n_dofs_per_quad = this->n_dofs_per_quad();
487 const unsigned int n_quad_dofs =
490 for (
unsigned int dof_index = 0; dof_index < this->n_dofs_per_cell();
497 const bool is_quad_dof =
499 (first_quad_index <= dof_index) &&
500 (dof_index < first_quad_index + n_quad_dofs));
510 double dof_sign = 1.0;
513 dof_sign = fe_data.dof_sign_change[dof_index];
521 const unsigned int face_index_from_dof_index =
522 (dof_index - first_quad_index) / (n_dofs_per_quad);
524 const unsigned int local_quad_dof_index = dof_index % n_dofs_per_quad;
527 if (adjust_quad_dof_sign_for_face_orientation(
528 local_quad_dof_index,
529 face_index_from_dof_index,
530 cell->face_orientation(face_index_from_dof_index),
531 cell->face_flip(face_index_from_dof_index),
532 cell->face_rotation(face_index_from_dof_index)))
536 const MappingKind mapping_kind = get_mapping_kind(dof_index);
538 const unsigned int first =
539 output_data.shape_function_to_row_table
540 [dof_index * this->n_components() +
541 this->get_nonzero_components(dof_index).first_selected_component()];
555 switch (mapping_kind)
559 for (
unsigned int k = 0; k < n_q_points; ++k)
560 for (
unsigned int d = 0;
d < dim; ++
d)
561 output_data.shape_values(
first + d, k) =
562 fe_data.shape_values[dof_index][k][
d];
575 for (
unsigned int k = 0; k < n_q_points; ++k)
576 for (
unsigned int d = 0;
d < dim; ++
d)
577 output_data.shape_values(
first + d, k) =
578 fe_data.transformed_shape_values[k][
d];
591 for (
unsigned int k = 0; k < n_q_points; ++k)
592 for (
unsigned int d = 0;
d < dim; ++
d)
593 output_data.shape_values(
first + d, k) =
594 dof_sign * fe_data.transformed_shape_values[k][
d];
606 for (
unsigned int k = 0; k < n_q_points; ++k)
607 for (
unsigned int d = 0;
d < dim; ++
d)
608 output_data.shape_values(
first + d, k) =
609 dof_sign * fe_data.transformed_shape_values[k][
d];
627 switch (mapping_kind)
636 for (
unsigned int k = 0; k < n_q_points; ++k)
637 for (
unsigned int d = 0;
d < dim; ++
d)
638 output_data.shape_gradients[
first + d][k] =
639 fe_data.transformed_shape_grads[k][d];
650 for (
unsigned int k = 0; k < n_q_points; ++k)
651 for (
unsigned int d = 0;
d < spacedim; ++
d)
652 for (
unsigned int n = 0; n < spacedim; ++n)
653 fe_data.transformed_shape_grads[k][d] -=
654 output_data.shape_values(
first + n, k) *
657 for (
unsigned int k = 0; k < n_q_points; ++k)
658 for (
unsigned int d = 0;
d < dim; ++
d)
659 output_data.shape_gradients[
first + d][k] =
660 fe_data.transformed_shape_grads[k][d];
666 for (
unsigned int k = 0; k < n_q_points; ++k)
667 fe_data.untransformed_shape_grads[k] =
668 fe_data.shape_grads[dof_index][k];
675 for (
unsigned int k = 0; k < n_q_points; ++k)
676 for (
unsigned int d = 0;
d < spacedim; ++
d)
677 for (
unsigned int n = 0; n < spacedim; ++n)
678 fe_data.transformed_shape_grads[k][d] +=
679 output_data.shape_values(
first + n, k) *
683 for (
unsigned int k = 0; k < n_q_points; ++k)
684 for (
unsigned int d = 0;
d < dim; ++
d)
685 output_data.shape_gradients[
first + d][k] =
686 fe_data.transformed_shape_grads[k][d];
693 for (
unsigned int k = 0; k < n_q_points; ++k)
694 fe_data.untransformed_shape_grads[k] =
695 fe_data.shape_grads[dof_index][k];
702 for (
unsigned int k = 0; k < n_q_points; ++k)
703 for (
unsigned int d = 0;
d < spacedim; ++
d)
704 for (
unsigned int n = 0; n < spacedim; ++n)
705 fe_data.transformed_shape_grads[k][d] +=
706 (output_data.shape_values(
first + n, k) *
709 (output_data.shape_values(
first + d, k) *
712 for (
unsigned int k = 0; k < n_q_points; ++k)
713 for (
unsigned int d = 0;
d < dim; ++
d)
714 output_data.shape_gradients[
first + d][k] =
715 dof_sign * fe_data.transformed_shape_grads[k][d];
732 for (
unsigned int k = 0; k < n_q_points; ++k)
733 fe_data.untransformed_shape_grads[k] =
734 fe_data.shape_grads[dof_index][k];
742 for (
unsigned int k = 0; k < n_q_points; ++k)
743 for (
unsigned int d = 0;
d < spacedim; ++
d)
744 for (
unsigned int n = 0; n < spacedim; ++n)
745 fe_data.transformed_shape_grads[k][d] -=
746 output_data.shape_values(
first + n, k) *
749 for (
unsigned int k = 0; k < n_q_points; ++k)
750 for (
unsigned int d = 0;
d < dim; ++
d)
751 output_data.shape_gradients[
first + d][k] =
752 dof_sign * fe_data.transformed_shape_grads[k][d];
770 switch (mapping_kind)
780 for (
unsigned int k = 0; k < n_q_points; ++k)
781 for (
unsigned int d = 0;
d < spacedim; ++
d)
782 for (
unsigned int n = 0; n < spacedim; ++n)
783 fe_data.transformed_shape_hessians[k][d] -=
784 output_data.shape_gradients[
first + d][k][n] *
787 for (
unsigned int k = 0; k < n_q_points; ++k)
788 for (
unsigned int d = 0;
d < dim; ++
d)
789 output_data.shape_hessians[
first + d][k] =
790 fe_data.transformed_shape_hessians[k][d];
796 for (
unsigned int k = 0; k < n_q_points; ++k)
797 fe_data.untransformed_shape_hessian_tensors[k] =
798 fe_data.shape_grad_grads[dof_index][k];
802 fe_data.untransformed_shape_hessian_tensors),
807 for (
unsigned int k = 0; k < n_q_points; ++k)
808 for (
unsigned int d = 0;
d < spacedim; ++
d)
809 for (
unsigned int n = 0; n < spacedim; ++n)
810 for (
unsigned int i = 0; i < spacedim; ++i)
811 for (
unsigned int j = 0; j < spacedim; ++j)
813 fe_data.transformed_shape_hessians[k][
d][i][j] -=
814 (output_data.shape_values(
first + n, k) *
818 (output_data.shape_gradients[
first + d][k][n] *
821 (output_data.shape_gradients[
first + n][k][i] *
824 (output_data.shape_gradients[
first + n][k][j] *
829 for (
unsigned int k = 0; k < n_q_points; ++k)
830 for (
unsigned int d = 0;
d < dim; ++
d)
831 output_data.shape_hessians[
first + d][k] =
832 fe_data.transformed_shape_hessians[k][d];
838 for (
unsigned int k = 0; k < n_q_points; ++k)
839 fe_data.untransformed_shape_hessian_tensors[k] =
840 fe_data.shape_grad_grads[dof_index][k];
844 fe_data.untransformed_shape_hessian_tensors),
849 for (
unsigned int k = 0; k < n_q_points; ++k)
850 for (
unsigned int d = 0;
d < spacedim; ++
d)
851 for (
unsigned int n = 0; n < spacedim; ++n)
852 for (
unsigned int i = 0; i < spacedim; ++i)
853 for (
unsigned int j = 0; j < spacedim; ++j)
855 fe_data.transformed_shape_hessians[k][
d][i][j] +=
856 (output_data.shape_values(
first + n, k) *
860 (output_data.shape_gradients[
first + n][k][i] *
863 (output_data.shape_gradients[
first + n][k][j] *
866 (output_data.shape_gradients[
first + d][k][n] *
869 for (
unsigned int m = 0; m < spacedim; ++m)
871 .transformed_shape_hessians[k][d][i][j] -=
873 .jacobian_pushed_forward_grads[k][d][i]
876 .jacobian_pushed_forward_grads[k][m][n]
878 output_data.shape_values(
first + n, k)) +
880 .jacobian_pushed_forward_grads[k][d][m]
883 .jacobian_pushed_forward_grads[k][m][i]
885 output_data.shape_values(
first + n, k));
888 for (
unsigned int k = 0; k < n_q_points; ++k)
889 for (
unsigned int d = 0;
d < dim; ++
d)
890 output_data.shape_hessians[
first + d][k] =
891 fe_data.transformed_shape_hessians[k][d];
898 for (
unsigned int k = 0; k < n_q_points; ++k)
899 fe_data.untransformed_shape_hessian_tensors[k] =
900 fe_data.shape_grad_grads[dof_index][k];
904 fe_data.untransformed_shape_hessian_tensors),
909 for (
unsigned int k = 0; k < n_q_points; ++k)
910 for (
unsigned int d = 0;
d < spacedim; ++
d)
911 for (
unsigned int n = 0; n < spacedim; ++n)
912 for (
unsigned int i = 0; i < spacedim; ++i)
913 for (
unsigned int j = 0; j < spacedim; ++j)
915 fe_data.transformed_shape_hessians[k][
d][i][j] +=
916 (output_data.shape_values(
first + n, k) *
920 (output_data.shape_gradients[
first + n][k][i] *
923 (output_data.shape_gradients[
first + n][k][j] *
926 (output_data.shape_gradients[
first + d][k][n] *
930 fe_data.transformed_shape_hessians[k][
d][i][j] -=
931 (output_data.shape_values(
first + d, k) *
935 (output_data.shape_gradients[
first + d][k][i] *
938 (output_data.shape_gradients[
first +
d][k][j] *
942 for (
unsigned int m = 0; m < spacedim; ++m)
945 .transformed_shape_hessians[k][
d][i][j] -=
952 output_data.shape_values(
first + n, k)) +
954 .jacobian_pushed_forward_grads[k][d][m]
957 .jacobian_pushed_forward_grads[k][m][i]
959 output_data.shape_values(
first + n, k));
962 .transformed_shape_hessians[k][
d][i][j] +=
969 output_data.shape_values(
first + d, k)) +
971 .jacobian_pushed_forward_grads[k][n][m]
974 .jacobian_pushed_forward_grads[k][m][i]
976 output_data.shape_values(
first + d, k));
980 for (
unsigned int k = 0; k < n_q_points; ++k)
981 for (
unsigned int d = 0;
d < dim; ++
d)
982 output_data.shape_hessians[
first + d][k] =
983 dof_sign * fe_data.transformed_shape_hessians[k][d];
990 for (
unsigned int k = 0; k < n_q_points; ++k)
991 fe_data.untransformed_shape_hessian_tensors[k] =
992 fe_data.shape_grad_grads[dof_index][k];
996 fe_data.untransformed_shape_hessian_tensors),
1001 for (
unsigned int k = 0; k < n_q_points; ++k)
1002 for (
unsigned int d = 0;
d < spacedim; ++
d)
1003 for (
unsigned int n = 0; n < spacedim; ++n)
1004 for (
unsigned int i = 0; i < spacedim; ++i)
1005 for (
unsigned int j = 0; j < spacedim; ++j)
1007 fe_data.transformed_shape_hessians[k][
d][i][j] -=
1008 (output_data.shape_values(
first + n, k) *
1012 (output_data.shape_gradients[
first + d][k][n] *
1015 (output_data.shape_gradients[
first + n][k][i] *
1018 (output_data.shape_gradients[
first + n][k][j] *
1023 for (
unsigned int k = 0; k < n_q_points; ++k)
1024 for (
unsigned int d = 0;
d < dim; ++
d)
1025 output_data.shape_hessians[
first + d][k] =
1026 dof_sign * fe_data.transformed_shape_hessians[k][d];
1050template <
int dim,
int spacedim>
1054 const unsigned int face_no,
1071 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
1073 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
1075 const unsigned int n_q_points = quadrature[0].
size();
1083 cell->face_orientation(face_no),
1084 cell->face_flip(face_no),
1085 cell->face_rotation(face_no),
1095 std::fill(fe_data.dof_sign_change.begin(),
1096 fe_data.dof_sign_change.end(),
1099 internal::FE_PolyTensor::get_dof_sign_change_nedelec(
1100 cell, *
this, this->mapping_kind, fe_data.dof_sign_change);
1107 internal::FE_PolyTensor::get_dof_sign_change_h_div(cell,
1110 fe_data.dof_sign_change);
1113 const unsigned int first_quad_index = this->get_first_quad_index();
1115 const unsigned int n_dofs_per_quad = this->n_dofs_per_quad();
1116 const unsigned int n_quad_dofs =
1119 for (
unsigned int dof_index = 0; dof_index < this->n_dofs_per_cell();
1126 const bool is_quad_dof =
1128 (first_quad_index <= dof_index) &&
1129 (dof_index < first_quad_index + n_quad_dofs));
1139 double dof_sign = 1.0;
1142 dof_sign = fe_data.dof_sign_change[dof_index];
1150 unsigned int face_index_from_dof_index =
1151 (dof_index - first_quad_index) / (n_dofs_per_quad);
1153 unsigned int local_quad_dof_index = dof_index % n_dofs_per_quad;
1156 if (adjust_quad_dof_sign_for_face_orientation(
1157 local_quad_dof_index,
1158 face_index_from_dof_index,
1159 cell->face_orientation(face_index_from_dof_index),
1160 cell->face_flip(face_index_from_dof_index),
1161 cell->face_rotation(face_index_from_dof_index)))
1165 const MappingKind mapping_kind = get_mapping_kind(dof_index);
1167 const unsigned int first =
1168 output_data.shape_function_to_row_table
1169 [dof_index * this->n_components() +
1170 this->get_nonzero_components(dof_index).first_selected_component()];
1174 switch (mapping_kind)
1178 for (
unsigned int k = 0; k < n_q_points; ++k)
1179 for (
unsigned int d = 0;
d < dim; ++
d)
1180 output_data.shape_values(
first + d, k) =
1181 fe_data.shape_values[dof_index][k + offset][
d];
1189 transformed_shape_values =
1199 transformed_shape_values);
1201 for (
unsigned int k = 0; k < n_q_points; ++k)
1202 for (
unsigned int d = 0;
d < dim; ++
d)
1203 output_data.shape_values(
first + d, k) =
1204 transformed_shape_values[k][
d];
1212 transformed_shape_values =
1222 transformed_shape_values);
1223 for (
unsigned int k = 0; k < n_q_points; ++k)
1224 for (
unsigned int d = 0;
d < dim; ++
d)
1225 output_data.shape_values(
first + d, k) =
1226 dof_sign * transformed_shape_values[k][
d];
1233 transformed_shape_values =
1243 transformed_shape_values);
1245 for (
unsigned int k = 0; k < n_q_points; ++k)
1246 for (
unsigned int d = 0;
d < dim; ++
d)
1247 output_data.shape_values(
first + d, k) =
1248 dof_sign * transformed_shape_values[k][
d];
1260 switch (mapping_kind)
1274 transformed_shape_grads);
1275 for (
unsigned int k = 0; k < n_q_points; ++k)
1276 for (
unsigned int d = 0;
d < dim; ++
d)
1277 output_data.shape_gradients[
first + d][k] =
1278 transformed_shape_grads[k][d];
1294 transformed_shape_grads);
1296 for (
unsigned int k = 0; k < n_q_points; ++k)
1297 for (
unsigned int d = 0;
d < spacedim; ++
d)
1298 for (
unsigned int n = 0; n < spacedim; ++n)
1299 transformed_shape_grads[k][d] -=
1300 output_data.shape_values(
first + n, k) *
1303 for (
unsigned int k = 0; k < n_q_points; ++k)
1304 for (
unsigned int d = 0;
d < dim; ++
d)
1305 output_data.shape_gradients[
first + d][k] =
1306 transformed_shape_grads[k][d];
1316 for (
unsigned int k = 0; k < n_q_points; ++k)
1317 fe_data.untransformed_shape_grads[k + offset] =
1318 fe_data.shape_grads[dof_index][k + offset];
1325 transformed_shape_grads);
1327 for (
unsigned int k = 0; k < n_q_points; ++k)
1328 for (
unsigned int d = 0;
d < spacedim; ++
d)
1329 for (
unsigned int n = 0; n < spacedim; ++n)
1330 transformed_shape_grads[k][d] +=
1331 output_data.shape_values(
first + n, k) *
1334 for (
unsigned int k = 0; k < n_q_points; ++k)
1335 for (
unsigned int d = 0;
d < dim; ++
d)
1336 output_data.shape_gradients[
first + d][k] =
1337 transformed_shape_grads[k][d];
1349 for (
unsigned int k = 0; k < n_q_points; ++k)
1350 fe_data.untransformed_shape_grads[k + offset] =
1351 fe_data.shape_grads[dof_index][k + offset];
1358 transformed_shape_grads);
1360 for (
unsigned int k = 0; k < n_q_points; ++k)
1361 for (
unsigned int d = 0;
d < spacedim; ++
d)
1362 for (
unsigned int n = 0; n < spacedim; ++n)
1363 transformed_shape_grads[k][d] +=
1364 (output_data.shape_values(
first + n, k) *
1367 (output_data.shape_values(
first + d, k) *
1370 for (
unsigned int k = 0; k < n_q_points; ++k)
1371 for (
unsigned int d = 0;
d < dim; ++
d)
1372 output_data.shape_gradients[
first + d][k] =
1373 dof_sign * transformed_shape_grads[k][d];
1390 for (
unsigned int k = 0; k < n_q_points; ++k)
1391 fe_data.untransformed_shape_grads[k + offset] =
1392 fe_data.shape_grads[dof_index][k + offset];
1404 transformed_shape_grads);
1406 for (
unsigned int k = 0; k < n_q_points; ++k)
1407 for (
unsigned int d = 0;
d < spacedim; ++
d)
1408 for (
unsigned int n = 0; n < spacedim; ++n)
1409 transformed_shape_grads[k][d] -=
1410 output_data.shape_values(
first + n, k) *
1413 for (
unsigned int k = 0; k < n_q_points; ++k)
1414 for (
unsigned int d = 0;
d < dim; ++
d)
1415 output_data.shape_gradients[
first + d][k] =
1416 dof_sign * transformed_shape_grads[k][d];
1428 switch (mapping_kind)
1433 transformed_shape_hessians =
1443 transformed_shape_hessians);
1445 for (
unsigned int k = 0; k < n_q_points; ++k)
1446 for (
unsigned int d = 0;
d < spacedim; ++
d)
1447 for (
unsigned int n = 0; n < spacedim; ++n)
1448 transformed_shape_hessians[k][d] -=
1449 output_data.shape_gradients[
first + d][k][n] *
1452 for (
unsigned int k = 0; k < n_q_points; ++k)
1453 for (
unsigned int d = 0;
d < dim; ++
d)
1454 output_data.shape_hessians[
first + d][k] =
1455 transformed_shape_hessians[k][d];
1461 for (
unsigned int k = 0; k < n_q_points; ++k)
1462 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1463 fe_data.shape_grad_grads[dof_index][k + offset];
1466 transformed_shape_hessians =
1476 transformed_shape_hessians);
1478 for (
unsigned int k = 0; k < n_q_points; ++k)
1479 for (
unsigned int d = 0;
d < spacedim; ++
d)
1480 for (
unsigned int n = 0; n < spacedim; ++n)
1481 for (
unsigned int i = 0; i < spacedim; ++i)
1482 for (
unsigned int j = 0; j < spacedim; ++j)
1484 transformed_shape_hessians[k][
d][i][j] -=
1485 (output_data.shape_values(
first + n, k) *
1489 (output_data.shape_gradients[
first + d][k][n] *
1492 (output_data.shape_gradients[
first + n][k][i] *
1495 (output_data.shape_gradients[
first + n][k][j] *
1500 for (
unsigned int k = 0; k < n_q_points; ++k)
1501 for (
unsigned int d = 0;
d < dim; ++
d)
1502 output_data.shape_hessians[
first + d][k] =
1503 transformed_shape_hessians[k][d];
1510 for (
unsigned int k = 0; k < n_q_points; ++k)
1511 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1512 fe_data.shape_grad_grads[dof_index][k + offset];
1515 transformed_shape_hessians =
1525 transformed_shape_hessians);
1527 for (
unsigned int k = 0; k < n_q_points; ++k)
1528 for (
unsigned int d = 0;
d < spacedim; ++
d)
1529 for (
unsigned int n = 0; n < spacedim; ++n)
1530 for (
unsigned int i = 0; i < spacedim; ++i)
1531 for (
unsigned int j = 0; j < spacedim; ++j)
1533 transformed_shape_hessians[k][
d][i][j] +=
1534 (output_data.shape_values(
first + n, k) *
1538 (output_data.shape_gradients[
first + n][k][i] *
1541 (output_data.shape_gradients[
first + n][k][j] *
1544 (output_data.shape_gradients[
first + d][k][n] *
1547 for (
unsigned int m = 0; m < spacedim; ++m)
1548 transformed_shape_hessians[k][d][i][j] -=
1550 .jacobian_pushed_forward_grads[k][d][i]
1553 .jacobian_pushed_forward_grads[k][m][n]
1555 output_data.shape_values(
first + n, k)) +
1557 .jacobian_pushed_forward_grads[k][d][m]
1560 .jacobian_pushed_forward_grads[k][m][i]
1562 output_data.shape_values(
first + n, k));
1565 for (
unsigned int k = 0; k < n_q_points; ++k)
1566 for (
unsigned int d = 0;
d < dim; ++
d)
1567 output_data.shape_hessians[
first + d][k] =
1568 transformed_shape_hessians[k][d];
1576 for (
unsigned int k = 0; k < n_q_points; ++k)
1577 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1578 fe_data.shape_grad_grads[dof_index][k + offset];
1581 transformed_shape_hessians =
1591 transformed_shape_hessians);
1593 for (
unsigned int k = 0; k < n_q_points; ++k)
1594 for (
unsigned int d = 0;
d < spacedim; ++
d)
1595 for (
unsigned int n = 0; n < spacedim; ++n)
1596 for (
unsigned int i = 0; i < spacedim; ++i)
1597 for (
unsigned int j = 0; j < spacedim; ++j)
1599 transformed_shape_hessians[k][
d][i][j] +=
1600 (output_data.shape_values(
first + n, k) *
1604 (output_data.shape_gradients[
first + n][k][i] *
1607 (output_data.shape_gradients[
first + n][k][j] *
1610 (output_data.shape_gradients[
first + d][k][n] *
1614 transformed_shape_hessians[k][
d][i][j] -=
1615 (output_data.shape_values(
first + d, k) *
1619 (output_data.shape_gradients[
first + d][k][i] *
1622 (output_data.shape_gradients[
first +
d][k][j] *
1626 for (
unsigned int m = 0; m < spacedim; ++m)
1628 transformed_shape_hessians[k][
d][i][j] -=
1635 output_data.shape_values(
first + n, k)) +
1637 .jacobian_pushed_forward_grads[k][d][m]
1640 .jacobian_pushed_forward_grads[k][m][i]
1642 output_data.shape_values(
first + n, k));
1644 transformed_shape_hessians[k][
d][i][j] +=
1651 output_data.shape_values(
first + d, k)) +
1653 .jacobian_pushed_forward_grads[k][n][m]
1656 .jacobian_pushed_forward_grads[k][m][i]
1658 output_data.shape_values(
first + d, k));
1662 for (
unsigned int k = 0; k < n_q_points; ++k)
1663 for (
unsigned int d = 0;
d < dim; ++
d)
1664 output_data.shape_hessians[
first + d][k] =
1665 dof_sign * transformed_shape_hessians[k][d];
1672 for (
unsigned int k = 0; k < n_q_points; ++k)
1673 fe_data.untransformed_shape_hessian_tensors[k + offset] =
1674 fe_data.shape_grad_grads[dof_index][k + offset];
1677 transformed_shape_hessians =
1687 transformed_shape_hessians);
1689 for (
unsigned int k = 0; k < n_q_points; ++k)
1690 for (
unsigned int d = 0;
d < spacedim; ++
d)
1691 for (
unsigned int n = 0; n < spacedim; ++n)
1692 for (
unsigned int i = 0; i < spacedim; ++i)
1693 for (
unsigned int j = 0; j < spacedim; ++j)
1695 transformed_shape_hessians[k][
d][i][j] -=
1696 (output_data.shape_values(
first + n, k) *
1700 (output_data.shape_gradients[
first + d][k][n] *
1703 (output_data.shape_gradients[
first + n][k][i] *
1706 (output_data.shape_gradients[
first + n][k][j] *
1711 for (
unsigned int k = 0; k < n_q_points; ++k)
1712 for (
unsigned int d = 0;
d < dim; ++
d)
1713 output_data.shape_hessians[
first + d][k] =
1714 dof_sign * transformed_shape_hessians[k][d];
1734template <
int dim,
int spacedim>
1738 const unsigned int face_no,
1739 const unsigned int sub_no,
1754 Assert(
dynamic_cast<const InternalData *
>(&fe_internal) !=
nullptr,
1756 const InternalData &fe_data =
static_cast<const InternalData &
>(fe_internal);
1758 const unsigned int n_q_points = quadrature.
size();
1767 cell->face_orientation(face_no),
1768 cell->face_flip(face_no),
1769 cell->face_rotation(face_no),
1771 cell->subface_case(face_no));
1780 std::fill(fe_data.dof_sign_change.begin(),
1781 fe_data.dof_sign_change.end(),
1784 internal::FE_PolyTensor::get_dof_sign_change_nedelec(
1785 cell, *
this, this->mapping_kind, fe_data.dof_sign_change);
1792 internal::FE_PolyTensor::get_dof_sign_change_h_div(cell,
1795 fe_data.dof_sign_change);
1798 const unsigned int first_quad_index = this->get_first_quad_index();
1800 const unsigned int n_dofs_per_quad = this->n_dofs_per_quad();
1801 const unsigned int n_quad_dofs =
1804 for (
unsigned int dof_index = 0; dof_index < this->n_dofs_per_cell();
1811 const bool is_quad_dof =
1813 (first_quad_index <= dof_index) &&
1814 (dof_index < first_quad_index + n_quad_dofs));
1824 double dof_sign = 1.0;
1827 dof_sign = fe_data.dof_sign_change[dof_index];
1835 unsigned int face_index_from_dof_index =
1836 (dof_index - first_quad_index) / (n_dofs_per_quad);
1838 unsigned int local_quad_dof_index = dof_index % n_dofs_per_quad;
1841 if (adjust_quad_dof_sign_for_face_orientation(
1842 local_quad_dof_index,
1843 face_index_from_dof_index,
1844 cell->face_orientation(face_index_from_dof_index),
1845 cell->face_flip(face_index_from_dof_index),
1846 cell->face_rotation(face_index_from_dof_index)))
1850 const MappingKind mapping_kind = get_mapping_kind(dof_index);
1852 const unsigned int first =
1853 output_data.shape_function_to_row_table
1854 [dof_index * this->n_components() +
1855 this->get_nonzero_components(dof_index).first_selected_component()];
1859 switch (mapping_kind)
1863 for (
unsigned int k = 0; k < n_q_points; ++k)
1864 for (
unsigned int d = 0;
d < dim; ++
d)
1865 output_data.shape_values(
first + d, k) =
1866 fe_data.shape_values[dof_index][k + offset][
d];
1874 transformed_shape_values =
1884 transformed_shape_values);
1886 for (
unsigned int k = 0; k < n_q_points; ++k)
1887 for (
unsigned int d = 0;
d < dim; ++
d)
1888 output_data.shape_values(
first + d, k) =
1889 transformed_shape_values[k][
d];
1898 transformed_shape_values =
1909 transformed_shape_values);
1910 for (
unsigned int k = 0; k < n_q_points; ++k)
1911 for (
unsigned int d = 0;
d < dim; ++
d)
1912 output_data.shape_values(
first + d, k) =
1913 dof_sign * transformed_shape_values[k][
d];
1920 transformed_shape_values =
1931 transformed_shape_values);
1933 for (
unsigned int k = 0; k < n_q_points; ++k)
1934 for (
unsigned int d = 0;
d < dim; ++
d)
1935 output_data.shape_values(
first + d, k) =
1936 dof_sign * transformed_shape_values[k][
d];
1952 switch (mapping_kind)
1962 transformed_shape_grads);
1963 for (
unsigned int k = 0; k < n_q_points; ++k)
1964 for (
unsigned int d = 0;
d < dim; ++
d)
1965 output_data.shape_gradients[
first + d][k] =
1966 transformed_shape_grads[k][d];
1978 transformed_shape_grads);
1980 for (
unsigned int k = 0; k < n_q_points; ++k)
1981 for (
unsigned int d = 0;
d < spacedim; ++
d)
1982 for (
unsigned int n = 0; n < spacedim; ++n)
1983 transformed_shape_grads[k][d] -=
1984 output_data.shape_values(
first + n, k) *
1987 for (
unsigned int k = 0; k < n_q_points; ++k)
1988 for (
unsigned int d = 0;
d < dim; ++
d)
1989 output_data.shape_gradients[
first + d][k] =
1990 transformed_shape_grads[k][d];
1997 for (
unsigned int k = 0; k < n_q_points; ++k)
1998 fe_data.untransformed_shape_grads[k + offset] =
1999 fe_data.shape_grads[dof_index][k + offset];
2007 transformed_shape_grads);
2009 for (
unsigned int k = 0; k < n_q_points; ++k)
2010 for (
unsigned int d = 0;
d < spacedim; ++
d)
2011 for (
unsigned int n = 0; n < spacedim; ++n)
2012 transformed_shape_grads[k][d] +=
2013 output_data.shape_values(
first + n, k) *
2016 for (
unsigned int k = 0; k < n_q_points; ++k)
2017 for (
unsigned int d = 0;
d < dim; ++
d)
2018 output_data.shape_gradients[
first + d][k] =
2019 transformed_shape_grads[k][d];
2027 for (
unsigned int k = 0; k < n_q_points; ++k)
2028 fe_data.untransformed_shape_grads[k + offset] =
2029 fe_data.shape_grads[dof_index][k + offset];
2037 transformed_shape_grads);
2039 for (
unsigned int k = 0; k < n_q_points; ++k)
2040 for (
unsigned int d = 0;
d < spacedim; ++
d)
2041 for (
unsigned int n = 0; n < spacedim; ++n)
2042 transformed_shape_grads[k][d] +=
2043 (output_data.shape_values(
first + n, k) *
2046 (output_data.shape_values(
first + d, k) *
2049 for (
unsigned int k = 0; k < n_q_points; ++k)
2050 for (
unsigned int d = 0;
d < dim; ++
d)
2051 output_data.shape_gradients[
first + d][k] =
2052 dof_sign * transformed_shape_grads[k][d];
2068 for (
unsigned int k = 0; k < n_q_points; ++k)
2069 fe_data.untransformed_shape_grads[k + offset] =
2070 fe_data.shape_grads[dof_index][k + offset];
2078 transformed_shape_grads);
2080 for (
unsigned int k = 0; k < n_q_points; ++k)
2081 for (
unsigned int d = 0;
d < spacedim; ++
d)
2082 for (
unsigned int n = 0; n < spacedim; ++n)
2083 transformed_shape_grads[k][d] -=
2084 output_data.shape_values(
first + n, k) *
2087 for (
unsigned int k = 0; k < n_q_points; ++k)
2088 for (
unsigned int d = 0;
d < dim; ++
d)
2089 output_data.shape_gradients[
first + d][k] =
2090 dof_sign * transformed_shape_grads[k][d];
2102 switch (mapping_kind)
2107 transformed_shape_hessians =
2117 transformed_shape_hessians);
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_hessians[k][d] -=
2123 output_data.shape_gradients[
first + d][k][n] *
2126 for (
unsigned int k = 0; k < n_q_points; ++k)
2127 for (
unsigned int d = 0;
d < dim; ++
d)
2128 output_data.shape_hessians[
first + d][k] =
2129 transformed_shape_hessians[k][d];
2135 for (
unsigned int k = 0; k < n_q_points; ++k)
2136 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2137 fe_data.shape_grad_grads[dof_index][k + offset];
2140 transformed_shape_hessians =
2150 transformed_shape_hessians);
2152 for (
unsigned int k = 0; k < n_q_points; ++k)
2153 for (
unsigned int d = 0;
d < spacedim; ++
d)
2154 for (
unsigned int n = 0; n < spacedim; ++n)
2155 for (
unsigned int i = 0; i < spacedim; ++i)
2156 for (
unsigned int j = 0; j < spacedim; ++j)
2158 transformed_shape_hessians[k][
d][i][j] -=
2159 (output_data.shape_values(
first + n, k) *
2163 (output_data.shape_gradients[
first + d][k][n] *
2166 (output_data.shape_gradients[
first + n][k][i] *
2169 (output_data.shape_gradients[
first + n][k][j] *
2174 for (
unsigned int k = 0; k < n_q_points; ++k)
2175 for (
unsigned int d = 0;
d < dim; ++
d)
2176 output_data.shape_hessians[
first + d][k] =
2177 transformed_shape_hessians[k][d];
2184 for (
unsigned int k = 0; k < n_q_points; ++k)
2185 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2186 fe_data.shape_grad_grads[dof_index][k + offset];
2189 transformed_shape_hessians =
2199 transformed_shape_hessians);
2201 for (
unsigned int k = 0; k < n_q_points; ++k)
2202 for (
unsigned int d = 0;
d < spacedim; ++
d)
2203 for (
unsigned int n = 0; n < spacedim; ++n)
2204 for (
unsigned int i = 0; i < spacedim; ++i)
2205 for (
unsigned int j = 0; j < spacedim; ++j)
2207 transformed_shape_hessians[k][
d][i][j] +=
2208 (output_data.shape_values(
first + n, k) *
2212 (output_data.shape_gradients[
first + n][k][i] *
2215 (output_data.shape_gradients[
first + n][k][j] *
2218 (output_data.shape_gradients[
first + d][k][n] *
2221 for (
unsigned int m = 0; m < spacedim; ++m)
2222 transformed_shape_hessians[k][d][i][j] -=
2224 .jacobian_pushed_forward_grads[k][d][i]
2227 .jacobian_pushed_forward_grads[k][m][n]
2229 output_data.shape_values(
first + n, k)) +
2231 .jacobian_pushed_forward_grads[k][d][m]
2234 .jacobian_pushed_forward_grads[k][m][i]
2236 output_data.shape_values(
first + n, k));
2239 for (
unsigned int k = 0; k < n_q_points; ++k)
2240 for (
unsigned int d = 0;
d < dim; ++
d)
2241 output_data.shape_hessians[
first + d][k] =
2242 transformed_shape_hessians[k][d];
2250 for (
unsigned int k = 0; k < n_q_points; ++k)
2251 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2252 fe_data.shape_grad_grads[dof_index][k + offset];
2255 transformed_shape_hessians =
2265 transformed_shape_hessians);
2267 for (
unsigned int k = 0; k < n_q_points; ++k)
2268 for (
unsigned int d = 0;
d < spacedim; ++
d)
2269 for (
unsigned int n = 0; n < spacedim; ++n)
2270 for (
unsigned int i = 0; i < spacedim; ++i)
2271 for (
unsigned int j = 0; j < spacedim; ++j)
2273 transformed_shape_hessians[k][
d][i][j] +=
2274 (output_data.shape_values(
first + n, k) *
2278 (output_data.shape_gradients[
first + n][k][i] *
2281 (output_data.shape_gradients[
first + n][k][j] *
2284 (output_data.shape_gradients[
first + d][k][n] *
2288 transformed_shape_hessians[k][
d][i][j] -=
2289 (output_data.shape_values(
first + d, k) *
2293 (output_data.shape_gradients[
first + d][k][i] *
2296 (output_data.shape_gradients[
first +
d][k][j] *
2299 for (
unsigned int m = 0; m < spacedim; ++m)
2301 transformed_shape_hessians[k][
d][i][j] -=
2308 output_data.shape_values(
first + n, k)) +
2310 .jacobian_pushed_forward_grads[k][d][m]
2313 .jacobian_pushed_forward_grads[k][m][i]
2315 output_data.shape_values(
first + n, k));
2317 transformed_shape_hessians[k][
d][i][j] +=
2324 output_data.shape_values(
first + d, k)) +
2326 .jacobian_pushed_forward_grads[k][n][m]
2329 .jacobian_pushed_forward_grads[k][m][i]
2331 output_data.shape_values(
first + d, k));
2335 for (
unsigned int k = 0; k < n_q_points; ++k)
2336 for (
unsigned int d = 0;
d < dim; ++
d)
2337 output_data.shape_hessians[
first + d][k] =
2338 dof_sign * transformed_shape_hessians[k][d];
2345 for (
unsigned int k = 0; k < n_q_points; ++k)
2346 fe_data.untransformed_shape_hessian_tensors[k + offset] =
2347 fe_data.shape_grad_grads[dof_index][k + offset];
2350 transformed_shape_hessians =
2360 transformed_shape_hessians);
2362 for (
unsigned int k = 0; k < n_q_points; ++k)
2363 for (
unsigned int d = 0;
d < spacedim; ++
d)
2364 for (
unsigned int n = 0; n < spacedim; ++n)
2365 for (
unsigned int i = 0; i < spacedim; ++i)
2366 for (
unsigned int j = 0; j < spacedim; ++j)
2368 transformed_shape_hessians[k][
d][i][j] -=
2369 (output_data.shape_values(
first + n, k) *
2373 (output_data.shape_gradients[
first + d][k][n] *
2376 (output_data.shape_gradients[
first + n][k][i] *
2379 (output_data.shape_gradients[
first + n][k][j] *
2384 for (
unsigned int k = 0; k < n_q_points; ++k)
2385 for (
unsigned int d = 0;
d < dim; ++
d)
2386 output_data.shape_hessians[
first + d][k] =
2387 dof_sign * transformed_shape_hessians[k][d];
2407template <
int dim,
int spacedim>
2414 for (
unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
2416 const MappingKind mapping_kind = get_mapping_kind(i);
2418 switch (mapping_kind)
2509#include "fe_poly_tensor.inst"
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, 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
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
Abstract base class for mapping classes.
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 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
#define DEAL_II_NOT_IMPLEMENTED()
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)
unsigned char combined_face_orientation(const bool face_orientation, const bool face_rotation, const bool face_flip)
const InputIterator OutputIterator out