deal.II version GIT relicensing-2165-gc91f007519 2024-11-20 01:40:00+00:00
|
#include <deal.II/matrix_free/fe_point_evaluation.h>
Public Types | |
using | number_type = Number |
using | ScalarNumber = typename internal::VectorizedArrayTrait< Number >::value_type |
using | VectorizedArrayType = typename ::internal::VectorizedArrayTrait< Number >::vectorized_value_type |
using | ETT = typename internal::FEPointEvaluation::EvaluatorTypeTraits< dim, spacedim, n_components, Number > |
using | value_type = typename ETT::value_type |
using | scalar_value_type = typename ETT::scalar_value_type |
using | vectorized_value_type = typename ETT::vectorized_value_type |
using | gradient_type = typename ETT::real_gradient_type |
using | interface_vectorized_unit_gradient_type = typename ETT::interface_vectorized_unit_gradient_type |
Public Member Functions | |
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) |
Number | get_divergence (const unsigned int point_index) const |
void | submit_divergence (const Number &value, const unsigned int point_index) |
Tensor< 1,(dim==2 ? 1 :dim), Number > | get_curl (const unsigned int point_index) const |
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 |
Point< spacedim, Number > | real_point (const unsigned int point_index) const |
Point< spacedim, Number > | quadrature_point (const unsigned int point_index) const |
Point< dim, Number > | unit_point (const unsigned int point_index) const |
scalar_value_type | integrate_value () const |
std_cxx20::ranges::iota_view< unsigned int, unsigned int > | quadrature_point_indices () const |
unsigned int | n_active_entries_per_quadrature_batch (unsigned int q) |
Static Public Attributes | |
static constexpr unsigned int | dimension = dim |
static constexpr unsigned int | n_components = n_components_ |
Protected Member Functions | |
FEPointEvaluationBase (const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const UpdateFlags update_flags, const unsigned int first_selected_component=0) | |
FEPointEvaluationBase (NonMatching::MappingInfo< dim, spacedim, Number > &mapping_info, const FiniteElement< dim, spacedim > &fe, const unsigned int first_selected_component=0, const bool is_interior=true) | |
FEPointEvaluationBase (FEPointEvaluationBase &other) noexcept | |
FEPointEvaluationBase (FEPointEvaluationBase &&other) noexcept | |
void | setup (const unsigned int first_selected_component) |
template<bool is_face, bool is_linear> | |
void | do_reinit () |
Static Protected Attributes | |
static constexpr std::size_t | n_lanes_user_interface |
static constexpr std::size_t | n_lanes_internal |
static constexpr std::size_t | stride |
Base class of FEPointEvaluation and FEFacePointEvaluation. This class needs usually not be called in user code and does not have any public constructor. The usage is through the class FEPointEvaluation/FEFacePointEvaluation instead.
Definition at line 623 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::number_type = Number |
Definition at line 629 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::ScalarNumber = typename internal::VectorizedArrayTrait<Number>::value_type |
Definition at line 631 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::VectorizedArrayType = typename ::internal::VectorizedArrayTrait< Number>::vectorized_value_type |
Definition at line 633 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::ETT = typename internal::FEPointEvaluation:: EvaluatorTypeTraits<dim, spacedim, n_components, Number> |
Definition at line 635 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::value_type = typename ETT::value_type |
Definition at line 637 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::scalar_value_type = typename ETT::scalar_value_type |
Definition at line 638 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::vectorized_value_type = typename ETT::vectorized_value_type |
Definition at line 639 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::gradient_type = typename ETT::real_gradient_type |
Definition at line 640 of file fe_point_evaluation.h.
using FEPointEvaluationBase< n_components_, dim, spacedim, Number >::interface_vectorized_unit_gradient_type = typename ETT::interface_vectorized_unit_gradient_type |
Definition at line 641 of file fe_point_evaluation.h.
|
protected |
Constructor.
mapping | The Mapping class describing the actual geometry of a cell passed to the evaluate() function. |
fe | The FiniteElement object that is used for the evaluation, which is typically the same on all cells to be evaluated. |
update_flags | Specify the quantities to be computed by the mapping during the call of reinit(). During evaluate() or integrate(), this data is queried to produce the desired result (e.g., the gradient of a finite element solution). |
first_selected_component | For multi-component FiniteElement objects, this parameter allows to select a range of n_components components starting from this parameter. |
Definition at line 1827 of file fe_point_evaluation.h.
|
protected |
Constructor to make the present class able to re-use the geometry data also used by other FEPointEvaluationBase
objects.
mapping_info | The MappingInfo class describes the geometry-related data for evaluating finite-element solutions. This object enables to construct such an object on the outside, possibly re-using it between several objects or between several calls to the same cell and unit points. |
fe | The FiniteElement object that is used for the evaluation, which is typically the same on all cells to be evaluated. |
first_selected_component | For multi-component FiniteElement objects, this parameter allows to select a range of n_components components starting from this parameter. |
is_interior | Defines if interior or exterior. Only makes sense for faces. |
Definition at line 1855 of file fe_point_evaluation.h.
|
protectednoexcept |
Copy constructor.
Definition at line 1880 of file fe_point_evaluation.h.
|
protectednoexcept |
Move constructor.
Definition at line 1918 of file fe_point_evaluation.h.
|
inline |
Return the value at quadrature point number point_index
after a call to FEPointEvaluation::evaluate() with EvaluationFlags::values set, or the value that has been stored there with a call to FEPointEvaluationBase::submit_value(). If the object is vector-valued, a vector-valued return argument is given.
Definition at line 2167 of file fe_point_evaluation.h.
|
inline |
Write a value to the field containing the values on points with component point_index. Access to the same field as through get_value(). If applied before the function FEPointEvaluation::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.
Definition at line 2206 of file fe_point_evaluation.h.
|
inline |
Return the gradient in real coordinates at the point with index point_index
after a call to FEPointEvaluation::evaluate() with EvaluationFlags::gradients set, or the gradient that has been stored there with a call to FEPointEvaluationBase::submit_gradient(). The gradient in real coordinates is obtained by taking the unit gradient (also accessible via get_unit_gradient()) and applying the inverse Jacobian of the mapping. If the object is vector-valued, a vector-valued return argument is given.
Definition at line 2181 of file fe_point_evaluation.h.
|
inline |
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 2218 of file fe_point_evaluation.h.
|
inline |
Return the divergence in real coordinates at the point with index point_index
after a call to FEPointEvaluation::evaluate() with EvaluationFlags::gradients set, or the divergence that has been stored there with a call to FEPointEvaluationBase::submit_divergence(). This functions only makes sense for a vector field with dim components.
Definition at line 2192 of file fe_point_evaluation.h.
|
inline |
Write a contribution that is tested by the divergence to the field containing the values on points with the given point_index
. Access to the same field as through get_divergence(). 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 2230 of file fe_point_evaluation.h.
Tensor< 1,(dim==2 ? 1 :dim), Number > FEPointEvaluationBase< n_components_, dim, spacedim, Number >::get_curl | ( | const unsigned int | point_index | ) | const |
Return the curl in real coordinates at the point with index point_index
after a call to FEPointEvaluation::evaluate() with EvaluationFlags::gradients set. This functions only makes sense for a vector field with dim components and dim > 1.
Definition at line 2247 of file fe_point_evaluation.h.
|
inline |
Return 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_jacobian
.
Definition at line 2276 of file fe_point_evaluation.h.
|
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 2294 of file fe_point_evaluation.h.
|
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 2313 of file fe_point_evaluation.h.
|
inline |
Return the position in real coordinates of the given point index among the points passed to reinit().
Definition at line 2328 of file fe_point_evaluation.h.
|
inline |
Return the position in real coordinates of the given point index among the points passed to reinit().
Definition at line 2339 of file fe_point_evaluation.h.
|
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 2354 of file fe_point_evaluation.h.
|
inline |
 * Take values collected at quadrature points via the submit_value() function, multiply by the Jacobian determinant and quadrature weights (JxW) and sum the values for all quadrature  * points on the cell. The result is a scalar, representing the integral  * of the function over the cell. 
Definition at line 2578 of file fe_point_evaluation.h.
|
inline |
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.
Definition at line 2370 of file fe_point_evaluation.h.
unsigned int FEPointEvaluationBase< n_components_, dim, spacedim, Number >::n_active_entries_per_quadrature_batch | ( | unsigned int | q | ) |
Returns how many lanes of a quadrature batch are active.
Definition at line 2593 of file fe_point_evaluation.h.
|
protected |
Common setup function for both constructors. Does the setup for both fast and slow path.
first_selected_component | For multi-component FiniteElement objects, this parameter allows to select a range of n_components components starting from this parameter. |
Definition at line 1951 of file fe_point_evaluation.h.
|
inlineprotected |
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 2038 of file fe_point_evaluation.h.
|
staticconstexpr |
Definition at line 626 of file fe_point_evaluation.h.
|
staticconstexpr |
Definition at line 627 of file fe_point_evaluation.h.
|
staticconstexprprotected |
Definition at line 858 of file fe_point_evaluation.h.
|
staticconstexprprotected |
Definition at line 860 of file fe_point_evaluation.h.
|
staticconstexprprotected |
Definition at line 862 of file fe_point_evaluation.h.
|
protected |
Number of quadrature batches of the current cell/face.
Definition at line 888 of file fe_point_evaluation.h.
|
protected |
Number of quadrature points/batches of the current cell/face.
Definition at line 893 of file fe_point_evaluation.h.
|
protected |
Number of quadrature points of the current cell/face.
Definition at line 898 of file fe_point_evaluation.h.
|
protected |
Pointer to the Mapping object passed to the constructor.
Definition at line 903 of file fe_point_evaluation.h.
|
protected |
Pointer to the FiniteElement object passed to the constructor.
Definition at line 908 of file fe_point_evaluation.h.
|
protected |
Description of the 1d polynomial basis for tensor product elements used for the fast path of this class using tensor product evaluators.
Definition at line 914 of file fe_point_evaluation.h.
|
protected |
Store whether the linear path should be used.
Definition at line 919 of file fe_point_evaluation.h.
|
protected |
Renumbering between the unknowns of unknowns implied by the FiniteElement class and a lexicographic numbering used for the tensorized code path.
Definition at line 925 of file fe_point_evaluation.h.
|
protected |
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 934 of file fe_point_evaluation.h.
|
protected |
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 942 of file fe_point_evaluation.h.
|
protected |
Temporary array for the face path (scalar).
Definition at line 947 of file fe_point_evaluation.h.
|
protected |
Temporary array to store the values at the points.
Definition at line 952 of file fe_point_evaluation.h.
|
protected |
Temporary array to store the gradients in real coordinates at the points.
Definition at line 957 of file fe_point_evaluation.h.
|
protected |
Pointer to first unit point batch of current cell/face from MappingInfo, set internally during do_reinit().
Definition at line 963 of file fe_point_evaluation.h.
|
protected |
Pointer to first unit point batch of current face from MappingInfo, set internally during do_reinit(). Needed for face path.
Definition at line 969 of file fe_point_evaluation.h.
|
protected |
Pointer to real point of first quadrature point of current cell/face from MappingInfo, set internally during do_reinit().
Definition at line 975 of file fe_point_evaluation.h.
|
protected |
Pointer to Jacobian of first quadrature point of current cell/face from MappingInfo, set internally during do_reinit().
Definition at line 981 of file fe_point_evaluation.h.
|
protected |
Pointer to inverse Jacobian of first quadrature point of current cell/face from MappingInfo, set internally during do_reinit().
Definition at line 987 of file fe_point_evaluation.h.
|
protected |
Pointer to normal vector of first quadrature point of current cell/face from MappingInfo, set internally during do_reinit().
Definition at line 993 of file fe_point_evaluation.h.
|
protected |
Pointer to Jacobian determinant times quadrature weight of first quadrature point of current cell/face from MappingInfo, set internally during do_reinit().
Definition at line 1000 of file fe_point_evaluation.h.
|
protected |
Cell type describing the geometry of the cell and compression of jacobians.
Definition at line 1005 of file fe_point_evaluation.h.
|
protected |
Number of unknowns per component, i.e., number of unique basis functions, for the chosen FiniteElement (or base element).
Definition at line 1011 of file fe_point_evaluation.h.
|
protected |
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 1018 of file fe_point_evaluation.h.
|
protected |
Scalar ShapeInfo object needed for face path.
Definition at line 1023 of file fe_point_evaluation.h.
|
protected |
The first selected component in the active base element.
Definition at line 1028 of file fe_point_evaluation.h.
|
protected |
For complicated FiniteElement objects this variable informs us about which unknowns actually carry degrees of freedom in the selected components.
Definition at line 1035 of file fe_point_evaluation.h.
|
protected |
The desired update flags for the evaluation.
Definition at line 1040 of file fe_point_evaluation.h.
|
protected |
The FEValues object underlying the slow evaluation path.
Definition at line 1045 of file fe_point_evaluation.h.
|
protected |
Pointer to mapping info on the fly computed during reinit.
Definition at line 1051 of file fe_point_evaluation.h.
|
protected |
Pointer to currently used mapping info (either on the fly or external precomputed).
Definition at line 1057 of file fe_point_evaluation.h.
|
protected |
The current cell index to access mapping data from mapping info.
Definition at line 1062 of file fe_point_evaluation.h.
|
protected |
The current face number to access mapping data from mapping info.
Definition at line 1067 of file fe_point_evaluation.h.
|
protected |
Bool indicating if fast path is chosen.
Definition at line 1072 of file fe_point_evaluation.h.
|
protected |
Bool indicating if class needs to call reinit() inside evaluate()/integrate() functions, which is the case when the present class does not own the MappingInfo object but shares evaluation points with another object.
Definition at line 1080 of file fe_point_evaluation.h.
|
protected |
Vector containing tensor product shape functions evaluated (during reinit()) at the vectorized unit points.
Definition at line 1086 of file fe_point_evaluation.h.
|
protected |
Vector containing tensor product shape functions evaluated (during reinit()) at the vectorized unit points on faces.
Definition at line 1092 of file fe_point_evaluation.h.
|
protected |
Definition at line 1094 of file fe_point_evaluation.h.