59 template <
int dim,
int spacedim,
typename VectorType>
64 , n_shape_functions(fe.n_dofs_per_cell())
66 , local_dof_indices(fe.n_dofs_per_cell())
67 , local_dof_values(fe.n_dofs_per_cell())
72 template <
int dim,
int spacedim,
typename VectorType>
83 template <
int dim,
int spacedim,
typename VectorType>
86 const unsigned int qpoint,
87 const unsigned int shape_nr)
89 AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
90 return shape_values[qpoint * n_shape_functions + shape_nr];
94 template <
int dim,
int spacedim,
typename VectorType>
97 const unsigned int qpoint,
98 const unsigned int shape_nr)
const
101 shape_derivatives.size());
102 return shape_derivatives[qpoint * n_shape_functions + shape_nr];
107 template <
int dim,
int spacedim,
typename VectorType>
110 const unsigned int qpoint,
111 const unsigned int shape_nr)
114 shape_derivatives.size());
115 return shape_derivatives[qpoint * n_shape_functions + shape_nr];
119 template <
int dim,
int spacedim,
typename VectorType>
122 const unsigned int qpoint,
123 const unsigned int shape_nr)
const
126 shape_second_derivatives.size());
127 return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
132 template <
int dim,
int spacedim,
typename VectorType>
135 const unsigned int qpoint,
136 const unsigned int shape_nr)
139 shape_second_derivatives.size());
140 return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
144 template <
int dim,
int spacedim,
typename VectorType>
147 const unsigned int qpoint,
148 const unsigned int shape_nr)
const
151 shape_third_derivatives.size());
152 return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
157 template <
int dim,
int spacedim,
typename VectorType>
160 const unsigned int qpoint,
161 const unsigned int shape_nr)
164 shape_third_derivatives.size());
165 return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
169 template <
int dim,
int spacedim,
typename VectorType>
172 const unsigned int qpoint,
173 const unsigned int shape_nr)
const
176 shape_fourth_derivatives.size());
177 return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
182 template <
int dim,
int spacedim,
typename VectorType>
185 const unsigned int qpoint,
186 const unsigned int shape_nr)
189 shape_fourth_derivatives.size());
190 return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
195 template <
int dim,
int spacedim,
typename VectorType>
203 , euler_dof_handler(&euler_dof_handler)
204 , fe_mask(mask.size() != 0u ?
207 euler_dof_handler.get_fe().get_nonzero_components(0).size(),
210 , fe_values(this->euler_dof_handler->get_fe(),
214 unsigned int size = 0;
215 for (
unsigned int i = 0; i < fe_mask.size(); ++i)
218 fe_to_real[i] = size++;
225 template <
int dim,
int spacedim,
typename VectorType>
228 const std::vector<VectorType> &euler_vector,
231 , uses_level_dofs(true)
232 , euler_dof_handler(&euler_dof_handler)
233 , fe_mask(mask.size() != 0u ?
236 euler_dof_handler.get_fe().get_nonzero_components(0).size(),
239 , fe_values(this->euler_dof_handler->get_fe(),
243 unsigned int size = 0;
252 ExcMessage(
"The underlying DoFHandler object did not call "
253 "distribute_mg_dofs(). In this case, the construction via "
254 "level vectors does not make sense."));
257 this->euler_vector.clear();
265 template <
int dim,
int spacedim,
typename VectorType>
271 , uses_level_dofs(true)
272 , euler_dof_handler(&euler_dof_handler)
273 , fe_mask(mask.size() != 0u ?
276 euler_dof_handler.get_fe().get_nonzero_components(0).size(),
279 , fe_values(this->euler_dof_handler->get_fe(),
283 unsigned int size = 0;
292 ExcMessage(
"The underlying DoFHandler object did not call "
293 "distribute_mg_dofs(). In this case, the construction via "
294 "level vectors does not make sense."));
297 this->euler_vector.
clear();
298 this->euler_vector.
resize(
307 template <
int dim,
int spacedim,
typename VectorType>
311 , uses_level_dofs(mapping.uses_level_dofs)
312 , euler_vector(mapping.euler_vector)
313 , euler_dof_handler(mapping.euler_dof_handler)
314 , fe_mask(mapping.fe_mask)
315 , fe_to_real(mapping.fe_to_real)
316 , fe_values(mapping.euler_dof_handler->get_fe(),
323 template <
int dim,
int spacedim,
typename VectorType>
324 inline const double &
326 const unsigned int qpoint,
327 const unsigned int shape_nr)
const
335 template <
int dim,
int spacedim,
typename VectorType>
344 template <
int dim,
int spacedim,
typename VectorType>
350 ExcMessage(
"The dimension of your mapping (" +
352 ") and the reference cell cell_type (" +
354 " ) do not agree."));
361 template <
int dim,
int spacedim,
typename VectorType>
362 boost::container::small_vector<Point<spacedim>,
389 const unsigned int dofs_per_cell =
391 std::vector<types::global_dof_index> dof_indices(dofs_per_cell);
393 dof_cell->get_mg_dof_indices(dof_indices);
395 dof_cell->get_dof_indices(dof_indices);
397 const VectorType &vector =
400 boost::container::small_vector<Point<spacedim>,
403 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
409 typename VectorType::value_type value =
412 for (
const unsigned int v : cell->vertex_indices())
424 template <
int dim,
int spacedim,
typename VectorType>
431 const unsigned int n_points = unit_points.size();
446 fe->shape_grad_grad(i, unit_points[
point]);
451 fe->shape_3rd_derivative(i, unit_points[
point]);
456 fe->shape_4th_derivative(i, unit_points[
point]);
461 template <
int dim,
int spacedim,
typename VectorType>
473 for (
unsigned int i = 0; i < 5; ++i)
516 template <
int dim,
int spacedim,
typename VectorType>
521 const unsigned int n_original_q_points,
529 const unsigned int n_q_points = q.
size();
536 if (data.update_each &
543 data.
covariant.resize(n_original_q_points);
551 if (data.update_each &
571 template <
int dim,
int spacedim,
typename VectorType>
576 const unsigned int n_original_q_points,
579 compute_data(update_flags, q, n_original_q_points, data);
593 for (
unsigned int i = 0; i < n_faces; ++i)
603 n_original_q_points);
615 template <
int dim,
int spacedim,
typename VectorType>
616 typename std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
621 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
631 template <
int dim,
int spacedim,
typename VectorType>
632 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
639 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
650 template <
int dim,
int spacedim,
typename VectorType>
651 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
656 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
670 namespace MappingFEFieldImplementation
680 template <
int dim,
int spacedim,
typename VectorType>
682 maybe_compute_q_points(
683 const typename ::QProjector<dim>::DataSetDescriptor data_set,
684 const typename ::MappingFEField<dim, spacedim, VectorType>::
699 const double *shape = &data.shape(
point + data_set, 0);
701 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
703 const unsigned int comp_k =
707 data.local_dof_values[k] * shape[k];
722 template <
int dim,
int spacedim,
typename VectorType>
724 maybe_update_Jacobians(
725 const typename ::QProjector<dim>::DataSetDescriptor data_set,
726 const typename ::MappingFEField<dim, spacedim, VectorType>::
737 const unsigned int n_q_points = data.contravariant.size();
744 &data.derivative(
point + data_set, 0);
748 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
750 const unsigned int comp_k =
754 data.local_dof_values[k] * data_derv[k];
758 for (
unsigned int i = 0; i < spacedim; ++i)
760 data.contravariant[
point][i] = result[i];
768 for (
unsigned int point = 0;
point < data.contravariant.size();
770 data.covariant[
point] =
771 (data.contravariant[
point]).covariant_form();
777 data.volume_elements.size());
778 for (
unsigned int point = 0;
point < data.contravariant.size();
780 data.volume_elements[
point] =
781 data.contravariant[
point].determinant();
791 template <
int dim,
int spacedim,
typename VectorType>
794 const typename ::QProjector<dim>::DataSetDescriptor data_set,
795 const typename ::MappingFEField<dim, spacedim, VectorType>::
805 const unsigned int n_q_points = jacobian_grads.size();
810 &data.second_derivative(
point + data_set, 0);
814 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
816 const unsigned int comp_k =
819 for (
unsigned int j = 0; j < dim; ++j)
820 for (
unsigned int l = 0;
l < dim; ++
l)
822 (
second[k][j][
l] * data.local_dof_values[k]);
827 for (
unsigned int i = 0; i < spacedim; ++i)
828 for (
unsigned int j = 0; j < dim; ++j)
829 for (
unsigned int l = 0;
l < dim; ++
l)
830 jacobian_grads[
point][i][j][
l] = result[i][j][
l];
841 template <
int dim,
int spacedim,
typename VectorType>
844 const typename ::QProjector<dim>::DataSetDescriptor data_set,
845 const typename ::MappingFEField<dim, spacedim, VectorType>::
855 const unsigned int n_q_points =
856 jacobian_pushed_forward_grads.size();
858 double tmp[spacedim][spacedim][spacedim];
862 &data.second_derivative(
point + data_set, 0);
866 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
868 const unsigned int comp_k =
871 for (
unsigned int j = 0; j < dim; ++j)
872 for (
unsigned int l = 0;
l < dim; ++
l)
874 (
second[k][j][
l] * data.local_dof_values[k]);
878 for (
unsigned int i = 0; i < spacedim; ++i)
879 for (
unsigned int j = 0; j < spacedim; ++j)
880 for (
unsigned int l = 0;
l < dim; ++
l)
883 result[i][0][
l] * data.covariant[
point][j][0];
884 for (
unsigned int jr = 1; jr < dim; ++jr)
887 result[i][jr][
l] * data.covariant[
point][j][jr];
892 for (
unsigned int i = 0; i < spacedim; ++i)
893 for (
unsigned int j = 0; j < spacedim; ++j)
894 for (
unsigned int l = 0;
l < spacedim; ++
l)
896 jacobian_pushed_forward_grads[
point][i][j][
l] =
897 tmp[i][j][0] * data.covariant[
point][
l][0];
898 for (
unsigned int lr = 1; lr < dim; ++lr)
900 jacobian_pushed_forward_grads[
point][i][j][
l] +=
901 tmp[i][j][lr] * data.covariant[
point][
l][lr];
914 template <
int dim,
int spacedim,
typename VectorType>
917 const typename ::QProjector<dim>::DataSetDescriptor data_set,
918 const typename ::MappingFEField<dim, spacedim, VectorType>::
928 const unsigned int n_q_points = jacobian_2nd_derivatives.size();
933 &data.third_derivative(
point + data_set, 0);
937 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
939 const unsigned int comp_k =
942 for (
unsigned int j = 0; j < dim; ++j)
943 for (
unsigned int l = 0;
l < dim; ++
l)
944 for (
unsigned int m = 0; m < dim; ++m)
946 (third[k][j][
l][m] * data.local_dof_values[k]);
951 for (
unsigned int i = 0; i < spacedim; ++i)
952 for (
unsigned int j = 0; j < dim; ++j)
953 for (
unsigned int l = 0;
l < dim; ++
l)
954 for (
unsigned int m = 0; m < dim; ++m)
955 jacobian_2nd_derivatives[
point][i][j][
l][m] =
968 template <
int dim,
int spacedim,
typename VectorType>
971 const typename ::QProjector<dim>::DataSetDescriptor data_set,
972 const typename ::MappingFEField<dim, spacedim, VectorType>::
978 &jacobian_pushed_forward_2nd_derivatives)
983 const unsigned int n_q_points =
984 jacobian_pushed_forward_2nd_derivatives.size();
986 double tmp[spacedim][spacedim][spacedim][spacedim];
990 &data.third_derivative(
point + data_set, 0);
994 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
996 const unsigned int comp_k =
999 for (
unsigned int j = 0; j < dim; ++j)
1000 for (
unsigned int l = 0;
l < dim; ++
l)
1001 for (
unsigned int m = 0; m < dim; ++m)
1003 (third[k][j][
l][m] * data.local_dof_values[k]);
1007 for (
unsigned int i = 0; i < spacedim; ++i)
1008 for (
unsigned int j = 0; j < spacedim; ++j)
1009 for (
unsigned int l = 0;
l < dim; ++
l)
1010 for (
unsigned int m = 0; m < dim; ++m)
1012 jacobian_pushed_forward_2nd_derivatives
1014 result[i][0][
l][m] * data.covariant[
point][j][0];
1015 for (
unsigned int jr = 1; jr < dim; ++jr)
1016 jacobian_pushed_forward_2nd_derivatives[
point][i][j]
1018 result[i][jr][
l][m] *
1019 data.covariant[
point][j][jr];
1023 for (
unsigned int i = 0; i < spacedim; ++i)
1024 for (
unsigned int j = 0; j < spacedim; ++j)
1025 for (
unsigned int l = 0;
l < spacedim; ++
l)
1026 for (
unsigned int m = 0; m < dim; ++m)
1029 jacobian_pushed_forward_2nd_derivatives[
point][i][j]
1031 data.covariant[
point][
l][0];
1032 for (
unsigned int lr = 1; lr < dim; ++lr)
1034 jacobian_pushed_forward_2nd_derivatives[
point][i]
1037 data.covariant[
point][
l][lr];
1041 for (
unsigned int i = 0; i < spacedim; ++i)
1042 for (
unsigned int j = 0; j < spacedim; ++j)
1043 for (
unsigned int l = 0;
l < spacedim; ++
l)
1044 for (
unsigned int m = 0; m < spacedim; ++m)
1046 jacobian_pushed_forward_2nd_derivatives
1048 tmp[i][j][
l][0] * data.covariant[
point][m][0];
1049 for (
unsigned int mr = 1; mr < dim; ++mr)
1050 jacobian_pushed_forward_2nd_derivatives[
point][i][j]
1052 tmp[i][j][
l][mr] * data.covariant[
point][m][mr];
1064 template <
int dim,
int spacedim,
typename VectorType>
1067 const typename ::QProjector<dim>::DataSetDescriptor data_set,
1068 const typename ::MappingFEField<dim, spacedim, VectorType>::
1075 const UpdateFlags update_flags = data.update_each;
1078 const unsigned int n_q_points = jacobian_3rd_derivatives.size();
1083 &data.fourth_derivative(
point + data_set, 0);
1087 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
1089 const unsigned int comp_k =
1092 for (
unsigned int j = 0; j < dim; ++j)
1093 for (
unsigned int l = 0;
l < dim; ++
l)
1094 for (
unsigned int m = 0; m < dim; ++m)
1095 for (
unsigned int n = 0; n < dim; ++n)
1097 (fourth[k][j][
l][m][n] *
1098 data.local_dof_values[k]);
1104 for (
unsigned int i = 0; i < spacedim; ++i)
1105 for (
unsigned int j = 0; j < dim; ++j)
1106 for (
unsigned int l = 0;
l < dim; ++
l)
1107 for (
unsigned int m = 0; m < dim; ++m)
1108 for (
unsigned int n = 0; n < dim; ++n)
1109 jacobian_3rd_derivatives[
point][i][j][
l][m][n] =
1110 result[i][j][
l][m][n];
1122 template <
int dim,
int spacedim,
typename VectorType>
1125 const typename ::QProjector<dim>::DataSetDescriptor data_set,
1126 const typename ::MappingFEField<dim, spacedim, VectorType>::
1132 &jacobian_pushed_forward_3rd_derivatives)
1134 const UpdateFlags update_flags = data.update_each;
1137 const unsigned int n_q_points =
1138 jacobian_pushed_forward_3rd_derivatives.size();
1140 double tmp[spacedim][spacedim][spacedim][spacedim][spacedim];
1144 &data.fourth_derivative(
point + data_set, 0);
1148 for (
unsigned int k = 0; k < data.n_shape_functions; ++k)
1150 const unsigned int comp_k =
1153 for (
unsigned int j = 0; j < dim; ++j)
1154 for (
unsigned int l = 0;
l < dim; ++
l)
1155 for (
unsigned int m = 0; m < dim; ++m)
1156 for (
unsigned int n = 0; n < dim; ++n)
1158 (fourth[k][j][
l][m][n] *
1159 data.local_dof_values[k]);
1163 for (
unsigned int i = 0; i < spacedim; ++i)
1164 for (
unsigned int j = 0; j < spacedim; ++j)
1165 for (
unsigned int l = 0;
l < dim; ++
l)
1166 for (
unsigned int m = 0; m < dim; ++m)
1167 for (
unsigned int n = 0; n < dim; ++n)
1169 tmp[i][j][
l][m][n] = result[i][0][
l][m][n] *
1170 data.covariant[
point][j][0];
1171 for (
unsigned int jr = 1; jr < dim; ++jr)
1172 tmp[i][j][
l][m][n] +=
1173 result[i][jr][
l][m][n] *
1174 data.covariant[
point][j][jr];
1178 for (
unsigned int i = 0; i < spacedim; ++i)
1179 for (
unsigned int j = 0; j < spacedim; ++j)
1180 for (
unsigned int l = 0;
l < spacedim; ++
l)
1181 for (
unsigned int m = 0; m < dim; ++m)
1182 for (
unsigned int n = 0; n < dim; ++n)
1184 jacobian_pushed_forward_3rd_derivatives
1186 tmp[i][j][0][m][n] *
1187 data.covariant[
point][
l][0];
1188 for (
unsigned int lr = 1; lr < dim; ++lr)
1189 jacobian_pushed_forward_3rd_derivatives[
point][i]
1192 tmp[i][j][lr][m][n] *
1193 data.covariant[
point][
l][lr];
1197 for (
unsigned int i = 0; i < spacedim; ++i)
1198 for (
unsigned int j = 0; j < spacedim; ++j)
1199 for (
unsigned int l = 0;
l < spacedim; ++
l)
1200 for (
unsigned int m = 0; m < spacedim; ++m)
1201 for (
unsigned int n = 0; n < dim; ++n)
1203 tmp[i][j][
l][m][n] =
1204 jacobian_pushed_forward_3rd_derivatives[
point][i]
1207 data.covariant[
point][m][0];
1208 for (
unsigned int mr = 1; mr < dim; ++mr)
1209 tmp[i][j][
l][m][n] +=
1210 jacobian_pushed_forward_3rd_derivatives[
point]
1213 data.covariant[
point][m][mr];
1217 for (
unsigned int i = 0; i < spacedim; ++i)
1218 for (
unsigned int j = 0; j < spacedim; ++j)
1219 for (
unsigned int l = 0;
l < spacedim; ++
l)
1220 for (
unsigned int m = 0; m < spacedim; ++m)
1221 for (
unsigned int n = 0; n < spacedim; ++n)
1223 jacobian_pushed_forward_3rd_derivatives
1225 tmp[i][j][
l][m][0] *
1226 data.covariant[
point][n][0];
1227 for (
unsigned int nr = 1; nr < dim; ++nr)
1228 jacobian_pushed_forward_3rd_derivatives[
point][i]
1231 tmp[i][j][
l][m][nr] *
1232 data.covariant[
point][n][nr];
1248 template <
int dim,
int spacedim,
typename VectorType>
1251 const ::Mapping<dim, spacedim> &mapping,
1252 const typename ::Triangulation<dim, spacedim>::cell_iterator
1254 const unsigned int face_no,
1255 const unsigned int subface_no,
1257 const typename ::MappingFEField<dim, spacedim, VectorType>::
1262 const UpdateFlags update_flags = data.update_each;
1266 const unsigned int n_q_points = output_data.
boundary_forms.size();
1275 for (
unsigned int d = 0;
d != dim - 1; ++
d)
1277 Assert(face_no + cell->n_faces() *
d <
1278 data.unit_tangentials.size(),
1281 data.aux[
d].size() <=
1282 data.unit_tangentials[face_no + cell->n_faces() *
d].size(),
1287 data.unit_tangentials[face_no + cell->n_faces() *
d]),
1295 if (dim == spacedim)
1297 for (
unsigned int i = 0; i < n_q_points; ++i)
1306 (face_no == 0 ? -1 : +1);
1310 cross_product_2d(data.aux[0][i]);
1314 cross_product_3d(data.aux[0][i], data.aux[1][i]);
1336 data.contravariant[
point].transpose()[0];
1338 (face_no == 0 ? -1. : +1.) *
1348 cross_product_3d(DX_t[0], DX_t[1]);
1349 cell_normal /= cell_normal.
norm();
1354 cross_product_3d(data.aux[0][
point], cell_normal);
1367 data.quadrature_weights[i + data_set];
1372 const double area_ratio =
1374 cell->subface_case(face_no), subface_no);
1393 template <
int dim,
int spacedim,
typename VectorType>
1396 const ::Mapping<dim, spacedim> &mapping,
1397 const typename ::Triangulation<dim, spacedim>::cell_iterator
1399 const unsigned int face_no,
1400 const unsigned int subface_no,
1401 const typename ::QProjector<dim>::DataSetDescriptor data_set,
1402 const typename ::MappingFEField<dim, spacedim, VectorType>::
1410 maybe_compute_q_points<dim, spacedim, VectorType>(
1418 maybe_update_Jacobians<dim, spacedim, VectorType>(
1421 const UpdateFlags update_flags = data.update_each;
1422 const unsigned int n_q_points = data.contravariant.size();
1431 data.covariant[
point].transpose();
1433 maybe_update_jacobian_grads<dim, spacedim, VectorType>(
1436 maybe_update_jacobian_pushed_forward_grads<dim, spacedim, VectorType>(
1444 maybe_update_jacobian_2nd_derivatives<dim, spacedim, VectorType>(
1462 maybe_update_jacobian_3rd_derivatives<dim, spacedim, VectorType>(
1480 maybe_compute_face_data<dim, spacedim, VectorType>(
1481 mapping, cell, face_no, subface_no, data_set, data, output_data);
1490 template <
int dim,
int spacedim,
typename VectorType>
1506 const unsigned int n_q_points = quadrature.
size();
1510 internal::MappingFEFieldImplementation::
1511 maybe_compute_q_points<dim, spacedim, VectorType>(
1519 internal::MappingFEFieldImplementation::
1520 maybe_update_Jacobians<dim, spacedim, VectorType>(
1527 const UpdateFlags update_flags = data.update_each;
1528 const std::vector<double> &weights = quadrature.
get_weights();
1545 if (dim == spacedim)
1554 Assert(det > 1
e-12 * Utilities::fixed_power<dim>(
1555 cell->diameter() / std::sqrt(
double(dim))),
1557 cell->center(), det,
point)));
1566 for (
unsigned int i = 0; i < spacedim; ++i)
1567 for (
unsigned int j = 0; j < dim; ++j)
1571 for (
unsigned int i = 0; i < dim; ++i)
1572 for (
unsigned int j = 0; j < dim; ++j)
1573 G[i][j] = DX_t[i] * DX_t[j];
1580 Assert(spacedim - dim == 1,
1581 ExcMessage(
"There is no cell normal in codim 2."));
1585 cross_product_2d(-DX_t[0]);
1593 cross_product_3d(DX_t[0], DX_t[dim - 1]);
1599 if (cell->direction_flag() ==
false)
1624 internal::MappingFEFieldImplementation::
1625 maybe_update_jacobian_grads<dim, spacedim, VectorType>(
1635 internal::MappingFEFieldImplementation::
1636 maybe_update_jacobian_pushed_forward_grads<dim, spacedim, VectorType>(
1645 internal::MappingFEFieldImplementation::
1646 maybe_update_jacobian_2nd_derivatives<dim, spacedim, VectorType>(
1655 internal::MappingFEFieldImplementation::
1656 maybe_update_jacobian_pushed_forward_2nd_derivatives<dim,
1667 internal::MappingFEFieldImplementation::
1668 maybe_update_jacobian_3rd_derivatives<dim, spacedim, VectorType>(
1678 internal::MappingFEFieldImplementation::
1679 maybe_update_jacobian_pushed_forward_3rd_derivatives<dim,
1694 template <
int dim,
int spacedim,
typename VectorType>
1698 const unsigned int face_no,
1714 internal::MappingFEFieldImplementation::
1715 do_fill_fe_face_values<dim, spacedim, VectorType>(
1722 cell->face_orientation(face_no),
1723 cell->face_flip(face_no),
1724 cell->face_rotation(face_no),
1725 quadrature[0].
size()),
1734 template <
int dim,
int spacedim,
typename VectorType>
1738 const unsigned int face_no,
1739 const unsigned int subface_no,
1763 cell->face_orientation(face_no),
1764 cell->face_flip(face_no),
1765 cell->face_rotation(face_no),
1767 cell->subface_case(face_no)),
1777 template <
int dim,
int spacedim,
typename VectorType>
1791 const unsigned int n_q_points = quadrature.
size();
1795 internal::MappingFEFieldImplementation::
1796 maybe_compute_q_points<dim, spacedim, VectorType>(
1804 internal::MappingFEFieldImplementation::
1805 maybe_update_Jacobians<dim, spacedim, VectorType>(
1812 const UpdateFlags update_flags = data.update_each;
1813 const std::vector<double> &weights = quadrature.
get_weights();
1834 Assert(det > 1
e-12 * Utilities::fixed_power<dim>(
1835 cell->diameter() / std::sqrt(
double(dim))),
1837 cell->center(), det,
point)));
1841 for (
unsigned int d = 0;
d < spacedim;
d++)
1849 normal /= normal.
norm();
1872 internal::MappingFEFieldImplementation::
1873 maybe_update_jacobian_grads<dim, spacedim, VectorType>(
1883 internal::MappingFEFieldImplementation::
1884 maybe_update_jacobian_pushed_forward_grads<dim, spacedim, VectorType>(
1893 internal::MappingFEFieldImplementation::
1894 maybe_update_jacobian_2nd_derivatives<dim, spacedim, VectorType>(
1904 internal::MappingFEFieldImplementation::
1905 maybe_update_jacobian_pushed_forward_2nd_derivatives<dim,
1916 internal::MappingFEFieldImplementation::
1917 maybe_update_jacobian_3rd_derivatives<dim, spacedim, VectorType>(
1927 internal::MappingFEFieldImplementation::
1928 maybe_update_jacobian_pushed_forward_3rd_derivatives<dim,
1942 namespace MappingFEFieldImplementation
1946 template <
int dim,
int spacedim,
int rank,
typename VectorType>
1957 MappingFEField<dim, spacedim, VectorType>::InternalData *
>(
1958 &mapping_data) !=
nullptr),
1960 const typename ::MappingFEField<dim, spacedim, VectorType>::
1961 InternalData &data =
static_cast<
1962 const typename ::MappingFEField<dim, spacedim, VectorType>::
1963 InternalData &
>(mapping_data);
1965 switch (mapping_kind)
1972 "update_contravariant_transformation"));
1974 for (
unsigned int i = 0; i < output.size(); ++i)
1986 "update_contravariant_transformation"));
1990 "update_volume_elements"));
1992 for (
unsigned int i = 0; i < output.size(); ++i)
1996 output[i] /= data.volume_elements[i];
2010 "update_contravariant_transformation"));
2012 for (
unsigned int i = 0; i < output.size(); ++i)
2024 template <
int dim,
int spacedim,
int rank,
typename VectorType>
2035 MappingFEField<dim, spacedim, VectorType>::InternalData *
>(
2036 &mapping_data) !=
nullptr),
2038 const typename ::MappingFEField<dim, spacedim, VectorType>::
2039 InternalData &data =
static_cast<
2040 const typename ::MappingFEField<dim, spacedim, VectorType>::
2041 InternalData &
>(mapping_data);
2043 switch (mapping_kind)
2050 "update_contravariant_transformation"));
2052 for (
unsigned int i = 0; i < output.size(); ++i)
2067 template <
int dim,
int spacedim,
typename VectorType>
2077 internal::MappingFEFieldImplementation::
2078 transform_fields<dim, spacedim, 1, VectorType>(input,
2086 template <
int dim,
int spacedim,
typename VectorType>
2096 internal::MappingFEFieldImplementation::
2097 transform_differential_forms<dim, spacedim, 1, VectorType>(input,
2105 template <
int dim,
int spacedim,
typename VectorType>
2123 template <
int dim,
int spacedim,
typename VectorType>
2136 switch (mapping_kind)
2142 "update_covariant_transformation"));
2144 for (
unsigned int q = 0; q < output.size(); ++q)
2145 for (
unsigned int i = 0; i < spacedim; ++i)
2146 for (
unsigned int j = 0; j < spacedim; ++j)
2147 for (
unsigned int k = 0; k < spacedim; ++k)
2149 output[q][i][j][k] = data.
covariant[q][j][0] *
2152 for (
unsigned int J = 0; J < dim; ++J)
2154 const unsigned int K0 = (0 == J) ? 1 : 0;
2155 for (
unsigned int K = K0;
K < dim; ++
K)
2156 output[q][i][j][k] += data.
covariant[q][j][J] *
2171 template <
int dim,
int spacedim,
typename VectorType>
2189 template <
int dim,
int spacedim,
typename VectorType>
2199 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> mdata(
2210 template <
int dim,
int spacedim,
typename VectorType>
2219 unsigned int comp_i =
2231 template <
int dim,
int spacedim,
typename VectorType>
2250 for (
unsigned int d = 0;
d < dim; ++
d)
2251 initial_p_unit[
d] = 0.5;
2254 initial_p_unit = cell->reference_cell().closest_point(initial_p_unit);
2261 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> mdata(
2262 get_data(update_flags, point_quadrature));
2275 template <
int dim,
int spacedim,
typename VectorType>
2283 const unsigned int n_shapes = mdata.
shape_values.size();
2303 const double eps = 1.e-12 * cell->diameter();
2304 const unsigned int newton_iteration_limit = 20;
2305 unsigned int newton_iteration = 0;
2314 const unsigned int comp_k =
2317 for (
unsigned int j = 0; j < dim; ++j)
2321 for (
unsigned int j = 0; j < dim; ++j)
2323 f[j] = DF[j] * p_minus_F;
2324 for (
unsigned int l = 0;
l < dim; ++
l)
2325 df[j][
l] = -DF[j] * DF[
l];
2331 double step_length = 1;
2344 for (
unsigned int i = 0; i < dim; ++i)
2345 p_unit_trial[i] -= step_length * delta[i];
2355 if (f_trial.
norm() < p_minus_F.
norm())
2357 p_real = p_real_trial;
2358 p_unit = p_unit_trial;
2359 p_minus_F = f_trial;
2362 else if (step_length > 0.05)
2369 if (newton_iteration > newton_iteration_limit)
2386 template <
int dim,
int spacedim,
typename VectorType>
2395 template <
int dim,
int spacedim,
typename VectorType>
2403 template <
int dim,
int spacedim,
typename VectorType>
2404 std::unique_ptr<Mapping<dim, spacedim>>
2407 return std::make_unique<MappingFEField<dim, spacedim, VectorType>>(*this);
2411 template <
int dim,
int spacedim,
typename VectorType>
2438 const VectorType &vector =
2448 #define SPLIT_INSTANTIATIONS_COUNT 2
2449 #ifndef SPLIT_INSTANTIATIONS_INDEX
2450 # define SPLIT_INSTANTIATIONS_INDEX 0
2452 #include "mapping_fe_field.inst"
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
unsigned int size() const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
void resize(const unsigned int new_minlevel, const unsigned int new_maxlevel, Args &&...args)
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
InternalData(const FiniteElement< dim, spacedim > &fe, const ComponentMask &mask)
std::vector< double > quadrature_weights
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
std::vector< std::vector< Tensor< 1, spacedim > > > aux
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
std::vector< Tensor< 3, dim > > shape_third_derivatives
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< types::global_dof_index > local_dof_indices
virtual std::size_t memory_consumption() const override
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< Tensor< 1, dim > > shape_derivatives
unsigned int n_shape_functions
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< double > local_dof_values
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< double > shape_values
std::vector< Tensor< 2, dim > > shape_second_derivatives
std::vector< double > volume_elements
virtual void compute_shapes_virtual(const std::vector< Point< dim >> &unit_points, typename MappingFEField< dim, spacedim, VectorType >::InternalData &data) const
virtual boost::container::small_vector< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
std::vector< unsigned int > fe_to_real
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
const ComponentMask fe_mask
void update_internal_dofs(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename MappingFEField< dim, spacedim, VectorType >::InternalData &data) const
void compute_face_data(const UpdateFlags update_flags, const Quadrature< dim > &q, const unsigned int n_original_q_points, InternalData &data) const
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
std::vector< SmartPointer< const VectorType, MappingFEField< dim, spacedim, VectorType > > > euler_vector
virtual bool preserves_vertex_locations() const override
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
virtual void fill_fe_immersed_surface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const NonMatching::ImmersedSurfaceQuadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
friend class MappingFEField
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
virtual CellSimilarity::Similarity fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
void compute_data(const UpdateFlags update_flags, const Quadrature< dim > &q, const unsigned int n_original_q_points, InternalData &data) const
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 override
ReferenceCell reference_cell
SmartPointer< const DoFHandler< dim, spacedim >, MappingFEField< dim, spacedim, VectorType > > euler_dof_handler
unsigned int get_degree() const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Point< spacedim > do_transform_unit_to_real_cell(const InternalData &mdata) const
Threads::Mutex fe_values_mutex
const bool uses_level_dofs
ComponentMask get_component_mask() const
Point< dim > do_transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p, const Point< dim > &initial_p_unit, InternalData &mdata) 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 typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
FEValues< dim, spacedim > fe_values
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const Quadrature< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Abstract base class for mapping classes.
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const =0
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
static DataSetDescriptor cell()
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)
const std::vector< Point< dim > > & get_points() const
const std::vector< double > & get_weights() const
unsigned int size() const
unsigned int n_faces() const
unsigned int get_dimension() const
constexpr numbers::NumberTraits< Number >::real_type norm_square() const
numbers::NumberTraits< Number >::real_type norm() const
Triangulation< dim, spacedim > & get_triangulation()
unsigned int n_vertices() const
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcInactiveCell()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_volume_elements
Determinant of the Jacobian.
@ update_contravariant_transformation
Contravariant transformation.
@ update_jacobian_pushed_forward_grads
@ update_jacobian_3rd_derivatives
@ update_values
Shape function values.
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_quadrature_points
Transformed quadrature points.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
@ mapping_covariant_gradient
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt K
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim >>> &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double >> &properties={})
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
std::string to_string(const number value, const unsigned int digits=numbers::invalid_unsigned_int)
void do_fill_fe_face_values(const ::MappingQ< dim, spacedim > &mapping, const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const typename QProjector< dim >::DataSetDescriptor data_set, const Quadrature< dim - 1 > &quadrature, const typename ::MappingQ< dim, spacedim >::InternalData &data, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber_lexicographic_to_hierarchic, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data)
void maybe_update_jacobian_3rd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, const ArrayView< const Point< dim >> &unit_points, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber_lexicographic_to_hierarchic, std::vector< DerivativeForm< 4, dim, spacedim >> &jacobian_3rd_derivatives)
void transform_differential_forms(const ArrayView< const DerivativeForm< rank, dim, spacedim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank+1, spacedim >> &output)
void maybe_update_jacobian_grads(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, const ArrayView< const Point< dim >> &unit_points, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber_lexicographic_to_hierarchic, std::vector< DerivativeForm< 2, dim, spacedim >> &jacobian_grads)
void maybe_update_jacobian_pushed_forward_grads(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, const ArrayView< const Point< dim >> &unit_points, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber_lexicographic_to_hierarchic, std::vector< Tensor< 3, spacedim >> &jacobian_pushed_forward_grads)
void maybe_update_jacobian_2nd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, const ArrayView< const Point< dim >> &unit_points, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber_lexicographic_to_hierarchic, std::vector< DerivativeForm< 3, dim, spacedim >> &jacobian_2nd_derivatives)
void maybe_update_jacobian_pushed_forward_2nd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, const ArrayView< const Point< dim >> &unit_points, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber_lexicographic_to_hierarchic, std::vector< Tensor< 4, spacedim >> &jacobian_pushed_forward_2nd_derivatives)
void maybe_update_jacobian_pushed_forward_3rd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, const ArrayView< const Point< dim >> &unit_points, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber_lexicographic_to_hierarchic, std::vector< Tensor< 5, spacedim >> &jacobian_pushed_forward_3rd_derivatives)
void maybe_compute_face_data(const ::MappingQ< dim, spacedim > &mapping, const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const unsigned int n_q_points, const std::vector< double > &weights, const typename ::MappingQ< dim, spacedim >::InternalData &data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data)
void transform_fields(const ArrayView< const Tensor< rank, dim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank, spacedim >> &output)
static const unsigned int invalid_unsigned_int
static double subface_ratio(const internal::SubfaceCase< dim > &subface_case, const unsigned int subface_no)
static VectorType::value_type get(const VectorType &V, const types::global_dof_index i)
constexpr DEAL_II_HOST SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &)
constexpr DEAL_II_HOST Number determinant(const SymmetricTensor< 2, dim, Number > &)