/**
- * This struct performs the evaluation of function values, gradients and
- * Hessians for tensor-product finite elements. The operation is used for
- * both the symmetric and non-symmetric case, which use different apply
- * functions 'values', 'gradients' in the individual coordinate
- * directions. The apply functions for values are provided through one of
- * the template classes EvaluatorTensorProduct which in turn are selected
- * from the MatrixFreeFunctions::ElementType template argument.
+ * This struct performs the evaluation of function values and gradients for
+ * tensor-product finite elements. The operation is used for both the
+ * symmetric and non-symmetric case, which use different apply functions
+ * 'values', 'gradients' in the individual coordinate directions. The apply
+ * functions for values are provided through one of the template classes
+ * EvaluatorTensorProduct which in turn are selected from the
+ * MatrixFreeFunctions::ElementType template argument.
*
* There are two specialized implementation classes
* FEEvaluationImplCollocation (for Gauss-Lobatto elements where the nodal
* identity) and FEEvaluationImplTransformToCollocation (which can be
* transformed to a collocation space and can then use the identity in these
* spaces), which both allow for shorter code.
+ *
+ * @note Hessians of the solution are handled in the general
+ * FEEvaluationImplSelector struct below, because they can be implemented
+ * with the only two code paths for all supported cases, including the
+ * specialized cases below.
*/
template <MatrixFreeFunctions::ElementType type,
int dim,
* @param values_out The array of size basis_size_2^dim where the results
* of the transformation are stored. It may alias with
* the values_in array.
- * @param basis_size_1_variable In case the template argument basis_size_1
- * is zero, the size of the first basis can alternatively be passed in as a
- * run time argument. The template argument takes precedence in case it is
- * nonzero for efficiency reasons.
- * @param basis_size_2_variable In case the template argument basis_size_1
- * is zero, the size of the second basis can alternatively be passed in as a
- * run time argument.
+ * @param basis_size_1_variable In case the template argument
+ * @p basis_size_1 is zero, the size of the first basis can alternatively
+ * be passed in as a run time argument. The template argument takes
+ * precedence in case it is nonzero for efficiency reasons.
+ * @param basis_size_2_variable In case the template argument
+ * @p basis_size_1 is zero, the size of the second basis can alternatively
+ * be passed in as a run time argument.
*/
template <typename Number, typename Number2>
#ifndef DEBUG
// basis_size_1==basis_size_2. The latter optimization increases
// optimization possibilities for the compiler but does only work for
// aliased pointers if the sizes are equal.
- constexpr int next_dim =
- (dim > 2 ||
- ((basis_size_1 == 0 || basis_size_2 > basis_size_1) && dim > 1)) ?
- dim - 1 :
- dim;
+ constexpr int next_dim = (dim == 1 || (dim == 2 && basis_size_1 > 0 &&
+ basis_size_1 == basis_size_2)) ?
+ dim :
+ dim - 1;
EvaluatorTensorProduct<variant,
dim,
* @param values_out The array of size basis_size_1^dim where the results
* of the transformation are stored. It may alias with
* the @p values_in array.
- * @param basis_size_1_variable In case the template argument basis_size_1
- * is zero, the size of the first basis can alternatively be passed in as a
- * run time argument. The template argument takes precedence in case it is
- * nonzero for efficiency reasons.
- * @param basis_size_2_variable In case the template argument basis_size_1
- * is zero, the size of the second basis can alternatively be passed in as a
- * run time argument.
+ * @param basis_size_1_variable In case the template argument
+ * @p basis_size_1 is zero, the size of the first basis can alternatively
+ * be passed in as a run time argument. The template argument takes
+ * precedence in case it is nonzero for efficiency reasons.
+ * @param basis_size_2_variable In case the template argument
+ * @p basis_size_1 is zero, the size of the second basis can alternatively
+ * be passed in as a run time argument.
*/
template <typename Number, typename Number2>
#ifndef DEBUG
/**
- * This struct performs the evaluation of function values, gradients and
- * Hessians for tensor-product finite elements. This a specialization for
- * elements where the nodal points coincide with the quadrature points like
- * FE_Q shape functions on Gauss-Lobatto elements integrated with
- * Gauss-Lobatto quadrature. The assumption of this class is that the shape
- * 'values' operation is identity, which allows us to write shorter code.
+ * This struct performs the evaluation of function values and gradients for
+ * tensor-product finite elements. This is a specialization for elements
+ * where the nodal points coincide with the quadrature points like FE_Q
+ * shape functions on Gauss-Lobatto elements integrated with Gauss-Lobatto
+ * quadrature. The assumption of this class is that the shape 'values'
+ * operation is identity, which allows us to write shorter code.
*
* In literature, this form of evaluation is often called spectral
* evaluation, spectral collocation or simply collocation, meaning the same
/**
- * This struct performs the evaluation of function values, gradients and
- * Hessians for tensor-product finite elements. This a specialization for
- * symmetric basis functions about the mid point 0.5 of the unit interval
- * with the same number of quadrature points as degrees of freedom. In that
- * case, we can first transform the basis to one that has the nodal points
- * in the quadrature points (i.e., the collocation space) and then perform
- * the evaluation of the first and second derivatives in this transformed
- * space, using the identity operation for the shape values.
+ * This struct performs the evaluation of function values and gradients for
+ * tensor-product finite elements. This is a specialization for symmetric
+ * basis functions about the mid point 0.5 of the unit interval with the
+ * same number of quadrature points as degrees of freedom. In that case, we
+ * can first transform the basis to one that has the nodal points in the
+ * quadrature points (i.e., the collocation space) and then perform the
+ * evaluation of the first and second derivatives in this transformed space,
+ * using the identity operation for the shape values.
*/
template <int dim, int fe_degree, int n_q_points_1d, typename Number>
struct FEEvaluationImplTransformToCollocation
}
/**
- * This function applies the tensor product kernel, corresponding to a
- * multiplication of 1d stripes, along the given @p direction of the tensor
- * data in the input array. This function allows the @p in and @p out
- * arrays to alias for the case n_rows == n_columns, i.e., it is safe to
- * perform the contraction in place where @p in and @p out point to the
- * same address. For the case n_rows != n_columns, the output is in general
- * not correct.
+ * This function applies the tensor product kernel with sum factorization,
+ * corresponding to a matrix-vector multiplication of 1d stripes, along
+ * the given @p direction of the tensor data in the input array. This
+ * function allows the @p in and @p out arrays to alias for the case
+ * n_rows == n_columns, i.e., it is safe to perform the contraction in
+ * place where @p in and @p out point to the same address. For the case
+ * `n_rows != n_columns`, the output is in general not correct.
*
* @tparam direction Direction that is evaluated
* @tparam contract_over_rows If true, the tensor contraction sums
* @tparam one_line If 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 @p false case.
+ * @tparam quantity Specify whether values, gradients or Hessians should
+ * be interpolated, allowing specialized algorithms
+ * for some class template parameters of `variant` to
+ * find the right path.
+ * @tparam extra_stride This parameter enables to place the result of the
+ * tensor product evaluation in the output array (if
+ * `contract_over_rows == true`) or input array (if
+ * `contract_over_rows == false`), which is used to
+ * group all components of a gradient adjacent in
+ * memory. If the stride is one, the data will form a
+ * contiguous range in memory.
*
* @param shape_data Transformation matrix with @p n_rows rows and
* @p n_columns columns, stored in row-major format
* @param in Pointer to the start of the input data vector
* @param out Pointer to the start of the output data vector
*/
- template <int direction,
- bool contract_over_rows,
- bool add,
- bool one_line = false,
- EvaluatorQuantity = EvaluatorQuantity::value,
- int extra_stride = 1>
+ template <int direction,
+ bool contract_over_rows,
+ bool add,
+ bool one_line = false,
+ EvaluatorQuantity quantity = EvaluatorQuantity::value,
+ int extra_stride = 1>
static void
apply(const Number2 *DEAL_II_RESTRICT shape_data,
const Number *in,
/**
- * This function applies the tensor product operation to produce face
- * values from cell values. The algorithm involved here can be interpreted
- * the first sweep in sum factorization, reducing the dimensionality of
- * the data set from dim-dimensional cell values to (dim-1)-dimensional
- * face values. This step is always done before we evaluate within the
- * face, as it reduces the dimensionality.
+ * This function applies the tensor product operation to produce face values
+ * from cell values. The algorithm involved here can be interpreted as the
+ * first sweep in sum factorization, reducing the dimensionality of the data
+ * set from dim-dimensional cell values to (dim-1)-dimensional face
+ * values. This step is always done before we evaluate within the face, as
+ * it reduces the length of the loops for the successive steps.
*
* @tparam n_rows_template The number of entries within the interpolation,
* typically equal to the polynomial degree plus one, if known
* the computed values overwrite the content in the output.
* @tparam max_derivative Sets 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, 2 up to second derivates. Note that all the
* derivatives access the data in @p shape_values passed to
* the constructor of the class.
*
- * @param shape_values address of the interpolation matrix
+ * @param shape_values Address of the interpolation matrix.
* @param n_blocks Number of interpolation layer used along the two other
- * dimensions tangential to the interpolation direction
+ * dimensions tangential to the interpolation direction.
* @param steps Increments in the input array from one step to the next,
* varied in conjunction with the @p stride variable.
* @param input Address of the input data vector.
+ /**
+ * This function performs the opposite operation to the interpolate_to_face
+ * function, done as the last step in sum factorization to embed face values
+ * and gradients back to values on all degrees of freedom of the cell.
+ */
template <int n_rows_template,
int stride_template,
bool contract_onto_face,