|
void | set_data_pointers (AlignedVector< VectorizedArrayType > *scratch_data, const unsigned int n_components) |
|
void | reinit_face (const internal::MatrixFreeFunctions::FaceToCellTopology< n_lanes > &face) |
|
|
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 |
|
|
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 > | get_normal_vector (const unsigned int q_point) const |
|
|
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 () |
|
|
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 |
|
|
AlignedVector< VectorizedArrayType > * | scratch_data_array |
|
const MatrixFree< dim, Number, VectorizedArrayType > * | matrix_free |
|
std::vector< types::global_dof_index > | local_dof_indices |
|
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 |
|
const MatrixFree< dim, Number, VectorizedArrayType > & | get_matrix_free () const |
|
| FEEvaluationBase (const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free, const unsigned int dof_no, const unsigned int first_selected_component, const unsigned int quad_no, const unsigned int fe_degree, const unsigned int n_q_points, const bool is_interior_face, const unsigned int active_fe_index, const unsigned int active_quad_index, const unsigned int face_type) |
|
| FEEvaluationBase (const Mapping< dim > &mapping, const FiniteElement< dim > &fe, const Quadrature< 1 > &quadrature, const UpdateFlags update_flags, const unsigned int first_selected_component, const FEEvaluationData< dim, VectorizedArrayType, is_face > *other) |
|
| FEEvaluationBase (const FEEvaluationBase &other) |
|
FEEvaluationBase & | operator= (const FEEvaluationBase &other) |
|
| ~FEEvaluationBase () |
|
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 |
|
|
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, VectorizedArrayType > * | quadrature_points |
|
const Tensor< 2, dim, VectorizedArrayType > * | jacobian |
|
const Tensor< 1, dim *(dim+1)/2, Tensor< 1, dim, VectorizedArrayType > > * | jacobian_gradients |
|
const VectorizedArrayType * | J_value |
|
const Tensor< 1, dim, VectorizedArrayType > * | normal_vectors |
|
const Tensor< 1, dim, VectorizedArrayType > * | normal_x_jacobian |
|
const ScalarNumber * | quadrature_weights |
|
ArrayView< VectorizedArrayType > | scratch_data |
|
VectorizedArrayType * | values_dofs |
|
VectorizedArrayType * | values_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_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, VectorizedArrayType > > | mapped_geometry |
|
VectorizedArrayType | read_cell_data (const AlignedVector< VectorizedArrayType > &array) const |
|
std::array< T, Number::size()> | read_cell_data (const AlignedVector< std::array< T, Number::size()>> &array) const |
|
void | set_cell_data (AlignedVector< VectorizedArrayType > &array, const VectorizedArrayType &value) const |
|
void | set_cell_data (AlignedVector< std::array< T, Number::size()>> &array, const std::array< T, Number::size()> &value) const |
|
VectorizedArrayType | read_face_data (const AlignedVector< VectorizedArrayType > &array) const |
|
std::array< T, Number::size()> | read_face_data (const AlignedVector< std::array< T, Number::size()>> &array) const |
|
void | set_face_data (AlignedVector< VectorizedArrayType > &array, const VectorizedArrayType &value) const |
|
void | set_face_data (AlignedVector< std::array< T, Number::size()>> &array, const std::array< T, Number::size()> &value) const |
|
template<int dim, int n_components_, typename Number, bool is_face, typename VectorizedArrayType>
class FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >
This is the base class for the FEEvaluation classes. This class needs usually not be called in user code and does not have any public constructor. The usage is through the class FEEvaluation instead. It implements a reinit method that is used to set pointers so that operations on quadrature points can be performed quickly, access functions to vectors for the FEEvaluationBase::read_dof_values(), FEEvaluationBase::set_dof_values(), and FEEvaluationBase::distribute_local_to_global() functions, as well as methods to access values and gradients of finite element functions. It also inherits the geometry access functions provided by the class FEEvaluationData.
This class has five template arguments:
- Template Parameters
-
dim | Dimension in which this class is to be used |
n_components | Number 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) |
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 94 of file fe_evaluation.h.
template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::FEEvaluationBase |
( |
const MatrixFree< dim, Number, VectorizedArrayType > & |
matrix_free, |
|
|
const unsigned int |
dof_no, |
|
|
const unsigned int |
first_selected_component, |
|
|
const unsigned int |
quad_no, |
|
|
const unsigned int |
fe_degree, |
|
|
const unsigned int |
n_q_points, |
|
|
const bool |
is_interior_face, |
|
|
const unsigned int |
active_fe_index, |
|
|
const unsigned int |
active_quad_index, |
|
|
const unsigned int |
face_type |
|
) |
| |
|
protected |
Constructor. Made protected to prevent users from directly using this class. 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
, dof_no
, first_selected_component
and quad_no
allow to select the appropriate components.
template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
FEEvaluationBase< dim, n_components_, Number, is_face, VectorizedArrayType >::FEEvaluationBase |
( |
const Mapping< dim > & |
mapping, |
|
|
const FiniteElement< dim > & |
fe, |
|
|
const Quadrature< 1 > & |
quadrature, |
|
|
const UpdateFlags |
update_flags, |
|
|
const unsigned int |
first_selected_component, |
|
|
const FEEvaluationData< dim, VectorizedArrayType, is_face > * |
other |
|
) |
| |
|
protected |
Constructor that comes with reduced functionality and works similar as FEValues. The arguments are similar to the ones passed to the constructor of FEValues, with the notable difference that FEEvaluation expects a one-dimensional quadrature formula, Quadrature<1>, instead of a dim
dimensional one. The finite element can be both scalar or vector valued, but this method always only selects a scalar base element at a time (with n_components
copies as specified by the class template argument). For vector-valued elements, the optional argument first_selected_component
allows to specify the index of the base element to be used for evaluation. Note that the internal data structures always assume that the base element is primitive, non-primitive are not supported currently.
As known from FEValues, a call to the reinit method with a Triangulation::cell_iterator is necessary to make the geometry and degrees of freedom of the current class known. If the iterator includes DoFHandler information (i.e., it is a DoFHandler::cell_iterator or similar), the initialization allows to also read from or write to vectors in the standard way for DoFHandler::active_cell_iterator types for one cell at a time. However, this approach is much slower than the path with MatrixFree with MPI since index translation has to be done. As only one cell at a time is used, this method does not vectorize over several elements (which is most efficient for vector operations), but only possibly within the element if the evaluate/integrate routines are combined inside user code (e.g. for computing cell matrices).
The optional FEEvaluationData object allows several FEEvaluation objects to share the geometry evaluation, i.e., the underlying mapping and quadrature points do only need to be evaluated once. This only works if the quadrature formulas are the same. Otherwise, a new evaluation object is created. Make sure to not pass an optional object around when you intend to use the FEEvaluation object in parallel with another one because otherwise the intended sharing may create race conditions.
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() |
|
) |
| |
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).
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() |
|
) |
| |
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).
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 |
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).
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 |
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).
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 |
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 that 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).
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 |
|
) |
| |
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 that 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).
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 |
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 that 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).
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 |
|
) |
| |
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 that 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).
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 |
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() * get_normal_vector() but will use a more efficient internal representation of data.
Note that 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).
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 |
|
) |
| |
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 that 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).
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 |
|
) |
| |
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 that 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).
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 |
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 that 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).
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 |
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 that 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).
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 |
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 that 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).
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 |
|
) |
| |
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(...,true)
is called, this specifies what is tested by all basis function gradients on the current cell and integrated over.
- Note
- Only available for 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().
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 |
|
) |
| |
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(...,true)
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 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().
template<int dim, int n_components_, typename Number , bool is_face, typename VectorizedArrayType >
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.
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 |
|
protected |
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.
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 |
|
protected |
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.
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 |
|
protected |
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.
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.
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.
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 750 of file fe_evaluation_data.h.