|
| FEPointEvaluation (const Mapping< dim > &mapping, const FiniteElement< dim > &fe, const UpdateFlags update_flags, const unsigned int first_selected_component=0) |
|
| FEPointEvaluation (NonMatching::MappingInfo< dim, spacedim, Number > &mapping_info, const FiniteElement< dim > &fe, const unsigned int first_selected_component=0) |
|
| FEPointEvaluation (FEPointEvaluation &other) noexcept |
|
| FEPointEvaluation (FEPointEvaluation &&other) noexcept |
|
| ~FEPointEvaluation () |
|
void | reinit (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< dim >> &unit_points) |
|
void | reinit () |
|
void | reinit (const unsigned int cell_index) |
|
void | reinit (const unsigned int cell_index, const unsigned int face_number) |
|
template<std::size_t stride_view> |
void | evaluate (const StridedArrayView< const ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &evaluation_flags) |
|
void | evaluate (const ArrayView< const ScalarNumber > &solution_values, const EvaluationFlags::EvaluationFlags &evaluation_flags) |
|
template<std::size_t stride_view> |
void | integrate (const StridedArrayView< ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) |
|
void | integrate (const ArrayView< ScalarNumber > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) |
|
template<std::size_t stride_view> |
void | test_and_sum (const StridedArrayView< ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) |
|
void | test_and_sum (const ArrayView< ScalarNumber > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) |
|
const value_type & | get_value (const unsigned int point_index) const |
|
void | submit_value (const value_type &value, const unsigned int point_index) |
|
const gradient_type & | get_gradient (const unsigned int point_index) const |
|
void | submit_gradient (const gradient_type &, const unsigned int point_index) |
|
DerivativeForm< 1, dim, spacedim, Number > | jacobian (const unsigned int point_index) const |
|
DerivativeForm< 1, spacedim, dim, Number > | inverse_jacobian (const unsigned int point_index) const |
|
Number | JxW (const unsigned int point_index) const |
|
Tensor< 1, spacedim, Number > | normal_vector (const unsigned int point_index) const |
|
Point< spacedim, Number > | real_point (const unsigned int point_index) const |
|
Point< dim, Number > | unit_point (const unsigned int point_index) const |
|
std_cxx20::ranges::iota_view< unsigned int, unsigned int > | quadrature_point_indices () const |
|
|
void | setup (const unsigned int first_selected_component) |
|
template<bool is_face, bool is_linear> |
void | do_reinit () |
|
template<bool is_face_path, bool is_linear, std::size_t stride_view> |
void | prepare_evaluate_fast (const StridedArrayView< const ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &evaluation_flags) |
|
template<bool is_face_path, bool is_linear, std::size_t stride_view> |
void | compute_evaluate_fast (const StridedArrayView< const ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &evaluation_flags, const unsigned int n_shapes, const unsigned int qb, vectorized_value_type &value, interface_vectorized_gradient_type &gradient) |
|
template<bool is_face_path, bool is_linear, std::size_t stride_view> |
void | evaluate_fast (const StridedArrayView< const ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &evaluation_flags) |
|
template<std::size_t stride_view> |
void | evaluate_slow (const StridedArrayView< const ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &evaluation_flags) |
|
template<bool is_face_path, bool is_linear> |
void | compute_integrate_fast (const EvaluationFlags::EvaluationFlags &integration_flags, const unsigned int n_shapes, const unsigned int qb, const vectorized_value_type value, const interface_vectorized_gradient_type gradient, vectorized_value_type *solution_values_vectorized_linear) |
|
template<bool is_face_path, bool is_linear, std::size_t stride_view> |
void | finish_integrate_fast (const StridedArrayView< ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags, vectorized_value_type *solution_values_vectorized_linear) |
|
template<bool do_JxW, bool is_face_path, bool is_linear, std::size_t stride_view> |
void | integrate_fast (const StridedArrayView< ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) |
|
template<bool do_JxW, std::size_t stride_view> |
void | integrate_slow (const StridedArrayView< ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) |
|
template<bool do_JxW, std::size_t stride_view> |
void | do_integrate (const StridedArrayView< ScalarNumber, stride_view > &solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) |
|
|
const unsigned int | n_q_batches |
|
const unsigned int | n_q_points |
|
const unsigned int | n_q_points_scalar |
|
SmartPointer< const Mapping< dim, spacedim > > | mapping |
|
SmartPointer< const FiniteElement< dim > > | fe |
|
std::vector< Polynomials::Polynomial< double > > | poly |
|
bool | use_linear_path |
|
std::vector< unsigned int > | renumber |
|
std::vector< scalar_value_type > | solution_renumbered |
|
AlignedVector< vectorized_value_type > | solution_renumbered_vectorized |
|
AlignedVector< ScalarNumber > | scratch_data_scalar |
|
std::vector< value_type > | values |
|
std::vector< gradient_type > | gradients |
|
const Point< dim, VectorizedArrayType > * | unit_point_ptr |
|
const Point< dim - 1, VectorizedArrayType > * | unit_point_faces_ptr |
|
const Point< spacedim, Number > * | real_point_ptr |
|
const DerivativeForm< 1, dim, spacedim, Number > * | jacobian_ptr |
|
const DerivativeForm< 1, spacedim, dim, Number > * | inverse_jacobian_ptr |
|
const Tensor< 1, spacedim, Number > * | normal_ptr |
|
const Number * | JxW_ptr |
|
internal::MatrixFreeFunctions::GeometryType | cell_type |
|
unsigned int | dofs_per_component |
|
unsigned int | dofs_per_component_face |
|
bool | use_face_path |
|
internal::MatrixFreeFunctions::ShapeInfo< ScalarNumber > | shape_info |
|
unsigned int | component_in_base_element |
|
std::vector< std::array< bool, n_components > > | nonzero_shape_function_component |
|
const UpdateFlags | update_flags |
|
std::shared_ptr< FEValues< dim, spacedim > > | fe_values |
|
std::unique_ptr< NonMatching::MappingInfo< dim, spacedim, Number > > | mapping_info_on_the_fly |
|
SmartPointer< NonMatching::MappingInfo< dim, spacedim, Number > > | mapping_info |
|
unsigned int | current_cell_index |
|
unsigned int | current_face_number |
|
bool | fast_path |
|
boost::signals2::connection | connection_is_reinitialized |
|
bool | is_reinitialized |
|
AlignedVector<::ndarray< VectorizedArrayType, 2, dim > > | shapes |
|
AlignedVector<::ndarray< VectorizedArrayType, 2, dim - 1 > > | shapes_faces |
|
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
class FEPointEvaluation< n_components_, dim, spacedim, Number >
This class provides an interface to the evaluation of interpolated solution values and gradients on cells on arbitrary reference point positions. These points can change from cell to cell, both with respect to their quantity as well to the location. The two typical use cases are evaluations on non-matching grids and particle simulations.
The use of this class is similar to FEValues or FEEvaluation: The class is first initialized to a cell by calling FEPointEvaluation::reinit(cell, unit_points)
, with the main difference to the other concepts that the underlying points in reference coordinates need to be passed along. Then, upon call to evaluate() or integrate(), the user can compute information at the give points. Eventually, the access functions get_value() or get_gradient() allow to query this information at a specific point index.
The functionality is similar to creating an FEValues object with a Quadrature object on the unit_points
on every cell separately and then calling FEValues::get_function_values or FEValues::get_function_gradients, and for some elements and mappings this is what actually happens internally. For specific combinations of Mapping and FiniteElement realizations, however, there is a much more efficient implementation that avoids the memory allocation and other expensive start-up cost of FEValues. Currently, the functionality is specialized for mappings derived from MappingQ and MappingCartesian and for finite elements with tensor product structure that work with the Matrix-free infrastructure module. In those cases, the cost implied by this class is similar (or sometimes even somewhat lower) than using FEValues::reinit(cell)
followed by FEValues::get_function_gradients
.
Definition at line 718 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
Set up the mapping information for the given cell, e.g., by computing the Jacobian of the mapping for the given points if gradients of the functions are requested.
- Parameters
-
[in] | cell | An iterator to the current cell |
[in] | unit_points | List of points in the reference locations of the current cell where the FiniteElement object should be evaluated/integrated in the evaluate() and integrate() functions. |
Definition at line 1657 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<std::size_t stride_view>
This function interpolates the finite element solution, represented by solution_values
, on the cell and unit_points
passed to reinit().
- Parameters
-
[in] | solution_values | This array is supposed to contain the unknown values on the element read out by FEEvaluation::read_dof_values(global_vector) . |
[in] | evaluation_flags | Flags specifying which quantities should be evaluated at the points. |
Definition at line 2233 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
This function interpolates the finite element solution, represented by solution_values
, on the cell and unit_points
passed to reinit().
- Parameters
-
[in] | solution_values | This array is supposed to contain the unknown values on the element as returned by cell->get_dof_values(global_vector, solution_values) . |
[in] | evaluation_flags | Flags specifying which quantities should be evaluated at the points. |
Definition at line 863 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<std::size_t stride_view>
This function multiplies the quantities passed in by previous submit_value() or submit_gradient() calls by the value or gradient of the test functions, and performs summation over all given points multiplied be the Jacobian determinant times the quadrature weight (JxW).
- Parameters
-
[out] | solution_values | This array will contain the result of the integral, which can be used during FEEvaluation::set_dof_values(global_vector) or FEEvaluation::distribute_local_to_global(global_vector) . Note that for multi-component systems where only some of the components are selected by the present class, the entries in solution_values not touched by this class will be set to zero. |
[in] | integration_flags | Flags specifying which quantities should be integrated at the points. |
Definition at line 2690 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
This function multiplies the quantities passed in by previous submit_value() or submit_gradient() calls by the value or gradient of the test functions, and performs summation over all given points multiplied be the Jacobian determinant times the quadrature weight (JxW).
- Parameters
-
[out] | solution_values | This array will contain the result of the integral, which can be used to during cell->set_dof_values(solution_values, global_vector) or cell->distribute_local_to_global(solution_values, global_vector) . Note that for multi-component systems where only some of the components are selected by the present class, the entries in solution_values not touched by this class will be set to zero. |
[in] | integration_flags | Flags specifying which quantities should be integrated at the points. |
Definition at line 914 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<std::size_t stride_view>
This function multiplies the quantities passed in by previous submit_value() or submit_gradient() calls by the value or gradient of the test functions, and performs summation over all given points. This is similar to the integration of a bilinear form in terms of the test function, with the difference that this formula does not include a JxW
factor (in contrast to the integrate function of this class). This allows the class to naturally embed point information (e.g. particles) into a finite element formulation.
- Parameters
-
[out] | solution_values | This array will contain the result of the integral, which can be used during FEEvaluation::set_dof_values(global_vector) or FEEvaluation::distribute_local_to_global(global_vector) . Note that for multi-component systems where only some of the components are selected by the present class, the entries in solution_values not touched by this class will be set to zero. |
[in] | integration_flags | Flags specifying which quantities should be integrated at the points. |
Definition at line 2702 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
This function multiplies the quantities passed in by previous submit_value() or submit_gradient() calls by the value or gradient of the test functions, and performs summation over all given points. This is similar to the integration of a bilinear form in terms of the test function, with the difference that this formula does not include a JxW
factor (in contrast to the integrate function of this class). This allows the class to naturally embed point information (e.g. particles) into a finite element formulation.
- Parameters
-
[out] | solution_values | This array will contain the result of the integral, which can be used during cell->set_dof_values(solution_values, global_vector) or cell->distribute_local_to_global(solution_values, global_vector) . Note that for multi-component systems where only some of the components are selected by the present class, the entries in solution_values not touched by this class will be set to zero. |
[in] | integration_flags | Flags specifying which quantities should be integrated at the points. |
Definition at line 973 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<int n_components_, int dim, int spacedim, typename Number >
Write a contribution that is tested by the gradient to the field containing the values on points with the given point_index
. Access to the same field as through get_gradient(). If applied before the function FEPointEvaluation::integrate(EvaluationFlags::gradients) is called, this specifies what is tested by all basis function gradients on the current cell and integrated over.
Definition at line 2749 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
DerivativeForm< 1, spacedim, dim, Number > FEPointEvaluation< n_components_, dim, spacedim, Number >::inverse_jacobian |
( |
const unsigned int |
point_index | ) |
const |
|
inline |
Return the inverse of the Jacobian of the transformation on the current cell with the given point index. Prerequisite: This class needs to be constructed with UpdateFlags containing update_inverse_jacobian
or update_gradients
.
Definition at line 2779 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
Number FEPointEvaluation< n_components_, dim, spacedim, Number >::JxW |
( |
const unsigned int |
point_index | ) |
const |
|
inline |
Return the Jacobian determinant multiplied by the quadrature weight. This class or the MappingInfo object passed to this function needs to be constructed with UpdateFlags containing update_JxW_values
.
Definition at line 2798 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
Point< dim, Number > FEPointEvaluation< n_components_, dim, spacedim, Number >::unit_point |
( |
const unsigned int |
point_index | ) |
const |
|
inline |
Return the position in unit/reference coordinates of the given point index, i.e., the respective point passed to the reinit() function.
Definition at line 2843 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<bool is_face, bool is_linear>
Shared functionality of all reinit()
functions. Resizes data fields and precomputes the shapes
vector, holding the evaluation of 1D basis functions of tensor product polynomials, if necessary.
Definition at line 1755 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<bool is_face_path, bool is_linear, std::size_t stride_view>
Resizes necessary data fields, reads in and renumbers solution values. Interpolates onto face if face path is selected.
Definition at line 1868 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<bool is_face_path, bool is_linear, std::size_t stride_view>
Addition across the lanes of VectorizedArray as accumulated by the compute_integrate_fast_function(), writing the sum into the result vector. Applies face contributions to cell contributions for face path.
Definition at line 2393 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim, typename Number >
template<bool do_JxW, bool is_face_path, bool is_linear, std::size_t stride_view>
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
Temporary array to store the solution_values
passed to the evaluate() function in a format compatible with the tensor product evaluators. For vector-valued setups, this array uses a Tensor<1, n_components, ScalarNumber>
type to collect the unknowns for a particular basis function.
Definition at line 1256 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
Temporary array to store a vectorized version of the solution_values
computed during integrate()
in a format compatible with the tensor product evaluators. For vector-valued setups, this array uses a Tensor<1, n_components, VectorizedArrayType>
format.
Definition at line 1264 of file fe_point_evaluation.h.
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
unsigned int FEPointEvaluation< n_components_, dim, spacedim, Number >::dofs_per_component |
|
private |
template<int n_components_, int dim, int spacedim = dim, typename Number = double>
unsigned int FEPointEvaluation< n_components_, dim, spacedim, Number >::dofs_per_component_face |
|
private |
Number of unknowns per component, i.e., number of unique basis functions, for a restriction to the face of the chosen FiniteElement (or base element). This means a (dim-1)-dimensional basis.
Definition at line 1340 of file fe_point_evaluation.h.