|
| FEEvaluationData (const ShapeInfoType &shape_info, const bool is_interior_face=true) |
|
| FEEvaluationData (const FEEvaluationData &other)=default |
|
FEEvaluationData & | operator= (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) |
|
|
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 |
|
|
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 () |
|
|
unsigned int | get_mapping_data_index_offset () const |
|
internal::MatrixFreeFunctions::GeometryType | get_cell_type () const |
|
const ShapeInfoType & | get_shape_info () const |
|
const internal::MatrixFreeFunctions::DoFInfo & | get_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 int > | quadrature_point_indices () const |
|
|
template<typename T > |
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 > |
T | read_face_data (const AlignedVector< T > &array) const |
|
template<typename T > |
void | set_face_data (AlignedVector< T > &array, const T &value) const |
|
|
const ShapeInfoType * | data |
|
const DoFInfo * | dof_info |
|
const MappingInfoStorageType * | mapping_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::QuadratureDescriptor * | descriptor |
|
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 ScalarNumber * | quadrature_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_lanes > | face_numbers |
|
std::array< std::uint8_t, n_lanes > | face_orientations |
|
unsigned int | subface_index |
|
internal::MatrixFreeFunctions::GeometryType | cell_type |
|
std::array< unsigned int, n_lanes > | cell_ids |
|
std::array< unsigned int, n_lanes > | face_ids |
|
std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly< dim, Number > > | mapped_geometry |
|
bool | divergence_is_requested |
|
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
-
dim | Dimension in which this class is to be used |
Number | Number format, usually double or float |
is_face | Whether 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) |
VectorizedArrayType | Type 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.
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\).
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
.
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
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.
template<int dim, typename Number , bool is_face>
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.
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.
template<int dim, typename Number , bool is_face>
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.
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.
template<int dim, typename Number , bool is_face>
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.