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
Classes | Public Types | Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Types | Friends | List of all members
FEEvaluationData< dim, Number, is_face > Class Template Reference

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

Inheritance diagram for FEEvaluationData< dim, Number, is_face >:
[legend]

Classes

struct  InitializationData
 

Public Types

using NumberType = Number
 
using shape_info_number_type = Number
 
using ScalarNumber = typename internal::VectorizedArrayTrait< Number >::value_type
 

Public Member Functions

 FEEvaluationData (const ShapeInfoType &shape_info, const bool is_interior_face=true)
 
 FEEvaluationData (const FEEvaluationData &other)=default
 
FEEvaluationDataoperator= (const FEEvaluationData &other)
 
virtual ~FEEvaluationData ()=default
 
void set_data_pointers (AlignedVector< Number > *scratch_data, const unsigned int n_components)
 
void reinit_face (const internal::MatrixFreeFunctions::FaceToCellTopology< n_lanes > &face)
 
1: Access to geometry data at quadrature points
ArrayView< Number > get_scratch_data () const
 
Number JxW (const unsigned int q_point) const
 
Point< dim, Number > quadrature_point (const unsigned int q) const
 
Tensor< 2, dim, Number > inverse_jacobian (const unsigned int q_point) const
 
Tensor< 1, dim, Number > normal_vector (const unsigned int q_point) const
 
Tensor< 1, dim, Number > get_normal_vector (const unsigned int q_point) const
 
2: Access to internal data arrays
const Number * begin_dof_values () const
 
Number * begin_dof_values ()
 
const Number * begin_values () const
 
Number * begin_values ()
 
const Number * begin_gradients () const
 
Number * begin_gradients ()
 
const Number * begin_hessians () const
 
Number * 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.
template<typename T >
read_cell_data (const AlignedVector< T > &array) const
 
template<typename T >
void set_cell_data (AlignedVector< T > &array, const T &value) const
 
template<typename T >
read_face_data (const AlignedVector< T > &array) const
 
template<typename T >
void set_face_data (AlignedVector< T > &array, const T &value) const
 

Static Public Attributes

static constexpr unsigned int dimension = dim
 
static constexpr unsigned int n_lanes
 

Protected Member Functions

 FEEvaluationData (const InitializationData &initialization_data, const bool is_interior_face, const unsigned int quad_no, const unsigned int first_selected_component)
 
 FEEvaluationData (const std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly< dim, Number > > &mapping_data, const unsigned int n_fe_components, const unsigned int first_selected_component)
 

Protected Attributes

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, Number > * quadrature_points
 
const Tensor< 2, dim, Number > * jacobian
 
const Tensor< 1, dim *(dim+1)/2, Tensor< 1, dim, Number > > * jacobian_gradients
 
const Tensor< 1, dim *(dim+1)/2, Tensor< 1, dim, Number > > * jacobian_gradients_non_inverse
 
const Number * J_value
 
const Tensor< 1, dim, Number > * normal_vectors
 
const Tensor< 1, dim, Number > * normal_x_jacobian
 
const ScalarNumberquadrature_weights
 
ArrayView< Number > scratch_data
 
Number * values_dofs
 
Number * values_quad
 
Number * values_from_gradients_quad
 
Number * gradients_quad
 
Number * gradients_from_hessians_quad
 
Number * 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, Number > > mapped_geometry
 
bool divergence_is_requested
 

Private Types

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

Friends

template<int , int , typename , bool , typename >
class FEEvaluationBase
 
template<int , int , int , int , typename , typename >
class FEEvaluation
 

Detailed Description

template<int dim, typename Number, bool is_face>
class FEEvaluationData< dim, Number, is_face >

This base class of the FEEvaluation and FEFaceEvaluation classes handles mapping-related information independent of the degrees of freedom and finite element in use. This class provides access functionality for user code. The main usage is through the class FEEvaluation. However, there is functionality to construct an object of type FEEvaluationData given suitable data pointers to internal data, which allows usage in some scenarios where no full MatrixFree object is available.

This class has four template arguments:

Template Parameters
dimDimension in which this class is to be used
NumberNumber format, usually double or float
is_faceWhether the class is used for a cell integrator (with quadrature dimension the same as the space dimension) or for a face integrator (with quadrature dimension one less)
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 114 of file fe_evaluation_data.h.

Member Typedef Documentation

◆ ShapeInfoType

template<int dim, typename Number , bool is_face>
using FEEvaluationData< dim, Number, is_face >::ShapeInfoType = internal::MatrixFreeFunctions::ShapeInfo<Number>
private

Definition at line 116 of file fe_evaluation_data.h.

◆ MappingInfoStorageType

template<int dim, typename Number , bool is_face>
using FEEvaluationData< dim, Number, is_face >::MappingInfoStorageType = internal::MatrixFreeFunctions:: MappingInfoStorage<(is_face ? dim - 1 : dim), dim, Number>
private

Definition at line 117 of file fe_evaluation_data.h.

◆ DoFInfo

template<int dim, typename Number , bool is_face>
using FEEvaluationData< dim, Number, is_face >::DoFInfo = internal::MatrixFreeFunctions::DoFInfo
private

Definition at line 119 of file fe_evaluation_data.h.

◆ NumberType

template<int dim, typename Number , bool is_face>
using FEEvaluationData< dim, Number, is_face >::NumberType = Number

Definition at line 124 of file fe_evaluation_data.h.

◆ shape_info_number_type

template<int dim, typename Number , bool is_face>
using FEEvaluationData< dim, Number, is_face >::shape_info_number_type = Number

Definition at line 125 of file fe_evaluation_data.h.

◆ ScalarNumber

template<int dim, typename Number , bool is_face>
using FEEvaluationData< dim, Number, is_face >::ScalarNumber = typename internal::VectorizedArrayTrait<Number>::value_type

Definition at line 126 of file fe_evaluation_data.h.

Constructor & Destructor Documentation

◆ FEEvaluationData() [1/4]

template<int dim, typename Number , bool is_face>
FEEvaluationData< dim, Number, is_face >::FEEvaluationData ( const ShapeInfoType shape_info,
const bool  is_interior_face = true 
)

Constructor, taking a single ShapeInfo object to inject the capability for evaluation and set up some data containers, compatible with the internal evaluation kernels of FEEvaluationImpl and friends. For actual use, select one of the derived classes FEEvaluation or FEFaceEvaluation with their respective arguments.

◆ FEEvaluationData() [2/4]

template<int dim, typename Number , bool is_face>
FEEvaluationData< dim, Number, is_face >::FEEvaluationData ( const FEEvaluationData< dim, Number, is_face > &  other)
default

Copy constructor.

◆ ~FEEvaluationData()

template<int dim, typename Number , bool is_face>
virtual FEEvaluationData< dim, Number, is_face >::~FEEvaluationData ( )
virtualdefault

Destructor.

◆ FEEvaluationData() [3/4]

template<int dim, typename Number , bool is_face>
FEEvaluationData< dim, Number, is_face >::FEEvaluationData ( const InitializationData initialization_data,
const bool  is_interior_face,
const unsigned int  quad_no,
const unsigned int  first_selected_component 
)
protected

Constructor filling information available from a MatrixFree class, collected in an internal data structure to reduce overhead of initialization. Used in derived classes when this class is initialized from a MatrixFree object.

◆ FEEvaluationData() [4/4]

template<int dim, typename Number , bool is_face>
FEEvaluationData< dim, Number, is_face >::FEEvaluationData ( const std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly< dim, Number > > &  mapping_data,
const unsigned int  n_fe_components,
const unsigned int  first_selected_component 
)
protected

Constructor that comes with reduced functionality and works similarly as FEValues.

Member Function Documentation

◆ operator=()

template<int dim, typename Number , bool is_face>
FEEvaluationData & FEEvaluationData< dim, Number, is_face >::operator= ( const FEEvaluationData< dim, Number, is_face > &  other)

Copy assignment operator.

◆ set_data_pointers()

template<int dim, typename Number , bool is_face>
void FEEvaluationData< dim, Number, is_face >::set_data_pointers ( AlignedVector< Number > *  scratch_data,
const unsigned int  n_components 
)

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()

template<int dim, typename Number , bool is_face>
void FEEvaluationData< dim, Number, is_face >::reinit_face ( const internal::MatrixFreeFunctions::FaceToCellTopology< n_lanes > &  face)

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

◆ get_scratch_data()

template<int dim, typename Number , bool is_face>
ArrayView< Number > FEEvaluationData< dim, Number, is_face >::get_scratch_data ( ) const

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()

template<int dim, typename Number , bool is_face>
Number FEEvaluationData< dim, Number, is_face >::JxW ( const unsigned int  q_point) const

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

◆ quadrature_point()

template<int dim, typename Number , bool is_face>
Point< dim, Number > FEEvaluationData< dim, Number, is_face >::quadrature_point ( const unsigned int  q) const

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

◆ inverse_jacobian()

template<int dim, typename Number , bool is_face>
Tensor< 2, dim, Number > FEEvaluationData< dim, Number, is_face >::inverse_jacobian ( const unsigned int  q_point) const

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()

template<int dim, typename Number , bool is_face>
Tensor< 1, dim, Number > FEEvaluationData< dim, Number, is_face >::normal_vector ( const unsigned int  q_point) const

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()

template<int dim, typename Number , bool is_face>
Tensor< 1, dim, Number > FEEvaluationData< dim, Number, is_face >::get_normal_vector ( const unsigned int  q_point) const

Same as normal_vector(const unsigned int q_point).

Warning
This function will be deprecated!

◆ begin_dof_values() [1/2]

template<int dim, typename Number , bool is_face>
const Number * FEEvaluationData< dim, Number, is_face >::begin_dof_values ( ) const

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]

template<int dim, typename Number , bool is_face>
Number * FEEvaluationData< dim, Number, is_face >::begin_dof_values ( )

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]

template<int dim, typename Number , bool is_face>
const Number * FEEvaluationData< dim, Number, is_face >::begin_values ( ) const

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]

template<int dim, typename Number , bool is_face>
Number * FEEvaluationData< dim, Number, is_face >::begin_values ( )

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]

template<int dim, typename Number , bool is_face>
const Number * FEEvaluationData< dim, Number, is_face >::begin_gradients ( ) const

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]

template<int dim, typename Number , bool is_face>
Number * FEEvaluationData< dim, Number, is_face >::begin_gradients ( )

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]

template<int dim, typename Number , bool is_face>
const Number * FEEvaluationData< dim, Number, is_face >::begin_hessians ( ) const

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]

template<int dim, typename Number , bool is_face>
Number * FEEvaluationData< dim, Number, is_face >::begin_hessians ( )

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()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_mapping_data_index_offset ( ) const

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()

template<int dim, typename Number , bool is_face>
internal::MatrixFreeFunctions::GeometryType FEEvaluationData< dim, Number, is_face >::get_cell_type ( ) const

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()

template<int dim, typename Number , bool is_face>
const ShapeInfoType & FEEvaluationData< dim, Number, is_face >::get_shape_info ( ) const

Return a reference to the ShapeInfo object currently in use.

◆ get_dof_info()

template<int dim, typename Number , bool is_face>
const internal::MatrixFreeFunctions::DoFInfo & FEEvaluationData< dim, Number, is_face >::get_dof_info ( ) const

Return a reference to the DoFInfo object currently in use.

◆ get_internal_dof_numbering()

template<int dim, typename Number , bool is_face>
const std::vector< unsigned int > & FEEvaluationData< dim, Number, is_face >::get_internal_dof_numbering ( ) const

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()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_quadrature_index ( ) const

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

◆ get_current_cell_index()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_current_cell_index ( ) const

Return index of the current cell or face.

◆ get_active_fe_index()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_active_fe_index ( ) const

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

◆ get_active_quadrature_index()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_active_quadrature_index ( ) const

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

◆ get_first_selected_component()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_first_selected_component ( ) const

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

◆ get_face_no()

template<int dim, typename Number , bool is_face>
std::uint8_t FEEvaluationData< dim, Number, is_face >::get_face_no ( const unsigned int  v = 0) const

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()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_subface_index ( ) const

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()

template<int dim, typename Number , bool is_face>
std::uint8_t FEEvaluationData< dim, Number, is_face >::get_face_orientation ( const unsigned int  v = 0) const

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()

template<int dim, typename Number , bool is_face>
internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex FEEvaluationData< dim, Number, is_face >::get_dof_access_index ( ) const

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()

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::is_interior_face ( ) const

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()

template<int dim, typename Number , bool is_face>
const std::array< unsigned int, n_lanes > & FEEvaluationData< dim, Number, is_face >::get_cell_ids ( ) const
inline

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()

template<int dim, typename Number , bool is_face>
const std::array< unsigned int, n_lanes > & FEEvaluationData< dim, Number, is_face >::get_face_ids ( ) const
inline

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()

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::get_cell_or_face_batch_id ( ) const
inline

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()

template<int dim, typename Number , bool is_face>
const std::array< unsigned int, n_lanes > & FEEvaluationData< dim, Number, is_face >::get_cell_or_face_ids ( ) const
inline

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()

template<int dim, typename Number , bool is_face>
std_cxx20::ranges::iota_view< unsigned int, unsigned int > FEEvaluationData< dim, Number, is_face >::quadrature_point_indices ( ) const

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()

template<int dim, typename Number , bool is_face>
template<typename T >
T FEEvaluationData< dim, Number, is_face >::read_cell_data ( const AlignedVector< T > &  array) const

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()

template<int dim, typename Number , bool is_face>
template<typename T >
void FEEvaluationData< dim, Number, is_face >::set_cell_data ( AlignedVector< T > &  array,
const T &  value 
) const

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()

template<int dim, typename Number , bool is_face>
template<typename T >
T FEEvaluationData< dim, Number, is_face >::read_face_data ( const AlignedVector< T > &  array) const

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()

template<int dim, typename Number , bool is_face>
template<typename T >
void FEEvaluationData< dim, Number, is_face >::set_face_data ( AlignedVector< T > &  array,
const T &  value 
) const

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.

Friends And Related Symbol Documentation

◆ FEEvaluationBase

template<int dim, typename Number , bool is_face>
template<int , int , typename , bool , typename >
friend class FEEvaluationBase
friend

Definition at line 1000 of file fe_evaluation_data.h.

◆ FEEvaluation

template<int dim, typename Number , bool is_face>
template<int , int , int , int , typename , typename >
friend class FEEvaluation
friend

Definition at line 1003 of file fe_evaluation_data.h.

Member Data Documentation

◆ dimension

template<int dim, typename Number , bool is_face>
constexpr unsigned int FEEvaluationData< dim, Number, is_face >::dimension = dim
staticconstexpr

Definition at line 122 of file fe_evaluation_data.h.

◆ n_lanes

template<int dim, typename Number , bool is_face>
constexpr unsigned int FEEvaluationData< dim, Number, is_face >::n_lanes
staticconstexpr
Initial value:

Definition at line 129 of file fe_evaluation_data.h.

◆ data

template<int dim, typename Number , bool is_face>
const ShapeInfoType* FEEvaluationData< dim, Number, is_face >::data
protected

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

template<int dim, typename Number , bool is_face>
const DoFInfo* FEEvaluationData< dim, Number, is_face >::dof_info
protected

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

template<int dim, typename Number , bool is_face>
const MappingInfoStorageType* FEEvaluationData< dim, Number, is_face >::mapping_data
protected

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

template<int dim, typename Number , bool is_face>
const unsigned int FEEvaluationData< dim, Number, is_face >::quad_no
protected

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

template<int dim, typename Number , bool is_face>
const unsigned int FEEvaluationData< dim, Number, is_face >::n_fe_components
protected

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

template<int dim, typename Number , bool is_face>
const unsigned int FEEvaluationData< dim, Number, is_face >::first_selected_component
protected

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

template<int dim, typename Number , bool is_face>
const unsigned int FEEvaluationData< dim, Number, is_face >::active_fe_index
protected

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

template<int dim, typename Number , bool is_face>
const unsigned int FEEvaluationData< dim, Number, is_face >::active_quad_index
protected

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

template<int dim, typename Number , bool is_face>
const MappingInfoStorageType::QuadratureDescriptor* FEEvaluationData< dim, Number, is_face >::descriptor
protected

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

template<int dim, typename Number , bool is_face>
const unsigned int FEEvaluationData< dim, Number, is_face >::n_quadrature_points
protected

The number of quadrature points in the current evaluation context.

Definition at line 717 of file fe_evaluation_data.h.

◆ quadrature_points

template<int dim, typename Number , bool is_face>
const Point<dim, Number>* FEEvaluationData< dim, Number, is_face >::quadrature_points
protected

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

template<int dim, typename Number , bool is_face>
const Tensor<2, dim, Number>* FEEvaluationData< dim, Number, is_face >::jacobian
protected

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

template<int dim, typename Number , bool is_face>
const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number> >* FEEvaluationData< dim, Number, is_face >::jacobian_gradients
protected

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

template<int dim, typename Number , bool is_face>
const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number> >* FEEvaluationData< dim, Number, is_face >::jacobian_gradients_non_inverse
protected

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

template<int dim, typename Number , bool is_face>
const Number* FEEvaluationData< dim, Number, is_face >::J_value
protected

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

template<int dim, typename Number , bool is_face>
const Tensor<1, dim, Number>* FEEvaluationData< dim, Number, is_face >::normal_vectors
protected

A pointer to the normal vectors at faces.

Definition at line 764 of file fe_evaluation_data.h.

◆ normal_x_jacobian

template<int dim, typename Number , bool is_face>
const Tensor<1, dim, Number>* FEEvaluationData< dim, Number, is_face >::normal_x_jacobian
protected

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

Definition at line 769 of file fe_evaluation_data.h.

◆ quadrature_weights

template<int dim, typename Number , bool is_face>
const ScalarNumber* FEEvaluationData< dim, Number, is_face >::quadrature_weights
protected

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

Definition at line 774 of file fe_evaluation_data.h.

◆ scratch_data

template<int dim, typename Number , bool is_face>
ArrayView<Number> FEEvaluationData< dim, Number, is_face >::scratch_data
mutableprotected

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

template<int dim, typename Number , bool is_face>
Number* FEEvaluationData< dim, Number, is_face >::values_dofs
protected

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

template<int dim, typename Number , bool is_face>
Number* FEEvaluationData< dim, Number, is_face >::values_quad
protected

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

template<int dim, typename Number , bool is_face>
Number* FEEvaluationData< dim, Number, is_face >::values_from_gradients_quad
protected

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

template<int dim, typename Number , bool is_face>
Number* FEEvaluationData< dim, Number, is_face >::gradients_quad
protected

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

template<int dim, typename Number , bool is_face>
Number* FEEvaluationData< dim, Number, is_face >::gradients_from_hessians_quad
protected

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

template<int dim, typename Number , bool is_face>
Number* FEEvaluationData< dim, Number, is_face >::hessians_quad
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::is_reinitialized
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::dof_values_initialized
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::values_quad_initialized
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::gradients_quad_initialized
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::hessians_quad_initialized
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::values_quad_submitted
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::gradients_quad_submitted
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::hessians_quad_submitted
protected

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

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::cell
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::interior_face
protected

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

template<int dim, typename Number , bool is_face>
internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex FEEvaluationData< dim, Number, is_face >::dof_access_index
protected

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

template<int dim, typename Number , bool is_face>
std::array<std::uint8_t, n_lanes> FEEvaluationData< dim, Number, is_face >::face_numbers
protected

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

template<int dim, typename Number , bool is_face>
std::array<std::uint8_t, n_lanes> FEEvaluationData< dim, Number, is_face >::face_orientations
protected

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

template<int dim, typename Number , bool is_face>
unsigned int FEEvaluationData< dim, Number, is_face >::subface_index
protected

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

template<int dim, typename Number , bool is_face>
internal::MatrixFreeFunctions::GeometryType FEEvaluationData< dim, Number, is_face >::cell_type
protected

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

template<int dim, typename Number , bool is_face>
std::array<unsigned int, n_lanes> FEEvaluationData< dim, Number, is_face >::cell_ids
protected

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

template<int dim, typename Number , bool is_face>
std::array<unsigned int, n_lanes> FEEvaluationData< dim, Number, is_face >::face_ids
protected

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

template<int dim, typename Number , bool is_face>
std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly<dim, Number> > FEEvaluationData< dim, Number, is_face >::mapped_geometry
protected

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

template<int dim, typename Number , bool is_face>
bool FEEvaluationData< dim, Number, is_face >::divergence_is_requested
protected

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 files: