Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Member Functions | Static Public Member Functions | Static Public Attributes | Private Attributes | List of all members
internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 > Struct Template Reference

#include <deal.II/matrix_free/tensor_product_kernels.h>

Public Member Functions

 EvaluatorTensorProduct ()
 
 EvaluatorTensorProduct (const AlignedVector< Number2 > &shape_values, const AlignedVector< Number2 > &shape_gradients, const AlignedVector< Number2 > &shape_hessians, const unsigned int dummy1=0, const unsigned int dummy2=0)
 
 EvaluatorTensorProduct (const Number2 *shape_values, const Number2 *shape_gradients, const Number2 *shape_hessians, const unsigned int dummy1=0, const unsigned int dummy2=0)
 
template<int direction, bool contract_over_rows, bool add>
void values (const Number in[], Number out[]) const
 
template<int direction, bool contract_over_rows, bool add>
void gradients (const Number in[], Number out[]) const
 
template<int direction, bool contract_over_rows, bool add>
void hessians (const Number in[], Number out[]) const
 
template<int direction, bool contract_over_rows, bool add>
void values_one_line (const Number in[], Number out[]) const
 
template<int direction, bool contract_over_rows, bool add>
void gradients_one_line (const Number in[], Number out[]) const
 
template<int direction, bool contract_over_rows, bool add>
void hessians_one_line (const Number in[], Number out[]) const
 
template<int face_direction, bool contract_onto_face, bool add, int max_derivative>
void apply_face (const Number *DEAL_II_RESTRICT in, Number *DEAL_II_RESTRICT out) const
 

Static Public Member Functions

template<int direction, bool contract_over_rows, bool add, bool one_line = false>
static void apply (const Number2 *DEAL_II_RESTRICT shape_data, const Number *in, Number *out)
 

Static Public Attributes

static constexpr unsigned int n_rows_of_product
 
static constexpr unsigned int n_columns_of_product
 

Private Attributes

const Number2 * shape_values
 
const Number2 * shape_gradients
 
const Number2 * shape_hessians
 

Detailed Description

template<int dim, int n_rows, int n_columns, typename Number, typename Number2>
struct internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >

Internal evaluator for shape function in arbitrary dimension using the tensor product form of the basis functions.

Template Parameters
dimSpace dimension in which this class is applied
n_rowsNumber of rows in the transformation matrix, which corresponds to the number of 1d shape functions in the usual tensor contraction setting
n_columnsNumber of columns in the transformation matrix, which corresponds to the number of 1d shape functions in the usual tensor contraction setting
NumberAbstract number type for input and output arrays
Number2Abstract number type for coefficient arrays (defaults to same type as the input/output arrays); must implement operator* with Number and produce Number as an output to be a valid type

Definition at line 182 of file tensor_product_kernels.h.

Constructor & Destructor Documentation

◆ EvaluatorTensorProduct() [1/3]

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::EvaluatorTensorProduct ( )
inline

Empty constructor. Does nothing. Be careful when using 'values' and related methods because they need to be filled with the other pointer

Definition at line 198 of file tensor_product_kernels.h.

◆ EvaluatorTensorProduct() [2/3]

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::EvaluatorTensorProduct ( const AlignedVector< Number2 > &  shape_values,
const AlignedVector< Number2 > &  shape_gradients,
const AlignedVector< Number2 > &  shape_hessians,
const unsigned int  dummy1 = 0,
const unsigned int  dummy2 = 0 
)
inline

Constructor, taking the data from ShapeInfo

Definition at line 207 of file tensor_product_kernels.h.

◆ EvaluatorTensorProduct() [3/3]

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::EvaluatorTensorProduct ( const Number2 *  shape_values,
const Number2 *  shape_gradients,
const Number2 *  shape_hessians,
const unsigned int  dummy1 = 0,
const unsigned int  dummy2 = 0 
)
inline

Constructor, taking the data from ShapeInfo via raw pointers

Definition at line 237 of file tensor_product_kernels.h.

Member Function Documentation

◆ values()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int direction, bool contract_over_rows, bool add>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::values ( const Number  in[],
Number  out[] 
) const
inline

Definition at line 252 of file tensor_product_kernels.h.

◆ gradients()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int direction, bool contract_over_rows, bool add>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::gradients ( const Number  in[],
Number  out[] 
) const
inline

Definition at line 259 of file tensor_product_kernels.h.

◆ hessians()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int direction, bool contract_over_rows, bool add>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::hessians ( const Number  in[],
Number  out[] 
) const
inline

Definition at line 266 of file tensor_product_kernels.h.

◆ values_one_line()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int direction, bool contract_over_rows, bool add>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::values_one_line ( const Number  in[],
Number  out[] 
) const
inline

Definition at line 273 of file tensor_product_kernels.h.

◆ gradients_one_line()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int direction, bool contract_over_rows, bool add>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::gradients_one_line ( const Number  in[],
Number  out[] 
) const
inline

Definition at line 281 of file tensor_product_kernels.h.

◆ hessians_one_line()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int direction, bool contract_over_rows, bool add>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::hessians_one_line ( const Number  in[],
Number  out[] 
) const
inline

Definition at line 289 of file tensor_product_kernels.h.

◆ apply()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int direction, bool contract_over_rows, bool add, bool one_line>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::apply ( const Number2 *DEAL_II_RESTRICT  shape_data,
const Number *  in,
Number *  out 
)
inlinestatic

This function applies the tensor product kernel, corresponding to a multiplication of 1d stripes, along the given direction of the tensor data in the input array. This function allows the in and out arrays to alias for the case n_rows == n_columns, i.e., it is safe to perform the contraction in place where in and out point to the same address. For the case n_rows != n_columns, the output is in general not correct.

Template Parameters
directionDirection that is evaluated
contract_over_rowsIf true, the tensor contraction sums over the rows in the given shape_data array, otherwise it sums over the columns
addIf true, the result is added to the output vector, else the computed values overwrite the content in the output
one_lineIf true, the kernel is only applied along a single 1d stripe within a dim-dimensional tensor, not the full n_rows^dim points as in the false case.
Parameters
shape_dataTransformation matrix with n_rows rows and n_columns columns, stored in row-major format
inPointer to the start of the input data vector
outPointer to the start of the output data vector

Definition at line 385 of file tensor_product_kernels.h.

◆ apply_face()

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
template<int face_direction, bool contract_onto_face, bool add, int max_derivative>
void internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::apply_face ( const Number *DEAL_II_RESTRICT  in,
Number *DEAL_II_RESTRICT  out 
) const
inline

This function applies the tensor product operation to produce face values from cell values. As opposed to the apply method, this method assumes that the directions orthogonal to the face have n_rows degrees of freedom per direction and not n_columns for those directions lower than the one currently applied. In other words, apply_face() must be called before calling any interpolation within the face.

Template Parameters
face_directionDirection of the normal vector (0=x, 1=y, etc)
contract_onto_faceIf true, the input vector is of size n_rows^dim and interpolation into n_rows^(dim-1) points is performed. This is a typical scenario in FEFaceEvaluation::evaluate() calls. If false, data from n_rows^(dim-1) points is expanded into the n_rows^dim points of the higher- dimensional data array. Derivatives in the case contract_onto_face==false are summed together
addIf true, the result is added to the output vector, else the computed values overwrite the content in the output
max_derivativeSets the number of derivatives that should be computed. 0 means only values, 1 means values and first derivatives, 2 second derivates. Note that all the derivatives access the data in shape_values passed to the constructor of the class
Parameters
inaddress of the input data vector
outaddress of the output data vector

Definition at line 468 of file tensor_product_kernels.h.

Member Data Documentation

◆ n_rows_of_product

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
constexpr unsigned int internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::n_rows_of_product
staticconstexpr
Initial value:
=
Utilities::pow(n_rows, dim)
constexpr T pow(const T base, const int iexp)
Definition utilities.h:447

Definition at line 189 of file tensor_product_kernels.h.

◆ n_columns_of_product

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
constexpr unsigned int internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::n_columns_of_product
staticconstexpr
Initial value:
=
Utilities::pow(n_columns, dim)

Definition at line 191 of file tensor_product_kernels.h.

◆ shape_values

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
const Number2* internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::shape_values
private

Definition at line 366 of file tensor_product_kernels.h.

◆ shape_gradients

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
const Number2* internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::shape_gradients
private

Definition at line 367 of file tensor_product_kernels.h.

◆ shape_hessians

template<int dim, int n_rows, int n_columns, typename Number , typename Number2 >
const Number2* internal::EvaluatorTensorProduct< evaluate_general, dim, n_rows, n_columns, Number, Number2 >::shape_hessians
private

Definition at line 368 of file tensor_product_kernels.h.


The documentation for this struct was generated from the following file: