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 Types | Public Member Functions | Static Public Member Functions | Public Attributes | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Types | List of all members
FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType > Class Template Reference

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

Inheritance diagram for FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >:
[legend]

Public Types

using BaseClass = FEEvaluationAccess< dim, n_components_, Number, true, VectorizedArrayType >
 
using number_type = Number
 
using value_type = typename BaseClass::value_type
 
using gradient_type = typename BaseClass::gradient_type
 
using hessian_type = Tensor< 1, n_components_, Tensor< 2, dim, VectorizedArrayType > >
 
using NumberType = VectorizedArrayType
 
using shape_info_number_type = VectorizedArrayType
 
using ScalarNumber = typename internal::VectorizedArrayTrait< VectorizedArrayType >::value_type
 

Public Member Functions

 FEFaceEvaluation (const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free, const bool is_interior_face=true, const unsigned int dof_no=0, const unsigned int quad_no=0, const unsigned int first_selected_component=0, const unsigned int active_fe_index=numbers::invalid_unsigned_int, const unsigned int active_quad_index=numbers::invalid_unsigned_int, const unsigned int face_type=numbers::invalid_unsigned_int)
 
 FEFaceEvaluation (const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free, const std::pair< unsigned int, unsigned int > &range, const bool is_interior_face=true, const unsigned int dof_no=0, const unsigned int quad_no=0, const unsigned int first_selected_component=0)
 
void reinit (const unsigned int face_batch_number)
 
void reinit (const unsigned int cell_batch_number, const unsigned int face_number)
 
void evaluate (const EvaluationFlags::EvaluationFlags evaluation_flag)
 
void evaluate (const bool evaluate_values, const bool evaluate_gradients)
 
void evaluate (const VectorizedArrayType *values_array, const EvaluationFlags::EvaluationFlags evaluation_flag)
 
void evaluate (const VectorizedArrayType *values_array, const bool evaluate_values, const bool evaluate_gradients)
 
template<typename VectorType >
void gather_evaluate (const VectorType &input_vector, const EvaluationFlags::EvaluationFlags evaluation_flag)
 
template<typename VectorType >
void gather_evaluate (const VectorType &input_vector, const bool evaluate_values, const bool evaluate_gradients)
 
void integrate (const EvaluationFlags::EvaluationFlags integration_flag)
 
void integrate (const bool integrate_values, const bool integrate_gradients)
 
void integrate (const EvaluationFlags::EvaluationFlags integration_flag, VectorizedArrayType *values_array)
 
void integrate (const bool integrate_values, const bool integrate_gradients, VectorizedArrayType *values_array)
 
template<typename VectorType >
void integrate_scatter (const EvaluationFlags::EvaluationFlags integration_flag, VectorType &output_vector)
 
template<typename VectorType >
void integrate_scatter (const bool integrate_values, const bool integrate_gradients, VectorType &output_vector)
 
std_cxx20::ranges::iota_view< unsigned int, unsigned intdof_indices () const
 
bool at_boundary () const
 
types::boundary_id boundary_id () const
 
const MatrixFree< dim, Number, VectorizedArrayType > & get_matrix_free () const
 
void set_data_pointers (AlignedVector< VectorizedArrayType > *scratch_data, const unsigned int n_components)
 
void reinit_face (const internal::MatrixFreeFunctions::FaceToCellTopology< n_lanes > &face)
 
1: Reading from and writing to vectors
template<typename VectorType >
void read_dof_values (const VectorType &src, const unsigned int first_index=0, const std::bitset< VectorizedArrayType::size()> &mask=std::bitset< VectorizedArrayType::size()>().flip())
 
template<typename VectorType >
void read_dof_values_plain (const VectorType &src, const unsigned int first_index=0, const std::bitset< VectorizedArrayType::size()> &mask=std::bitset< VectorizedArrayType::size()>().flip())
 
template<typename VectorType >
void distribute_local_to_global (VectorType &dst, const unsigned int first_index=0, const std::bitset< VectorizedArrayType::size()> &mask=std::bitset< VectorizedArrayType::size()>().flip()) const
 
template<typename VectorType >
void set_dof_values (VectorType &dst, const unsigned int first_index=0, const std::bitset< VectorizedArrayType::size()> &mask=std::bitset< VectorizedArrayType::size()>().flip()) const
 
template<typename VectorType >
void set_dof_values_plain (VectorType &dst, const unsigned int first_index=0, const std::bitset< VectorizedArrayType::size()> &mask=std::bitset< VectorizedArrayType::size()>().flip()) const
 
2: Access to data at quadrature points or the gather vector data
value_type get_dof_value (const unsigned int dof) const
 
void submit_dof_value (const value_type val_in, const unsigned int dof)
 
value_type get_value (const unsigned int q_point) const
 
void submit_value (const value_type val_in, const unsigned int q_point)
 
gradient_type get_gradient (const unsigned int q_point) const
 
value_type get_normal_derivative (const unsigned int q_point) const
 
void submit_gradient (const gradient_type grad_in, const unsigned int q_point)
 
void submit_normal_derivative (const value_type grad_in, const unsigned int q_point)
 
void submit_hessian (const hessian_type hessian_in, const unsigned int q_point)
 
Tensor< 1, n_components_, Tensor< 2, dim, VectorizedArrayType > > get_hessian (const unsigned int q_point) const
 
gradient_type get_hessian_diagonal (const unsigned int q_point) const
 
value_type get_laplacian (const unsigned int q_point) const
 
VectorizedArrayType get_divergence (const unsigned int q_point) const
 
SymmetricTensor< 2, dim, VectorizedArrayType > get_symmetric_gradient (const unsigned int q_point) const
 
Tensor< 1,(dim==2 ? 1 :dim), VectorizedArrayType > get_curl (const unsigned int q_point) const
 
void submit_divergence (const VectorizedArrayType div_in, const unsigned int q_point)
 
void submit_symmetric_gradient (const SymmetricTensor< 2, dim, VectorizedArrayType > grad_in, const unsigned int q_point)
 
void submit_curl (const Tensor< 1, dim==2 ? 1 :dim, VectorizedArrayType > curl_in, const unsigned int q_point)
 
value_type integrate_value () const
 
1: Access to geometry data at quadrature points
ArrayView< VectorizedArrayType > get_scratch_data () const
 
VectorizedArrayType JxW (const unsigned int q_point) const
 
Point< dim, VectorizedArrayType > quadrature_point (const unsigned int q) const
 
Tensor< 2, dim, VectorizedArrayType > inverse_jacobian (const unsigned int q_point) const
 
Tensor< 1, dim, VectorizedArrayType > normal_vector (const unsigned int q_point) const
 
Tensor< 1, dim, VectorizedArrayType > get_normal_vector (const unsigned int q_point) const
 
2: Access to internal data arrays
const VectorizedArrayType * begin_dof_values () const
 
VectorizedArrayType * begin_dof_values ()
 
const VectorizedArrayType * begin_values () const
 
VectorizedArrayType * begin_values ()
 
const VectorizedArrayType * begin_gradients () const
 
VectorizedArrayType * begin_gradients ()
 
const VectorizedArrayType * begin_hessians () const
 
VectorizedArrayType * begin_hessians ()
 
3: Information about the current cell this class operates on
unsigned int get_mapping_data_index_offset () const
 
internal::MatrixFreeFunctions::GeometryType get_cell_type () const
 
const ShapeInfoTypeget_shape_info () const
 
const internal::MatrixFreeFunctions::DoFInfoget_dof_info () const
 
const std::vector< unsigned int > & get_internal_dof_numbering () const
 
unsigned int get_quadrature_index () const
 
unsigned int get_current_cell_index () const
 
unsigned int get_active_fe_index () const
 
unsigned int get_active_quadrature_index () const
 
unsigned int get_first_selected_component () const
 
std::uint8_t get_face_no (const unsigned int v=0) const
 
unsigned int get_subface_index () const
 
std::uint8_t get_face_orientation (const unsigned int v=0) const
 
internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex get_dof_access_index () const
 
bool is_interior_face () const
 
const std::array< unsigned int, n_lanes > & get_cell_ids () const
 
const std::array< unsigned int, n_lanes > & get_face_ids () const
 
unsigned int get_cell_or_face_batch_id () const
 
const std::array< unsigned int, n_lanes > & get_cell_or_face_ids () const
 
std_cxx20::ranges::iota_view< unsigned int, unsigned intquadrature_point_indices () const
 
3: Functions to access cell- and face-data vectors.
read_cell_data (const AlignedVector< T > &array) const
 
void set_cell_data (AlignedVector< T > &array, const T &value) const
 
read_face_data (const AlignedVector< T > &array) const
 
void set_face_data (AlignedVector< T > &array, const T &value) const
 

Static Public Member Functions

static bool fast_evaluation_supported (const unsigned int given_degree, const unsigned int given_n_q_points_1d)
 

Public Attributes

const unsigned int dofs_per_component
 
const unsigned int dofs_per_cell
 
const unsigned int n_q_points
 

Static Public Attributes

static constexpr unsigned int dimension = dim
 
static constexpr unsigned int n_components = n_components_
 
static constexpr unsigned int static_n_q_points
 
static constexpr unsigned int static_n_q_points_cell
 
static constexpr unsigned int static_dofs_per_component
 
static constexpr unsigned int tensor_dofs_per_cell
 
static constexpr unsigned int static_dofs_per_cell
 
static constexpr unsigned int n_lanes
 

Protected Member Functions

template<typename VectorType , typename VectorOperation >
void read_write_operation (const VectorOperation &operation, const std::array< VectorType *, n_components_ > &vectors, const std::array< const std::vector< ArrayView< const typename VectorType::value_type > > *, n_components_ > &vectors_sm, const std::bitset< VectorizedArrayType::size()> &mask, const bool apply_constraints=true) const
 
template<typename VectorType , typename VectorOperation >
void read_write_operation_contiguous (const VectorOperation &operation, const std::array< VectorType *, n_components_ > &vectors, const std::array< const std::vector< ArrayView< const typename VectorType::value_type > > *, n_components_ > &vectors_sm, const std::bitset< VectorizedArrayType::size()> &mask) const
 
template<typename VectorType , typename VectorOperation >
void read_write_operation_global (const VectorOperation &operation, const std::array< VectorType *, n_components_ > &vectors) const
 
void apply_hanging_node_constraints (const bool transpose) const
 

Protected Attributes

AlignedVector< VectorizedArrayType > * scratch_data_array
 
const MatrixFree< dim, Number, VectorizedArrayType > * matrix_free
 
std::vector< types::global_dof_indexlocal_dof_indices
 
const ShapeInfoTypedata
 
const DoFInfodof_info
 
const MappingInfoStorageTypemapping_data
 
const unsigned int quad_no
 
const unsigned int n_fe_components
 
const unsigned int first_selected_component
 
const unsigned int active_fe_index
 
const unsigned int active_quad_index
 
const MappingInfoStorageType::QuadratureDescriptordescriptor
 
const unsigned int n_quadrature_points
 
const Point< dim, VectorizedArrayType > * quadrature_points
 
const Tensor< 2, dim, VectorizedArrayType > * jacobian
 
const Tensor< 1, dim *(dim+1)/2, Tensor< 1, dim, VectorizedArrayType > > * jacobian_gradients
 
const Tensor< 1, dim *(dim+1)/2, Tensor< 1, dim, VectorizedArrayType > > * jacobian_gradients_non_inverse
 
const VectorizedArrayType * J_value
 
const Tensor< 1, dim, VectorizedArrayType > * normal_vectors
 
const Tensor< 1, dim, VectorizedArrayType > * normal_x_jacobian
 
const ScalarNumberquadrature_weights
 
ArrayView< VectorizedArrayType > scratch_data
 
VectorizedArrayType * values_dofs
 
VectorizedArrayType * values_quad
 
VectorizedArrayType * values_from_gradients_quad
 
VectorizedArrayType * gradients_quad
 
VectorizedArrayType * gradients_from_hessians_quad
 
VectorizedArrayType * hessians_quad
 
bool is_reinitialized
 
bool dof_values_initialized
 
bool values_quad_initialized
 
bool gradients_quad_initialized
 
bool hessians_quad_initialized
 
bool values_quad_submitted
 
bool gradients_quad_submitted
 
bool hessians_quad_submitted
 
unsigned int cell
 
bool interior_face
 
internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex dof_access_index
 
std::array< std::uint8_t, n_lanesface_numbers
 
std::array< std::uint8_t, n_lanesface_orientations
 
unsigned int subface_index
 
internal::MatrixFreeFunctions::GeometryType cell_type
 
std::array< unsigned int, n_lanescell_ids
 
std::array< unsigned int, n_lanesface_ids
 
std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly< dim, VectorizedArrayType > > mapped_geometry
 
bool divergence_is_requested
 

Private Types

using ShapeInfoType = internal::MatrixFreeFunctions::ShapeInfo< VectorizedArrayType >
 
using MappingInfoStorageType = internal::MatrixFreeFunctions::MappingInfoStorage<(is_face ? dim - 1 :dim), dim, VectorizedArrayType >
 
using DoFInfo = internal::MatrixFreeFunctions::DoFInfo
 

Detailed Description

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
class FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >

The class that provides all functions necessary to evaluate functions at quadrature points and face integrations. The design of the class is similar to FEEvaluation and most of the interfaces are shared with that class, in particular most access functions that come from the common base classes FEEvaluationAccess and FEEvaluationBase. Furthermore, the relation of this class to FEEvaluation is similar to the relation between FEValues and FEFaceValues.

Template Parameters
dimDimension in which this class is to be used
fe_degreeDegree of the tensor product finite element with fe_degree+1 degrees of freedom per coordinate direction. If set to -1, the degree of the underlying element will be used, which acts as a run time constant rather than a compile time constant that slows down the execution.
n_q_points_1dNumber of points in the quadrature formula in 1d, usually chosen as fe_degree+1
n_componentsNumber of vector components when solving a system of PDEs. If the same operation is applied to several components of a PDE (e.g. a vector Laplace equation), they can be applied simultaneously with one call (and often more efficiently)
NumberNumber format, usually double or float
VectorizedArrayTypeType of array to be woked on in a vectorized fashion, defaults to VectorizedArray<Number>
Note
Currently only VectorizedArray<Number, width> is supported as VectorizedArrayType.

Definition at line 2422 of file fe_evaluation.h.

Member Typedef Documentation

◆ BaseClass

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
using FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::BaseClass = FEEvaluationAccess<dim, n_components_, Number, true, VectorizedArrayType>

An alias to the base class.

Definition at line 2436 of file fe_evaluation.h.

◆ number_type

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
using FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::number_type = Number

A underlying number type specified as template argument.

Definition at line 2442 of file fe_evaluation.h.

◆ value_type

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
using FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::value_type = typename BaseClass::value_type

The type of function values, e.g. VectorizedArrayType for n_components=1 or Tensor<1,dim,VectorizedArrayType > for n_components=dim.

Definition at line 2449 of file fe_evaluation.h.

◆ gradient_type

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
using FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::gradient_type = typename BaseClass::gradient_type

The type of gradients, e.g. Tensor<1,dim,VectorizedArrayType> for n_components=1 or Tensor<2,dim,VectorizedArrayType > for n_components=dim.

Definition at line 2456 of file fe_evaluation.h.

◆ hessian_type

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
using FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::hessian_type = Tensor<1, n_components_, Tensor<2, dim, VectorizedArrayType> >
inherited

Definition at line 102 of file fe_evaluation.h.

◆ ShapeInfoType

using FEEvaluationData< dim, VectorizedArrayType , is_face >::ShapeInfoType = internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType >
privateinherited

Definition at line 116 of file fe_evaluation_data.h.

◆ MappingInfoStorageType

using FEEvaluationData< dim, VectorizedArrayType , is_face >::MappingInfoStorageType = internal::MatrixFreeFunctions:: MappingInfoStorage<(is_face ? dim - 1 : dim), dim, VectorizedArrayType >
privateinherited

Definition at line 117 of file fe_evaluation_data.h.

◆ DoFInfo

using FEEvaluationData< dim, VectorizedArrayType , is_face >::DoFInfo = internal::MatrixFreeFunctions::DoFInfo
privateinherited

Definition at line 119 of file fe_evaluation_data.h.

◆ NumberType

using FEEvaluationData< dim, VectorizedArrayType , is_face >::NumberType = VectorizedArrayType
inherited

Definition at line 124 of file fe_evaluation_data.h.

◆ shape_info_number_type

using FEEvaluationData< dim, VectorizedArrayType , is_face >::shape_info_number_type = VectorizedArrayType
inherited

Definition at line 125 of file fe_evaluation_data.h.

◆ ScalarNumber

using FEEvaluationData< dim, VectorizedArrayType , is_face >::ScalarNumber = typename internal::VectorizedArrayTrait<VectorizedArrayType >::value_type
inherited

Definition at line 126 of file fe_evaluation_data.h.

Constructor & Destructor Documentation

◆ FEFaceEvaluation() [1/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::FEFaceEvaluation ( const MatrixFree< dim, Number, VectorizedArrayType > &  matrix_free,
const bool  is_interior_face = true,
const unsigned int  dof_no = 0,
const unsigned int  quad_no = 0,
const unsigned int  first_selected_component = 0,
const unsigned int  active_fe_index = numbers::invalid_unsigned_int,
const unsigned int  active_quad_index = numbers::invalid_unsigned_int,
const unsigned int  face_type = numbers::invalid_unsigned_int 
)

Constructor. Takes all data stored in MatrixFree. If applied to problems with more than one finite element or more than one quadrature formula selected during construction of matrix_free, the appropriate component can be selected by the optional arguments.

Parameters
matrix_freeData object that contains all data
is_interior_faceThis selects which of the two cells of an internal face the current evaluator will be based upon. The interior face is the main face along which the normal vectors are oriented. The exterior face coming from the other side provides the same normal vector as the interior side, so if the outer normal vector to that side is desired, it must be multiplied by -1.
dof_noIf matrix_free was set up with multiple DoFHandler objects, this parameter selects to which DoFHandler/AffineConstraints pair the given evaluator should be attached to.
quad_noIf matrix_free was set up with multiple Quadrature objects, this parameter selects the appropriate number of the quadrature formula.
first_selected_componentIf the dof_handler selected by dof_no uses an FESystem consisting of more than one base element, this parameter selects the number of the base element in FESystem. Note that this does not directly relate to the component of the respective element due to the possibility for a multiplicity in the element.
active_fe_indexIf matrix_free was set up with DoFHandler objects with hp::FECollections, this parameter selects to which DoFHandler/AffineConstraints pair the given evaluator should be attached to.
face_typeIn the case of a face, indicate its reference-cell type (0 for line or quadrilateral 1 for triangle).
active_quad_indexIf matrix_free was set up with hp::Collection objects, this parameter selects the appropriate number of the quadrature formula.

◆ FEFaceEvaluation() [2/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::FEFaceEvaluation ( const MatrixFree< dim, Number, VectorizedArrayType > &  matrix_free,
const std::pair< unsigned int, unsigned int > &  range,
const bool  is_interior_face = true,
const unsigned int  dof_no = 0,
const unsigned int  quad_no = 0,
const unsigned int  first_selected_component = 0 
)

Constructor. Takes all data stored in MatrixFree for a given face range, which allows to automatically identify the active_fe_index and active_quad_index in case of a p-adaptive strategy.

The rest of the arguments are the same as in the constructor above.

Member Function Documentation

◆ reinit() [1/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::reinit ( const unsigned int  face_batch_number)

Initializes the operation pointer to the current face. This method is the default choice for face integration as the data stored in MappingInfo is stored according to this numbering. Unlike the reinit functions taking a cell iterator as argument below and the FEValues::reinit() methods, where the information related to a particular cell is generated in the reinit call, this function is very cheap since all data is pre-computed in matrix_free, and only a few indices and pointers have to be set appropriately.

◆ reinit() [2/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::reinit ( const unsigned int  cell_batch_number,
const unsigned int  face_number 
)

As opposed to the reinit() method from the base class, this reinit() method initializes for a given number of cells and a face number. This method is less efficient than the other reinit() method taking a numbering of the faces because it needs to copy the data associated with the faces to the cells in this call.

◆ fast_evaluation_supported()

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
static bool FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::fast_evaluation_supported ( const unsigned int  given_degree,
const unsigned int  given_n_q_points_1d 
)
static

Check if fast evaluation/integration is supported.

◆ evaluate() [1/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::evaluate ( const EvaluationFlags::EvaluationFlags  evaluation_flag)

Evaluates the function values, the gradients, and the Laplacians of the FE function given at the DoF values stored in the internal data field dof_values (that is usually filled by the read_dof_values() method) at the quadrature points on the unit cell. The function arguments specify which parts shall actually be computed. Needs to be called before the functions get_value(), get_gradient() or get_normal_derivative() give useful information (unless these values have been set manually by accessing the internal data pointers).

◆ evaluate() [2/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::evaluate ( const bool  evaluate_values,
const bool  evaluate_gradients 
)
Deprecated:
Please use the evaluate() function with the EvaluationFlags argument.

◆ evaluate() [3/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::evaluate ( const VectorizedArrayType *  values_array,
const EvaluationFlags::EvaluationFlags  evaluation_flag 
)

Evaluates the function values, the gradients, and the Laplacians of the FE function given at the DoF values in the input array values_array at the quadrature points on the unit cell. If multiple components are involved in the current FEEvaluation object, the sorting in values_array is such that all degrees of freedom for the first component come first, then all degrees of freedom for the second, and so on. The function arguments specify which parts shall actually be computed. Needs to be called before the functions get_value(), get_gradient(), or get_normal_derivative() give useful information (unless these values have been set manually).

◆ evaluate() [4/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::evaluate ( const VectorizedArrayType *  values_array,
const bool  evaluate_values,
const bool  evaluate_gradients 
)
Deprecated:
Please use the evaluate() function with the EvaluationFlags argument.

◆ gather_evaluate() [1/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
template<typename VectorType >
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::gather_evaluate ( const VectorType &  input_vector,
const EvaluationFlags::EvaluationFlags  evaluation_flag 
)

Reads from the input vector and evaluates the function values, the gradients, and the Laplacians of the FE function at the quadrature points on the unit cell. The function arguments specify which parts shall actually be computed. Needs to be called before the functions get_value(), get_gradient(), or get_normal_derivative() give useful information.

This call is equivalent to calling read_dof_values() followed by evaluate(), but might internally use some additional optimizations.

◆ gather_evaluate() [2/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
template<typename VectorType >
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::gather_evaluate ( const VectorType &  input_vector,
const bool  evaluate_values,
const bool  evaluate_gradients 
)
Deprecated:
Please use the gather_evaluate() function with the EvaluationFlags argument.

◆ integrate() [1/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::integrate ( const EvaluationFlags::EvaluationFlags  integration_flag)

This function takes the values and/or gradients that are stored on quadrature points, tests them by all the basis functions/gradients on the cell and performs the cell integration. The two function arguments integrate_val and integrate_grad are used to enable/disable some of values or gradients. The result is written into the internal data field dof_values (that is usually written into the result vector by the distribute_local_to_global() or set_dof_values() methods).

◆ integrate() [2/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::integrate ( const bool  integrate_values,
const bool  integrate_gradients 
)
Deprecated:
Please use the integrate() function with the EvaluationFlags argument.

◆ integrate() [3/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::integrate ( const EvaluationFlags::EvaluationFlags  integration_flag,
VectorizedArrayType *  values_array 
)

This function takes the values and/or gradients that are stored on quadrature points, tests them by all the basis functions/gradients on the cell and performs the cell integration. The two function arguments integrate_val and integrate_grad are used to enable/disable some of values or gradients. As opposed to the other integrate() method, this call stores the result of the testing in the given array values_array.

◆ integrate() [4/4]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::integrate ( const bool  integrate_values,
const bool  integrate_gradients,
VectorizedArrayType *  values_array 
)
Deprecated:
Please use the integrate() function with the EvaluationFlags argument.

◆ integrate_scatter() [1/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
template<typename VectorType >
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::integrate_scatter ( const EvaluationFlags::EvaluationFlags  integration_flag,
VectorType &  output_vector 
)

This function takes the values and/or gradients that are stored on quadrature points, tests them by all the basis functions/gradients on the cell and performs the cell integration. The two function arguments integrate_val and integrate_grad are used to enable/disable some of values or gradients.

This call is equivalent to calling integrate() followed by distribute_local_to_global(), but might internally use some additional optimizations.

◆ integrate_scatter() [2/2]

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
template<typename VectorType >
void FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::integrate_scatter ( const bool  integrate_values,
const bool  integrate_gradients,
VectorType &  output_vector 
)
Deprecated:
Please use the integrate_scatter() function with the EvaluationFlags argument.

◆ dof_indices()

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
std_cxx20::ranges::iota_view< unsigned int, unsigned int > FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::dof_indices ( ) const

Return an object that can be thought of as an array containing all indices from zero to dofs_per_cell. This allows to write code using range-based for loops.

◆ at_boundary()

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::at_boundary ( ) const

Return whether the face associated to this FEFaceEvaluation object is at the boundary.

◆ boundary_id()

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
types::boundary_id FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::boundary_id ( ) const

Return the boundary indicator of the face associated to this FEFaceEvaluation object.

If the return value is the special value numbers::internal_face_boundary_id, then the face is in the interior of the domain.

Note
Alternatively to this function, you can use MatrixFree::get_boundary_id() to get the same information if no FEFaceEvaluation object has been set up.

◆ read_dof_values()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::read_dof_values ( const VectorType &  src,
const unsigned int  first_index = 0,
const std::bitset< VectorizedArrayType::size()> &  mask = std::bitset< VectorizedArrayType::size()>().flip() 
)
inherited

For the vector src, read out the values on the degrees of freedom of the current cell, and store them internally. Similar functionality as the function DoFAccessor::get_interpolated_dof_values when no constraints are present, but it also includes constraints from hanging nodes, so one can see it as a similar function to AffineConstraints::read_dof_values as well. Note that if vectorization is enabled, the DoF values for several cells are set.

If some constraints on the vector are inhomogeneous, use the function read_dof_values_plain instead and provide the vector with useful data also in constrained positions by calling AffineConstraints::distribute. When accessing vector entries during the solution of linear systems, the temporary solution should always have homogeneous constraints and this method is the correct one.

If the given vector template class is a block vector (determined through the template function 'IsBlockVector<VectorType>::value', which checks for vectors derived from BlockVectorBase) or an std::vector<VectorType> or std::vector<VectorType *>, this function reads n_components blocks from the block vector starting at the index first_index. For non-block vectors, first_index is ignored.

Note
If this class was constructed without a MatrixFree object and the information is acquired on the fly through a DoFHandler<dim>::cell_iterator, only one single cell is used by this class and this function extracts the values of the underlying components on the given cell. This call is slower than the ones done through a MatrixFree object and lead to a structure that does not effectively use vectorization in the evaluate routines based on these values (instead, VectorizedArray::size() same copies are worked on).

◆ read_dof_values_plain()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::read_dof_values_plain ( const VectorType &  src,
const unsigned int  first_index = 0,
const std::bitset< VectorizedArrayType::size()> &  mask = std::bitset< VectorizedArrayType::size()>().flip() 
)
inherited

For the vector src, read out the values on the degrees of freedom of the current cell, and store them internally. Similar functionality as the function DoFAccessor::get_interpolated_dof_values. As opposed to the read_dof_values function, this function reads out the plain entries from vectors, without taking stored constraints into account. This way of access is appropriate when the constraints have been distributed on the vector by a call to AffineConstraints::distribute previously. This function is also necessary when inhomogeneous constraints are to be used, as MatrixFree can only handle homogeneous constraints. Note that if vectorization is enabled, the DoF values for several cells are set.

If the given vector template class is a block vector (determined through the template function 'IsBlockVector<VectorType>::value', which checks for vectors derived from BlockVectorBase) or an std::vector<VectorType> or std::vector<VectorType *>, this function reads n_components blocks from the block vector starting at the index first_index. For non-block vectors, first_index is ignored.

Note
If this class was constructed without a MatrixFree object and the information is acquired on the fly through a DoFHandler<dim>::cell_iterator, only one single cell is used by this class and this function extracts the values of the underlying components on the given cell. This call is slower than the ones done through a MatrixFree object and lead to a structure that does not effectively use vectorization in the evaluate routines based on these values (instead, VectorizedArray::size() same copies are worked on).

◆ distribute_local_to_global()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::distribute_local_to_global ( VectorType &  dst,
const unsigned int  first_index = 0,
const std::bitset< VectorizedArrayType::size()> &  mask = std::bitset< VectorizedArrayType::size()>().flip() 
) const
inherited

Takes the values stored internally on dof values of the current cell and sums them into the vector dst. The function also applies constraints during the write operation. The functionality is hence similar to the function AffineConstraints::distribute_local_to_global. If vectorization is enabled, the DoF values for several cells are used.

If the given vector template class is a block vector (determined through the template function 'IsBlockVector<VectorType>::value', which checks for vectors derived from BlockVectorBase) or an std::vector<VectorType> or std::vector<VectorType *>, this function writes to n_components blocks of the block vector starting at the index first_index. For non-block vectors, first_index is ignored.

The mask can be used to suppress the write access for some of the cells contained in the current cell vectorization batch, e.g. in case of local time stepping, where some cells are excluded from a call. A value of true in the bitset means that the respective lane index will be processed, whereas a value of false skips this index. The default setting is a bitset that contains all ones, which will write the accumulated integrals to all cells in the batch.

Note
If this class was constructed without a MatrixFree object and the information is acquired on the fly through a DoFHandler<dim>::cell_iterator, only one single cell is used by this class and this function extracts the values of the underlying components on the given cell. This call is slower than the ones done through a MatrixFree object and lead to a structure that does not effectively use vectorization in the evaluate routines based on these values (instead, VectorizedArray::size() same copies are worked on).

◆ set_dof_values()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::set_dof_values ( VectorType &  dst,
const unsigned int  first_index = 0,
const std::bitset< VectorizedArrayType::size()> &  mask = std::bitset< VectorizedArrayType::size()>().flip() 
) const
inherited

Takes the values stored internally on dof values of the current cell and writes them into the vector dst. The function skips the degrees of freedom which are constrained. As opposed to the distribute_local_to_global method, the old values at the position given by the current cell are overwritten. Thus, if a degree of freedom is associated to more than one cell (as usual in continuous finite elements), the values will be overwritten and only the value written last is retained. Please note that in a parallel context this function might also touch degrees of freedom owned by other MPI processes, so that a subsequent update or accumulation of ghost values as done by MatrixFree::loop() might invalidate the degrees of freedom set by this function.

If the given vector template class is a block vector (determined through the template function 'IsBlockVector<VectorType>::value', which checks for vectors derived from BlockVectorBase) or an std::vector<VectorType> or std::vector<VectorType *>, this function writes to n_components blocks of the block vector starting at the index first_index. For non-block vectors, first_index is ignored.

The mask can be used to suppress the write access for some of the cells contained in the current cell vectorization batch, e.g. in case of local time stepping, where some cells are excluded from a call. A value of true in the bitset means that the respective lane index will be processed, whereas a value of false skips this index. The default setting is a bitset that contains all ones, which will write the accumulated integrals to all cells in the batch.

Note
If this class was constructed without a MatrixFree object and the information is acquired on the fly through a DoFHandler<dim>::cell_iterator, only one single cell is used by this class and this function extracts the values of the underlying components on the given cell. This call is slower than the ones done through a MatrixFree object and lead to a structure that does not effectively use vectorization in the evaluate routines based on these values (instead, VectorizedArray::size() same copies are worked on).

◆ set_dof_values_plain()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::set_dof_values_plain ( VectorType &  dst,
const unsigned int  first_index = 0,
const std::bitset< VectorizedArrayType::size()> &  mask = std::bitset< VectorizedArrayType::size()>().flip() 
) const
inherited

Same as set_dof_values(), but without resolving constraints.

◆ get_dof_value()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
value_type FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_dof_value ( const unsigned int  dof) const
inherited

Return the value stored for the local degree of freedom with index dof. If the object is vector-valued, a vector-valued return argument is given. Thus, the argument dof can at most run until dofs_per_component rather than dofs_per_cell since the different components of a vector-valued FE are return together. Note that when vectorization is enabled, values from several cells are grouped together. If set_dof_values was called last, the value corresponds to the one set there. If integrate was called last, it instead corresponds to the value of the integrated function with the test function of the given index.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ submit_dof_value()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_dof_value ( const value_type  val_in,
const unsigned int  dof 
)
inherited

Write a value to the field containing the degrees of freedom with component dof. Writes to the same field as is accessed through get_dof_value. Therefore, the original data that was read from a vector is overwritten as soon as a value is submitted.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ get_value()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
value_type FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_value ( const unsigned int  q_point) const
inherited

Return the value of a finite element function at quadrature point number q_point after a call to FEEvaluation::evaluate() with EvaluationFlags::values set, or the value that has been stored there with a call to FEEvaluationBase::submit_value(). If the object is vector-valued, a vector-valued return argument is given. Note that when vectorization is enabled, values from several cells are grouped together.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ submit_value()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_value ( const value_type  val_in,
const unsigned int  q_point 
)
inherited

Write a value to the field containing the values on quadrature points with component q_point. Access to the same field as through get_value(). If applied before the function FEEvaluation::integrate() with EvaluationFlags::values set is called, this specifies the value which is tested by all basis function on the current cell and integrated over.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ get_gradient()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
gradient_type FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_gradient ( const unsigned int  q_point) const
inherited

Return the gradient of a finite element function at quadrature point number q_point after a call to FEEvaluation::evaluate() with EvaluationFlags::gradients, or the value that has been stored there with a call to FEEvaluationBase::submit_gradient().

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ get_normal_derivative()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
value_type FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_normal_derivative ( const unsigned int  q_point) const
inherited

Return the derivative of a finite element function at quadrature point number q_point after a call to FEEvaluation::evaluate(EvaluationFlags::gradients) the direction normal to the face: \(\boldsymbol \nabla u(\mathbf x_q) \cdot \mathbf n(\mathbf x_q)\)

This call is equivalent to calling get_gradient() * normal_vector() but will use a more efficient internal representation of data.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ submit_gradient()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_gradient ( const gradient_type  grad_in,
const unsigned int  q_point 
)
inherited

Write a contribution that is tested by the gradient to the field containing the values on quadrature points with component q_point. Access to the same field as through get_gradient(). If applied before the function FEEvaluation::integrate(EvaluationFlags::gradients) is called, this specifies what is tested by all basis function gradients on the current cell and integrated over.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ submit_normal_derivative()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_normal_derivative ( const value_type  grad_in,
const unsigned int  q_point 
)
inherited

Write a contribution that is tested by the gradient to the field containing the values on quadrature points with component q_point. Access to the same field as through get_gradient() or get_normal_derivative(). If applied before the function FEEvaluation::integrate(EvaluationFlags::gradients) is called, this specifies what is tested by all basis function gradients on the current cell and integrated over.

Note
This operation writes the data to the same field as submit_gradient(). As a consequence, only one of these two can be used. Usually, the contribution of a potential call to this function must be added into the contribution for submit_gradient().
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ submit_hessian()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_hessian ( const hessian_type  hessian_in,
const unsigned int  q_point 
)
inherited

Write a contribution that is tested by the Hessian to the field containing the values at quadrature points with component q_point. Access to the same field as through get_hessian(). If applied before the function FEEvaluation::integrate(EvaluationFlags::hessians) is called, this specifies what is tested by the Hessians of all basis functions on the current cell and integrated over.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ get_hessian()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
Tensor< 1, n_components_, Tensor< 2, dim, VectorizedArrayType > > FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_hessian ( const unsigned int  q_point) const
inherited

Return the Hessian of a finite element function at quadrature point number q_point after a call to FEEvaluation::evaluate(EvaluationFlags::hessians). If only the diagonal or even the trace of the Hessian, the Laplacian, is needed, use the other functions below.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ get_hessian_diagonal()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
gradient_type FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_hessian_diagonal ( const unsigned int  q_point) const
inherited

Return the diagonal of the Hessian of a finite element function at quadrature point number q_point after a call to FEEvaluation::evaluate(EvaluationFlags::hessians).

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ get_laplacian()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
value_type FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_laplacian ( const unsigned int  q_point) const
inherited

Return the Laplacian (i.e., the trace of the Hessian) of a finite element function at quadrature point number q_point after a call to FEEvaluation::evaluate(EvaluationFlags::hessians). Compared to the case when computing the full Hessian, some operations can be saved when only the Laplacian is requested.

Note
The derived class FEEvaluationAccess overloads this operation with specializations for the scalar case (n_components == 1) and for the vector-valued case (n_components == dim).

◆ get_divergence()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
VectorizedArrayType FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_divergence ( const unsigned int  q_point) const
inherited

Return the divergence of a vector-valued finite element at quadrature point number q_point after a call to evaluate(EvaluationFlags::gradients).

Note
Only available for the vector-valued case (n_components == dim).

◆ get_symmetric_gradient()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
SymmetricTensor< 2, dim, VectorizedArrayType > FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_symmetric_gradient ( const unsigned int  q_point) const
inherited

Return the symmetric gradient of a vector-valued finite element at quadrature point number q_point after a call to evaluate(EvaluationFlags::gradients). It corresponds to 0.5 (grad+gradT).

Note
Only available for the vector-valued case (n_components == dim).

◆ get_curl()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
Tensor< 1,(dim==2 ? 1 :dim), VectorizedArrayType > FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_curl ( const unsigned int  q_point) const
inherited

Return the curl of the vector field, \(\nabla \times v\) after a call to evaluate(EvaluationFlags::gradients).

Note
Only available for the vector-valued case (n_components == dim).

◆ submit_divergence()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_divergence ( const VectorizedArrayType  div_in,
const unsigned int  q_point 
)
inherited

Write a contribution that is tested by the divergence to the field containing the values on quadrature points with component q_point. Access to the same field as through get_gradient. If applied before the function integrate(EvaluationFlags::gradients) is called, this specifies what is tested by all basis function gradients on the current cell and integrated over.

Note
Only available for the vector-valued case (n_components == dim).
This operation writes the data to the same field as submit_gradient(). As a consequence, only one of these two can be used. Usually, the contribution of a potential call to this function must be added into the diagonal of the contribution for submit_gradient().

◆ submit_symmetric_gradient()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_symmetric_gradient ( const SymmetricTensor< 2, dim, VectorizedArrayType >  grad_in,
const unsigned int  q_point 
)
inherited

Write a contribution that is tested by the symmetric gradient to the field containing the values on quadrature points with component q_point. Access to the same field as through get_symmetric_gradient. If applied before the function integrate(EvaluationFlags::gradients) is called, this specifies the symmetric gradient which is tested by all basis function symmetric gradients on the current cell and integrated over.

Note
Only available for the vector-valued case (n_components == dim).
This operation writes the data to the same field as submit_gradient(). As a consequence, only one of these two can be used. Usually, the contribution of a potential call to this function must be added to the respective entries of the rank-2 tensor for submit_gradient().

◆ submit_curl()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::submit_curl ( const Tensor< 1, dim==2 ? 1 :dim, VectorizedArrayType >  curl_in,
const unsigned int  q_point 
)
inherited

Write the components of a curl containing the values on quadrature point q_point. Access to the same data field as through get_gradient.

Note
Only available for the vector-valued case (n_components == dim).
This operation writes the data to the same field as submit_gradient(). As a consequence, only one of these two can be used. Usually, the contribution of a potential call to this function must be added to the respective entries of the rank-2 tensor for submit_gradient().

◆ integrate_value()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
value_type FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::integrate_value ( ) const
inherited

Takes values at quadrature points, multiplies by the Jacobian determinant and quadrature weights (JxW) and sums the values for all quadrature points on the cell. The result is a scalar, representing the integral over the function over the cell. If a vector-element is used, the resulting components are still separated. Moreover, if vectorization is enabled, the integral values of several cells are contained in the slots of the returned VectorizedArray field.

Note
In case the FEEvaluation object is initialized with a batch of cells where not all lanes in the SIMD vector VectorizedArray are representing actual data, this method performs computations on dummy data (that is copied from the last valid lane) that will not make sense. Thus, the user needs to make sure that it is not used in any computation explicitly, like when summing the results of several cells.

◆ get_matrix_free()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
const MatrixFree< dim, Number, VectorizedArrayType > & FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::get_matrix_free ( ) const
inherited

Return the underlying MatrixFree object.

◆ read_write_operation()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType , typename VectorOperation >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::read_write_operation ( const VectorOperation operation,
const std::array< VectorType *, n_components_ > &  vectors,
const std::array< const std::vector< ArrayView< const typename VectorType::value_type > > *, n_components_ > &  vectors_sm,
const std::bitset< VectorizedArrayType::size()> &  mask,
const bool  apply_constraints = true 
) const
protectedinherited

A unified function to read from and write into vectors based on the given template operation. It can perform the operation for read_dof_values, distribute_local_to_global, and set_dof_values. It performs the operation for several vectors at a time.

◆ read_write_operation_contiguous()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType , typename VectorOperation >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::read_write_operation_contiguous ( const VectorOperation operation,
const std::array< VectorType *, n_components_ > &  vectors,
const std::array< const std::vector< ArrayView< const typename VectorType::value_type > > *, n_components_ > &  vectors_sm,
const std::bitset< VectorizedArrayType::size()> &  mask 
) const
protectedinherited

A unified function to read from and write into vectors based on the given template operation for DG-type schemes where all degrees of freedom on cells are contiguous. It can perform the operation for read_dof_values(), distribute_local_to_global(), and set_dof_values() for several vectors at a time, depending on n_components.

◆ read_write_operation_global()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
template<typename VectorType , typename VectorOperation >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::read_write_operation_global ( const VectorOperation operation,
const std::array< VectorType *, n_components_ > &  vectors 
) const
protectedinherited

A unified function to read from and write into vectors based on the given template operation for the case when we do not have an underlying MatrixFree object. It can perform the operation for read_dof_values, distribute_local_to_global, and set_dof_values. It performs the operation for several vectors at a time, depending on n_components.

◆ apply_hanging_node_constraints()

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
void FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::apply_hanging_node_constraints ( const bool  transpose) const
protectedinherited

Apply hanging-node constraints.

◆ set_data_pointers()

void FEEvaluationData< dim, VectorizedArrayType , is_face >::set_data_pointers ( AlignedVector< VectorizedArrayType > *  scratch_data,
const unsigned int  n_components 
)
inherited

Sets the pointers for values, gradients, hessians to the central scratch_data_array inside the given scratch array, for a given number of components as provided by one of the derived classes.

◆ reinit_face()

void FEEvaluationData< dim, VectorizedArrayType , is_face >::reinit_face ( const internal::MatrixFreeFunctions::FaceToCellTopology< n_lanes > &  face)
inherited

Initialize the indices related to the given face description. Used during the reinit() call of the FEFaceEvaluation class.

◆ get_scratch_data()

ArrayView< VectorizedArrayType > FEEvaluationData< dim, VectorizedArrayType , is_face >::get_scratch_data ( ) const
inherited

Return an ArrayView to internal memory for temporary use. Note that some of this memory is overwritten during evaluate() and integrate() calls so do not assume it to be stable over those calls. The maximum size you can write into is 3*dofs_per_cell+2*n_q_points.

◆ JxW()

VectorizedArrayType FEEvaluationData< dim, VectorizedArrayType , is_face >::JxW ( const unsigned int  q_point) const
inherited

Return the determinant of the Jacobian from the unit to the real cell times the quadrature weight.

◆ quadrature_point()

Point< dim, VectorizedArrayType > FEEvaluationData< dim, VectorizedArrayType , is_face >::quadrature_point ( const unsigned int  q) const
inherited

Returns the q-th quadrature point on the face in real coordinates stored in MappingInfo.

◆ inverse_jacobian()

Tensor< 2, dim, VectorizedArrayType > FEEvaluationData< dim, VectorizedArrayType , is_face >::inverse_jacobian ( const unsigned int  q_point) const
inherited

Return the inverse and transposed version \(J^{-\mathrm T}\) of the Jacobian of the mapping between the unit to the real cell defined as \(J_{ij} = d x_i / d\hat x_j\). The \((i,j)\) entry of the returned tensor contains \(d\hat x_j/dx_i\), i.e., columns refer to reference space coordinates and rows to real cell coordinates. Thus, the returned tensor represents a covariant transformation, which is used in the FEEvaluationBase::get_gradient() function to transform the unit cell gradients to gradients on the real cell by a multiplication \(J^{-\mathrm T} \hat{\nabla} u_h\).

◆ normal_vector()

Tensor< 1, dim, VectorizedArrayType > FEEvaluationData< dim, VectorizedArrayType , is_face >::normal_vector ( const unsigned int  q_point) const
inherited

Return the unit normal vector on a face. Note that both sides of a face use the same orientation of the normal vector: For the faces enumerated as interior in FaceToCellTopology and selected with the is_interior_face=true flag of the constructor, this corresponds to the outer normal vector, whereas for faces enumerated as exterior in FaceToCellTopology and selected with the is_interior_face=false flag of the constructor, the normal points into the element as a consequence of the single normal vector.

Note
Only implemented in case is_face == true.

◆ get_normal_vector()

Tensor< 1, dim, VectorizedArrayType > FEEvaluationData< dim, VectorizedArrayType , is_face >::get_normal_vector ( const unsigned int  q_point) const
inherited

Same as normal_vector(const unsigned int q_point).

Warning
This function will be deprecated!

◆ begin_dof_values() [1/2]

const VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_dof_values ( ) const
inherited

Return a read-only pointer to the first field of the dof values. This is the data field the read_dof_values() functions write into. First come the dof values for the first component, then all values for the second component, and so on. This is related to the internal data structures used in this class. In general, it is safer to use the get_dof_value() function instead.

◆ begin_dof_values() [2/2]

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_dof_values ( )
inherited

Return a read and write pointer to the first field of the dof values. This is the data field the read_dof_values() functions write into. First come the dof values for the first component, then all values for the second component, and so on. This is related to the internal data structures used in this class. In general, it is safer to use the get_dof_value() function instead.

◆ begin_values() [1/2]

const VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_values ( ) const
inherited

Return a read-only pointer to the first field of function values on quadrature points. First come the function values on all quadrature points for the first component, then all values for the second component, and so on. This is related to the internal data structures used in this class. The raw data after a call to evaluate only contains unit cell operations, so possible transformations, quadrature weights etc. must be applied manually. In general, it is safer to use the get_value() function instead, which does all the transformation internally.

◆ begin_values() [2/2]

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_values ( )
inherited

Return a read and write pointer to the first field of function values on quadrature points. First come the function values on all quadrature points for the first component, then all values for the second component, and so on. This is related to the internal data structures used in this class. The raw data after a call to evaluate only contains unit cell operations, so possible transformations, quadrature weights etc. must be applied manually. In general, it is safer to use the get_value() function instead, which does all the transformation internally.

◆ begin_gradients() [1/2]

const VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_gradients ( ) const
inherited

Return a read-only pointer to the first field of function gradients on quadrature points. First comes the x-component of the gradient for the first component on all quadrature points, then the y-component, and so on. Next comes the x-component of the second component, and so on. This is related to the internal data structures used in this class. The raw data after a call to evaluate only contains unit cell operations, so possible transformations, quadrature weights etc. must be applied manually. In general, it is safer to use the get_gradient() function instead, which does all the transformation internally.

◆ begin_gradients() [2/2]

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_gradients ( )
inherited

Return a read and write pointer to the first field of function gradients on quadrature points. First comes the x-component of the gradient for the first component on all quadrature points, then the y-component, and so on. Next comes the x-component of the second component, and so on. This is related to the internal data structures used in this class. The raw data after a call to evaluate only contains unit cell operations, so possible transformations, quadrature weights etc. must be applied manually. In general, it is safer to use the get_gradient() function instead, which does all the transformation internally.

◆ begin_hessians() [1/2]

const VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_hessians ( ) const
inherited

Return a read-only pointer to the first field of function hessians on quadrature points. First comes the xx-component of the hessian for the first component on all quadrature points, then the yy-component, zz-component in (3d), then the xy-component, and so on. Next comes the xx- component of the second component, and so on. This is related to the internal data structures used in this class. The raw data after a call to evaluate only contains unit cell operations, so possible transformations, quadrature weights etc. must be applied manually. In general, it is safer to use the get_laplacian() or get_hessian() functions instead, which does all the transformation internally.

◆ begin_hessians() [2/2]

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::begin_hessians ( )
inherited

Return a read and write pointer to the first field of function hessians on quadrature points. First comes the xx-component of the hessian for the first component on all quadrature points, then the yy-component, zz-component in (3d), then the xy-component, and so on. Next comes the xx-component of the second component, and so on. This is related to the internal data structures used in this class. The raw data after a call to evaluate only contains unit cell operations, so possible transformations, quadrature weights etc. must be applied manually. In general, it is safer to use the get_laplacian() or get_hessian() functions instead, which does all the transformation internally.

◆ get_mapping_data_index_offset()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_mapping_data_index_offset ( ) const
inherited

Return the index offset within the geometry fields for the cell the reinit() function has been called for. This index can be used to access an index into a field that has the same compression behavior as the Jacobian of the geometry, e.g., to store an effective coefficient tensors that combines a coefficient with the geometry for lower memory transfer as the available data fields.

◆ get_cell_type()

internal::MatrixFreeFunctions::GeometryType FEEvaluationData< dim, VectorizedArrayType , is_face >::get_cell_type ( ) const
inherited

Return the type of the cell the reinit() function has been called for. Valid values are cartesian for Cartesian cells (which allows for considerable data compression), affine for cells with affine mappings, and general for general cells without any compressed storage applied.

◆ get_shape_info()

const ShapeInfoType & FEEvaluationData< dim, VectorizedArrayType , is_face >::get_shape_info ( ) const
inherited

Return a reference to the ShapeInfo object currently in use.

◆ get_dof_info()

const internal::MatrixFreeFunctions::DoFInfo & FEEvaluationData< dim, VectorizedArrayType , is_face >::get_dof_info ( ) const
inherited

Return a reference to the DoFInfo object currently in use.

◆ get_internal_dof_numbering()

const std::vector< unsigned int > & FEEvaluationData< dim, VectorizedArrayType , is_face >::get_internal_dof_numbering ( ) const
inherited

Return the numbering of local degrees of freedom within the evaluation routines of FEEvaluation in terms of the standard numbering on finite elements.

◆ get_quadrature_index()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_quadrature_index ( ) const
inherited

Return the number of the quadrature formula of the present cell.

◆ get_current_cell_index()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_current_cell_index ( ) const
inherited

Return index of the current cell or face.

◆ get_active_fe_index()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_active_fe_index ( ) const
inherited

Return the active FE index for this class for efficient indexing in the hp-case.

◆ get_active_quadrature_index()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_active_quadrature_index ( ) const
inherited

Return the active quadrature index for this class for efficient indexing in the hp-case.

◆ get_first_selected_component()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_first_selected_component ( ) const
inherited

Return the first selected component in a multi-component system.

◆ get_face_no()

std::uint8_t FEEvaluationData< dim, VectorizedArrayType , is_face >::get_face_no ( const unsigned int  v = 0) const
inherited

If is_face is true, this function returns the face number of the given lane within a cell for the face this object was initialized to. On cells where the face makes no sense, an exception is thrown.

Note
Only available for dof_access_index == dof_access_cell and is_interior_face == false.
This function depends on the internal representation of data, which is not stable between releases of deal.II, and is hence mostly for internal use.

◆ get_subface_index()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_subface_index ( ) const
inherited

If is_face is true, this function returns the index of a subface along a face in case this object was initialized to a face with hanging nodes. On faces On cells where the face makes no sense, an exception is thrown.

Note
This function depends on the internal representation of data, which is not stable between releases of deal.II, and is hence mostly for internal use.

◆ get_face_orientation()

std::uint8_t FEEvaluationData< dim, VectorizedArrayType , is_face >::get_face_orientation ( const unsigned int  v = 0) const
inherited

If is_face is true, this function returns for a given lane the orientation index within an array of orientations as stored in ShapeInfo for unknowns and quadrature points.

Note
Only available for dof_access_index == dof_access_cell and is_interior_face == false.

◆ get_dof_access_index()

internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex FEEvaluationData< dim, VectorizedArrayType , is_face >::get_dof_access_index ( ) const
inherited

Return the current index in the access to compressed indices.

Note
This function depends on the internal representation of data, which is not stable between releases of deal.II, and is hence mostly for internal use.

◆ is_interior_face()

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::is_interior_face ( ) const
inherited

Return whether this object was set up to work on what is called "interior" faces in the evaluation context.

Note
This function depends on the internal representation of data, which is not stable between releases of deal.II, and is hence mostly for internal use.

◆ get_cell_ids()

const std::array< unsigned int, n_lanes > & FEEvaluationData< dim, VectorizedArrayType , is_face >::get_cell_ids ( ) const
inlineinherited

Return the id of the cells this FEEvaluation or FEFaceEvaluation is associated with.

Definition at line 495 of file fe_evaluation_data.h.

◆ get_face_ids()

const std::array< unsigned int, n_lanes > & FEEvaluationData< dim, VectorizedArrayType , is_face >::get_face_ids ( ) const
inlineinherited

Return the id of the faces this FEFaceEvaluation is associated with.

Definition at line 509 of file fe_evaluation_data.h.

◆ get_cell_or_face_batch_id()

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::get_cell_or_face_batch_id ( ) const
inlineinherited

Return the id of the cell/face batch this FEEvaluation/FEFaceEvaluation is associated with.

Definition at line 523 of file fe_evaluation_data.h.

◆ get_cell_or_face_ids()

const std::array< unsigned int, n_lanes > & FEEvaluationData< dim, VectorizedArrayType , is_face >::get_cell_or_face_ids ( ) const
inlineinherited

Return the id of the cells/faces this FEEvaluation/FEFaceEvaluation is associated with.

Definition at line 538 of file fe_evaluation_data.h.

◆ quadrature_point_indices()

std_cxx20::ranges::iota_view< unsigned int, unsigned int > FEEvaluationData< dim, VectorizedArrayType , is_face >::quadrature_point_indices ( ) const
inherited

Return an object that can be thought of as an array containing all indices from zero to n_quadrature_points. This allows to write code using range-based for loops.

◆ read_cell_data()

T FEEvaluationData< dim, VectorizedArrayType , is_face >::read_cell_data ( const AlignedVector< T > &  array) const
inherited

Provides a unified interface to access data in a vector of VectorizedArray fields of length MatrixFree::n_cell_batches() + MatrixFree::n_ghost_cell_batches() for cell data. It is implemented both for cells and faces (access data to the associated cell).

The underlying type can be both the Number type as given by the class template (i.e., VectorizedArrayType) as well as std::array with as many entries as n_lanes.

◆ set_cell_data()

void FEEvaluationData< dim, VectorizedArrayType , is_face >::set_cell_data ( AlignedVector< T > &  array,
const T &  value 
) const
inherited

Provides a unified interface to set data in a vector of VectorizedArray fields of length MatrixFree::n_cell_batches() + MatrixFree::n_ghost_cell_batches() for cell data. It is implemented both for cells and faces (access data to the associated cell).

◆ read_face_data()

T FEEvaluationData< dim, VectorizedArrayType , is_face >::read_face_data ( const AlignedVector< T > &  array) const
inherited

Provides a unified interface to access data in a vector of VectorizedArray fields of length MatrixFree::n_inner_face_batches() + MatrixFree::n_boundary_face_batches() + MatrixFree::n_ghost_inner_face_batches() for face data.

Note
Only implemented for faces.

◆ set_face_data()

void FEEvaluationData< dim, VectorizedArrayType , is_face >::set_face_data ( AlignedVector< T > &  array,
const T &  value 
) const
inherited

Provides a unified interface to set data in a vector of VectorizedArray fields of length MatrixFree::n_inner_face_batches() + MatrixFree::n_boundary_face_batches() + MatrixFree::n_ghost_inner_face_batches() for face data.

Note
Only implemented for faces.

Member Data Documentation

◆ dimension

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
constexpr unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::dimension = dim
staticconstexpr

The dimension given as template argument.

Definition at line 2461 of file fe_evaluation.h.

◆ n_components

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
constexpr unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::n_components = n_components_
staticconstexpr

The number of solution components of the evaluator given as template argument.

Definition at line 2467 of file fe_evaluation.h.

◆ static_n_q_points

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
constexpr unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::static_n_q_points
staticconstexpr
Initial value:
=
Utilities::pow(n_q_points_1d, dim - 1)
constexpr T pow(const T base, const int iexp)
Definition utilities.h:447

The static number of quadrature points determined from the given template argument n_q_points_1d taken to the power of dim-1.

Note
The actual number of quadrature points, n_q_points, can be different if fe_degree=-1 is given and run-time loop lengths are used rather than compile time ones.

Definition at line 2478 of file fe_evaluation.h.

◆ static_n_q_points_cell

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
constexpr unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::static_n_q_points_cell
staticconstexpr
Initial value:
=
Utilities::pow(n_q_points_1d, dim)

The static number of quadrature points on a cell with the same quadrature formula.

Note
This value is only present for simpler comparison with the cell quadrature, as the actual number of points is given to a face by the static_n_q_points variable.

Definition at line 2489 of file fe_evaluation.h.

◆ static_dofs_per_component

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
constexpr unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::static_dofs_per_component
staticconstexpr
Initial value:
=
Utilities::pow(fe_degree + 1, dim)

The static number of degrees of freedom of a scalar component determined from the given template argument fe_degree.

Note
The actual number of degrees of freedom dofs_per_component can be different if fe_degree=-1 is given.

Definition at line 2500 of file fe_evaluation.h.

◆ tensor_dofs_per_cell

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
constexpr unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::tensor_dofs_per_cell
staticconstexpr
Initial value:
=
static constexpr unsigned int n_components
static constexpr unsigned int static_dofs_per_component

The static number of degrees of freedom of all components determined from the given template argument fe_degree.

Note
The actual number of degrees of freedom dofs_per_cell can be different if fe_degree=-1 is given.

Definition at line 2511 of file fe_evaluation.h.

◆ static_dofs_per_cell

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
constexpr unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::static_dofs_per_cell
staticconstexpr
Initial value:

The static number of degrees of freedom of all components determined from the given template argument fe_degree.

Note
The actual number of degrees of freedom dofs_per_cell can be different if fe_degree=-1 is given.

Definition at line 2522 of file fe_evaluation.h.

◆ dofs_per_component

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
const unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::dofs_per_component

The number of degrees of freedom of a single component on the cell for the underlying evaluation object. Usually close to static_dofs_per_component, but the number depends on the actual element selected and is thus not static.

Definition at line 2788 of file fe_evaluation.h.

◆ dofs_per_cell

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
const unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::dofs_per_cell

The number of degrees of freedom on the cell accumulated over all components in the current evaluation object. Usually close to static_dofs_per_cell = static_dofs_per_component*n_components, but the number depends on the actual element selected and is thus not static.

Definition at line 2796 of file fe_evaluation.h.

◆ n_q_points

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
const unsigned int FEFaceEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number, VectorizedArrayType >::n_q_points

The number of quadrature points in use. If the number of quadrature points in 1d is given as a template, this number is simply the dim-1-th power of that value. If the element degree is set to -1 (dynamic selection of element degree), the static value of quadrature points is inaccurate and this value must be used instead.

Definition at line 2805 of file fe_evaluation.h.

◆ scratch_data_array

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
AlignedVector<VectorizedArrayType>* FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::scratch_data_array
protectedinherited

This is the general array for all data fields.

Definition at line 732 of file fe_evaluation.h.

◆ matrix_free

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
const MatrixFree<dim, Number, VectorizedArrayType>* FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::matrix_free
protectedinherited

A pointer to the underlying data.

Definition at line 737 of file fe_evaluation.h.

◆ local_dof_indices

template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
std::vector<types::global_dof_index> FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::local_dof_indices
mutableprotectedinherited

A temporary data structure necessary to read degrees of freedom when no MatrixFree object was given at initialization.

Definition at line 743 of file fe_evaluation.h.

◆ n_lanes

constexpr unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::n_lanes
staticconstexprinherited

Definition at line 129 of file fe_evaluation_data.h.

◆ data

const ShapeInfoType* FEEvaluationData< dim, VectorizedArrayType , is_face >::data
protectedinherited

A pointer to the unit cell shape data, i.e., values, gradients and Hessians in 1d at the quadrature points that constitute the tensor product. Also contained in matrix_info, but it simplifies code if we store a reference to it.

Definition at line 660 of file fe_evaluation_data.h.

◆ dof_info

const DoFInfo* FEEvaluationData< dim, VectorizedArrayType , is_face >::dof_info
protectedinherited

A pointer to the underlying DoF indices and constraint description for the component specified at construction. Also contained in matrix_info, but it simplifies code if we store a reference to it.

Definition at line 667 of file fe_evaluation_data.h.

◆ mapping_data

const MappingInfoStorageType* FEEvaluationData< dim, VectorizedArrayType , is_face >::mapping_data
protectedinherited

A pointer to the underlying transformation data from unit to real cells for the given quadrature formula specified at construction. Also contained in matrix_info, but it simplifies code if we store a reference to it.

Definition at line 675 of file fe_evaluation_data.h.

◆ quad_no

const unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::quad_no
protectedinherited

The number of the quadrature formula of the present cell among all quadrature formulas available in the MatrixFree objects pass to derived classes.

Definition at line 682 of file fe_evaluation_data.h.

◆ n_fe_components

const unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::n_fe_components
protectedinherited

Stores the number of components in the finite element as detected in the MatrixFree storage class for comparison with the template argument.

Definition at line 688 of file fe_evaluation_data.h.

◆ first_selected_component

const unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::first_selected_component
protectedinherited

For a FiniteElement with more than one base element, select at which component this data structure should start.

Definition at line 694 of file fe_evaluation_data.h.

◆ active_fe_index

const unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::active_fe_index
protectedinherited

The active FE index for this class for efficient indexing in the hp-case.

Definition at line 699 of file fe_evaluation_data.h.

◆ active_quad_index

const unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::active_quad_index
protectedinherited

The active quadrature index for this class for efficient indexing in the hp-case.

Definition at line 705 of file fe_evaluation_data.h.

◆ descriptor

const MappingInfoStorageType::QuadratureDescriptor* FEEvaluationData< dim, VectorizedArrayType , is_face >::descriptor
protectedinherited

A pointer to the underlying quadrature formula specified at construction. Also contained in mapping_data, but it simplifies code if we store a reference to it.

Definition at line 712 of file fe_evaluation_data.h.

◆ n_quadrature_points

const unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::n_quadrature_points
protectedinherited

The number of quadrature points in the current evaluation context.

Definition at line 717 of file fe_evaluation_data.h.

◆ quadrature_points

const Point<dim, VectorizedArrayType >* FEEvaluationData< dim, VectorizedArrayType , is_face >::quadrature_points
protectedinherited

A pointer to the quadrature-point information of the present cell. Only set to a useful value if on a non-Cartesian cell.

Definition at line 723 of file fe_evaluation_data.h.

◆ jacobian

const Tensor<2, dim, VectorizedArrayType >* FEEvaluationData< dim, VectorizedArrayType , is_face >::jacobian
protectedinherited

A pointer to the inverse transpose Jacobian information of the present cell. Only the first inverse transpose Jacobian (q_point = 0) is set for Cartesian/affine cells, and the actual Jacobian is stored at index 1 instead. For faces on hypercube elements, the derivatives are reorder s.t. the derivative orthogonal to the face is stored last, i.e for dim = 3 and face_no = 0 or 1, the derivatives are ordered as [dy, dz, dx], face_no = 2 or 3: [dz, dx, dy], and face_no = 5 or 6: [dx, dy, dz]. If the Jacobian also is stored, the components are instead reordered in the same way. Filled from MappingInfoStorage.jacobians in include/deal.II/matrix_free/mapping_info.templates.h

Definition at line 737 of file fe_evaluation_data.h.

◆ jacobian_gradients

const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, VectorizedArrayType > >* FEEvaluationData< dim, VectorizedArrayType , is_face >::jacobian_gradients
protectedinherited

A pointer to the gradients of the inverse Jacobian transformation of the present cell.

Definition at line 744 of file fe_evaluation_data.h.

◆ jacobian_gradients_non_inverse

const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, VectorizedArrayType > >* FEEvaluationData< dim, VectorizedArrayType , is_face >::jacobian_gradients_non_inverse
protectedinherited

A pointer to the gradients of the Jacobian transformation of the present cell.

Definition at line 751 of file fe_evaluation_data.h.

◆ J_value

const VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::J_value
protectedinherited

A pointer to the Jacobian determinant of the present cell. If on a Cartesian cell or on a cell with constant Jacobian, this is just the Jacobian determinant, otherwise the Jacobian determinant times the quadrature weight.

Definition at line 759 of file fe_evaluation_data.h.

◆ normal_vectors

const Tensor<1, dim, VectorizedArrayType >* FEEvaluationData< dim, VectorizedArrayType , is_face >::normal_vectors
protectedinherited

A pointer to the normal vectors at faces.

Definition at line 764 of file fe_evaluation_data.h.

◆ normal_x_jacobian

const Tensor<1, dim, VectorizedArrayType >* FEEvaluationData< dim, VectorizedArrayType , is_face >::normal_x_jacobian
protectedinherited

A pointer to the normal vectors times the jacobian at faces.

Definition at line 769 of file fe_evaluation_data.h.

◆ quadrature_weights

const ScalarNumber* FEEvaluationData< dim, VectorizedArrayType , is_face >::quadrature_weights
protectedinherited

A pointer to the quadrature weights of the underlying quadrature formula.

Definition at line 774 of file fe_evaluation_data.h.

◆ scratch_data

ArrayView<VectorizedArrayType > FEEvaluationData< dim, VectorizedArrayType , is_face >::scratch_data
mutableprotectedinherited

This is the user-visible part of FEEvaluationBase::scratch_data_array, only showing the part that can be consumed by various users. It is set during the allocation of the internal data structures in FEEvaluationBase.

Definition at line 782 of file fe_evaluation_data.h.

◆ values_dofs

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::values_dofs
protectedinherited

This field stores the values for local degrees of freedom (e.g. after reading out from a vector but before applying unit cell transformations or before distributing them into a result vector). The methods get_dof_value() and submit_dof_value() read from or write to this field.

The values of this array are stored in the start section of scratch_data_array. Due to its access as a thread local memory, the memory can get reused between different calls.

Definition at line 794 of file fe_evaluation_data.h.

◆ values_quad

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::values_quad
protectedinherited

This field stores the values of the finite element function on quadrature points after applying unit cell transformations or before integrating. The methods get_value() and submit_value() access this field.

The values of this array are stored in the start section of scratch_data_array. Due to its access as a thread local memory, the memory can get reused between different calls.

Definition at line 805 of file fe_evaluation_data.h.

◆ values_from_gradients_quad

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::values_from_gradients_quad
protectedinherited

This field stores the values of the finite element function on quadrature points after applying unit cell transformations or before integrating. This field is accessed when performing the contravariant Piola transform for gradients on general cells. This is done by the functions get_gradient() and submit_gradient() when using a H(div)- conforming finite element such as FE_RaviartThomasNodal.

The values of this array are stored in the start section of scratch_data_array. Due to its access as a thread local memory, the memory can get reused between different calls.

Definition at line 819 of file fe_evaluation_data.h.

◆ gradients_quad

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::gradients_quad
protectedinherited

This field stores the gradients of the finite element function on quadrature points after applying unit cell transformations or before integrating. The methods get_gradient() and submit_gradient() (as well as some specializations like get_symmetric_gradient() or get_divergence()) access this field.

The values of this array are stored in the start section of scratch_data_array. Due to its access as a thread local memory, the memory can get reused between different calls.

Definition at line 832 of file fe_evaluation_data.h.

◆ gradients_from_hessians_quad

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::gradients_from_hessians_quad
protectedinherited

This field stores the gradients of the finite element function on quadrature points after applying unit cell transformations or before integrating. The methods get_hessian() and submit_hessian() (as well as some specializations like get_hessian_diagonal() or get_laplacian()) access this field for general cell/face types.

The values of this array are stored in the start section of scratch_data_array. Due to its access as a thread local memory, the memory can get reused between different calls.

Definition at line 845 of file fe_evaluation_data.h.

◆ hessians_quad

VectorizedArrayType * FEEvaluationData< dim, VectorizedArrayType , is_face >::hessians_quad
protectedinherited

This field stores the Hessians of the finite element function on quadrature points after applying unit cell transformations. The methods get_hessian(), get_laplacian(), get_hessian_diagonal() access this field.

The values of this array are stored in the start section of scratch_data_array. Due to its access as a thread local memory, the memory can get reused between different calls.

Definition at line 856 of file fe_evaluation_data.h.

◆ is_reinitialized

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::is_reinitialized
protectedinherited

Debug information to track whether the reinit() function of derived classes was called.

Definition at line 862 of file fe_evaluation_data.h.

◆ dof_values_initialized

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::dof_values_initialized
protectedinherited

Debug information to track whether dof values have been initialized before accessed. Used to control exceptions when uninitialized data is used.

Definition at line 869 of file fe_evaluation_data.h.

◆ values_quad_initialized

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::values_quad_initialized
protectedinherited

Debug information to track whether values on quadrature points have been initialized before accessed. Used to control exceptions when uninitialized data is used.

Definition at line 876 of file fe_evaluation_data.h.

◆ gradients_quad_initialized

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::gradients_quad_initialized
protectedinherited

Debug information to track whether gradients on quadrature points have been initialized before accessed. Used to control exceptions when uninitialized data is used.

Definition at line 883 of file fe_evaluation_data.h.

◆ hessians_quad_initialized

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::hessians_quad_initialized
protectedinherited

Debug information to track whether Hessians on quadrature points have been initialized before accessed. Used to control exceptions when uninitialized data is used.

Definition at line 890 of file fe_evaluation_data.h.

◆ values_quad_submitted

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::values_quad_submitted
protectedinherited

Debug information to track whether values on quadrature points have been submitted for integration before the integration is actually stared. Used to control exceptions when uninitialized data is used.

Definition at line 897 of file fe_evaluation_data.h.

◆ gradients_quad_submitted

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::gradients_quad_submitted
protectedinherited

Debug information to track whether gradients on quadrature points have been submitted for integration before the integration is actually stared. Used to control exceptions when uninitialized data is used.

Definition at line 904 of file fe_evaluation_data.h.

◆ hessians_quad_submitted

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::hessians_quad_submitted
protectedinherited

Debug information to track whether hessians on quadrature points have been submitted for integration before the integration is actually stared. Used to control exceptions when uninitialized data is used.

Definition at line 911 of file fe_evaluation_data.h.

◆ cell

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::cell
protectedinherited

After a call to reinit(), stores the number of the cell we are currently working with.

Definition at line 917 of file fe_evaluation_data.h.

◆ interior_face

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::interior_face
protectedinherited

Flag holding information whether a face is an interior or exterior face according to the defined direction of the normal. For cells it defines if the dof values should be read from the actual cell corresponding to the interior face or the neighboring cell corresponding to the exterior face.

Definition at line 925 of file fe_evaluation_data.h.

◆ dof_access_index

internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex FEEvaluationData< dim, VectorizedArrayType , is_face >::dof_access_index
protectedinherited

Stores the index an FEFaceEvaluation object is currently pointing into (interior face, exterior face, data associated with cell).

Definition at line 931 of file fe_evaluation_data.h.

◆ face_numbers

std::array<std::uint8_t, n_lanes> FEEvaluationData< dim, VectorizedArrayType , is_face >::face_numbers
protectedinherited

Stores the (non-vectorized) number of faces within cells in case of ECL for the exterior cells where the single number face_no is not enough to cover all cases.

Note
Only available for dof_access_index == dof_access_cell and is_interior_face == false.

Definition at line 941 of file fe_evaluation_data.h.

◆ face_orientations

std::array<std::uint8_t, n_lanes> FEEvaluationData< dim, VectorizedArrayType , is_face >::face_orientations
protectedinherited

Store the orientation of the neighbor's faces with respect to the current cell for the case of exterior faces on ECL with possibly different orientations behind different cells.

Note
Only available for dof_access_index == dof_access_cell and is_interior_face == false.

Definition at line 951 of file fe_evaluation_data.h.

◆ subface_index

unsigned int FEEvaluationData< dim, VectorizedArrayType , is_face >::subface_index
protectedinherited

Stores the subface index of the given face. Usually, this variable takes the value numbers::invalid_unsigned_int to indicate integration over the full face, but in case the current physical face has a neighbor that is more refined, it is a subface and must scale the entries in ShapeInfo appropriately.

Definition at line 960 of file fe_evaluation_data.h.

◆ cell_type

internal::MatrixFreeFunctions::GeometryType FEEvaluationData< dim, VectorizedArrayType , is_face >::cell_type
protectedinherited

Stores the type of the cell we are currently working with after a call to reinit(). Valid values are cartesian, affine and general, which have different implications on how the Jacobian transformations are stored internally in MappingInfo.

Definition at line 968 of file fe_evaluation_data.h.

◆ cell_ids

std::array<unsigned int, n_lanes> FEEvaluationData< dim, VectorizedArrayType , is_face >::cell_ids
protectedinherited

Stores the (non-vectorized) id of the cells or faces this object is initialized to. Relevant for ECL.

Definition at line 974 of file fe_evaluation_data.h.

◆ face_ids

std::array<unsigned int, n_lanes> FEEvaluationData< dim, VectorizedArrayType , is_face >::face_ids
protectedinherited

Stores the (non-vectorized) id of the cells or faces this object is initialized to. Relevant for ECL.

Definition at line 980 of file fe_evaluation_data.h.

◆ mapped_geometry

std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly<dim, VectorizedArrayType > > FEEvaluationData< dim, VectorizedArrayType , is_face >::mapped_geometry
protectedinherited

Geometry data that can be generated FEValues on the fly with the respective constructor, as an alternative to the entry point with MatrixFree.

Definition at line 989 of file fe_evaluation_data.h.

◆ divergence_is_requested

bool FEEvaluationData< dim, VectorizedArrayType , is_face >::divergence_is_requested
protectedinherited

Bool indicating if the divergence is requested. Used internally in the case of the Piola transform.

Definition at line 995 of file fe_evaluation_data.h.


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