deal.II version GIT relicensing-2901-g19332422bd 2025-03-23 19:50:00+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Static Public Attributes | Private Attributes | List of all members
Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number > Class Template Reference

#include <deal.II/matrix_free/portable_fe_evaluation.h>

Public Types

using value_type = std::conditional_t<(n_components_==1), Number, Tensor< 1, n_components_, Number > >
 
using gradient_type = std::conditional_t< n_components_==1, Tensor< 1, dim, Number >, std::conditional_t< n_components_==dim, Tensor< 2, dim, Number >, Tensor< 1, n_components_, Tensor< 1, dim, Number > > > >
 
using data_type = typename MatrixFree< dim, Number >::Data
 

Public Member Functions

 FEEvaluation (const data_type *data)
 
int get_current_cell_index ()
 
const data_typeget_matrix_free_data ()
 
void read_dof_values (const Number *src)
 
void distribute_local_to_global (Number *dst) const
 
void evaluate (const EvaluationFlags::EvaluationFlags evaluate_flag)
 
void integrate (const EvaluationFlags::EvaluationFlags integration_flag)
 
value_type get_value (int q_point) const
 
value_type get_dof_value (int q_point) const
 
void submit_value (const value_type &val_in, int q_point)
 
void submit_dof_value (const value_type &val_in, int q_point)
 
gradient_type get_gradient (int q_point) const
 
void submit_gradient (const gradient_type &grad_in, int q_point)
 
template<typename Functor >
void apply_for_each_quad_point (const Functor &func)
 

Static Public Attributes

static constexpr unsigned int dimension = dim
 
static constexpr unsigned int n_components = n_components_
 
static constexpr unsigned int n_q_points
 
static constexpr unsigned int tensor_dofs_per_component
 
static constexpr unsigned int tensor_dofs_per_cell
 

Private Attributes

const MatrixFree< dim, Number >::Data * data
 
const MatrixFree< dim, Number >::PrecomputedData * precomputed_data
 
SharedData< dim, Number > * shared_data
 
int cell_id
 

Detailed Description

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
class Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >

This class provides all the functions necessary to evaluate functions at quadrature points and cell integrations. In functionality, this class is similar to FEValues<dim>.

This class has five template arguments:

Template Parameters
dimDimension in which this class is to be used
fe_degreeDegree of the tensor prodict finite element with fe_degree+1 degrees of freedom per coordinate direction
n_q_points_1dNumber of points in the quadrature formular in 1d, defaults to fe_degree+1
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). Defaults to 1
NumberNumber format, double or float. Defaults to double.

Definition at line 68 of file portable_fe_evaluation.h.

Member Typedef Documentation

◆ value_type

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
using Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::value_type = std::conditional_t<(n_components_ == 1), Number, Tensor<1, n_components_, Number> >

An alias for scalar quantities.

Definition at line 74 of file portable_fe_evaluation.h.

◆ gradient_type

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
using Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::gradient_type = std::conditional_t< n_components_ == 1, Tensor<1, dim, Number>, std::conditional_t<n_components_ == dim, Tensor<2, dim, Number>, Tensor<1, n_components_, Tensor<1, dim, Number> >> >

An alias for vectorial quantities.

Definition at line 81 of file portable_fe_evaluation.h.

◆ data_type

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
using Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::data_type = typename MatrixFree<dim, Number>::Data

An alias to kernel specific information.

Definition at line 91 of file portable_fe_evaluation.h.

Constructor & Destructor Documentation

◆ FEEvaluation()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::FEEvaluation ( const data_type data)
explicit

Constructor.

Definition at line 257 of file portable_fe_evaluation.h.

Member Function Documentation

◆ get_current_cell_index()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
int FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::get_current_cell_index ( )

Return the index of the current cell.

Definition at line 273 of file portable_fe_evaluation.h.

◆ get_matrix_free_data()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
const FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::data_type * FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::get_matrix_free_data ( )

Return a pointer to the MatrixFree<dim, Number>::Data object on device that contains necessary constraint, dof index, and shape function information for evaluation used in the matrix-free kernels.

Definition at line 291 of file portable_fe_evaluation.h.

◆ read_dof_values()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::read_dof_values ( const Number *  src)

For the vector src, read out the values on the degrees of freedom of the current cell, and store them internally. Similar functionality as the function DoFAccessor::get_interpolated_dof_values when no constraints are present, but it also includes constraints from hanging nodes, so once can see it as a similar function to AffineConstraints::read_dof_values() as well.

Definition at line 305 of file portable_fe_evaluation.h.

◆ distribute_local_to_global()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::distribute_local_to_global ( Number *  dst) const

Take the value stored internally on dof values of the current cell and sum them into the vector dst. The function also applies constraints during the write operation. The functionality is hence similar to the function AffineConstraints::distribute_local_to_global.

Definition at line 337 of file portable_fe_evaluation.h.

◆ evaluate()

template<int dim, int fe_degree, int n_q_points_1d, int n_components, typename Number >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components, Number >::evaluate ( const EvaluationFlags::EvaluationFlags  evaluate_flag)

Evaluate the function values and the gradients of the FE function given at the DoF values in the input vector at the quadrature points on the unit cell. The function argument evaluate_flag specifies which parts shall actually be computed. This function needs to be called before the functions get_value() or get_gradient() give useful information.

Definition at line 381 of file portable_fe_evaluation.h.

◆ integrate()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::integrate ( const EvaluationFlags::EvaluationFlags  integration_flag)

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 as specified by the integration_flag argument.

Definition at line 426 of file portable_fe_evaluation.h.

◆ get_value()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::value_type FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::get_value ( int  q_point) const

Same as above, except that the quadrature point is computed from thread id.

Definition at line 475 of file portable_fe_evaluation.h.

◆ get_dof_value()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::value_type FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::get_dof_value ( int  q_point) const

Same as above, except that the local dof index is computed from the thread id.

Definition at line 504 of file portable_fe_evaluation.h.

◆ submit_value()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::submit_value ( const value_type val_in,
int  q_point 
)

Same as above, except that the quadrature point is computed from the thread id.

Definition at line 529 of file portable_fe_evaluation.h.

◆ submit_dof_value()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::submit_dof_value ( const value_type val_in,
int  q_point 
)

Same as above, except that the local dof index is computed from the thread id.

Definition at line 554 of file portable_fe_evaluation.h.

◆ get_gradient()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::gradient_type FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::get_gradient ( int  q_point) const

Same as above, except that the quadrature point is computed from the thread id.

Definition at line 581 of file portable_fe_evaluation.h.

◆ submit_gradient()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::submit_gradient ( const gradient_type grad_in,
int  q_point 
)

Same as above, except that the quadrature point is computed from the thread id.

Definition at line 624 of file portable_fe_evaluation.h.

◆ apply_for_each_quad_point()

template<int dim, int fe_degree, int n_q_points_1d, int n_components_, typename Number >
template<typename Functor >
void FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::apply_for_each_quad_point ( const Functor &  func)

Same as above, except that the functor func only takes a single input argument (fe_eval) and computes the quadrature point from the thread id.

func needs to define

Definition at line 666 of file portable_fe_evaluation.h.

Member Data Documentation

◆ dimension

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
constexpr unsigned int Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::dimension = dim
staticconstexpr

Dimension.

Definition at line 96 of file portable_fe_evaluation.h.

◆ n_components

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
constexpr unsigned int Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::n_components = n_components_
staticconstexpr

Number of components.

Definition at line 101 of file portable_fe_evaluation.h.

◆ n_q_points

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
constexpr unsigned int Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::n_q_points
staticconstexpr
Initial value:
=
Utilities::pow(n_q_points_1d, dim)
constexpr T pow(const T base, const int iexp)
Definition utilities.h:967

Number of quadrature points per cell.

Definition at line 106 of file portable_fe_evaluation.h.

◆ tensor_dofs_per_component

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
constexpr unsigned int Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::tensor_dofs_per_component
staticconstexpr
Initial value:
=
Utilities::pow(fe_degree + 1, dim)

Number of tensor degrees of freedom of a scalar component determined from the given template argument fe_degree.

Definition at line 113 of file portable_fe_evaluation.h.

◆ tensor_dofs_per_cell

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
constexpr unsigned int Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::tensor_dofs_per_cell
staticconstexpr
Initial value:
=
static constexpr unsigned int n_components
static constexpr unsigned int tensor_dofs_per_component

Number of tensor degrees of freedom of all component determined from the given template argument fe_degree.

Definition at line 120 of file portable_fe_evaluation.h.

◆ data

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
const MatrixFree<dim,Number>::Data* Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::data
private

Definition at line 243 of file portable_fe_evaluation.h.

◆ precomputed_data

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
const MatrixFree<dim,Number>::PrecomputedData* Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::precomputed_data
private

Definition at line 244 of file portable_fe_evaluation.h.

◆ shared_data

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
SharedData<dim, Number>* Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::shared_data
private

Definition at line 245 of file portable_fe_evaluation.h.

◆ cell_id

template<int dim, int fe_degree, int n_q_points_1d = fe_degree + 1, int n_components_ = 1, typename Number = double>
int Portable::FEEvaluation< dim, fe_degree, n_q_points_1d, n_components_, Number >::cell_id
private

Definition at line 246 of file portable_fe_evaluation.h.


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