Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Types | Public Member Functions | Static Public Attributes | Protected Member Functions | List of all members
FEEvaluationGL< dim, fe_degree, n_components_, Number > Class Template Reference

#include <fe_evaluation.h>

Inheritance diagram for FEEvaluationGL< dim, fe_degree, n_components_, Number >:
[legend]

Public Types

typedef FEEvaluation< dim,
fe_degree, fe_degree+1,
n_components_, Number > 
BaseClass
 
typedef Number number_type
 
typedef BaseClass::value_type value_type
 
typedef BaseClass::gradient_type gradient_type
 
- Public Types inherited from FEEvaluation< dim, fe_degree, fe_degree+1, n_components_, Number >
typedef FEEvaluationGeneral
< dim, fe_degree,
n_q_points_1d, n_components_,
Number > 
BaseClass
 
typedef Number number_type
 
typedef BaseClass::value_type value_type
 
typedef BaseClass::gradient_type gradient_type
 
- Public Types inherited from FEEvaluationGeneral< dim, fe_degree, n_q_points_1d, n_components_, Number >
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
 
- Public Types inherited from FEEvaluationAccess< dim, Utilities::fixed_int_power< fe_degree+1, dim >::value, Utilities::fixed_int_power< n_q_points_1d, dim >::value, n_components_, Number >
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
 
- Public Types inherited from FEEvaluationBase< dim, dofs_per_cell_, n_q_points_, n_components_, Number >
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

 FEEvaluationGL (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_lapl=false)
 
void integrate (const bool integrate_val, const bool integrate_grad)
 
- Public Member Functions inherited from FEEvaluation< dim, fe_degree, fe_degree+1, n_components_, Number >
 FEEvaluation (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)
 
- Public Member Functions inherited from 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)
 
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
 
- Public Member Functions inherited from FEEvaluationBase< dim, dofs_per_cell_, n_q_points_, n_components_, Number >
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 Public Attributes inherited from FEEvaluation< dim, fe_degree, fe_degree+1, n_components_, Number >
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 Public Attributes inherited from FEEvaluationGeneral< dim, fe_degree, n_q_points_1d, n_components_, Number >
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 Public Attributes inherited from FEEvaluationAccess< dim, Utilities::fixed_int_power< fe_degree+1, dim >::value, Utilities::fixed_int_power< n_q_points_1d, dim >::value, n_components_, Number >
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 Public Attributes inherited from FEEvaluationBase< dim, dofs_per_cell_, n_q_points_, n_components_, Number >
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_gradients (const VectorizedArray< Number > in[], VectorizedArray< Number > out[])
 
- Protected Member Functions inherited from FEEvaluation< dim, fe_degree, fe_degree+1, n_components_, Number >
void apply_values (const VectorizedArray< Number > in[], VectorizedArray< Number > out[])
 
void apply_gradients (const VectorizedArray< Number > in[], VectorizedArray< Number > out[])
 
void apply_hessians (const VectorizedArray< Number > in[], VectorizedArray< Number > out[])
 
- Protected Member Functions inherited from FEEvaluationGeneral< dim, fe_degree, n_q_points_1d, n_components_, Number >
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[])
 
- Protected Member Functions inherited from FEEvaluationAccess< dim, Utilities::fixed_int_power< fe_degree+1, dim >::value, Utilities::fixed_int_power< n_q_points_1d, dim >::value, n_components_, Number >
 FEEvaluationAccess (const MatrixFree< dim, Number > &matrix_free, const unsigned int fe_no=0, const unsigned int quad_no=0)
 
- Protected Member Functions inherited from FEEvaluationBase< dim, dofs_per_cell_, n_q_points_, n_components_, Number >
 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[])
 

Additional Inherited Members

- Protected Attributes inherited from FEEvaluation< dim, fe_degree, fe_degree+1, n_components_, Number >
VectorizedArray< Number > shape_val_evenodd [fe_degree+1][(n_q_points_1d+1)/2]
 
VectorizedArray< Number > shape_gra_evenodd [fe_degree+1][(n_q_points_1d+1)/2]
 
VectorizedArray< Number > shape_hes_evenodd [fe_degree+1][(n_q_points_1d+1)/2]
 
- Protected Attributes inherited from FEEvaluationBase< dim, dofs_per_cell_, n_q_points_, n_components_, Number >
VectorizedArray< Number > values_dofs [n_components][dofs_per_cell >0?dofs_per_cell:1]
 
VectorizedArray< Number > values_quad [n_components][n_q_points >0?n_q_points:1]
 
VectorizedArray< Number > gradients_quad [n_components][dim][n_q_points >0?n_q_points:1]
 
VectorizedArray< Number > hessians_quad [n_components][(dim *(dim+1))/2][n_q_points >0?n_q_points:1]
 
const unsigned int quad_no
 
const unsigned int n_fe_components
 
const unsigned int active_fe_index
 
const unsigned int active_quad_index
 
const MatrixFree< dim, Number > & matrix_info
 
const
internal::MatrixFreeFunctions::DoFInfo & 
dof_info
 
const
internal::MatrixFreeFunctions::MappingInfo
< dim, Number > & 
mapping_info
 
const
internal::MatrixFreeFunctions::ShapeInfo
< Number > & 
data
 
const Tensor< 1, dim,
VectorizedArray< Number > > * 
cartesian_data
 
const Tensor< 2, dim,
VectorizedArray< Number > > * 
jacobian
 
const VectorizedArray< Number > * J_value
 
const VectorizedArray< Number > * quadrature_weights
 
const Point< dim,
VectorizedArray< Number > > * 
quadrature_points
 
const Tensor< 2, dim,
VectorizedArray< Number > > * 
jacobian_grad
 
const Tensor< 1,(dim >1?dim
*(dim-1)/2:1), Tensor< 1, dim,
VectorizedArray< Number > > > * 
jacobian_grad_upper
 
unsigned int cell
 
internal::MatrixFreeFunctions::CellType cell_type
 
unsigned int cell_data_number
 
bool dof_values_initialized
 
bool values_quad_initialized
 
bool gradients_quad_initialized
 
bool hessians_quad_initialized
 
bool values_quad_submitted
 
bool gradients_quad_submitted
 

Detailed Description

template<int dim, int fe_degree, int n_components_ = 1, typename Number = double>
class FEEvaluationGL< dim, fe_degree, n_components_, Number >

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, depending on the polynomial order).

This class is a specialization of FEEvaluation for elements where quadrature formula and support points are chosen so that a orthonormal 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 (QGaussLobatto). In that case, application of values is trivial (as the transformation is the identity matrix), and application of gradients is considerably simpler (since all value applications in directions other than the gradient direction are again identity operations).

This class has four template arguments:

Parameters
dimDimension in which this class is to be used
fe_degreeDegree of the tensor product finite element with fe_degree+1 degrees of freedom per coordinate direction. The quadrature formula is tied to the choice of the element by setting n_q_points_1d = fe_degree+1, which gives a diagonal mass matrix
n_componentsNumber 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)
NumberNumber format, usually double or float
Author
Katharina Kormann and Martin Kronbichler, 2010, 2011

Definition at line 1360 of file fe_evaluation.h.

Constructor & Destructor Documentation

template<int dim, int fe_degree, int n_components_ = 1, typename Number = double>
FEEvaluationGL< dim, fe_degree, n_components_, Number >::FEEvaluationGL ( 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.

Member Function Documentation

template<int dim, int fe_degree, int n_components_ = 1, typename Number = double>
void FEEvaluationGL< dim, fe_degree, n_components_, Number >::evaluate ( const bool  evaluate_val,
const bool  evaluate_grad,
const bool  evaluate_lapl = false 
)

Evaluates the function values, the gradients, and the Hessians of the FE function given at the DoF values in the input vector at the quadrature points of the unit cell. The function arguments specify which parts shall actually be computed. Needs to be called before the functions get_value(), get_gradient() or get_laplacian give useful information (unless these values have been set manually).

template<int dim, int fe_degree, int n_components_ = 1, typename Number = double>
void FEEvaluationGL< dim, fe_degree, 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.

template<int dim, int fe_degree, int n_components_ = 1, typename Number = double>
template<int direction, bool dof_to_quad, bool add>
void FEEvaluationGL< dim, fe_degree, n_components_, Number >::apply_gradients ( const VectorizedArray< Number >  in[],
VectorizedArray< Number >  out[] 
)
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 some previous results or not.


The documentation for this class was generated from the following file: