15#ifndef dealii_mapping_q_internal_h
16#define dealii_mapping_q_internal_h
63 template <
int spacedim>
71 return Point<1>((p[0] - vertices[0][0]) /
72 (vertices[1][0] - vertices[0][0]));
77 template <
int spacedim>
89 const long double x = p[0];
90 const long double y = p[1];
92 const long double x0 = vertices[0][0];
93 const long double x1 = vertices[1][0];
94 const long double x2 = vertices[2][0];
95 const long double x3 = vertices[3][0];
97 const long double y0 = vertices[0][1];
98 const long double y1 = vertices[1][1];
99 const long double y2 = vertices[2][1];
100 const long double y3 = vertices[3][1];
141 const long double eta =
156 return {
xi,
static_cast<double>(
eta)};
160 const long double max_y =
169 return {
xi,
static_cast<double>(
eta)};
176 spacedim>::ExcTransformationFailed()));
182 return {std::numeric_limits<double>::quiet_NaN(),
183 std::numeric_limits<double>::quiet_NaN()};
188 template <
int spacedim>
197 return {std::numeric_limits<double>::quiet_NaN(),
198 std::numeric_limits<double>::quiet_NaN(),
199 std::numeric_limits<double>::quiet_NaN()};
210 namespace MappingQImplementation
217 std::vector<Point<dim>>
223 std::vector<Point<dim>> points(
renumbering.size());
224 const unsigned int n1 = line_support_points.size();
225 for (
unsigned int q2 = 0,
q = 0;
q2 < (dim > 2 ?
n1 : 1); ++
q2)
226 for (
unsigned int q1 = 0;
q1 < (dim > 1 ?
n1 : 1); ++
q1)
227 for (
unsigned int q0 = 0;
q0 <
n1; ++
q0, ++
q)
247 inline ::Table<2, double>
254 if (polynomial_degree == 1)
257 const unsigned int M = polynomial_degree - 1;
264 for (
unsigned int i = 0; i <
M; ++i)
265 for (
unsigned int j = 0;
j <
M; ++
j)
268 gl.point((i + 1) * (polynomial_degree + 1) + (
j + 1));
270 for (
unsigned int v = 0; v < 4; ++v)
281 for (
unsigned int unit_point = 0; unit_point <
n_inner_2d; ++unit_point)
282 Assert(std::fabs(std::accumulate(
loqvs[unit_point].begin(),
283 loqvs[unit_point].end(),
285 1) < 1e-13 * polynomial_degree,
299 inline ::Table<2, double>
306 if (polynomial_degree == 1)
309 const unsigned int M = polynomial_degree - 1;
311 const unsigned int n_inner = Utilities::fixed_power<3>(
M);
312 const unsigned int n_outer = 8 + 12 *
M + 6 *
M *
M;
317 for (
unsigned int i = 0; i <
M; ++i)
318 for (
unsigned int j = 0;
j <
M; ++
j)
319 for (
unsigned int k = 0;
k <
M; ++
k)
321 const Point<3> &p =
gl.point((i + 1) * (
M + 2) * (
M + 2) +
322 (
j + 1) * (
M + 2) + (
k + 1));
326 for (
unsigned int v = 0; v < 8; ++v)
335 for (
unsigned int l = 0; l < 4; ++l)
344 for (
unsigned int l = 0; l < 4; ++l)
353 for (
unsigned int l = 0; l < 4; ++l)
372 for (
unsigned int unit_point = 0; unit_point <
n_inner; ++unit_point)
373 Assert(std::fabs(std::accumulate(
lohvs[unit_point].begin(),
374 lohvs[unit_point].end(),
376 1) < 1e-13 * polynomial_degree,
388 inline std::vector<::Table<2, double>>
390 const unsigned int polynomial_degree,
391 const unsigned int dim)
394 std::vector<::Table<2, double>> output(dim);
395 if (polynomial_degree <= 1)
400 output[0].reinit(polynomial_degree - 1,
402 for (
unsigned int q = 0;
q < polynomial_degree - 1; ++
q)
423 inline ::Table<2, double>
427 if (polynomial_degree <= 1)
428 return ::Table<2, double>();
431 const std::vector<unsigned int>
h2l =
432 FETools::hierarchic_to_lexicographic_numbering<dim>(polynomial_degree);
437 for (
unsigned int q = 0;
q < output.size(0); ++
q)
454 template <
int dim,
int spacedim>
457 const typename ::MappingQ<dim, spacedim>::InternalData &
data)
460 data.mapping_support_points.size());
464 for (
unsigned int i = 0; i <
data.mapping_support_points.size(); ++i)
476 template <
int dim,
int spacedim,
typename Number>
483 const std::vector<unsigned int> &renumber,
487 deallog <<
"Start MappingQ::do_transform_real_to_unit_cell for real "
488 <<
"point [ " << p <<
" ] " << std::endl;
505 polynomials_1d, points,
p_unit, polynomials_1d.size() == 2, renumber);
552 const double eps = 1.e-11;
567 <<
" for unit point guess " <<
p_unit << std::endl;
571 for (
unsigned int d = 0; d < spacedim; ++d)
572 for (
unsigned int e = 0; e < dim; ++e)
583 Number(std::numeric_limits<double>::min())) ==
584 Number(std::numeric_limits<double>::min())))
599 polynomials_1d.size() == 2,
617 deallog <<
" delta=" << delta << std::endl;
628 for (
unsigned int i = 0; i < dim; ++i)
637 polynomials_1d.size() == 2,
698 polynomials_1d.size() == 2,
728 <<
" ] and iteration error "
746 const std::vector<unsigned int> &renumber)
748 const unsigned int spacedim = dim + 1;
755 const double eps = 1.e-12;
758 unsigned int loop = 0;
768 polynomials_1d.size() == 2,
774 polynomials_1d, points,
p_unit, renumber);
777 for (
unsigned int j = 0;
j < dim; ++
j)
780 for (
unsigned int l = 0; l < dim; ++l)
821 template <
int dim,
int spacedim>
829 (spacedim == 1 ? 3 : (spacedim == 2 ? 6 : 10));
865 const auto affine = GridTools::affine_cell_approximation<dim>(
868 affine.first.covariant_form().transpose();
875 for (
unsigned int d = 0; d < spacedim; ++d)
876 for (
unsigned int e = 0; e < dim; ++e)
884 std::array<double, n_functions> shape_values;
890 shape_values[0] = 1.;
894 for (
unsigned int d = 0; d < spacedim; ++d)
896 for (
unsigned int d = 0, c = 0; d < spacedim; ++d)
897 for (
unsigned int e = 0; e <= d; ++e, ++c)
907 matrix[i][
j] += shape_values[i] * shape_values[
j];
920 for (
unsigned int j = 0;
j < i; ++
j)
923 for (
unsigned int k = 0;
k <
j; ++
k)
929 ExcMessage(
"Matrix of normal equations not positive "
942 for (
unsigned int j = 0;
j < i; ++
j)
960 for (
unsigned int i = dim + 1; i <
n_functions; ++i)
977 template <
typename Number>
982 for (
unsigned int d = 0; d < dim; ++d)
990 for (
unsigned int d = 0; d < spacedim; ++d)
993 for (
unsigned int d = 0; d < spacedim; ++d)
999 for (
unsigned int d = 0, c = 0; d < spacedim; ++d)
1000 for (
unsigned int e = 0; e <= d; ++e, ++c)
1011 const Number distance_to_unit_cell =
result.distance_square(
1016 for (
unsigned int d = 0; d < dim; ++d)
1018 distance_to_unit_cell,
1061 template <
int dim,
int spacedim>
1065 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1073 using VectorizedArrayType =
1074 typename ::MappingQ<dim,
1075 spacedim>::InternalData::VectorizedArrayType;
1077 const unsigned int n_q_points =
data.shape_info.n_q_points;
1078 constexpr unsigned int n_lanes = VectorizedArrayType::size();
1079 constexpr unsigned int n_comp = 1 + (spacedim - 1) / n_lanes;
1080 constexpr unsigned int n_hessians = (dim * (dim + 1)) / 2;
1087 jacobians.resize(n_q_points);
1089 inverse_jacobians.resize(n_q_points);
1106 n_q_points == quadrature_points.size(),
1109 data.n_shape_functions > 0,
1112 n_q_points == jacobian_grads.size(),
1118 data.shape_info.element_type ==
1121 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1122 quadrature_points[
q] =
1123 data.mapping_support_points[
data.shape_info
1124 .lexicographic_numbering[
q]];
1138 eval.begin_dof_values()[i] = VectorizedArrayType();
1141 data.shape_info.lexicographic_numbering;
1143 for (
unsigned int d = 0; d < spacedim; ++d)
1145 const unsigned int in_comp = d % n_lanes;
1146 const unsigned int out_comp = d / n_lanes;
1161 for (
unsigned int i = 0; i < n_q_points; ++i)
1162 for (
unsigned int in_comp = 0;
1173 for (
unsigned int point = 0; point < n_q_points; ++point)
1174 for (
unsigned int in_comp = 0;
1177 for (
unsigned int j = 0;
j < dim; ++
j)
1180 eval.begin_gradients()[(
out_comp * n_q_points + point) *
1188 for (
unsigned int point = 0; point < n_q_points; ++point)
1189 data.volume_elements[point] = jacobians[point].determinant();
1197 for (
unsigned int point = 0; point < n_q_points; ++point)
1198 inverse_jacobians[point] =
1199 jacobians[point].covariant_form().transpose();
1205 {0, 0}, {1, 1}, {2, 2}, {0, 1}, {0, 2}, {1, 2}};
1210 for (
unsigned int point = 0; point < n_q_points; ++point)
1211 for (
unsigned int j = 0;
j < n_hessians; ++
j)
1212 for (
unsigned int in_comp = 0;
1227 const double value =
1228 eval.begin_hessians()[(
out_comp * n_q_points + point) *
1243 template <
int dim,
int spacedim>
1247 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1250 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
1257 data.mapping_support_points);
1259 const unsigned int n_points = unit_points.size();
1267 jacobians.resize(n_points);
1269 inverse_jacobians.resize(n_points);
1282 for (
unsigned int i = 0; i < n_points; i += n_lanes)
1283 if (n_points - i > 1)
1286 for (
unsigned int j = 0;
j < n_lanes; ++
j)
1287 if (i +
j < n_points)
1288 for (
unsigned int d = 0; d < dim; ++d)
1289 p_vec[d][
j] = unit_points[i +
j][d];
1291 for (
unsigned int d = 0; d < dim; ++d)
1292 p_vec[d][
j] = unit_points[i][d];
1304 polynomials_1d.size() == 2,
1305 renumber_lexicographic_to_hierarchic);
1309 for (
unsigned int d = 0; d < spacedim; ++d)
1310 for (
unsigned int e = 0; e < dim; ++e)
1311 derivative[d][e] =
result.second[e][d];
1318 polynomials_1d.size() == 2,
1319 renumber_lexicographic_to_hierarchic);
1322 for (
unsigned int j = 0;
j < n_lanes && i +
j < n_points; ++
j)
1323 for (
unsigned int d = 0; d < spacedim; ++d)
1324 quadrature_points[i +
j][d] =
value[d][
j];
1330 for (
unsigned int j = 0;
j < n_lanes && i +
j < n_points; ++
j)
1331 for (
unsigned int d = 0; d < spacedim; ++d)
1332 for (
unsigned int e = 0; e < dim; ++e)
1333 jacobians[i +
j][d][e] = derivative[d][e][
j];
1337 const auto determinant = derivative.determinant();
1338 for (
unsigned int j = 0;
j < n_lanes && i +
j < n_points; ++
j)
1344 const auto covariant = derivative.covariant_form();
1345 for (
unsigned int j = 0;
j < n_lanes && i +
j < n_points; ++
j)
1346 for (
unsigned int d = 0; d < dim; ++d)
1347 for (
unsigned int e = 0; e < spacedim; ++e)
1348 inverse_jacobians[i +
j][d][e] = covariant[e][d][
j];
1362 polynomials_1d.size() == 2,
1363 renumber_lexicographic_to_hierarchic);
1367 for (
unsigned int d = 0; d < spacedim; ++d)
1368 for (
unsigned int e = 0; e < dim; ++e)
1369 derivative[d][e] =
result.second[e][d];
1376 polynomials_1d.size() == 2,
1377 renumber_lexicographic_to_hierarchic);
1380 quadrature_points[i] =
value;
1386 data.volume_elements[i] = derivative.determinant();
1389 jacobians[i] = derivative;
1392 inverse_jacobians[i] = jacobians[i].covariant_form().transpose();
1404 template <
int dim,
int spacedim>
1408 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1411 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
1417 data.mapping_support_points);
1418 const unsigned int n_q_points = jacobian_grads.size();
1421 for (
unsigned int point = 0; point < n_q_points; ++point)
1428 renumber_lexicographic_to_hierarchic);
1430 for (
unsigned int i = 0; i < spacedim; ++i)
1431 for (
unsigned int j = 0;
j < dim; ++
j)
1432 for (
unsigned int l = 0; l < dim; ++l)
1433 jacobian_grads[point][i][
j][l] =
second[
j][l][i];
1446 template <
int dim,
int spacedim>
1450 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1453 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
1459 data.mapping_support_points);
1460 const unsigned int n_q_points = jacobian_pushed_forward_grads.size();
1464 double tmp[spacedim][spacedim][spacedim];
1465 for (
unsigned int point = 0; point < n_q_points; ++point)
1472 renumber_lexicographic_to_hierarchic);
1474 data.output_data->inverse_jacobians[point].transpose();
1477 for (
unsigned int i = 0; i < spacedim; ++i)
1478 for (
unsigned int j = 0;
j < spacedim; ++
j)
1479 for (
unsigned int l = 0; l < dim; ++l)
1481 tmp[i][
j][l] =
second[0][l][i] * covariant[
j][0];
1482 for (
unsigned int jr = 1;
jr < dim; ++
jr)
1490 for (
unsigned int i = 0; i < spacedim; ++i)
1491 for (
unsigned int j = 0;
j < spacedim; ++
j)
1492 for (
unsigned int l = 0; l < spacedim; ++l)
1494 jacobian_pushed_forward_grads[point][i][
j][l] =
1495 tmp[i][
j][0] * covariant[l][0];
1496 for (
unsigned int lr = 1;
lr < dim; ++
lr)
1498 jacobian_pushed_forward_grads[point][i][
j][l] +=
1499 tmp[i][
j][
lr] * covariant[l][
lr];
1509 template <
int dim,
int spacedim,
int length_tensor>
1516 for (
unsigned int i = 0; i < spacedim; ++i)
1519 result[i][0][0][0] = compressed[0][i];
1522 for (
unsigned int d = 0; d < 2; ++d)
1523 for (
unsigned int e = 0; e < 2; ++e)
1524 for (
unsigned int f = 0; f < 2; ++f)
1525 result[i][d][e][f] = compressed[d + e + f][i];
1533 for (
unsigned int d = 0; d < 2; ++d)
1534 for (
unsigned int e = 0; e < 2; ++e)
1536 result[i][d][e][2] = compressed[4 + d + e][i];
1537 result[i][d][2][e] = compressed[4 + d + e][i];
1538 result[i][2][d][e] = compressed[4 + d + e][i];
1540 for (
unsigned int d = 0; d < 2; ++d)
1542 result[i][d][2][2] = compressed[7 + d][i];
1543 result[i][2][d][2] = compressed[7 + d][i];
1544 result[i][2][2][d] = compressed[7 + d][i];
1546 result[i][2][2][2] = compressed[9][i];
1561 template <
int dim,
int spacedim>
1565 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1568 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
1574 data.mapping_support_points);
1575 const unsigned int n_q_points = jacobian_2nd_derivatives.size();
1579 for (
unsigned int point = 0; point < n_q_points; ++point)
1582 internal::evaluate_tensor_product_higher_derivatives<3>(
1586 renumber_lexicographic_to_hierarchic));
1601 template <
int dim,
int spacedim>
1605 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1608 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
1614 data.mapping_support_points);
1615 const unsigned int n_q_points =
1616 jacobian_pushed_forward_2nd_derivatives.size();
1622 for (
unsigned int point = 0; point < n_q_points; ++point)
1626 internal::evaluate_tensor_product_higher_derivatives<3>(
1630 renumber_lexicographic_to_hierarchic));
1632 data.output_data->inverse_jacobians[point].transpose();
1635 for (
unsigned int i = 0; i < spacedim; ++i)
1636 for (
unsigned int j = 0;
j < spacedim; ++
j)
1637 for (
unsigned int l = 0; l < dim; ++l)
1638 for (
unsigned int m = 0; m < dim; ++m)
1641 third[i][0][l][m] * covariant[
j][0];
1642 for (
unsigned int jr = 1;
jr < dim; ++
jr)
1648 for (
unsigned int i = 0; i < spacedim; ++i)
1649 for (
unsigned int j = 0;
j < spacedim; ++
j)
1650 for (
unsigned int l = 0; l < spacedim; ++l)
1651 for (
unsigned int m = 0; m < dim; ++m)
1654 tmp[i][
j][0][m] * covariant[l][0];
1655 for (
unsigned int lr = 1;
lr < dim; ++
lr)
1657 tmp[i][
j][
lr][m] * covariant[l][
lr];
1661 for (
unsigned int i = 0; i < spacedim; ++i)
1662 for (
unsigned int j = 0;
j < spacedim; ++
j)
1663 for (
unsigned int l = 0; l < spacedim; ++l)
1664 for (
unsigned int m = 0; m < spacedim; ++m)
1666 jacobian_pushed_forward_2nd_derivatives
1667 [point][i][
j][l][m] =
1668 tmp2[i][
j][l][0] * covariant[m][0];
1669 for (
unsigned int mr = 1;
mr < dim; ++
mr)
1670 jacobian_pushed_forward_2nd_derivatives[point][i]
1682 template <
int dim,
int spacedim,
int length_tensor>
1689 for (
unsigned int i = 0; i < spacedim; ++i)
1692 result[i][0][0][0][0] = compressed[0][i];
1695 for (
unsigned int d = 0; d < 2; ++d)
1696 for (
unsigned int e = 0; e < 2; ++e)
1697 for (
unsigned int f = 0; f < 2; ++f)
1698 for (
unsigned int g = 0;
g < 2; ++
g)
1699 result[i][d][e][f][
g] = compressed[d + e + f +
g][i];
1707 for (
unsigned int d = 0; d < 2; ++d)
1708 for (
unsigned int e = 0; e < 2; ++e)
1709 for (
unsigned int f = 0; f < 2; ++f)
1711 result[i][d][e][f][2] = compressed[5 + d + e + f][i];
1712 result[i][d][e][2][f] = compressed[5 + d + e + f][i];
1713 result[i][d][2][e][f] = compressed[5 + d + e + f][i];
1714 result[i][2][d][e][f] = compressed[5 + d + e + f][i];
1716 for (
unsigned int d = 0; d < 2; ++d)
1717 for (
unsigned int e = 0; e < 2; ++e)
1719 result[i][d][e][2][2] = compressed[9 + d + e][i];
1720 result[i][d][2][e][2] = compressed[9 + d + e][i];
1721 result[i][d][2][2][e] = compressed[9 + d + e][i];
1722 result[i][2][d][e][2] = compressed[9 + d + e][i];
1723 result[i][2][d][2][e] = compressed[9 + d + e][i];
1724 result[i][2][2][d][e] = compressed[9 + d + e][i];
1726 for (
unsigned int d = 0; d < 2; ++d)
1728 result[i][d][2][2][2] = compressed[12 + d][i];
1729 result[i][2][d][2][2] = compressed[12 + d][i];
1730 result[i][2][2][d][2] = compressed[12 + d][i];
1731 result[i][2][2][2][d] = compressed[12 + d][i];
1733 result[i][2][2][2][2] = compressed[14][i];
1748 template <
int dim,
int spacedim>
1752 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1755 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
1761 data.mapping_support_points);
1762 const unsigned int n_q_points = jacobian_3rd_derivatives.size();
1766 for (
unsigned int point = 0; point < n_q_points; ++point)
1769 internal::evaluate_tensor_product_higher_derivatives<4>(
1773 renumber_lexicographic_to_hierarchic));
1788 template <
int dim,
int spacedim>
1792 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1795 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
1801 data.mapping_support_points);
1802 const unsigned int n_q_points =
1803 jacobian_pushed_forward_3rd_derivatives.size();
1812 for (
unsigned int point = 0; point < n_q_points; ++point)
1816 internal::evaluate_tensor_product_higher_derivatives<4>(
1820 renumber_lexicographic_to_hierarchic));
1823 data.output_data->inverse_jacobians[point].transpose();
1826 for (
unsigned int i = 0; i < spacedim; ++i)
1827 for (
unsigned int j = 0;
j < spacedim; ++
j)
1828 for (
unsigned int l = 0; l < dim; ++l)
1829 for (
unsigned int m = 0; m < dim; ++m)
1830 for (
unsigned int n = 0; n < dim; ++n)
1832 tmp[i][
j][l][m][n] =
1833 fourth[i][0][l][m][n] * covariant[
j][0];
1834 for (
unsigned int jr = 1;
jr < dim; ++
jr)
1835 tmp[i][
j][l][m][n] +=
1840 for (
unsigned int i = 0; i < spacedim; ++i)
1841 for (
unsigned int j = 0;
j < spacedim; ++
j)
1842 for (
unsigned int l = 0; l < spacedim; ++l)
1843 for (
unsigned int m = 0; m < dim; ++m)
1844 for (
unsigned int n = 0; n < dim; ++n)
1846 tmp2[i][
j][l][m][n] =
1847 tmp[i][
j][0][m][n] * covariant[l][0];
1848 for (
unsigned int lr = 1;
lr < dim; ++
lr)
1849 tmp2[i][
j][l][m][n] +=
1850 tmp[i][
j][
lr][m][n] * covariant[l][
lr];
1854 for (
unsigned int i = 0; i < spacedim; ++i)
1855 for (
unsigned int j = 0;
j < spacedim; ++
j)
1856 for (
unsigned int l = 0; l < spacedim; ++l)
1857 for (
unsigned int m = 0; m < spacedim; ++m)
1858 for (
unsigned int n = 0; n < dim; ++n)
1860 tmp[i][
j][l][m][n] =
1861 tmp2[i][
j][l][0][n] * covariant[m][0];
1862 for (
unsigned int mr = 1;
mr < dim; ++
mr)
1863 tmp[i][
j][l][m][n] +=
1864 tmp2[i][
j][l][
mr][n] * covariant[m][
mr];
1868 for (
unsigned int i = 0; i < spacedim; ++i)
1869 for (
unsigned int j = 0;
j < spacedim; ++
j)
1870 for (
unsigned int l = 0; l < spacedim; ++l)
1871 for (
unsigned int m = 0; m < spacedim; ++m)
1872 for (
unsigned int n = 0; n < spacedim; ++n)
1874 jacobian_pushed_forward_3rd_derivatives
1875 [point][i][
j][l][m][n] =
1876 tmp[i][
j][l][m][0] * covariant[n][0];
1877 for (
unsigned int nr = 1;
nr < dim; ++
nr)
1878 jacobian_pushed_forward_3rd_derivatives[point]
1881 tmp[i][
j][l][m][
nr] * covariant[n][
nr];
1899 template <
int dim,
int spacedim>
1903 const typename ::Triangulation<dim, spacedim>::cell_iterator &cell,
1906 const unsigned int n_q_points,
1907 const std::vector<double> &weights,
1908 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
1932 for (
unsigned int d = 0; d != dim - 1; ++d)
1952 if (dim == spacedim)
1954 for (
unsigned int i = 0; i < n_q_points; ++i)
1989 for (
unsigned int point = 0; point < n_q_points; ++point)
1992 data.output_data->jacobians[point];
1997 contravariant.transpose()[0];
2006 contravariant.transpose();
2022 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
2036 for (
unsigned int i = 0; i < output_data.
normal_vectors.size(); ++i)
2051 template <
int dim,
int spacedim>
2055 const typename ::Triangulation<dim, spacedim>::cell_iterator &cell,
2060 const typename ::MappingQ<dim, spacedim>::InternalData &
data,
2062 const std::vector<unsigned int> &renumber_lexicographic_to_hierarchic,
2068 if (dim > 1 &&
data.tensor_product_quadrature)
2086 renumber_lexicographic_to_hierarchic,
2095 renumber_lexicographic_to_hierarchic,
2103 renumber_lexicographic_to_hierarchic,
2110 renumber_lexicographic_to_hierarchic,
2117 renumber_lexicographic_to_hierarchic,
2124 renumber_lexicographic_to_hierarchic,
2131 renumber_lexicographic_to_hierarchic,
2139 quadrature.get_weights(),
2149 template <
int dim,
int spacedim,
int rank>
2159 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2160 &mapping_data) !=
nullptr),
2162 const typename ::MappingQ<dim, spacedim>::InternalData &
data =
2164 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2167 switch (mapping_kind)
2173 "update_contravariant_transformation"));
2175 for (
unsigned int i = 0; i < output.size(); ++i)
2186 "update_contravariant_transformation"));
2189 "update_volume_elements"));
2194 for (
unsigned int i = 0; i < output.size(); ++i)
2199 output[i] /=
data.volume_elements[i];
2210 "update_covariant_transformation"));
2212 for (
unsigned int i = 0; i < output.size(); ++i)
2214 data.output_data->inverse_jacobians[i].transpose(), input[i]);
2229 template <
int dim,
int spacedim,
int rank>
2239 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2240 &mapping_data) !=
nullptr),
2242 const typename ::MappingQ<dim, spacedim>::InternalData &
data =
2244 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2247 switch (mapping_kind)
2253 "update_covariant_transformation"));
2256 "update_contravariant_transformation"));
2259 for (
unsigned int i = 0; i < output.size(); ++i)
2265 data.output_data->inverse_jacobians[i].transpose(),
2276 "update_covariant_transformation"));
2279 for (
unsigned int i = 0; i < output.size(); ++i)
2282 data.output_data->inverse_jacobians[i].transpose();
2295 "update_covariant_transformation"));
2298 "update_contravariant_transformation"));
2301 "update_volume_elements"));
2304 for (
unsigned int i = 0; i < output.size(); ++i)
2307 data.output_data->inverse_jacobians[i].transpose();
2315 output[i] /=
data.volume_elements[i];
2331 template <
int dim,
int spacedim>
2341 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2342 &mapping_data) !=
nullptr),
2344 const typename ::MappingQ<dim, spacedim>::InternalData &
data =
2346 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2349 switch (mapping_kind)
2355 "update_covariant_transformation"));
2358 "update_contravariant_transformation"));
2360 for (
unsigned int q = 0;
q < output.size(); ++
q)
2363 data.output_data->inverse_jacobians[
q].transpose();
2365 data.output_data->jacobians[
q];
2367 for (
unsigned int i = 0; i < spacedim; ++i)
2369 double tmp1[dim][dim];
2370 for (
unsigned int J = 0; J < dim; ++J)
2371 for (
unsigned int K = 0; K < dim; ++K)
2374 contravariant[i][0] * input[
q][0][J][K];
2375 for (
unsigned int I = 1; I < dim; ++I)
2377 contravariant[i][I] * input[
q][I][J][K];
2379 for (
unsigned int j = 0;
j < spacedim; ++
j)
2382 for (
unsigned int K = 0; K < dim; ++K)
2385 for (
unsigned int J = 1; J < dim; ++J)
2386 tmp2[K] += covariant[
j][J] *
tmp1[J][K];
2388 for (
unsigned int k = 0;
k < spacedim; ++
k)
2390 output[
q][i][
j][
k] = covariant[
k][0] *
tmp2[0];
2391 for (
unsigned int K = 1; K < dim; ++K)
2392 output[
q][i][
j][
k] += covariant[
k][K] *
tmp2[K];
2404 "update_covariant_transformation"));
2406 for (
unsigned int q = 0;
q < output.size(); ++
q)
2409 data.output_data->inverse_jacobians[
q].transpose();
2411 for (
unsigned int i = 0; i < spacedim; ++i)
2413 double tmp1[dim][dim];
2414 for (
unsigned int J = 0; J < dim; ++J)
2415 for (
unsigned int K = 0; K < dim; ++K)
2417 tmp1[J][K] = covariant[i][0] * input[
q][0][J][K];
2418 for (
unsigned int I = 1; I < dim; ++I)
2419 tmp1[J][K] += covariant[i][I] * input[
q][I][J][K];
2421 for (
unsigned int j = 0;
j < spacedim; ++
j)
2424 for (
unsigned int K = 0; K < dim; ++K)
2427 for (
unsigned int J = 1; J < dim; ++J)
2428 tmp2[K] += covariant[
j][J] *
tmp1[J][K];
2430 for (
unsigned int k = 0;
k < spacedim; ++
k)
2432 output[
q][i][
j][
k] = covariant[
k][0] *
tmp2[0];
2433 for (
unsigned int K = 1; K < dim; ++K)
2434 output[
q][i][
j][
k] += covariant[
k][K] *
tmp2[K];
2447 "update_covariant_transformation"));
2450 "update_contravariant_transformation"));
2453 "update_volume_elements"));
2455 for (
unsigned int q = 0;
q < output.size(); ++
q)
2458 data.output_data->inverse_jacobians[
q].transpose();
2460 data.output_data->jacobians[
q];
2461 for (
unsigned int i = 0; i < spacedim; ++i)
2464 for (
unsigned int I = 0; I < dim; ++I)
2466 contravariant[i][I] * (1. /
data.volume_elements[
q]);
2467 double tmp1[dim][dim];
2468 for (
unsigned int J = 0; J < dim; ++J)
2469 for (
unsigned int K = 0; K < dim; ++K)
2471 tmp1[J][K] = factor[0] * input[
q][0][J][K];
2472 for (
unsigned int I = 1; I < dim; ++I)
2473 tmp1[J][K] += factor[I] * input[
q][I][J][K];
2475 for (
unsigned int j = 0;
j < spacedim; ++
j)
2478 for (
unsigned int K = 0; K < dim; ++K)
2481 for (
unsigned int J = 1; J < dim; ++J)
2482 tmp2[K] += covariant[
j][J] *
tmp1[J][K];
2484 for (
unsigned int k = 0;
k < spacedim; ++
k)
2486 output[
q][i][
j][
k] = covariant[
k][0] *
tmp2[0];
2487 for (
unsigned int K = 1; K < dim; ++K)
2488 output[
q][i][
j][
k] += covariant[
k][K] *
tmp2[K];
2508 template <
int dim,
int spacedim,
int rank>
2518 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2519 &mapping_data) !=
nullptr),
2521 const typename ::MappingQ<dim, spacedim>::InternalData &
data =
2523 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2526 switch (mapping_kind)
2532 "update_covariant_transformation"));
2534 for (
unsigned int i = 0; i < output.size(); ++i)
2536 data.output_data->inverse_jacobians[i].transpose(), input[i]);
ArrayView< std::remove_reference_t< typename std::iterator_traits< Iterator >::reference >, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Abstract base class for mapping classes.
numbers::NumberTraits< Number >::real_type norm() const
constexpr numbers::NumberTraits< Number >::real_type norm_square() const
static constexpr std::size_t size()
const Point< spacedim > normalization_shift
const double normalization_length
InverseQuadraticApproximation(const InverseQuadraticApproximation &)=default
static constexpr unsigned int n_functions
Point< dim, Number > compute(const Point< spacedim, Number > &p) const
InverseQuadraticApproximation(const ArrayView< const Point< spacedim > > &real_support_points, const std::vector< Point< dim > > &unit_support_points)
std::array< Point< dim >, n_functions > coefficients
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ update_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_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ 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
@ mapping_contravariant_hessian
@ mapping_covariant_hessian
@ mapping_contravariant_gradient
@ tensor_symmetric_collocation
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
std::vector< index_type > data
EvaluationFlags
The EvaluationFlags enum.
constexpr T pow(const T base, const int iexp)
Point< 1 > transform_real_to_unit_cell(const std::array< Point< spacedim >, GeometryInfo< 1 >::vertices_per_cell > &vertices, const Point< spacedim > &p)
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)
Point< dim, Number > do_transform_real_to_unit_cell_internal(const Point< spacedim, Number > &p, const Point< dim, Number > &initial_p_unit, const ArrayView< const Point< spacedim > > &points, const std::vector< Polynomials::Polynomial< double > > &polynomials_1d, const std::vector< unsigned int > &renumber, const bool print_iterations_to_deallog=false)
DerivativeForm< 3, dim, spacedim > expand_3rd_derivatives(const Tensor< 1, length_tensor, Tensor< 1, spacedim > > &compressed)
inline ::Table< 2, double > compute_support_point_weights_cell(const unsigned int polynomial_degree)
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 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)
inline ::Table< 2, double > compute_support_point_weights_on_hex(const unsigned int polynomial_degree)
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)
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 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)
inline ::Table< 2, double > compute_support_point_weights_on_quad(const unsigned int polynomial_degree)
void maybe_update_q_points_Jacobians_generic(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< Point< spacedim > > &quadrature_points, std::vector< DerivativeForm< 1, dim, spacedim > > &jacobians, std::vector< DerivativeForm< 1, spacedim, dim > > &inverse_jacobians)
std::vector< Point< dim > > unit_support_points(const std::vector< Point< 1 > > &line_support_points, const std::vector< unsigned int > &renumbering)
void transform_gradients(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)
void maybe_update_q_points_Jacobians_and_grads_tensor(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< Point< spacedim > > &quadrature_points, std::vector< DerivativeForm< 1, dim, spacedim > > &jacobians, std::vector< DerivativeForm< 1, spacedim, dim > > &inverse_jacobians, std::vector< DerivativeForm< 2, dim, spacedim > > &jacobian_grads)
void transform_hessians(const ArrayView< const Tensor< 3, dim > > &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< 3, spacedim > > &output)
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)
std::vector<::Table< 2, double > > compute_support_point_weights_perimeter_to_interior(const unsigned int polynomial_degree, const unsigned int dim)
DerivativeForm< 4, dim, spacedim > expand_4th_derivatives(const Tensor< 1, length_tensor, Tensor< 1, spacedim > > &compressed)
Point< spacedim > compute_mapped_location_of_point(const typename ::MappingQ< dim, spacedim >::InternalData &data)
Point< dim > do_transform_real_to_unit_cell_internal_codim1(const Point< dim+1 > &p, const Point< dim > &initial_p_unit, const ArrayView< const Point< dim+1 > > &points, const std::vector< Polynomials::Polynomial< double > > &polynomials_1d, const std::vector< unsigned int > &renumber)
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 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)
SymmetricTensor< 2, dim, typename ProductTypeNoPoint< Number, Number2 >::type > evaluate_tensor_product_hessian(const std::vector< Polynomials::Polynomial< double > > &poly, const ArrayView< const Number > &values, const Point< dim, Number2 > &p, const std::vector< unsigned int > &renumber={})
std::pair< typename ProductTypeNoPoint< Number, Number2 >::type, Tensor< 1, dim, typename ProductTypeNoPoint< Number, Number2 >::type > > evaluate_tensor_product_value_and_gradient(const std::vector< Polynomials::Polynomial< double > > &poly, const ArrayView< const Number > &values, const Point< dim, Number2 > &p, const bool d_linear=false, const std::vector< unsigned int > &renumber={})
ProductTypeNoPoint< Number, Number2 >::type evaluate_tensor_product_value(const std::vector< Polynomials::Polynomial< double > > &poly, const ArrayView< const Number > &values, const Point< dim, Number2 > &p, const bool d_linear=false, const std::vector< unsigned int > &renumber={})
constexpr unsigned int invalid_unsigned_int
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
static double subface_ratio(const internal::SubfaceCase< dim > &subface_case, const unsigned int subface_no)
static Point< dim, Number > project_to_unit_cell(const Point< dim, Number > &p)
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices()
static void evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number *values_dofs, FEEvaluationData< dim, Number, false > &fe_eval)
DEAL_II_HOST constexpr Number determinant(const SymmetricTensor< 2, dim, Number > &)
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &)