16#ifndef dealii_matrix_free_evaluation_kernels_face_h
17#define dealii_matrix_free_evaluation_kernels_face_h
64 const unsigned int subface_index,
65 const unsigned int direction)
81 const unsigned int index =
82 direction == 0 ? subface_index % 2 : subface_index / 2;
93 const unsigned int n_components,
98 Number *gradients_quad,
99 Number *hessians_quad,
100 Number *scratch_data,
101 const unsigned int subface_index)
106 const std::size_t n_dofs = fe_degree > -1 ?
109 const std::size_t n_q_points =
118 for (
unsigned int c = 0; c < n_components; ++c)
133 values_quad[0] = values_dofs[0];
140 values_dofs += 3 * n_dofs;
141 values_quad += n_q_points;
144 for (
unsigned int c = 0; c < n_components; ++c)
164 values_quad, gradients_quad);
166 values_quad, gradients_quad + 1);
180 scratch_data, gradients_quad + 1);
203 values_quad[0] = values_dofs[0];
204 gradients_quad[0] = values_dofs[1];
209 values_dofs += 3 * n_dofs;
210 values_quad += n_q_points;
211 gradients_quad += dim * n_q_points;
217 for (
unsigned int c = 0; c < n_components; ++c)
272 values_dofs + 2 * n_dofs, hessians_quad + n_q_points);
275 values_dofs + n_dofs, hessians_quad + 2 * n_q_points);
278 hessians_quad[0] = values_dofs[2];
283 values_dofs += 3 * n_dofs;
284 hessians_quad += dim * (dim + 1) / 2 * n_q_points;
291 const unsigned int n_components,
296 Number *gradients_quad,
297 Number *hessians_quad,
298 Number *scratch_data,
299 const unsigned int subface_index)
304 const std::size_t n_dofs =
308 const std::size_t n_q_points =
317 for (
unsigned int c = 0; c < n_components; ++c)
332 values_dofs[0] = values_quad[0];
337 values_dofs += 3 * n_dofs;
338 values_quad += n_q_points;
341 for (
unsigned int c = 0; c < n_components; ++c)
350 values_dofs + n_dofs);
363 gradients_quad + 1, values_quad);
366 gradients_quad + 1, values_quad);
368 gradients_quad, values_quad);
381 gradients_quad + 1, scratch_data);
385 gradients_quad + 1, scratch_data);
409 values_dofs[0] = values_quad[0];
410 values_dofs[1] = gradients_quad[0];
415 values_dofs += 3 * n_dofs;
416 values_quad += n_q_points;
417 gradients_quad += dim * n_q_points;
423 for (
unsigned int c = 0; c < n_components; ++c)
479 values_dofs + n_dofs);
494 hessians_quad + n_q_points, values_dofs + 2 * n_dofs);
498 hessians_quad + 2 * n_q_points, values_dofs + n_dofs);
501 hessians_quad + 2 * n_q_points, values_dofs + n_dofs);
504 values_dofs[2] = hessians_quad[0];
513 values_dofs += 3 * n_dofs;
514 hessians_quad += dim * (dim + 1) / 2 * n_q_points;
522 template <
int dim,
int fe_degree,
int n_q_po
ints_1d,
typename Number>
532 template <
bool do_
integrate>
541 Number *scratch_data,
542 const unsigned int subface_index,
547 const int degree = fe_degree != -1 ? fe_degree :
shape_data[0].fe_degree;
564 (fe_degree > 0 ? fe_degree : 0),
586 std::array<unsigned int, dim> components;
590 for (
const unsigned int comp : components)
594 std::array<int, 2> n_blocks{
599 if constexpr (dim == 3)
608 shape_data[0].shape_gradients_collocation_eo.data(),
631 n_blocks[0] * n_blocks[1],
651 n_blocks[0] * n_blocks[1],
670 scratch_data, gradients + 2);
719 n_blocks[0] * n_blocks[1]);
739 n_blocks[0] * n_blocks[1]);
753 gradients + 2, scratch_data);
756 values_dofs + n_blocks[0] * n_blocks[1]);
781 values_dofs, gradients);
783 values_dofs +
n_rows_n, gradients + 1);
795 values_dofs, gradients);
797 values_dofs +
n_rows_t, gradients + 1);
816 gradients, values_dofs);
819 gradients, values_dofs);
821 gradients + 1, values_dofs +
n_rows_n);
836 gradients, values_dofs);
839 gradients, values_dofs);
841 gradients + 1, values_dofs +
n_rows_t);
854 template <
int dim,
int fe_degree,
typename Number>
860 template <
bool do_evaluate,
bool add_
into_output>
869 Assert(
static_cast<unsigned int>(fe_degree) ==
870 shape_info.
data.front().fe_degree ||
875 n_components, input, output, flags,
face_no, shape_info);
883 shape_info.
data.front().fe_degree + 1,
884 shape_info.
data.front().shape_data_on_face,
892 template <
bool do_evaluate,
bool add_
into_output>
895 const unsigned int n_components,
902 Assert(
static_cast<unsigned int>(fe_degree + 1) ==
903 shape_info.
data.front().n_q_points_1d ||
913 shape_info.
data.front().quadrature.size(),
914 shape_info.
data.front().quadrature_data_on_face,
920 template <
bool do_evaluate,
bool add_
into_output,
int face_direction = 0>
929 const unsigned int dofs_per_component_on_cell,
930 const unsigned int dofs_per_component_on_face)
936 const int n_rows = fe_degree != -1 ? fe_degree + 1 :
n_points_1d;
938 const std::array<int, 2> n_blocks{
939 {(dim > 1 ? n_rows : 1), (dim > 2 ? n_rows : 1)}};
940 std::array<int, 2> steps;
942 steps = {{n_rows, 0}};
947 steps = {{n_rows * n_rows, -n_rows * n_rows * n_rows + 1}};
951 for (
unsigned int c = 0; c < n_components; ++c)
991 input += dofs_per_component_on_cell;
992 output += dofs_per_component_on_face;
996 output += dofs_per_component_on_cell;
997 input += dofs_per_component_on_face;
1013 dofs_per_component_on_cell,
1014 dofs_per_component_on_face);
1018 template <
bool do_evaluate,
1024 const unsigned int n_components,
1025 const Number *input,
1049 constexpr int stride2 = fe_degree * (fe_degree + 1);
1052 fe_degree != -1 ? fe_degree : shape_info.
data[0].fe_degree;
1056 std::array<int, 3>
strides{{1, 1, 1}};
1069 std::array<int, 2> steps, n_blocks;
1092 .shape_data_on_face[
face_no % 2]
1104 output += 3 * n_blocks[0] * n_blocks[1];
1109 input += 3 * n_blocks[0] * n_blocks[1];
1114 steps = {{degree, 0}};
1126 .shape_data_on_face[
face_no % 2]
1135 if constexpr (dim > 2)
1140 output += 3 * n_blocks[0] * n_blocks[1];
1145 input += 3 * n_blocks[0] * n_blocks[1];
1149 steps = {{degree, 0}};
1166 .shape_data_on_face[
face_no % 2]
1183 n_components, input, output, flag,
face_no, shape_info);
1193 n_components, input, output, flag,
face_no, shape_info);
1201 n_components, input, output, flag,
face_no, shape_info);
1210 template <
typename VectorizedArrayType,
typename Number2>
1214 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
1222 template <
typename Number, std::
size_t w
idth>
1232 template <
typename VectorizedArrayType,
typename Number2>
1235 const unsigned int *indices,
1236 VectorizedArrayType &dst)
1238 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
1246 template <
typename Number, std::
size_t w
idth>
1249 const unsigned int *indices,
1258 template <
typename VectorizedArrayType,
typename Number2>
1262 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
1270 template <
typename Number, std::
size_t w
idth>
1282 template <
typename VectorizedArrayType,
typename Number2>
1285 const unsigned int *indices,
1288 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
1289 dst_ptr[indices[v]] += src[v];
1296 template <
typename Number, std::
size_t w
idth>
1299 const unsigned int *indices,
1302#if DEAL_II_VECTORIZATION_WIDTH_IN_BITS < 512
1303 for (
unsigned int v = 0; v < width; ++v)
1304 dst_ptr[indices[v]] += src[v];
1308 (tmp + src).scatter(indices,
dst_ptr);
1314 template <
typename Number>
1317 const unsigned int n_components,
1319 const unsigned int *orientation,
1320 const bool integrate,
1321 const std::size_t n_q_points,
1323 Number *values_quad,
1324 Number *gradients_quad,
1325 Number *hessians_quad)
1327 for (
unsigned int c = 0; c < n_components; ++c)
1332 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1333 tmp_values[
q] = values_quad[c * n_q_points + orientation[
q]];
1335 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1336 tmp_values[orientation[
q]] = values_quad[c * n_q_points +
q];
1337 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1341 for (
unsigned int d = 0; d < dim; ++d)
1344 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1346 gradients_quad[(c * n_q_points + orientation[
q]) * dim + d];
1348 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1350 gradients_quad[(c * n_q_points +
q) * dim + d];
1351 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1352 gradients_quad[(c * n_q_points +
q) * dim + d] =
tmp_values[
q];
1356 const unsigned int hdim = (dim * (dim + 1)) / 2;
1357 for (
unsigned int d = 0; d <
hdim; ++d)
1360 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1364 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1366 hessians_quad[(c *
hdim + d) * n_q_points +
q];
1367 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1368 hessians_quad[(c *
hdim + d) * n_q_points +
q] =
1377 template <
typename Number,
typename VectorizedArrayType>
1380 const unsigned int dim,
1381 const unsigned int n_components,
1382 const unsigned int v,
1384 const unsigned int *orientation,
1385 const bool integrate,
1386 const std::size_t n_q_points,
1388 VectorizedArrayType *values_quad,
1389 VectorizedArrayType *gradients_quad =
nullptr,
1390 VectorizedArrayType *hessians_quad =
nullptr)
1392 for (
unsigned int c = 0; c < n_components; ++c)
1397 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1398 tmp_values[
q] = values_quad[c * n_q_points + orientation[
q]][v];
1400 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1401 tmp_values[orientation[
q]] = values_quad[c * n_q_points +
q][v];
1402 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1406 for (
unsigned int d = 0; d < dim; ++d)
1410 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1412 gradients_quad[(c * n_q_points + orientation[
q]) * dim + d]
1415 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1417 gradients_quad[(c * n_q_points +
q) * dim + d][v];
1418 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1419 gradients_quad[(c * n_q_points +
q) * dim + d][v] =
1425 const unsigned int hdim = (dim * (dim + 1)) / 2;
1426 for (
unsigned int d = 0; d <
hdim; ++d)
1429 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1433 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1435 hessians_quad[(c *
hdim + d) * n_q_points +
q][v];
1436 for (
unsigned int q = 0;
q < n_q_points; ++
q)
1437 hessians_quad[(c *
hdim + d) * n_q_points +
q][v] =
1446 template <
int dim,
typename Number>
1452 const Number *values_dofs,
1455 const auto &shape_info = fe_eval.get_shape_info();
1456 const auto &
shape_data = shape_info.data.front();
1460 Assert((fe_eval.get_dof_access_index() ==
1462 fe_eval.is_interior_face() ==
false) ==
false,
1465 const unsigned int face_no = fe_eval.get_face_no();
1466 const unsigned int face_orientation = fe_eval.get_face_orientation();
1467 const std::size_t n_dofs = shape_info.dofs_per_component_on_cell;
1468 const std::size_t n_q_points = shape_info.n_q_points_faces[
face_no];
1475 const auto *
const shape_values =
1481 Eval
eval(shape_values,
nullptr,
nullptr, n_dofs, n_q_points);
1482 for (
unsigned int c = 0; c < n_components; ++c)
1497 std::array<const Number2 *, dim> shape_gradients;
1498 for (
unsigned int d = 0; d < dim; ++d)
1499 shape_gradients[d] =
1502 for (
unsigned int c = 0; c < n_components; ++c)
1504 for (
unsigned int d = 0; d < dim; ++d)
1507 nullptr, shape_gradients[d],
nullptr, n_dofs, n_q_points);
1523 template <
int fe_degree>
1530 const Number *values_dofs,
1534 Number *scratch_data)
1536 const auto &shape_info = fe_eval.get_shape_info();
1540 const auto &
shape_data = shape_info.data.front();
1548 for (
unsigned int v = 0; v < Number::size(); ++v)
1555 for (
unsigned int i = 0; i < 3 * dofs_per_face; ++i)
1566 fe_eval.get_face_no(v));
1568 for (
unsigned int i = 0; i < 3 * dofs_per_face; ++i)
1569 temp[i][v] = scratch_data[i][v];
1579 fe_eval.get_face_no());
1583 template <
int fe_degree,
int n_q_po
ints_1d>
1592 Number *scratch_data)
1594 const auto &shape_info = fe_eval.get_shape_info();
1595 const auto &
shape_data = shape_info.data.front();
1597 const unsigned int subface_index = fe_eval.get_subface_index();
1599 fe_degree > -1 ? n_q_points_1d : 0;
1609 fe_eval.get_shape_info().data,
1611 fe_eval.begin_values(),
1612 fe_eval.begin_gradients(),
1615 fe_eval.get_face_no() / 2);
1617 else if (fe_degree > -1 &&
1628 fe_eval.begin_values(),
1631 fe_eval.begin_hessians(),
1643 fe_eval.begin_values(),
1646 fe_eval.begin_hessians(),
1656 const unsigned int n_components,
1662 const auto &shape_info = fe_eval.get_shape_info();
1666 for (
unsigned int v = 0; v < Number::size(); ++v)
1674 if (fe_eval.get_face_orientation(v) != 0)
1680 &shape_info.face_orientations_quad(
1681 fe_eval.get_face_orientation(v), 0),
1683 shape_info.n_q_points_face,
1685 fe_eval.begin_values(),
1686 fe_eval.begin_gradients(),
1687 fe_eval.begin_hessians());
1690 else if (fe_eval.get_face_orientation() != 0)
1695 &shape_info.face_orientations_quad(fe_eval.get_face_orientation(), 0),
1697 shape_info.n_q_points_face,
1699 fe_eval.begin_values(),
1700 fe_eval.begin_gradients(),
1701 fe_eval.begin_hessians());
1706 template <
int fe_degree,
int n_q_po
ints_1d>
1710 const Number *values_dofs,
1713 const auto &shape_info = fe_eval.get_shape_info();
1714 const auto &
shape_data = shape_info.data.front();
1723 Number *
temp = fe_eval.get_scratch_data().begin();
1727 if (fe_eval.get_dof_access_index() ==
1729 fe_eval.is_interior_face() ==
false)
1730 for (
unsigned int v = 0; v < Number::size(); ++v)
1732 fe_eval.get_face_no(v) != fe_eval.get_face_no(0))
1753 template <
int fe_degree,
int n_q_po
ints_1d>
1755 run(
const unsigned int n_components,
1757 const Number *values_dofs,
1760 const auto &shape_info = fe_eval.get_shape_info();
1777 template <
int dim,
typename Number>
1780 template <
int fe_degree>
1782 run(
const unsigned int n_components,
1784 const Number *values_dofs,
1787 const auto &shape_info = fe_eval.get_shape_info();
1788 const auto &
shape_data = shape_info.data.front();
1797 Number *
temp = fe_eval.get_scratch_data().begin();
1801 if (fe_eval.get_dof_access_index() ==
1803 fe_eval.is_interior_face() ==
false)
1804 for (
unsigned int v = 0; v < Number::size(); ++v)
1806 fe_eval.get_face_no(v) != fe_eval.get_face_no(0))
1824 template <
int dim,
typename Number>
1827 template <
int fe_degree,
int n_q_po
ints_1d>
1829 run(
const unsigned int n_components,
1833 const auto &shape_info = fe_eval.get_shape_info();
1834 const auto &
shape_data = shape_info.data.front();
1843 Number *
temp = fe_eval.get_scratch_data().begin();
1856 template <
int dim,
typename Number>
1861 const unsigned int n_components,
1863 Number *values_dofs,
1867 const auto &shape_info = fe_eval.get_shape_info();
1868 const auto &
shape_data = shape_info.data.front();
1872 Assert((fe_eval.get_dof_access_index() ==
1874 fe_eval.is_interior_face() ==
false) ==
false,
1877 const unsigned int face_no = fe_eval.get_face_no();
1878 const unsigned int face_orientation = fe_eval.get_face_orientation();
1879 const std::size_t n_dofs = shape_info.dofs_per_component_on_cell;
1880 const std::size_t n_q_points = shape_info.n_q_points_faces[
face_no];
1887 const auto *
const shape_values =
1893 Eval
eval(shape_values,
nullptr,
nullptr, n_dofs, n_q_points);
1894 for (
unsigned int c = 0; c < n_components; ++c)
1912 std::array<const Number2 *, dim> shape_gradients;
1913 for (
unsigned int d = 0; d < dim; ++d)
1914 shape_gradients[d] =
1917 for (
unsigned int c = 0; c < n_components; ++c)
1919 for (
unsigned int d = 0; d < dim; ++d)
1922 nullptr, shape_gradients[d],
nullptr, n_dofs, n_q_points);
1948 const unsigned int n_components,
1954 const auto &shape_info = fe_eval.get_shape_info();
1958 for (
unsigned int v = 0; v < Number::size(); ++v)
1966 if (fe_eval.get_face_orientation(v) != 0)
1972 &fe_eval.get_shape_info().face_orientations_quad(
1973 fe_eval.get_face_orientation(v), 0),
1975 shape_info.n_q_points_face,
1977 fe_eval.begin_values(),
1978 fe_eval.begin_gradients(),
1979 fe_eval.begin_hessians());
1982 else if (fe_eval.get_face_orientation() != 0)
1987 &fe_eval.get_shape_info().face_orientations_quad(
1988 fe_eval.get_face_orientation(), 0),
1990 shape_info.n_q_points_face,
1992 fe_eval.begin_values(),
1993 fe_eval.begin_gradients(),
1994 fe_eval.begin_hessians());
1997 template <
int fe_degree,
int n_q_po
ints_1d>
2006 Number *scratch_data)
2008 const auto &shape_info = fe_eval.get_shape_info();
2009 const auto &
shape_data = shape_info.data.front();
2012 fe_degree > -1 ? n_q_points_1d : 0;
2013 const unsigned int subface_index = fe_eval.get_subface_index();
2023 fe_eval.get_shape_info().data,
2025 fe_eval.begin_values(),
2026 fe_eval.begin_gradients(),
2029 fe_eval.get_face_no() / 2);
2031 else if (fe_degree > -1 &&
2032 fe_eval.get_subface_index() >=
2044 fe_eval.begin_values(),
2045 fe_eval.begin_gradients(),
2046 fe_eval.begin_hessians(),
2059 fe_eval.begin_values(),
2060 fe_eval.begin_gradients(),
2061 fe_eval.begin_hessians(),
2066 template <
int fe_degree>
2073 Number *values_dofs,
2077 Number *scratch_data,
2080 const auto &shape_info = fe_eval.get_shape_info();
2081 const auto &
shape_data = shape_info.data.front();
2091 for (
unsigned int v = 0; v < Number::size(); ++v)
2105 fe_eval.get_face_no(v));
2108 for (
unsigned int i = 0; i < 3 * dofs_per_face; ++i)
2109 values_dofs[i][v] += scratch_data[i][v];
2111 for (
unsigned int i = 0; i < 3 * dofs_per_face; ++i)
2112 values_dofs[i][v] = scratch_data[i][v];
2124 fe_eval.get_face_no());
2132 fe_eval.get_face_no());
2136 template <
int fe_degree,
int n_q_po
ints_1d>
2140 Number *values_dofs,
2144 const auto &shape_info = fe_eval.get_shape_info();
2145 const auto &
shape_data = shape_info.data.front();
2152 Number *
temp = fe_eval.get_scratch_data().begin();
2157 if (fe_eval.get_dof_access_index() ==
2159 fe_eval.is_interior_face() ==
false)
2162 std::all_of(fe_eval.get_cell_ids().begin() + 1,
2163 fe_eval.get_cell_ids().end(),
2164 [&](
const auto &v) {
2165 return v == fe_eval.get_cell_ids()[0] ||
2166 v == numbers::invalid_unsigned_int;
2188 template <
int fe_degree,
int n_q_po
ints_1d>
2190 run(
const unsigned int n_components,
2192 Number *values_dofs,
2196 const auto &shape_info = fe_eval.get_shape_info();
2215 template <
int dim,
typename Number>
2218 template <
int fe_degree>
2220 run(
const unsigned int n_components,
2222 Number *values_dofs,
2226 const auto &shape_info = fe_eval.get_shape_info();
2227 const auto &
shape_data = shape_info.data.front();
2234 Number *
temp = fe_eval.get_scratch_data().begin();
2239 if (fe_eval.get_dof_access_index() ==
2241 fe_eval.is_interior_face() ==
false)
2244 std::all_of(fe_eval.get_cell_ids().begin() + 1,
2245 fe_eval.get_cell_ids().end(),
2246 [&](
const auto &v) {
2247 return v == fe_eval.get_cell_ids()[0] ||
2248 v == numbers::invalid_unsigned_int;
2267 template <
int dim,
typename Number>
2270 template <
int fe_degree,
int n_q_po
ints_1d>
2272 run(
const unsigned int n_components,
2277 const auto &shape_info = fe_eval.get_shape_info();
2278 const auto &
shape_data = shape_info.data.front();
2285 Number *
temp = fe_eval.get_scratch_data().begin();
2298 template <
int n_face_orientations,
2305 const unsigned int n_components,
2310 typename Processor::VectorizedArrayType_ *
temp1)
2312 constexpr int dim = Processor::dim_;
2313 constexpr int fe_degree = Processor::fe_degree_;
2314 using VectorizedArrayType =
typename Processor::VectorizedArrayType_;
2315 constexpr int n_lanes = VectorizedArrayType::size();
2317 using Number =
typename Processor::Number_;
2318 using Number2_ =
typename Processor::Number2_;
2320 const auto &
shape_data = fe_eval.get_shape_info().data.front();
2321 constexpr bool integrate = Processor::do_integrate;
2322 const unsigned int face_no = fe_eval.get_face_no();
2323 const auto &dof_info = fe_eval.get_dof_info();
2324 const unsigned int cell = fe_eval.get_cell_or_face_batch_id();
2326 fe_eval.get_dof_access_index();
2328 dof_info.index_storage_variants[dof_access_index].size());
2329 constexpr unsigned int dofs_per_face =
2331 const unsigned int subface_index = fe_eval.get_subface_index();
2334 dof_info.n_vectorization_lanes_filled[dof_access_index][cell];
2337 if (n_face_orientations == n_lanes)
2338 for (
unsigned int v = 1; v < n_lanes; ++v)
2339 if (fe_eval.get_face_no(v) != fe_eval.get_face_no(0) ||
2340 fe_eval.get_face_orientation(v) != fe_eval.get_face_orientation(0))
2347 std::array<const unsigned int *, n_face_orientations> orientation = {};
2350 fe_eval.is_interior_face() == 0)
2351 for (
unsigned int v = 0; v < n_lanes; ++v)
2360 fe_eval.get_face_orientation(v) != 0)
2379 orientation[v] = &fe_eval.get_shape_info().face_orientations_dofs(
2380 fe_eval.get_face_orientation(v), 0);
2383 else if (dim == 3 && fe_eval.get_face_orientation() != 0)
2401 for (
unsigned int v = 0; v < n_face_orientations; ++v)
2402 orientation[v] = &fe_eval.get_shape_info().face_orientations_dofs(
2403 fe_eval.get_face_orientation(), 0);
2411 .shape_data_on_face[0][fe_degree + (integrate ? (2 -
face_no % 2) :
2419 if (n_face_orientations == 1)
2421 &fe_eval.get_shape_info().face_to_cell_index_hermite(
face_no, 0);
2424 for (
unsigned int v = 0; v < n_lanes; ++v)
2429 const auto face_no = fe_eval.get_face_no(v);
2438 &fe_eval.get_shape_info().face_to_cell_index_hermite(
face_no,
2447 if (
shape_data.nodal_at_cell_boundaries ==
true)
2449 if (n_face_orientations == 1)
2451 &fe_eval.get_shape_info().face_to_cell_index_nodal(
face_no, 0);
2454 for (
unsigned int v = 0; v < n_lanes; ++v)
2459 const auto face_no = fe_eval.get_face_no(v);
2462 &fe_eval.get_shape_info().face_to_cell_index_nodal(
face_no,
2469 const auto reorientate = [&](
const unsigned int v,
const unsigned int i) {
2477 fe_eval.get_cell_ids()[0] :
2479 const unsigned int *dof_indices =
2480 &dof_info.dof_indices_contiguous[dof_access_index][
cell_index];
2482 for (
unsigned int comp = 0;
comp < n_components; ++
comp)
2485 dof_info.component_dof_indices_offset
2486 [fe_eval.get_active_fe_index()]
2487 [fe_eval.get_first_selected_component()] +
2491 if (n_face_orientations == 1 &&
2492 dof_info.index_storage_variants[dof_access_index][cell] ==
2494 interleaved_contiguous)
2497 dof_info.n_vectorization_lanes_filled[dof_access_index][cell],
2499 Number2_ *vector_ptr =
2504 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2512 temp1[i_ + dofs_per_face],
2513 vector_ptr +
ind1 * n_lanes,
2514 vector_ptr +
ind2 * n_lanes,
2520 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2527 vector_ptr +
ind * n_lanes);
2533 else if (n_face_orientations == 1 &&
2534 dof_info.index_storage_variants[dof_access_index][cell] ==
2536 interleaved_contiguous_strided)
2539 dof_info.n_vectorization_lanes_filled[dof_access_index][cell],
2544 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2549 const unsigned int ind1 =
2551 const unsigned int ind2 =
2553 proc.hermite_grad_vectorized_indexed(
2555 temp1[i_ + dofs_per_face],
2565 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2571 proc.value_vectorized_indexed(
temp1[i_],
2579 else if (n_face_orientations == 1 &&
2580 dof_info.index_storage_variants[dof_access_index][cell] ==
2582 interleaved_contiguous_mixed_strides)
2585 &dof_info.dof_indices_interleave_strides[dof_access_index]
2587 unsigned int indices[n_lanes];
2588 for (
unsigned int v = 0; v < n_lanes; ++v)
2591 dof_info.n_vectorization_lanes_filled[dof_access_index][cell];
2596 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2601 unsigned int ind1[n_lanes];
2603 for (
unsigned int v = 0; v < n_lanes; ++v)
2604 ind1[v] = indices[v] +
2606 unsigned int ind2[n_lanes];
2608 for (
unsigned int v = 0; v < n_lanes; ++v)
2613 proc.hermite_grad_vectorized_indexed(
2615 temp1[i_ + dofs_per_face],
2624 if (integrate ==
false)
2625 for (
unsigned int i = 0; i < 2 * dofs_per_face; ++i)
2626 temp1[i] = VectorizedArrayType();
2629 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2631 const unsigned int i_ =
2635 temp1[i_ + dofs_per_face][v],
2639 [n_face_orientations == 1 ? 0 : v][2 * i] *
2654 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2657 unsigned int ind[n_lanes];
2659 for (
unsigned int v = 0; v < n_lanes; ++v)
2663 proc.value_vectorized_indexed(
temp1[i_],
2669 if (integrate ==
false)
2670 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2671 temp1[i] = VectorizedArrayType();
2674 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2688 else if (n_face_orientations > 1 ||
2689 dof_info.index_storage_variants[dof_access_index][cell] ==
2704 if (n_face_orientations == 1)
2715 .dof_indices_contiguous_sm[dof_access_index]
2716 [cell * n_lanes + v];
2722 else if (n_face_orientations == n_lanes)
2724 const auto &cells = fe_eval.get_cell_ids();
2725 for (
unsigned int v = 0; v < n_lanes; ++v)
2733 .dof_indices_contiguous[dof_access_index]
2740 .dof_indices_contiguous_sm[dof_access_index]
2753 else if (n_face_orientations == n_lanes)
2755 for (
unsigned int v = 0; v < n_lanes; ++v)
2757 dof_info.dof_indices_contiguous[dof_access_index]
2758 [fe_eval.get_cell_ids()[v]];
2765 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2768 const unsigned int ind2 =
2772 proc.hermite_grad_vectorized_indexed(
2774 temp1[i_ + dofs_per_face],
2781 else if (n_face_orientations == 1)
2782 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2785 const unsigned int ind2 =
2791 temp1[i_ + dofs_per_face][v],
2796 if (integrate ==
false)
2800 temp1[i + dofs_per_face][v] = 0.0;
2806 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2807 temp1[i] =
temp1[i + dofs_per_face] = Number();
2810 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2822 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2827 proc.value_vectorized_indexed(
temp1[i_],
2831 else if (n_face_orientations == 1)
2832 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2840 if (integrate ==
false)
2847 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2848 temp1[i] = Number();
2851 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2863 temp1 += 3 * dofs_per_face;
2869 template <
int dim,
typename Number2,
typename VectorizedArrayType>
2872 using Number =
typename VectorizedArrayType::value_type;
2874 template <
int fe_degree,
int n_q_po
ints_1d>
2876 run(
const unsigned int n_components,
2883 Assert(fe_eval.get_shape_info().element_type <=
2887 const unsigned int dofs_per_face =
Utilities::pow(fe_degree + 1, dim - 1);
2889 VectorizedArrayType *
temp = fe_eval.get_scratch_data().begin();
2890 VectorizedArrayType *scratch_data =
2891 temp + 3 * n_components * dofs_per_face;
2895 if (fe_eval.get_dof_access_index() ==
2897 fe_eval.is_interior_face() ==
false)
2904 const unsigned int subface_index = fe_eval.get_subface_index();
2911 VectorizedArrayType>
::
2914 fe_eval.get_shape_info().data.front(),
2916 fe_eval.begin_values(),
2917 fe_eval.begin_gradients(),
2918 fe_eval.begin_hessians(),
2926 VectorizedArrayType>
::
2929 fe_eval.get_shape_info().data.front(),
2931 fe_eval.begin_values(),
2932 fe_eval.begin_gradients(),
2933 fe_eval.begin_hessians(),
2940 if (fe_eval.get_dof_access_index() ==
2942 fe_eval.is_interior_face() ==
false)
2944 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
2949 if (fe_eval.get_cell_ids()[v] ==
2953 if (fe_eval.get_face_orientation(v) != 0)
2959 &fe_eval.get_shape_info().face_orientations_quad(
2960 fe_eval.get_face_orientation(v), 0),
2964 fe_eval.begin_values(),
2965 fe_eval.begin_gradients(),
2966 fe_eval.begin_hessians());
2969 else if (fe_eval.get_face_orientation() != 0)
2974 &fe_eval.get_shape_info().face_orientations_quad(
2975 fe_eval.get_face_orientation(), 0),
2979 fe_eval.begin_values(),
2980 fe_eval.begin_gradients(),
2981 fe_eval.begin_hessians());
2987 template <
typename Number3>
2991 const Number2 *vector_ptr,
2994 const unsigned int fe_degree = shape_info.
data.front().fe_degree;
2995 if (fe_degree < 1 || !shape_info.
data.front().nodal_at_cell_boundaries ||
2998 shape_info.
data.front().element_type !=
3001 vector_ptr ==
nullptr ||
3002 shape_info.
data.front().element_type >
3012 template <
int fe_degree>
3022 template <
typename T0,
typename T1,
typename T2>
3035 template <
typename T1,
typename T2>
3042 template <
typename T0,
typename T1,
typename T2,
typename T3>
3057 template <
typename T0,
typename T1,
typename T2>
3064 template <
typename T0,
typename T1,
typename T2>
3077 template <
typename T1,
typename T2>
3089 template <
int dim,
typename Number2,
typename VectorizedArrayType>
3092 using Number =
typename VectorizedArrayType::value_type;
3094 template <
int fe_degree,
int n_q_po
ints_1d>
3096 run(
const unsigned int n_components,
3103 Assert(fe_eval.get_shape_info().element_type <=
3107 const unsigned int dofs_per_face =
Utilities::pow(fe_degree + 1, dim - 1);
3109 VectorizedArrayType *
temp = fe_eval.get_scratch_data().begin();
3110 VectorizedArrayType *scratch_data =
3111 temp + 3 * n_components * dofs_per_face;
3113 const unsigned int subface_index = fe_eval.get_subface_index();
3118 if (fe_eval.get_dof_access_index() ==
3120 fe_eval.is_interior_face() ==
false)
3121 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
3129 if (fe_eval.get_face_orientation(v) != 0)
3135 &fe_eval.get_shape_info().face_orientations_quad(
3136 fe_eval.get_face_orientation(v), 0),
3140 fe_eval.begin_values(),
3141 fe_eval.begin_gradients(),
3142 fe_eval.begin_hessians());
3144 else if (fe_eval.get_face_orientation() != 0)
3149 &fe_eval.get_shape_info().face_orientations_quad(
3150 fe_eval.get_face_orientation(), 0),
3154 fe_eval.begin_values(),
3155 fe_eval.begin_gradients(),
3156 fe_eval.begin_hessians());
3159 if (fe_degree > -1 && fe_eval.get_subface_index() >=
3165 VectorizedArrayType>::
3166 integrate_in_face(n_components,
3168 fe_eval.get_shape_info().data.front(),
3170 fe_eval.begin_values(),
3171 fe_eval.begin_gradients(),
3172 fe_eval.begin_hessians(),
3180 VectorizedArrayType>
::
3183 fe_eval.get_shape_info().data.front(),
3185 fe_eval.begin_values(),
3186 fe_eval.begin_gradients(),
3187 fe_eval.begin_hessians(),
3193 if (fe_eval.get_dof_access_index() ==
3195 fe_eval.is_interior_face() ==
false)
3206 template <
int fe_degree>
3216 template <
typename T0,
typename T1,
typename T2,
typename T3,
typename T4>
3231 template <
typename T0,
typename T1>
3239 template <
typename T0,
typename T1,
typename T2,
typename T3>
3256 template <
typename T0,
typename T1,
typename T2>
3264 template <
typename T0,
typename T1,
typename T2>
3279 template <
typename T0,
typename T1>
Number shape_info_number_type
#define DEAL_II_ALWAYS_INLINE
#define DEAL_II_OPENMP_SIMD_PRAGMA
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
#define AssertThrow(cond, exc)
@ tensor_symmetric_no_collocation
@ tensor_symmetric_hermite
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
EvaluationFlags
The EvaluationFlags enum.
constexpr T fixed_power(const T t)
constexpr T pow(const T base, const int iexp)
void do_vectorized_add(const VectorizedArrayType src, Number2 *dst_ptr)
constexpr bool use_collocation_evaluation(const unsigned int fe_degree, const unsigned int n_q_points_1d)
void adjust_for_face_orientation_per_lane(const unsigned int dim, const unsigned int n_components, const unsigned int v, const EvaluationFlags::EvaluationFlags flag, const unsigned int *orientation, const bool integrate, const std::size_t n_q_points, Number *tmp_values, VectorizedArrayType *values_quad, VectorizedArrayType *gradients_quad=nullptr, VectorizedArrayType *hessians_quad=nullptr)
void fe_face_evaluation_process_and_io(Processor &proc, const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, typename Processor::Number2_ *global_vector_ptr, const std::vector< ArrayView< const typename Processor::Number2_ > > *sm_ptr, const EvaluationData &fe_eval, typename Processor::VectorizedArrayType_ *temp1)
void do_vectorized_scatter_add(const VectorizedArrayType src, const unsigned int *indices, Number2 *dst_ptr)
void do_vectorized_gather(const Number2 *src_ptr, const unsigned int *indices, VectorizedArrayType &dst)
void do_vectorized_read(const Number2 *src_ptr, VectorizedArrayType &dst)
std::enable_if_t< contract_onto_face, void > interpolate_to_face(const Number2 *shape_values, const std::array< int, 2 > &n_blocks, const std::array< int, 2 > &steps, const Number *input, Number *DEAL_II_RESTRICT output, const int n_rows_runtime=0, const int stride_runtime=1)
void adjust_for_face_orientation(const unsigned int dim, const unsigned int n_components, const EvaluationFlags::EvaluationFlags flag, const unsigned int *orientation, const bool integrate, const std::size_t n_q_points, Number *tmp_values, Number *values_quad, Number *gradients_quad, Number *hessians_quad)
static const unsigned int invalid_unsigned_int
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval, const bool sum_into_values)
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, FEEvaluationData< dim, Number, true > &fe_eval)
static void evaluate_in_face(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, FEEvaluationData< dim, Number, true > &fe_eval, Number *temp, Number *scratch_data)
static void project_to_face(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval, const bool use_vectorization, Number *temp, Number *scratch_data)
static bool evaluate_tensor_none(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval)
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval)
static void adjust_quadrature_for_face_orientation(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, FEEvaluationData< dim, Number, true > &fe_eval, const bool use_vectorization, Number *temp)
static bool evaluate_tensor(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval)
void hermite_grad(T0 &temp_1, T0 &temp_2, const T1 &src_ptr_1, const T1 &src_ptr_2, const T2 &grad_weight)
static const int fe_degree_
static const bool do_integrate
void value_vectorized(T1 &temp, const T2 src_ptr)
void value(T1 &temp, const T2 &src_ptr)
VectorizedArrayType VectorizedArrayType_
void value_vectorized_indexed(T0 &temp, const T1 src_ptr, const T2 &indices)
void hermite_grad_vectorized(T0 &temp_1, T0 &temp_2, const T1 src_ptr_1, const T1 src_ptr_2, const T2 &grad_weight)
void hermite_grad_vectorized_indexed(T0 &temp_1, T0 &temp_2, const T1 src_ptr_1, const T1 src_ptr_2, const T2 &grad_weight, const T3 &indices_1, const T3 &indices_2)
typename VectorizedArrayType::value_type Number
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number2 *src_ptr, const std::vector< ArrayView< const Number2 > > *sm_ptr, FEEvaluationData< dim, VectorizedArrayType, true > &fe_eval)
static bool supports(const EvaluationFlags::EvaluationFlags evaluation_flag, const MatrixFreeFunctions::ShapeInfo< Number3 > &shape_info, const Number2 *vector_ptr, MatrixFreeFunctions::DoFInfo::IndexStorageVariants storage)
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, FEEvaluationData< dim, Number, true > &fe_eval)
static const int fe_degree_
void value(const T0 &temp, T1 &dst_ptr)
void hermite_grad_vectorized(const T0 &temp_1, const T1 &temp_2, T2 dst_ptr_1, T3 dst_ptr_2, const T4 &grad_weight)
void hermite_grad_vectorized_indexed(const T0 &temp_1, const T0 &temp_2, T1 dst_ptr_1, T1 dst_ptr_2, const T2 &grad_weight, const T3 &indices_1, const T3 &indices_2)
static const bool do_integrate
VectorizedArrayType VectorizedArrayType_
void hermite_grad(const T0 &temp_1, const T0 &temp_2, T1 &dst_ptr_1, T1 &dst_ptr_2, const T2 &grad_weight)
void value_vectorized(const T0 &temp, T1 dst_ptr)
void value_vectorized_indexed(const T0 &temp, T1 dst_ptr, const T2 &indices)
typename VectorizedArrayType::value_type Number
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number2 *dst_ptr, const std::vector< ArrayView< const Number2 > > *sm_ptr, FEEvaluationData< dim, VectorizedArrayType, true > &fe_eval)
static bool integrate_tensor_none(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval, const bool sum_into_values)
static bool integrate_tensor(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval, const bool sum_into_values)
static void collect_from_face(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval, const bool use_vectorization, const Number *temp, Number *scratch_data, const bool sum_into_values)
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval, const bool sum_into_values)
static void integrate_in_face(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, FEEvaluationData< dim, Number, true > &fe_eval, Number *temp, Number *scratch_data)
static void adjust_quadrature_for_face_orientation(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, FEEvaluationData< dim, Number, true > &fe_eval, const bool use_vectorization, Number *temp)
static bool run(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number *values_dofs, FEEvaluationData< dim, Number, true > &fe_eval)
static void evaluate_or_integrate_in_face(const EvaluationFlags::EvaluationFlags evaluation_flag, const std::vector< MatrixFreeFunctions::UnivariateShapeData< Number2 > > &shape_data, Number *values_dofs_in, Number *values, Number *gradients, Number *scratch_data, const unsigned int subface_index, const unsigned int face_direction)
typename FEEvaluationData< dim, Number, true >::shape_info_number_type Number2
EvaluatorTensorProduct< symmetric_evaluate ? evaluate_evenodd :evaluate_general, dim - 1, fe_degree+1, n_q_points_1d, Number, Number2 > Eval
static void evaluate_in_face(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const MatrixFreeFunctions::UnivariateShapeData< Number2 > &data, Number *values_dofs, Number *values_quad, Number *gradients_quad, Number *hessians_quad, Number *scratch_data, const unsigned int subface_index)
static Eval create_evaluator_tensor_product(const MatrixFreeFunctions::UnivariateShapeData< Number2 > &data, const unsigned int subface_index, const unsigned int direction)
static void integrate_in_face(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, const MatrixFreeFunctions::UnivariateShapeData< Number2 > &data, Number *values_dofs, Number *values_quad, Number *gradients_quad, Number *hessians_quad, Number *scratch_data, const unsigned int subface_index)
typename FEEvaluationData< dim, Number, true >::shape_info_number_type Number2
typename FEEvaluationData< dim, Number, true >::shape_info_number_type Number2
static void interpolate_quadrature(const unsigned int n_components, const EvaluationFlags::EvaluationFlags flags, const MatrixFreeFunctions::ShapeInfo< Number2 > &shape_info, const Number *input, Number *output, const unsigned int face_no)
static void interpolate_generic(const unsigned int n_components, const Number *input, Number *output, const EvaluationFlags::EvaluationFlags flag, const unsigned int face_no, const unsigned int n_points_1d, const std::array< AlignedVector< Number2 >, 2 > &shape_data, const unsigned int dofs_per_component_on_cell, const unsigned int dofs_per_component_on_face)
static void interpolate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags flags, const MatrixFreeFunctions::ShapeInfo< Number2 > &shape_info, const Number *input, Number *output, const unsigned int face_no)
static void interpolate_raviart_thomas(const unsigned int n_components, const Number *input, Number *output, const EvaluationFlags::EvaluationFlags flag, const unsigned int face_no, const MatrixFreeFunctions::ShapeInfo< Number2 > &shape_info)
unsigned int n_q_points_face
unsigned int dofs_per_component_on_cell
std::vector< UnivariateShapeData< Number > > data
unsigned int dofs_per_component_on_face
AlignedVector< Number > shape_values
AlignedVector< Number > shape_values_eo
AlignedVector< Number > shape_hessians_eo
AlignedVector< Number > shape_gradients_collocation_eo
unsigned int n_q_points_1d
AlignedVector< Number > shape_gradients_eo
AlignedVector< Number > shape_hessians
AlignedVector< Number > shape_gradients
std::array< AlignedVector< Number >, 2 > hessians_within_subface
std::array< AlignedVector< Number >, 2 > values_within_subface
std::array< AlignedVector< Number >, 2 > gradients_within_subface