![]() |
Reference documentation for deal.II version 8.1.0
|
#include <fe_evaluation.h>
Public Types | |
typedef FEEvaluationAccess < dim, Utilities::fixed_int_power < fe_degree+1, dim >::value, Utilities::fixed_int_power < n_q_points_1d, dim >::value, n_components_, Number > | BaseClass |
typedef Number | number_type |
typedef BaseClass::value_type | value_type |
typedef BaseClass::gradient_type | gradient_type |
![]() | |
typedef Number | number_type |
typedef Tensor < 1, n_components_, VectorizedArray< Number > > | value_type |
typedef Tensor < 1, n_components_, Tensor < 1, dim, VectorizedArray < Number > > > | gradient_type |
typedef FEEvaluationBase< dim, dofs_per_cell_, n_q_points_, n_components_, Number > | BaseClass |
![]() | |
typedef Number | number_type |
typedef Tensor < 1, n_components_, VectorizedArray< Number > > | value_type |
typedef Tensor < 1, n_components_, Tensor < 1, dim, VectorizedArray < Number > > > | gradient_type |
Public Member Functions | |
FEEvaluationGeneral (const MatrixFree< dim, Number > &matrix_free, const unsigned int fe_no=0, const unsigned int quad_no=0) | |
void | evaluate (const bool evaluate_val, const bool evaluate_grad, const bool evaluate_hess=false) |
void | integrate (const bool integrate_val, const bool integrate_grad) |
Point< dim, VectorizedArray < Number > > | quadrature_point (const unsigned int q_point) const |
![]() | |
void | reinit (const unsigned int cell) |
unsigned int | get_cell_data_number () const |
internal::MatrixFreeFunctions::CellType | get_cell_type () const |
template<typename VectorType > | |
void | read_dof_values (const VectorType &src) |
template<typename VectorType > | |
void | read_dof_values (const std::vector< VectorType > &src, const unsigned int first_index=0) |
template<typename VectorType > | |
void | read_dof_values (const std::vector< VectorType * > &src, const unsigned int first_index=0) |
template<typename VectorType > | |
void | read_dof_values_plain (const VectorType &src) |
template<typename VectorType > | |
void | read_dof_values_plain (const std::vector< VectorType > &src, const unsigned int first_index=0) |
template<typename VectorType > | |
void | read_dof_values_plain (const std::vector< VectorType * > &src, const unsigned int first_index=0) |
template<typename VectorType > | |
void | distribute_local_to_global (VectorType &dst) const |
template<typename VectorType > | |
void | distribute_local_to_global (std::vector< VectorType > &dst, const unsigned int first_index=0) const |
template<typename VectorType > | |
void | distribute_local_to_global (std::vector< VectorType * > &dst, const unsigned int first_index=0) const |
template<typename VectorType > | |
void | set_dof_values (VectorType &dst) const |
template<typename VectorType > | |
void | set_dof_values (std::vector< VectorType > &dst, const unsigned int first_index=0) const |
template<typename VectorType > | |
void | set_dof_values (std::vector< VectorType * > &dst, const unsigned int first_index=0) const |
value_type | get_dof_value (const unsigned int dof) const |
void | submit_dof_value (const value_type val_in, const unsigned int dof) |
value_type | get_value (const unsigned int q_point) const |
void | submit_value (const value_type val_in, const unsigned int q_point) |
gradient_type | get_gradient (const unsigned int q_point) const |
void | submit_gradient (const gradient_type grad_in, const unsigned int q_point) |
Tensor< 1, n_components_, Tensor< 2, dim, VectorizedArray< Number > > > | get_hessian (const unsigned int q_point) const |
gradient_type | get_hessian_diagonal (const unsigned int q_point) const |
value_type | get_laplacian (const unsigned int q_point) const |
value_type | integrate_value () const |
const VectorizedArray< Number > * | begin_dof_values () const |
VectorizedArray< Number > * | begin_dof_values () |
const VectorizedArray< Number > * | begin_values () const |
VectorizedArray< Number > * | begin_values () |
const VectorizedArray< Number > * | begin_gradients () const |
VectorizedArray< Number > * | begin_gradients () |
const VectorizedArray< Number > * | begin_hessians () const |
VectorizedArray< Number > * | begin_hessians () |
Static Public Attributes | |
static const unsigned int | dimension = dim |
static const unsigned int | n_components = n_components_ |
static const unsigned int | dofs_per_cell = BaseClass::dofs_per_cell |
static const unsigned int | n_q_points = BaseClass::n_q_points |
![]() | |
static const unsigned int | dimension |
static const unsigned int | n_components |
static const unsigned int | dofs_per_cell |
static const unsigned int | n_q_points |
![]() | |
static const unsigned int | dimension = dim |
static const unsigned int | n_components = n_components_ |
static const unsigned int | dofs_per_cell = dofs_per_cell_ |
static const unsigned int | n_q_points = n_q_points_ |
Protected Member Functions | |
template<int direction, bool dof_to_quad, bool add> | |
void | apply_values (const VectorizedArray< Number > in[], VectorizedArray< Number > out[]) |
template<int direction, bool dof_to_quad, bool add> | |
void | apply_gradients (const VectorizedArray< Number > in[], VectorizedArray< Number > out[]) |
template<int direction, bool dof_to_quad, bool add> | |
void | apply_hessians (const VectorizedArray< Number > in[], VectorizedArray< Number > out[]) |
![]() | |
FEEvaluationAccess (const MatrixFree< dim, Number > &matrix_free, const unsigned int fe_no=0, const unsigned int quad_no=0) | |
![]() | |
FEEvaluationBase (const MatrixFree< dim, Number > &matrix_free, const unsigned int fe_no=0, const unsigned int quad_no=0) | |
template<typename VectorType , typename VectorOperation > | |
void | read_write_operation (const VectorOperation &operation, VectorType *vectors[]) const |
template<typename VectorType > | |
void | read_dof_values_plain (const VectorType *src_data[]) |
Friends | |
template<typename FEEval > | |
void | internal::do_evaluate (FEEval &, const bool, const bool, const bool) |
template<typename FEEval > | |
void | internal::do_integrate (FEEval &, const bool, const bool) |
The class that provides all functions necessary to evaluate functions at quadrature points and cell integrations. In functionality, this class is similar to FEValues<dim>, however, it includes a lot of specialized functions that make it much faster (between 5 and 500 times as fast, depending on the polynomial order). Access to the data fields is provided through functionality in the class FEEvaluationAccess.
This class is designed for general local finite element operations based on tensor products of 1D polynomials and quadrature points. Often, there are some symmetries or zeros in the unit cell data that allow for a more efficient operator application. FEEvaluation is specialized to standard FE_Q/FE_DGQ elements and quadrature points symmetric around 0.5 (like Gauss quadrature), and hence the most common situation. FEEvaluationGL is a specialization for elements where quadrature formula and support points are chosen so that a orthogonal relation between the values holds. This is the case for FE_Q based on Gauss-Lobatto support points with Gauss-Lobatto quadrature formula of the same order.
This class has five template arguments:
dim | Dimension in which this class is to be used |
fe_degree | Degree of the tensor product finite element with fe_degree+1 degrees of freedom per coordinate direction |
n_q_points_1d | Number of points in the quadrature formula in 1D, usually chosen as fe_degree+1 |
n_components | Number of vector components when solving a system of PDEs. If the same operation is applied to several components of a PDE (e.g. a vector Laplace equation), they can be applied simultaneously with one call (and often more efficiently) |
Number | Number format, usually double or float |
Definition at line 1074 of file fe_evaluation.h.
FEEvaluationGeneral< dim, fe_degree, n_q_points_1d, n_components_, Number >::FEEvaluationGeneral | ( | const MatrixFree< dim, Number > & | matrix_free, |
const unsigned int | fe_no = 0 , |
||
const unsigned int | quad_no = 0 |
||
) |
Constructor. Takes all data stored in MatrixFree. If applied to problems with more than one finite element or more than one quadrature formula selected during construction of matrix_free
, fe_no
and quad_no
allow to select the appropriate components.
void FEEvaluationGeneral< dim, fe_degree, n_q_points_1d, n_components_, Number >::evaluate | ( | const bool | evaluate_val, |
const bool | evaluate_grad, | ||
const bool | evaluate_hess = false |
||
) |
Evaluates the function values, the gradients, and the Laplacians of the FE function given at the DoF values in the input vector at the quadrature points. The function arguments specify which parts shall actually be computed. Needs to be called before the functions get_value()
, get_gradient()
or get_laplacian
return useful information.
void FEEvaluationGeneral< dim, fe_degree, n_q_points_1d, n_components_, Number >::integrate | ( | const bool | integrate_val, |
const bool | integrate_grad | ||
) |
This function takes the values and/or gradients that are stored on quadrature points, tests them by all the basis functions/gradients on the cell and performs the cell integration. The two function arguments integrate_val
and integrate_grad
are used to enable/disable some of values or gradients.
Point<dim,VectorizedArray<Number> > FEEvaluationGeneral< dim, fe_degree, n_q_points_1d, n_components_, Number >::quadrature_point | ( | const unsigned int | q_point | ) | const |
Returns the q-th quadrature point stored in MappingInfo.
|
protected |
Internal function that applies the function values of the tensor product in a given coordinate direction (first template argument), from polynomials to values on quadrature points (second flag set to true) or in an integration loop from values on quadrature points to values tested by different test function (second flag set to false), and if the result is to be added to previous content in the data fields or not.
|
protected |
Internal function that applies the gradient operation of the tensor product in a given coordinate direction (first template argument), from polynomials to values on quadrature points (second flag set to true) or in an integration loop from values on quadrature points to values tested by different test function (second flag set to false), and if the result is to be added to previous content in the data fields or not.
|
protected |
Internal function that applies the second derivative operation (Hessian) of the tensor product in a given coordinate direction (first template argument), from polynomials to values on quadrature points (second flag set to true) or in an integration loop from values on quadrature points to values tested by different test function (second flag set to false), and if the result is to be added to previous content in the data fields or not.
|
friend |
Friend declaration.