Reference documentation for deal.II version 9.0.0
|
#include <deal.II/fe/fe_values.h>
Classes | |
class | CellIterator |
class | TriaCellIterator |
Public Member Functions | |
FEValuesBase (const unsigned int n_q_points, const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe) | |
~FEValuesBase () | |
ShapeAccess Access to shape function values. These fields are filled by the finite element. | |
const double & | shape_value (const unsigned int function_no, const unsigned int point_no) const |
double | shape_value_component (const unsigned int function_no, const unsigned int point_no, const unsigned int component) const |
const Tensor< 1, spacedim > & | shape_grad (const unsigned int function_no, const unsigned int quadrature_point) const |
Tensor< 1, spacedim > | shape_grad_component (const unsigned int function_no, const unsigned int point_no, const unsigned int component) const |
const Tensor< 2, spacedim > & | shape_hessian (const unsigned int function_no, const unsigned int point_no) const |
Tensor< 2, spacedim > | shape_hessian_component (const unsigned int function_no, const unsigned int point_no, const unsigned int component) const |
const Tensor< 3, spacedim > & | shape_3rd_derivative (const unsigned int function_no, const unsigned int point_no) const |
Tensor< 3, spacedim > | shape_3rd_derivative_component (const unsigned int function_no, const unsigned int point_no, const unsigned int component) const |
Access to values of global finite element fields | |
template<class InputVector > | |
void | get_function_values (const InputVector &fe_function, std::vector< typename InputVector::value_type > &values) const |
template<class InputVector > | |
void | get_function_values (const InputVector &fe_function, std::vector< Vector< typename InputVector::value_type > > &values) const |
template<class InputVector > | |
void | get_function_values (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< typename InputVector::value_type > &values) const |
template<class InputVector > | |
void | get_function_values (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< Vector< typename InputVector::value_type > > &values) const |
template<class InputVector > | |
void | get_function_values (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, VectorSlice< std::vector< std::vector< typename InputVector::value_type > > > values, const bool quadrature_points_fastest) const |
Access to derivatives of global finite element fields | |
template<class InputVector > | |
void | get_function_gradients (const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > &gradients) const |
template<class InputVector > | |
void | get_function_gradients (const InputVector &fe_function, std::vector< std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > > &gradients) const |
template<class InputVector > | |
void | get_function_gradients (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > &gradients) const |
template<class InputVector > | |
void | get_function_gradients (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, VectorSlice< std::vector< std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > > > gradients, bool quadrature_points_fastest=false) const |
Access to second derivatives (Hessian matrices and Laplacians) of global finite element fields | |
template<class InputVector > | |
void | get_function_hessians (const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const |
template<class InputVector > | |
void | get_function_hessians (const InputVector &fe_function, std::vector< std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > > &hessians, bool quadrature_points_fastest=false) const |
template<class InputVector > | |
void | get_function_hessians (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const |
template<class InputVector > | |
void | get_function_hessians (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, VectorSlice< std::vector< std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > > > hessians, bool quadrature_points_fastest=false) const |
template<class InputVector > | |
void | get_function_laplacians (const InputVector &fe_function, std::vector< typename InputVector::value_type > &laplacians) const |
template<class InputVector > | |
void | get_function_laplacians (const InputVector &fe_function, std::vector< Vector< typename InputVector::value_type > > &laplacians) const |
template<class InputVector > | |
void | get_function_laplacians (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< typename InputVector::value_type > &laplacians) const |
template<class InputVector > | |
void | get_function_laplacians (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< Vector< typename InputVector::value_type > > &laplacians) const |
template<class InputVector > | |
void | get_function_laplacians (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< std::vector< typename InputVector::value_type > > &laplacians, bool quadrature_points_fastest=false) const |
Access to third derivatives of global finite element fields | |
template<class InputVector > | |
void | get_function_third_derivatives (const InputVector &fe_function, std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > &third_derivatives) const |
template<class InputVector > | |
void | get_function_third_derivatives (const InputVector &fe_function, std::vector< std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > > &third_derivatives, bool quadrature_points_fastest=false) const |
template<class InputVector > | |
void | get_function_third_derivatives (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > &third_derivatives) const |
template<class InputVector > | |
void | get_function_third_derivatives (const InputVector &fe_function, const VectorSlice< const std::vector< types::global_dof_index > > &indices, VectorSlice< std::vector< std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > > > third_derivatives, bool quadrature_points_fastest=false) const |
Geometry of the cell | |
const Point< spacedim > & | quadrature_point (const unsigned int q) const |
const std::vector< Point< spacedim > > & | get_quadrature_points () const |
double | JxW (const unsigned int quadrature_point) const |
const std::vector< double > & | get_JxW_values () const |
const DerivativeForm< 1, dim, spacedim > & | jacobian (const unsigned int quadrature_point) const |
const std::vector< DerivativeForm< 1, dim, spacedim > > & | get_jacobians () const |
const DerivativeForm< 2, dim, spacedim > & | jacobian_grad (const unsigned int quadrature_point) const |
const std::vector< DerivativeForm< 2, dim, spacedim > > & | get_jacobian_grads () const |
const Tensor< 3, spacedim > & | jacobian_pushed_forward_grad (const unsigned int quadrature_point) const |
const std::vector< Tensor< 3, spacedim > > & | get_jacobian_pushed_forward_grads () const |
const DerivativeForm< 3, dim, spacedim > & | jacobian_2nd_derivative (const unsigned int quadrature_point) const |
const std::vector< DerivativeForm< 3, dim, spacedim > > & | get_jacobian_2nd_derivatives () const |
const Tensor< 4, spacedim > & | jacobian_pushed_forward_2nd_derivative (const unsigned int quadrature_point) const |
const std::vector< Tensor< 4, spacedim > > & | get_jacobian_pushed_forward_2nd_derivatives () const |
const DerivativeForm< 4, dim, spacedim > & | jacobian_3rd_derivative (const unsigned int quadrature_point) const |
const std::vector< DerivativeForm< 4, dim, spacedim > > & | get_jacobian_3rd_derivatives () const |
const Tensor< 5, spacedim > & | jacobian_pushed_forward_3rd_derivative (const unsigned int quadrature_point) const |
const std::vector< Tensor< 5, spacedim > > & | get_jacobian_pushed_forward_3rd_derivatives () const |
const DerivativeForm< 1, spacedim, dim > & | inverse_jacobian (const unsigned int quadrature_point) const |
const std::vector< DerivativeForm< 1, spacedim, dim > > & | get_inverse_jacobians () const |
const Tensor< 1, spacedim > & | normal_vector (const unsigned int i) const |
const std::vector< Tensor< 1, spacedim > > & | get_all_normal_vectors () const |
const std::vector< Tensor< 1, spacedim > > & | get_normal_vectors () const |
Extractors Methods to extract individual components | |
const FEValuesViews::Scalar< dim, spacedim > & | operator[] (const FEValuesExtractors::Scalar &scalar) const |
const FEValuesViews::Vector< dim, spacedim > & | operator[] (const FEValuesExtractors::Vector &vector) const |
const FEValuesViews::SymmetricTensor< 2, dim, spacedim > & | operator[] (const FEValuesExtractors::SymmetricTensor< 2 > &tensor) const |
const FEValuesViews::Tensor< 2, dim, spacedim > & | operator[] (const FEValuesExtractors::Tensor< 2 > &tensor) const |
Access to the raw data | |
const Mapping< dim, spacedim > & | get_mapping () const |
const FiniteElement< dim, spacedim > & | get_fe () const |
UpdateFlags | get_update_flags () const |
const Triangulation< dim, spacedim >::cell_iterator | get_cell () const |
CellSimilarity::Similarity | get_cell_similarity () const |
std::size_t | memory_consumption () const |
Public Member Functions inherited from Subscriptor | |
Subscriptor () | |
Subscriptor (const Subscriptor &) | |
Subscriptor (Subscriptor &&) noexcept | |
virtual | ~Subscriptor () |
Subscriptor & | operator= (const Subscriptor &) |
Subscriptor & | operator= (Subscriptor &&) noexcept |
void | subscribe (const char *identifier=nullptr) const |
void | unsubscribe (const char *identifier=nullptr) const |
unsigned int | n_subscriptions () const |
void | list_subscribers () const |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
Static Public Member Functions | |
static ::ExceptionBase & | ExcAccessToUninitializedField (char *arg1) |
static ::ExceptionBase & | ExcFEDontMatch () |
static ::ExceptionBase & | ExcShapeFunctionNotPrimitive (int arg1) |
static ::ExceptionBase & | ExcFENotPrimitive () |
Static Public Member Functions inherited from Subscriptor | |
static ::ExceptionBase & | ExcInUse (int arg1, std::string arg2, std::string arg3) |
static ::ExceptionBase & | ExcNoSubscriber (std::string arg1, std::string arg2) |
Public Attributes | |
const unsigned int | n_quadrature_points |
const unsigned int | dofs_per_cell |
Static Public Attributes | |
static const unsigned int | dimension = dim |
static const unsigned int | space_dimension = spacedim |
Protected Member Functions | |
void | invalidate_present_cell () |
void | maybe_invalidate_previous_present_cell (const typename Triangulation< dim, spacedim >::cell_iterator &cell) |
UpdateFlags | compute_update_flags (const UpdateFlags update_flags) const |
void | check_cell_similarity (const typename Triangulation< dim, spacedim >::cell_iterator &cell) |
Protected Attributes | |
std::unique_ptr< const CellIteratorBase > | present_cell |
boost::signals2::connection | tria_listener_refinement |
boost::signals2::connection | tria_listener_mesh_transform |
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > | mapping |
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > | mapping_data |
::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > | mapping_output |
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > | fe |
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > | fe_data |
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > | finite_element_output |
UpdateFlags | update_flags |
CellSimilarity::Similarity | cell_similarity |
Private Member Functions | |
FEValuesBase (const FEValuesBase &) | |
FEValuesBase & | operator= (const FEValuesBase &) |
Private Attributes | |
::internal::FEValuesViews::Cache< dim, spacedim > | fe_values_views_cache |
Friends | |
template<int , int > | |
class | FEValuesViews::Scalar |
template<int , int > | |
class | FEValuesViews::Vector |
FEValues, FEFaceValues and FESubfaceValues objects are interfaces to finite element and mapping classes on the one hand side, to cells and quadrature rules on the other side. They allow to evaluate values or derivatives of shape functions at the quadrature points of a quadrature formula when projected by a mapping from the unit cell onto a cell in real space. The reason for this abstraction is possible optimization: Depending on the type of finite element and mapping, some values can be computed once on the unit cell. Others must be computed on each cell, but maybe computation of several values at the same time offers ways for optimization. Since this interplay may be complex and depends on the actual finite element, it cannot be left to the applications programmer.
FEValues, FEFaceValues and FESubfaceValues provide only data handling: computations are left to objects of type Mapping and FiniteElement. These provide functions get_*_data
and fill_*_values
which are called by the constructor and reinit
functions of FEValues*
, respectively.
Usually, an object of FEValues*
is used in integration loops over all cells of a triangulation (or faces of cells). To take full advantage of the optimization features, it should be constructed before the loop so that information that does not depend on the location and shape of cells can be computed once and for all (this includes, for example, the values of shape functions at quadrature points for the most common elements: we can evaluate them on the unit cell and they will be the same when mapped to the real cell). Then, in the loop over all cells, it must be re-initialized for each grid cell to compute that part of the information that changes depending on the actual cell (for example, the gradient of shape functions equals the gradient on the unit cell – which can be computed once and for all – times the Jacobian matrix of the mapping between unit and real cell, which needs to be recomputed for each cell).
A typical piece of code, adding up local contributions to the Laplace matrix looks like this:
The individual functions used here are described below. Note that by design, the order of quadrature points used inside the FEValues object is the same as defined by the quadrature formula passed to the constructor of the FEValues object above.
The functions of this class fall into different categories:
shape_value(), shape_grad(), etc: return one of the values of this object at a time. These functions are inlined, so this is the suggested access to all finite element values. There should be no loss in performance with an optimizing compiler. If the finite element is vector valued, then these functions return the only non-zero component of the requested shape function. However, some finite elements have shape functions that have more than one non-zero component (we call them non-"primitive"), and in this case this set of functions will throw an exception since they cannot generate a useful result. Rather, use the next set of functions.
shape_value_component(), shape_grad_component(), etc: This is the same set of functions as above, except that for vector valued finite elements they return only one vector component. This is useful for elements of which shape functions have more than one non-zero component, since then the above functions cannot be used, and you have to walk over all (or only the non- zero) components of the shape function using this set of functions.
get_function_values(), get_function_gradients(), etc.: Compute a finite element function or its derivative in quadrature points.
The mechanisms by which this class work are discussed on the page on Update flags and about the How Mapping, FiniteElement, and FEValues work together.
FEValuesBase< dim, spacedim >::FEValuesBase | ( | const unsigned int | n_q_points, |
const unsigned int | dofs_per_cell, | ||
const UpdateFlags | update_flags, | ||
const Mapping< dim, spacedim > & | mapping, | ||
const FiniteElement< dim, spacedim > & | fe | ||
) |
Constructor. Set up the array sizes with n_q_points
quadrature points, dofs_per_cell
trial functions per cell and with the given pattern to update the fields when the reinit
function of the derived classes is called. The fields themselves are not set up, this must happen in the constructor of the derived class.
Definition at line 2733 of file fe_values.cc.
FEValuesBase< dim, spacedim >::~FEValuesBase | ( | ) |
Destructor.
Definition at line 2756 of file fe_values.cc.
|
private |
Copy constructor. Since objects of this class are not copyable, we make it private, and also do not implement it.
const double& FEValuesBase< dim, spacedim >::shape_value | ( | const unsigned int | function_no, |
const unsigned int | point_no | ||
) | const |
Value of a shape function at a quadrature point on the cell, face or subface selected the last time the reinit
function of the derived class was called.
If the shape function is vector-valued, then this returns the only non- zero component. If the shape function has more than one non-zero component (i.e. it is not primitive), then throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_value_component() function.
function_no | Number of the shape function to be evaluated. Note that this number runs from zero to dofs_per_cell, even in the case of an FEFaceValues or FESubfaceValues object. |
point_no | Number of the quadrature point at which function is to be evaluated |
update_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. double FEValuesBase< dim, spacedim >::shape_value_component | ( | const unsigned int | function_no, |
const unsigned int | point_no, | ||
const unsigned int | component | ||
) | const |
Compute one vector component of the value of a shape function at a quadrature point. If the finite element is scalar, then only component zero is allowed and the return value equals that of the shape_value() function. If the finite element is vector valued but all shape functions are primitive (i.e. they are non-zero in only one component), then the value returned by shape_value() equals that of this function for exactly one component. This function is therefore only of greater interest if the shape function is not primitive, but then it is necessary since the other function cannot be used.
function_no | Number of the shape function to be evaluated. |
point_no | Number of the quadrature point at which function is to be evaluated. |
component | vector component to be evaluated. |
update_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const Tensor<1,spacedim>& FEValuesBase< dim, spacedim >::shape_grad | ( | const unsigned int | function_no, |
const unsigned int | quadrature_point | ||
) | const |
Compute the gradient of the function_no
th shape function at the quadrature_point
th quadrature point with respect to real cell coordinates. If you want to get the derivative in one of the coordinate directions, use the appropriate function of the Tensor class to extract one component of the Tensor returned by this function. Since only a reference to the gradient's value is returned, there should be no major performance drawback.
If the shape function is vector-valued, then this returns the only non- zero component. If the shape function has more than one non-zero component (i.e. it is not primitive), then it will throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_grad_component() function.
The same holds for the arguments of this function as for the shape_value() function.
function_no | Number of the shape function to be evaluated. |
quadrature_point | Number of the quadrature point at which function is to be evaluated. |
update_gradients
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Tensor<1,spacedim> FEValuesBase< dim, spacedim >::shape_grad_component | ( | const unsigned int | function_no, |
const unsigned int | point_no, | ||
const unsigned int | component | ||
) | const |
Return one vector component of the gradient of a shape function at a quadrature point. If the finite element is scalar, then only component zero is allowed and the return value equals that of the shape_grad() function. If the finite element is vector valued but all shape functions are primitive (i.e. they are non-zero in only one component), then the value returned by shape_grad() equals that of this function for exactly one component. This function is therefore only of greater interest if the shape function is not primitive, but then it is necessary since the other function cannot be used.
The same holds for the arguments of this function as for the shape_value_component() function.
update_gradients
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const Tensor<2,spacedim>& FEValuesBase< dim, spacedim >::shape_hessian | ( | const unsigned int | function_no, |
const unsigned int | point_no | ||
) | const |
Second derivatives of the function_no
th shape function at the point_no
th quadrature point with respect to real cell coordinates. If you want to get the derivatives in one of the coordinate directions, use the appropriate function of the Tensor class to extract one component. Since only a reference to the hessian values is returned, there should be no major performance drawback.
If the shape function is vector-valued, then this returns the only non- zero component. If the shape function has more than one non-zero component (i.e. it is not primitive), then throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_hessian_component() function.
The same holds for the arguments of this function as for the shape_value() function.
update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Tensor<2,spacedim> FEValuesBase< dim, spacedim >::shape_hessian_component | ( | const unsigned int | function_no, |
const unsigned int | point_no, | ||
const unsigned int | component | ||
) | const |
Return one vector component of the hessian of a shape function at a quadrature point. If the finite element is scalar, then only component zero is allowed and the return value equals that of the shape_hessian() function. If the finite element is vector valued but all shape functions are primitive (i.e. they are non-zero in only one component), then the value returned by shape_hessian() equals that of this function for exactly one component. This function is therefore only of greater interest if the shape function is not primitive, but then it is necessary since the other function cannot be used.
The same holds for the arguments of this function as for the shape_value_component() function.
update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const Tensor<3,spacedim>& FEValuesBase< dim, spacedim >::shape_3rd_derivative | ( | const unsigned int | function_no, |
const unsigned int | point_no | ||
) | const |
Third derivatives of the function_no
th shape function at the point_no
th quadrature point with respect to real cell coordinates. If you want to get the 3rd derivatives in one of the coordinate directions, use the appropriate function of the Tensor class to extract one component. Since only a reference to the 3rd derivative values is returned, there should be no major performance drawback.
If the shape function is vector-valued, then this returns the only non- zero component. If the shape function has more than one non-zero component (i.e. it is not primitive), then throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_3rdderivative_component() function.
The same holds for the arguments of this function as for the shape_value() function.
update_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Tensor<3,spacedim> FEValuesBase< dim, spacedim >::shape_3rd_derivative_component | ( | const unsigned int | function_no, |
const unsigned int | point_no, | ||
const unsigned int | component | ||
) | const |
Return one vector component of the third derivative of a shape function at a quadrature point. If the finite element is scalar, then only component zero is allowed and the return value equals that of the shape_3rdderivative() function. If the finite element is vector valued but all shape functions are primitive (i.e. they are non-zero in only one component), then the value returned by shape_3rdderivative() equals that of this function for exactly one component. This function is therefore only of greater interest if the shape function is not primitive, but then it is necessary since the other function cannot be used.
The same holds for the arguments of this function as for the shape_value_component() function.
update_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. void FEValuesBase< dim, spacedim >::get_function_values | ( | const InputVector & | fe_function, |
std::vector< typename InputVector::value_type > & | values | ||
) | const |
Return the values of a finite element function restricted to the current cell, face or subface selected the last time the reinit
function of the derived class was called, at the quadrature points.
If the present cell is not active then values are interpolated to the current cell and point values are computed from that.
This function may only be used if the finite element in use is a scalar one, i.e. has only one vector component. To get values of multi- component elements, there is another get_function_values() below, returning a vector of vectors of results.
[in] | fe_function | A vector of values that describes (globally) the finite element function that this function should evaluate at the quadrature points of the current cell. |
[out] | values | The values of the function specified by fe_function at the quadrature points of the current cell. The object is assume to already have the correct size. The data type stored by this output vector must be what you get when you multiply the values of shape function times the type used to store the values of the unknowns \(U_j\) of your finite element vector \(U\) (represented by the fe_function argument). This happens to be equal to the type of the elements of the solution vector. |
values[q]
will contain the value of the field described by fe_function at the \(q\)th quadrature point.update_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3197 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_values | ( | const InputVector & | fe_function, |
std::vector< Vector< typename InputVector::value_type > > & | values | ||
) | const |
This function does the same as the other get_function_values(), but applied to multi-component (vector-valued) elements. The meaning of the arguments is as explained there.
values[q]
is a vector of values of the field described by fe_function at the \(q\)th quadrature point. The size of the vector accessed by values[q]
equals the number of components of the finite element, i.e. values[q](c)
returns the value of the \(c\)th vector component at the \(q\)th quadrature point.update_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3242 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_values | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< typename InputVector::value_type > & | values | ||
) | const |
Generate function values from an arbitrary vector.
This function offers the possibility to extract function values in quadrature points from vectors not corresponding to a whole discretization.
The vector indices
corresponds to the degrees of freedom on a single cell. Its length may even be a multiple of the number of dofs per cell. Then, the vectors in value
should allow for the same multiple of the components of the finite element.
You may want to use this function, if you want to access just a single block from a BlockVector, if you have a multi-level vector or if you already have a local representation of your finite element data.
update_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3221 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_values | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< Vector< typename InputVector::value_type > > & | values | ||
) | const |
Generate vector function values from an arbitrary vector.
This function offers the possibility to extract function values in quadrature points from vectors not corresponding to a whole discretization.
The vector indices
corresponds to the degrees of freedom on a single cell. Its length may even be a multiple of the number of dofs per cell. Then, the vectors in value
should allow for the same multiple of the components of the finite element.
You may want to use this function, if you want to access just a single block from a BlockVector, if you have a multi-level vector or if you already have a local representation of your finite element data.
Since this function allows for fairly general combinations of argument sizes, be aware that the checks on the arguments may not detect errors.
update_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3268 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_values | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
VectorSlice< std::vector< std::vector< typename InputVector::value_type > > > | values, | ||
const bool | quadrature_points_fastest | ||
) | const |
Generate vector function values from an arbitrary vector.
This function offers the possibility to extract function values in quadrature points from vectors not corresponding to a whole discretization.
The vector indices
corresponds to the degrees of freedom on a single cell. Its length may even be a multiple of the number of dofs per cell. Then, the vectors in value
should allow for the same multiple of the components of the finite element.
Depending on the value of the last argument, the outer vector of values
has either the length of the quadrature rule (quadrature_points_fastest == false
) or the length of components to be filled quadrature_points_fastest == true
. If p
is the current quadrature point number and i
is the vector component of the solution desired, the access to values
is values[p][i]
if quadrature_points_fastest == false
, and values[i][p]
otherwise.
You may want to use this function, if you want to access just a single block from a BlockVector, if you have a multi-level vector or if you already have a local representation of your finite element data.
Since this function allows for fairly general combinations of argument sizes, be aware that the checks on the arguments may not detect errors.
update_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3297 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_gradients | ( | const InputVector & | fe_function, |
std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > & | gradients | ||
) | const |
Compute the gradients of a finite element at the quadrature points of a cell. This function is the equivalent of the corresponding get_function_values() function (see there for more information) but evaluates the finite element field's gradient instead of its value.
This function may only be used if the finite element in use is a scalar one, i.e. has only one vector component. There is a corresponding function of the same name for vector-valued finite elements.
[in] | fe_function | A vector of values that describes (globally) the finite element function that this function should evaluate at the quadrature points of the current cell. |
[out] | gradients | The gradients of the function specified by fe_function at the quadrature points of the current cell. The gradients are computed in real space (as opposed to on the unit cell). The object is assume to already have the correct size. The data type stored by this output vector must be what you get when you multiply the gradients of shape function times the type used to store the values of the unknowns \(U_j\) of your finite element vector \(U\) (represented by the fe_function argument). |
gradients[q]
will contain the gradient of the field described by fe_function at the \(q\)th quadrature point. gradients[q][d]
represents the derivative in coordinate direction \(d\) at quadrature point \(q\).update_gradients
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3329 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_gradients | ( | const InputVector & | fe_function, |
std::vector< std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > > & | gradients | ||
) | const |
This function does the same as the other get_function_gradients(), but applied to multi-component (vector-valued) elements. The meaning of the arguments is as explained there.
gradients[q]
is a vector of gradients of the field described by fe_function at the \(q\)th quadrature point. The size of the vector accessed by gradients[q]
equals the number of components of the finite element, i.e. gradients[q][c]
returns the gradient of the \(c\)th vector component at the \(q\)th quadrature point. Consequently, gradients[q][c][d]
is the derivative in coordinate direction \(d\) of the \(c\)th vector component of the vector field at quadrature point \(q\) of the current cell.update_gradients
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3375 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_gradients | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > & | gradients | ||
) | const |
Function gradient access with more flexibility. See get_function_values() with corresponding arguments.
update_gradients
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3352 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_gradients | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
VectorSlice< std::vector< std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > > > | gradients, | ||
bool | quadrature_points_fastest = false |
||
) | const |
Function gradient access with more flexibility. See get_function_values() with corresponding arguments.
update_gradients
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3400 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_hessians | ( | const InputVector & | fe_function, |
std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > & | hessians | ||
) | const |
Compute the tensor of second derivatives of a finite element at the quadrature points of a cell. This function is the equivalent of the corresponding get_function_values() function (see there for more information) but evaluates the finite element field's second derivatives instead of its value.
This function may only be used if the finite element in use is a scalar one, i.e. has only one vector component. There is a corresponding function of the same name for vector-valued finite elements.
[in] | fe_function | A vector of values that describes (globally) the finite element function that this function should evaluate at the quadrature points of the current cell. |
[out] | hessians | The Hessians of the function specified by fe_function at the quadrature points of the current cell. The Hessians are computed in real space (as opposed to on the unit cell). The object is assume to already have the correct size. The data type stored by this output vector must be what you get when you multiply the Hessians of shape function times the type used to store the values of the unknowns \(U_j\) of your finite element vector \(U\) (represented by the fe_function argument). |
hessians[q]
will contain the Hessian of the field described by fe_function at the \(q\)th quadrature point. hessians[q][i][j]
represents the \((i,j)\)th component of the matrix of second derivatives at quadrature point \(q\).update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3432 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_hessians | ( | const InputVector & | fe_function, |
std::vector< std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > > & | hessians, | ||
bool | quadrature_points_fastest = false |
||
) | const |
This function does the same as the other get_function_hessians(), but applied to multi-component (vector-valued) elements. The meaning of the arguments is as explained there.
hessians[q]
is a vector of Hessians of the field described by fe_function at the \(q\)th quadrature point. The size of the vector accessed by hessians[q]
equals the number of components of the finite element, i.e. hessians[q][c]
returns the Hessian of the \(c\)th vector component at the \(q\)th quadrature point. Consequently, hessians[q][c][i][j]
is the \((i,j)\)th component of the matrix of second derivatives of the \(c\)th vector component of the vector field at quadrature point \(q\) of the current cell.update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3478 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_hessians | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > & | hessians | ||
) | const |
Access to the second derivatives of a function with more flexibility. See get_function_values() with corresponding arguments.
Definition at line 3454 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_hessians | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
VectorSlice< std::vector< std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > > > | hessians, | ||
bool | quadrature_points_fastest = false |
||
) | const |
Access to the second derivatives of a function with more flexibility. See get_function_values() with corresponding arguments.
update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3504 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_laplacians | ( | const InputVector & | fe_function, |
std::vector< typename InputVector::value_type > & | laplacians | ||
) | const |
Compute the (scalar) Laplacian (i.e. the trace of the tensor of second derivatives) of a finite element at the quadrature points of a cell. This function is the equivalent of the corresponding get_function_values() function (see there for more information) but evaluates the finite element field's second derivatives instead of its value.
This function may only be used if the finite element in use is a scalar one, i.e. has only one vector component. There is a corresponding function of the same name for vector-valued finite elements.
[in] | fe_function | A vector of values that describes (globally) the finite element function that this function should evaluate at the quadrature points of the current cell. |
[out] | laplacians | The Laplacians of the function specified by fe_function at the quadrature points of the current cell. The Laplacians are computed in real space (as opposed to on the unit cell). The object is assume to already have the correct size. The data type stored by this output vector must be what you get when you multiply the Laplacians of shape function times the type used to store the values of the unknowns \(U_j\) of your finite element vector \(U\) (represented by the fe_function argument). This happens to be equal to the type of the elements of the input vector. |
laplacians[q]
will contain the Laplacian of the field described by fe_function at the \(q\)th quadrature point.laplacians[q]=trace(hessians[q])
, where hessians
would be the output of the get_function_hessians() function.update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3532 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_laplacians | ( | const InputVector & | fe_function, |
std::vector< Vector< typename InputVector::value_type > > & | laplacians | ||
) | const |
This function does the same as the other get_function_laplacians(), but applied to multi-component (vector-valued) elements. The meaning of the arguments is as explained there.
laplacians[q]
is a vector of Laplacians of the field described by fe_function at the \(q\)th quadrature point. The size of the vector accessed by laplacians[q]
equals the number of components of the finite element, i.e. laplacians[q][c]
returns the Laplacian of the \(c\)th vector component at the \(q\)th quadrature point.laplacians[q][c]=trace(hessians[q][c])
, where hessians
would be the output of the get_function_hessians() function.update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3577 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_laplacians | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< typename InputVector::value_type > & | laplacians | ||
) | const |
Access to the second derivatives of a function with more flexibility. See get_function_values() with corresponding arguments.
update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3555 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_laplacians | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< Vector< typename InputVector::value_type > > & | laplacians | ||
) | const |
Access to the second derivatives of a function with more flexibility. See get_function_values() with corresponding arguments.
update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3600 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_laplacians | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< std::vector< typename InputVector::value_type > > & | laplacians, | ||
bool | quadrature_points_fastest = false |
||
) | const |
Access to the second derivatives of a function with more flexibility. See get_function_values() with corresponding arguments.
update_hessians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3626 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_third_derivatives | ( | const InputVector & | fe_function, |
std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > & | third_derivatives | ||
) | const |
Compute the tensor of third derivatives of a finite element at the quadrature points of a cell. This function is the equivalent of the corresponding get_function_values() function (see there for more information) but evaluates the finite element field's third derivatives instead of its value.
This function may only be used if the finite element in use is a scalar one, i.e. has only one vector component. There is a corresponding function of the same name for vector-valued finite elements.
[in] | fe_function | A vector of values that describes (globally) the finite element function that this function should evaluate at the quadrature points of the current cell. |
[out] | third_derivatives | The third derivatives of the function specified by fe_function at the quadrature points of the current cell. The third derivatives are computed in real space (as opposed to on the unit cell). The object is assumed to already have the correct size. The data type stored by this output vector must be what you get when you multiply the third derivatives of shape function times the type used to store the values of the unknowns \(U_j\) of your finite element vector \(U\) (represented by the fe_function argument). |
third_derivatives[q]
will contain the third derivatives of the field described by fe_function at the \(q\)th quadrature point. third_derivatives[q][i][j][k]
represents the \((i,j,k)\)th component of the 3rd order tensor of third derivatives at quadrature point \(q\).update_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3653 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_third_derivatives | ( | const InputVector & | fe_function, |
std::vector< std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > > & | third_derivatives, | ||
bool | quadrature_points_fastest = false |
||
) | const |
This function does the same as the other get_function_third_derivatives(), but applied to multi-component (vector- valued) elements. The meaning of the arguments is as explained there.
third_derivatives[q]
is a vector of third derivatives of the field described by fe_function at the \(q\)th quadrature point. The size of the vector accessed by third_derivatives[q]
equals the number of components of the finite element, i.e. third_derivatives[q][c]
returns the third derivative of the \(c\)th vector component at the \(q\)th quadrature point. Consequently, third_derivatives[q][c][i][j][k]
is the \((i,j,k)\)th component of the tensor of third derivatives of the \(c\)th vector component of the vector field at quadrature point \(q\) of the current cell.update_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3699 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_third_derivatives | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > & | third_derivatives | ||
) | const |
Access to the third derivatives of a function with more flexibility. See get_function_values() with corresponding arguments.
Definition at line 3675 of file fe_values.cc.
void FEValuesBase< dim, spacedim >::get_function_third_derivatives | ( | const InputVector & | fe_function, |
const VectorSlice< const std::vector< types::global_dof_index > > & | indices, | ||
VectorSlice< std::vector< std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > > > | third_derivatives, | ||
bool | quadrature_points_fastest = false |
||
) | const |
Access to the third derivatives of a function with more flexibility. See get_function_values() with corresponding arguments.
update_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3725 of file fe_values.cc.
const Point<spacedim>& FEValuesBase< dim, spacedim >::quadrature_point | ( | const unsigned int | q | ) | const |
Position of the q
th quadrature point in real space.
update_quadrature_points
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<Point<spacedim> >& FEValuesBase< dim, spacedim >::get_quadrature_points | ( | ) | const |
Return a reference to the vector of quadrature points in real space.
update_quadrature_points
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. double FEValuesBase< dim, spacedim >::JxW | ( | const unsigned int | quadrature_point | ) | const |
Mapped quadrature weight. If this object refers to a volume evaluation (i.e. the derived class is of type FEValues), then this is the Jacobi determinant times the weight of the *i
th unit quadrature point.
For surface evaluations (i.e. classes FEFaceValues or FESubfaceValues), it is the mapped surface element times the weight of the quadrature point.
You can think of the quantity returned by this function as the volume or surface element \(dx, ds\) in the integral that we implement here by quadrature.
update_JxW_values
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<double>& FEValuesBase< dim, spacedim >::get_JxW_values | ( | ) | const |
Return a reference to the array holding the values returned by JxW().
const DerivativeForm<1,dim,spacedim>& FEValuesBase< dim, spacedim >::jacobian | ( | const unsigned int | quadrature_point | ) | const |
Return the Jacobian of the transformation at the specified quadrature point, i.e. \(J_{ij}=dx_i/d\hat x_j\)
update_jacobians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<DerivativeForm<1,dim,spacedim> >& FEValuesBase< dim, spacedim >::get_jacobians | ( | ) | const |
Return a reference to the array holding the values returned by jacobian().
update_jacobians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const DerivativeForm<2,dim,spacedim>& FEValuesBase< dim, spacedim >::jacobian_grad | ( | const unsigned int | quadrature_point | ) | const |
Return the second derivative of the transformation from unit to real cell, i.e. the first derivative of the Jacobian, at the specified quadrature point, i.e. \(G_{ijk}=dJ_{jk}/d\hat x_i\).
update_jacobian_grads
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<DerivativeForm<2,dim,spacedim> >& FEValuesBase< dim, spacedim >::get_jacobian_grads | ( | ) | const |
Return a reference to the array holding the values returned by jacobian_grads().
update_jacobian_grads
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const Tensor<3,spacedim>& FEValuesBase< dim, spacedim >::jacobian_pushed_forward_grad | ( | const unsigned int | quadrature_point | ) | const |
Return the second derivative of the transformation from unit to real cell, i.e. the first derivative of the Jacobian, at the specified quadrature point, pushed forward to the real cell coordinates, i.e. \(G_{ijk}=dJ_{iJ}/d\hat x_K (J_{jJ})^{-1} (J_{kK})^{-1}\).
update_jacobian_pushed_forward_grads
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<Tensor<3,spacedim> >& FEValuesBase< dim, spacedim >::get_jacobian_pushed_forward_grads | ( | ) | const |
Return a reference to the array holding the values returned by jacobian_pushed_forward_grads().
update_jacobian_pushed_forward_grads
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const DerivativeForm<3,dim,spacedim>& FEValuesBase< dim, spacedim >::jacobian_2nd_derivative | ( | const unsigned int | quadrature_point | ) | const |
Return the third derivative of the transformation from unit to real cell, i.e. the second derivative of the Jacobian, at the specified quadrature point, i.e. \(G_{ijkl}=\frac{d^2J_{ij}}{d\hat x_k d\hat x_l}\).
update_jacobian_2nd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<DerivativeForm<3,dim,spacedim> >& FEValuesBase< dim, spacedim >::get_jacobian_2nd_derivatives | ( | ) | const |
Return a reference to the array holding the values returned by jacobian_2nd_derivatives().
update_jacobian_2nd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const Tensor<4,spacedim>& FEValuesBase< dim, spacedim >::jacobian_pushed_forward_2nd_derivative | ( | const unsigned int | quadrature_point | ) | const |
Return the third derivative of the transformation from unit to real cell, i.e. the second derivative of the Jacobian, at the specified quadrature point, pushed forward to the real cell coordinates, i.e. \(G_{ijkl}=\frac{d^2J_{iJ}}{d\hat x_K d\hat x_L} (J_{jJ})^{-1} (J_{kK})^{-1}(J_{lL})^{-1}\).
update_jacobian_pushed_forward_2nd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<Tensor<4,spacedim> >& FEValuesBase< dim, spacedim >::get_jacobian_pushed_forward_2nd_derivatives | ( | ) | const |
Return a reference to the array holding the values returned by jacobian_pushed_forward_2nd_derivatives().
update_jacobian_pushed_forward_2nd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const DerivativeForm<4,dim,spacedim>& FEValuesBase< dim, spacedim >::jacobian_3rd_derivative | ( | const unsigned int | quadrature_point | ) | const |
Return the fourth derivative of the transformation from unit to real cell, i.e. the third derivative of the Jacobian, at the specified quadrature point, i.e. \(G_{ijklm}=\frac{d^2J_{ij}}{d\hat x_k d\hat x_l d\hat x_m}\).
update_jacobian_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<DerivativeForm<4,dim,spacedim> >& FEValuesBase< dim, spacedim >::get_jacobian_3rd_derivatives | ( | ) | const |
Return a reference to the array holding the values returned by jacobian_3rd_derivatives().
update_jacobian_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const Tensor<5,spacedim>& FEValuesBase< dim, spacedim >::jacobian_pushed_forward_3rd_derivative | ( | const unsigned int | quadrature_point | ) | const |
Return the fourth derivative of the transformation from unit to real cell, i.e. the third derivative of the Jacobian, at the specified quadrature point, pushed forward to the real cell coordinates, i.e. \(G_{ijklm}=\frac{d^3J_{iJ}}{d\hat x_K d\hat x_L d\hat x_M} (J_{jJ})^{-1} (J_{kK})^{-1} (J_{lL})^{-1} (J_{mM})^{-1}\).
update_jacobian_pushed_forward_3rd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<Tensor<5,spacedim> >& FEValuesBase< dim, spacedim >::get_jacobian_pushed_forward_3rd_derivatives | ( | ) | const |
Return a reference to the array holding the values returned by jacobian_pushed_forward_3rd_derivatives().
update_jacobian_pushed_forward_2nd_derivatives
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const DerivativeForm<1,spacedim,dim>& FEValuesBase< dim, spacedim >::inverse_jacobian | ( | const unsigned int | quadrature_point | ) | const |
Return the inverse Jacobian of the transformation at the specified quadrature point, i.e. \(J_{ij}=d\hat x_i/dx_j\)
update_inverse_jacobians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector<DerivativeForm<1,spacedim,dim> >& FEValuesBase< dim, spacedim >::get_inverse_jacobians | ( | ) | const |
Return a reference to the array holding the values returned by inverse_jacobian().
update_inverse_jacobians
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const Tensor<1,spacedim>& FEValuesBase< dim, spacedim >::normal_vector | ( | const unsigned int | i | ) | const |
For a face, return the outward normal vector to the cell at the i
th quadrature point.
For a cell of codimension one, return the normal vector. There are of course two normal directions to a manifold in that case, and this function returns the "up" direction as induced by the numbering of the vertices.
The length of the vector is normalized to one.
update_normal_vectors
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. const std::vector< Tensor< 1, spacedim > > & FEValuesBase< dim, spacedim >::get_all_normal_vectors | ( | ) | const |
Return the normal vectors at the quadrature points. For a face, these are the outward normal vectors to the cell. For a cell of codimension one, the orientation is given by the numbering of vertices.
update_normal_vectors
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information.Definition at line 3762 of file fe_values.cc.
const std::vector< Tensor< 1, spacedim > > & FEValuesBase< dim, spacedim >::get_normal_vectors | ( | ) | const |
Return the normal vectors at the quadrature points. For a face, these are the outward normal vectors to the cell. For a cell of codimension one, the orientation is given by the numbering of vertices.
update_normal_vectors
flag must be an element of the list of UpdateFlags that you passed to the constructor of this object. See The interplay of UpdateFlags, Mapping, and FiniteElement in FEValues for more information. Definition at line 3773 of file fe_values.cc.
const FEValuesViews::Scalar<dim,spacedim>& FEValuesBase< dim, spacedim >::operator[] | ( | const FEValuesExtractors::Scalar & | scalar | ) | const |
Create a view of the current FEValues object that represents a particular scalar component of the possibly vector-valued finite element. The concept of views is explained in the documentation of the namespace FEValuesViews and in particular in the Handling vector valued problems module.
const FEValuesViews::Vector<dim,spacedim>& FEValuesBase< dim, spacedim >::operator[] | ( | const FEValuesExtractors::Vector & | vector | ) | const |
Create a view of the current FEValues object that represents a set of dim
scalar components (i.e. a vector) of the vector-valued finite element. The concept of views is explained in the documentation of the namespace FEValuesViews and in particular in the Handling vector valued problems module.
const FEValuesViews::SymmetricTensor<2,dim,spacedim>& FEValuesBase< dim, spacedim >::operator[] | ( | const FEValuesExtractors::SymmetricTensor< 2 > & | tensor | ) | const |
Create a view of the current FEValues object that represents a set of (dim*dim + dim)/2
scalar components (i.e. a symmetric 2nd order tensor) of the vector-valued finite element. The concept of views is explained in the documentation of the namespace FEValuesViews and in particular in the Handling vector valued problems module.
const FEValuesViews::Tensor<2,dim,spacedim>& FEValuesBase< dim, spacedim >::operator[] | ( | const FEValuesExtractors::Tensor< 2 > & | tensor | ) | const |
Create a view of the current FEValues object that represents a set of (dim*dim)
scalar components (i.e. a 2nd order tensor) of the vector-valued finite element. The concept of views is explained in the documentation of the namespace FEValuesViews and in particular in the Handling vector valued problems module.
const Mapping<dim,spacedim>& FEValuesBase< dim, spacedim >::get_mapping | ( | ) | const |
Constant reference to the selected mapping object.
const FiniteElement<dim,spacedim>& FEValuesBase< dim, spacedim >::get_fe | ( | ) | const |
Constant reference to the selected finite element object.
UpdateFlags FEValuesBase< dim, spacedim >::get_update_flags | ( | ) | const |
Return the update flags set for this object.
const Triangulation< dim, spacedim >::cell_iterator FEValuesBase< dim, spacedim >::get_cell | ( | ) | const |
Return a triangulation iterator to the current cell.
Definition at line 3753 of file fe_values.cc.
CellSimilarity::Similarity FEValuesBase< dim, spacedim >::get_cell_similarity | ( | ) | const |
Return the relation of the current cell to the previous cell. This allows re-use of some cell data (like local matrices for equations with constant coefficients) if the result is CellSimilarity::translation
.
Definition at line 3943 of file fe_values.cc.
std::size_t FEValuesBase< dim, spacedim >::memory_consumption | ( | ) | const |
Determine an estimate for the memory consumption (in bytes) of this object.
Definition at line 3785 of file fe_values.cc.
|
protected |
A function that is connected to the triangulation in order to reset the stored 'present_cell' iterator to an invalid one whenever the triangulation is changed and the iterator consequently becomes invalid.
Definition at line 3824 of file fe_values.cc.
|
protected |
This function is called by the various reinit() functions in derived classes. Given the cell indicated by the argument, test whether we have to throw away the previously stored present_cell argument because it would require us to compare cells from different triangulations. In checking all this, also make sure that we have tria_listener connected to the triangulation to which we will set present_cell right after calling this function.
Definition at line 3843 of file fe_values.cc.
|
protected |
Initialize some update flags. Called from the initialize
functions of derived classes, which are in turn called from their constructors.
Basically, this function finds out using the finite element and mapping object already stored which flags need to be set to compute everything the user wants, as expressed through the flags passed as argument.
Definition at line 3805 of file fe_values.cc.
|
inlineprotected |
A function that checks whether the new cell is similar to the one previously used. Then, a significant amount of the data can be reused, e.g. the derivatives of the basis functions in real space, shape_grad.
Definition at line 3889 of file fe_values.cc.
|
private |
Copy operator. Since objects of this class are not copyable, we make it private, and also do not implement it.
Make the view classes friends of this class, since they access internal data.
Definition at line 3122 of file fe_values.h.
|
static |
Dimension in which this object operates.
Definition at line 1830 of file fe_values.h.
|
static |
Dimension of the space in which this object operates.
Definition at line 1835 of file fe_values.h.
const unsigned int FEValuesBase< dim, spacedim >::n_quadrature_points |
Number of quadrature points.
Definition at line 1840 of file fe_values.h.
const unsigned int FEValuesBase< dim, spacedim >::dofs_per_cell |
Number of shape functions per cell. If we use this base class to evaluate a finite element on faces of cells, this is still the number of degrees of freedom per cell, not per face.
Definition at line 1847 of file fe_values.h.
|
protected |
Store the cell selected last time the reinit() function was called. This is necessary for the get_function_*
functions as well as the functions of same name in the extractor classes.
Definition at line 2985 of file fe_values.h.
|
protected |
A signal connection we use to ensure we get informed whenever the triangulation changes by refinement. We need to know about that because it invalidates all cell iterators and, as part of that, the 'present_cell' iterator we keep around between subsequent calls to reinit() in order to compute the cell similarity.
Definition at line 3001 of file fe_values.h.
|
protected |
A signal connection we use to ensure we get informed whenever the triangulation changes by mesh transformations. We need to know about that because it invalidates all cell iterators and, as part of that, the 'present_cell' iterator we keep around between subsequent calls to reinit() in order to compute the cell similarity.
Definition at line 3010 of file fe_values.h.
|
protected |
A pointer to the mapping object associated with this FEValues object.
Definition at line 3034 of file fe_values.h.
|
protected |
A pointer to the internal data object of mapping, obtained from Mapping::get_data(), Mapping::get_face_data(), or Mapping::get_subface_data().
Definition at line 3041 of file fe_values.h.
|
protected |
An object into which the Mapping::fill_fe_values() and similar functions place their output.
Definition at line 3047 of file fe_values.h.
|
protected |
A pointer to the finite element object associated with this FEValues object.
Definition at line 3054 of file fe_values.h.
|
protected |
A pointer to the internal data object of finite element, obtained from FiniteElement::get_data(), Mapping::get_face_data(), or FiniteElement::get_subface_data().
Definition at line 3061 of file fe_values.h.
|
protected |
An object into which the FiniteElement::fill_fe_values() and similar functions place their output.
Definition at line 3067 of file fe_values.h.
|
protected |
Original update flags handed to the constructor of FEValues.
Definition at line 3073 of file fe_values.h.
|
protected |
An enum variable that can store different states of the current cell in comparison to the previously visited cell. If wanted, additional states can be checked here and used in one of the methods used during reinit.
Definition at line 3090 of file fe_values.h.
|
private |
A cache for all possible FEValuesViews objects.
Definition at line 3116 of file fe_values.h.