16#ifndef dealii_matrix_free_tensor_product_point_kernels_h
17#define dealii_matrix_free_tensor_product_point_kernels_h
40 template <
typename Number,
typename Number2>
46 template <
int dim,
typename Number,
typename Number2>
58 template <
int dim,
typename Number>
64 const unsigned int derivative = 1)
69 std::array<Number, dim> point;
70 for (
unsigned int d = 0; d < dim; ++d)
73 poly[i].values_of_array(point, derivative, shapes[i].
data());
81 template <
typename Number>
100 template <
int dim,
typename Number>
114 std::array<Number, 2 * dim> point,
result;
115 for (
unsigned int d = 0; d < dim; ++d)
117 for (
unsigned int d = 0; d < dim; ++d)
118 point[dim + d] =
p1[d];
121 poly[i].values_of_array(point, 0, &
result);
122 for (
unsigned int j = 0;
j < 2; ++
j)
123 for (
unsigned int d = 0; d < dim; ++d)
144 std::array<typename ProductTypeNoPoint<Number, Number2>::type,
147 const std::vector<unsigned int> &renumber,
152 static_assert(0 <= dim && dim <= 3,
"Only dim=0,1,2,3 implemented");
153 static_assert(1 <= n_values && n_values <= 2,
154 "Only n_values=1,2 implemented");
163 values + stride * (length > 0 ?
168 std::array<Number3, 2 + n_values>
result = {};
181 shapes[
i0][1][0] * values[renumber[i] * stride];
184 shapes[
i0][0][0] * values[renumber[i] * stride];
187 shapes[
i0][0][0] *
values_2[renumber[i] * stride];
236 inline std::array<typename ProductTypeNoPoint<Number, Number2>::type,
241 const Number *values,
242 const std::vector<unsigned int> &renumber = {})
244 static_assert(0 <= dim && dim <= 3,
"Only dim=0,1,2,3 implemented");
245 static_assert(1 <= n_values && n_values <= 2,
246 "Only n_values=1,2 implemented");
250 std::array<Number3, dim + n_values>
result = {};
370 inline std::array<typename ProductTypeNoPoint<Number, Number2>::type,
373 const Number *values,
376 static_assert(0 <= dim && dim <= 3,
"Only dim=0,1,2,3 implemented");
377 static_assert(1 <= n_values && n_values <= 2,
378 "Only n_values=1,2 implemented");
386 std::array<Number3, dim + n_values>
result;
392 result[1] = values[1 * stride];
402 p[0] * (values[3 * stride] - values[2 * stride]);
422 p[0] * (values[5 * stride] - values[4 * stride]);
425 p[0] * (values[7 * stride] - values[6 * stride]);
495 template <
int dim,
typename Number,
typename Number2>
504 const std::vector<unsigned int> &renumber = {})
508 std::array<Number3, dim + 1>
result;
517 std::array<::ndarray<Number2, 2, dim>, 200> shapes;
522 shapes.data(), poly.size(),
values.data(), renumber);
524 return std::make_pair(
result[dim],
543 const std::vector<unsigned int> &renumber,
560 value += shapes[
i0][0][0] * values[renumber[i] * stride];
563 value += shapes[
i0][0][0] * values[i * stride];
584 const Number *values,
585 const std::vector<unsigned int> &renumber = {})
587 static_assert(dim >= 0 && dim <= 3,
"Only dim=0,1,2,3 implemented");
669 template <
int dim,
typename Number,
typename Number2,
int str
ide = 1>
674 static_assert(dim >= 0 && dim <= 3,
"Only dim=0,1,2,3 implemented");
685 return Number3(values[0]) + p[0] *
Number3(values[stride] - values[0]);
718 template <
int dim,
typename Number,
typename Number2>
725 const std::vector<unsigned int> &renumber = {})
735 std::array<::ndarray<Number2, 2, dim>, 200> shapes;
737 std::array<Number2, dim>
point;
738 for (
unsigned int d = 0;
d < dim; ++
d)
741 poly[i].values_of_array(point, 0, shapes[i].
data());
754 template <
int derivative_order,
typename Number,
typename Number2>
760 const std::vector<unsigned int> &renumber = {})
766 Assert(renumber.empty() || renumber.size() == values.size(),
769 std::array<Number2, derivative_order + 1> shapes;
771 if (renumber.empty())
792 template <
int derivative_order,
typename Number,
typename Number2>
800 const std::vector<unsigned int> &renumber = {})
803 constexpr int dim = 2;
807 Assert(renumber.empty() || renumber.size() == values.size(),
813 std::array<Number2, dim> point;
814 for (
unsigned int d = 0; d < dim; ++d)
823 if (renumber.empty())
830 result_x[d] += shapes[
i0][d][0] * values[renumber[i]];
844 template <
int derivative_order,
typename Number,
typename Number2>
852 const std::vector<unsigned int> &renumber = {})
855 constexpr int dim = 3;
861 Assert(renumber.empty() || renumber.size() == values.size(),
867 std::array<Number2, dim> point;
868 for (
unsigned int d = 0; d < dim; ++d)
881 if (renumber.empty())
888 result_x[d] += shapes[
i0][d][0] * values[renumber[i]];
910 template <
int dim,
typename Number,
typename Number2>
916 const std::vector<unsigned int> &renumber = {})
918 static_assert(dim >= 1 && dim <= 3,
"Only dim=1,2,3 implemented");
931 for (
unsigned int d = 0, c = 0; d < 2; ++d)
932 for (
unsigned int e = d; e < 2; ++e, ++c)
936 for (
unsigned int d = 0; d < 2; ++d)
968 static_assert(0 <= dim && dim <= 3,
"Only dim=0,1,2,3 implemented");
969 static_assert(1 <= n_values && n_values <= 2,
970 "Only n_values=1,2 implemented");
976 constexpr unsigned int array_size = length > 0 ? length : 1;
984 for (
int i1 = 0;
i1 < (dim > 1 ? length : 1); ++
i1)
1002 for (
int i0 = 0;
i0 < length; ++
i0)
1018 i += (dim > 1 ? length * length : length);
1077 const Number2 *
value,
1081 static_assert(0 <= dim && dim <= 3,
"Only dim=0,1,2,3 implemented");
1082 static_assert(1 <= n_values && n_values <= 2,
1083 "Only n_values=1,2 implemented");
1090 values[0] +=
value[0];
1092 values[0] =
value[0];
1096 values[1] +=
value[1];
1098 values[1] =
value[1];
1113 (dim > 1 ?
gradient[1] : Number2());
1159 const Number2 *
value,
1164 static_assert(0 <= dim && dim <= 3,
"Only dim=0,1,2,3 implemented");
1165 static_assert(1 <= n_values && n_values <= 2,
1166 "Only n_values=1,2 implemented");
1173 values[0] +=
value[0];
1175 values[0] =
value[0];
1179 values[1] +=
value[1];
1181 values[1] =
value[1];
1189 values[0] +=
value[0] - difference;
1190 values[1] += difference;
1194 values[0] =
value[0] - difference;
1195 values[1] = difference;
1325 const Number2 *
value,
1372 template <
int dim,
int length,
typename Number2,
typename Number,
bool add>
1387 constexpr unsigned int array_size = length > 0 ? length : 1;
1391 for (
int i1 = 0;
i1 < (dim > 1 ? length : 1); ++
i1)
1397 for (
int i0 = 0;
i0 < length; ++
i0)
1405 i += (dim > 1 ? length * length : length);
1432 template <
int dim,
typename Number,
typename Number2,
bool add>
1437 const Number2 &
value,
1440 static_assert(dim >= 0 && dim <= 3,
"Only dim=0,1,2,3 implemented");
1489 template <
int dim,
typename Number,
typename Number2,
bool add>
1495 static_assert(dim >= 0 && dim <= 3,
"Only dim=0,1,2,3 implemented");
1506 const auto x0 = Number(1.) - p[0],
x1 = p[0];
1521 const auto x0 = Number(1.) - p[0],
x1 = p[0],
y0 = Number(1.) - p[1],
1544 const auto x0 = Number(1.) - p[0],
x1 = p[0],
y0 = Number(1.) - p[1],
1545 y1 = p[1],
z0 = Number(1.) - p[2],
z1 = p[2];
1587 template <
bool is_linear,
int dim,
typename Number,
typename Number2>
1591 const Number2 &
value,
#define DEAL_II_ALWAYS_INLINE
#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 & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define DEAL_II_ASSERT_UNREACHABLE()
std::vector< index_type > data
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
constexpr T pow(const T base, const int iexp)
std::array< typename ProductTypeNoPoint< Number, Number2 >::type, dim+n_values > evaluate_tensor_product_value_and_gradient_shapes(const ::ndarray< Number2, 2, dim > *shapes, const int n_shapes, const Number *values, const std::vector< unsigned int > &renumber={})
std::array< typename ProductTypeNoPoint< Number, Number2 >::type, 2+n_values > do_interpolate_xy(const Number *values, const std::vector< unsigned int > &renumber, const ::ndarray< Number2, 2, dim > *shapes, const int n_shapes_runtime, int &i)
void integrate_add_tensor_product_value_linear(const Number2 &value, Number2 *values, const Point< dim, Number > &p)
void integrate_tensor_product_value_and_gradient(const ::ndarray< Number, 2, dim > *shapes, const unsigned int n_shapes, const Number2 *value, const Tensor< 1, dim, Number2 > &gradient, Number2 *values, const Point< dim, Number > &p, const bool do_add)
void do_apply_test_functions_xy_value(Number2 *values, const ::ndarray< Number, 2, dim > *shapes, const Number2 &test_value, const int n_shapes_runtime, int &i)
void compute_values_of_array_in_pairs(::ndarray< Number, 2, dim > *shapes, const std::vector< Polynomials::Polynomial< double > > &poly, const Point< dim, Number > &p0, const Point< dim, Number > &p1)
void integrate_add_tensor_product_value_and_gradient_linear(const Number2 *value, const Tensor< 1, dim, Number2 > &gradient, Number2 *values, const Point< dim, Number > &p)
void integrate_add_tensor_product_value_shapes(const ::ndarray< Number, 2, dim > *shapes, const int n_shapes, const Number2 &value, Number2 *values)
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={})
ProductTypeNoPoint< Number, Number2 >::type evaluate_tensor_product_value_shapes(const ::ndarray< Number2, 2, dim > *shapes, const int n_shapes, const Number *values, const std::vector< unsigned int > &renumber={})
void compute_values_of_array(::ndarray< Number, 2, dim > *shapes, const std::vector< Polynomials::Polynomial< double > > &poly, const Point< dim, Number > &p, const unsigned int derivative=1)
void integrate_add_tensor_product_value_and_gradient_shapes(const ::ndarray< Number, 2, dim > *shapes, const int n_shapes, const Number2 *value, const Tensor< 1, dim, Number2 > &gradient, Number2 *values)
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={})
void do_apply_test_functions_xy(Number2 *values, const ::ndarray< Number, 2, dim > *shapes, const std::array< Number2, 2+n_values > &test_grads_value, const int n_shapes_runtime, int &i)
ProductTypeNoPoint< Number, Number2 >::type evaluate_tensor_product_value_linear(const Number *values, const Point< dim, Number2 > &p)
std::array< typename ProductTypeNoPoint< Number, Number2 >::type, dim+n_values > evaluate_tensor_product_value_and_gradient_linear(const Number *values, const Point< dim, Number2 > &p)
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={})
Tensor< 1, 1, typename ProductTypeNoPoint< Number, Number2 >::type > evaluate_tensor_product_higher_derivatives(const std::vector< Polynomials::Polynomial< double > > &poly, const ArrayView< const Number > &values, const Point< 1, Number2 > &p, const std::vector< unsigned int > &renumber={})
void integrate_tensor_product_value(const ::ndarray< Number, 2, dim > *shapes, const unsigned int n_shapes, const Number2 &value, Number2 *values, const Point< dim, Number > &p, const bool do_add)
ProductTypeNoPoint< Number, Number2 >::type do_interpolate_xy_value(const Number *values, const std::vector< unsigned int > &renumber, const ::ndarray< Number2, 2, dim > *shapes, const int n_shapes_runtime, int &i)
typename internal::ProductTypeImpl< std::decay_t< T >, std::decay_t< U > >::type type
typename ProductType< Tensor< 1, dim, Number >, Number2 >::type type
typename ProductType< Number, Number2 >::type type