hyper.deal
Loading...
Searching...
No Matches
hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType > Class Template Reference

#include <fe_evaluation_cell.h>

Inheritance diagram for hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >:
hyperdeal::FEEvaluationBase< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType > hyperdeal::FEEvaluationInverse< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >

Public Types

using PARENT
 
- Public Types inherited from hyperdeal::FEEvaluationBase< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >
using NUMBER_ = Number
 
using VEC_NUMBER_ = VectorizedArrayType
 
using MF = MatrixFree<dim_x, dim_v, Number, VectorizedArrayType>
 
using VectorizedArrayTypeX = typename MF::VectorizedArrayTypeX
 
using VectorizedArrayTypeV = typename MF::VectorizedArrayTypeV
 

Public Member Functions

 FEEvaluation (const MatrixFree< dim_x, dim_v, Number, VectorizedArrayType > &matrix_free, const unsigned int dof_no_x, const unsigned int dof_no_v, const unsigned int quad_no_x, const unsigned int quad_no_v)
 
void reinit (typename MatrixFree< dim_x, dim_v, Number, VectorizedArrayType >::ID cell_index) override
 
void read_dof_values (const dealii::LinearAlgebra::distributed::Vector< Number > &src, VectorizedArrayType *data) const
 
void read_dof_values (const dealii::LinearAlgebra::distributed::Vector< Number > &src)
 
void set_dof_values (dealii::LinearAlgebra::distributed::Vector< Number > &dst) const
 
unsigned int n_vectorization_lanes_filled () const
 
VectorizedArrayType submit_inplace (const VectorizedArrayType data, const unsigned int q) const
 
VectorizedArrayType submit_inplace (const VectorizedArrayType data, unsigned int, unsigned int qx, unsigned int qv) const
 
DEAL_II_ALWAYS_INLINE VectorizedArrayType JxW (unsigned int qx, unsigned int qv) const
 
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, dim_x, VectorizedArrayType > get_gradient_x (const VectorizedArrayType *__restrict grad_in, const unsigned int q, const unsigned int qx, const unsigned int qv) const
 
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, dim_v, VectorizedArrayType > get_gradient_v (const VectorizedArrayType *__restrict grad_in, const unsigned int q, const unsigned int qx, const unsigned int qv) const
 
template<bool do_add>
DEAL_II_ALWAYS_INLINE void submit_value (VectorizedArrayType *__restrict data_ptr_out, const VectorizedArrayType values_in, const unsigned int q, const unsigned int qx, const unsigned int qv) const
 
DEAL_II_ALWAYS_INLINE void submit_gradient_x (VectorizedArrayType *__restrict data_ptr_out, const VectorizedArrayType *__restrict grad_in, const unsigned int q, const unsigned int qx, const unsigned int qv) const
 
DEAL_II_ALWAYS_INLINE void submit_gradient_v (VectorizedArrayType *__restrict data_ptr_out, const VectorizedArrayType *__restrict grad_in, const unsigned int q, const unsigned int qx, const unsigned int qv) const
 
dealii::Point< dim, VectorizedArrayType > get_quadrature_point (const unsigned int q) const
 
dealii::Point< dim, VectorizedArrayType > get_quadrature_point (const unsigned int qx, const unsigned int qv) const
 
dealii::Point< dim_v, VectorizedArrayType > get_quadrature_point_v (const unsigned int qv) const
 
- Public Member Functions inherited from hyperdeal::FEEvaluationBase< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >
 FEEvaluationBase (const MatrixFree< dim_x, dim_v, Number, VectorizedArrayType > &matrix_free)
 
virtual ~FEEvaluationBase ()=default
 
VectorizedArrayType * get_data_ptr ()
 
const dealii::AlignedVector< Number > * get_shape_values () const
 
const dealii::AlignedVector< Number > * get_shape_gradients () const
 

Static Public Attributes

static const unsigned int dim = dim_x + dim_v
 
static const unsigned int static_dofs_per_cell_x = dealii::Utilities::pow(degree + 1, dim_x)
 
static const unsigned int static_dofs_per_cell_v = dealii::Utilities::pow(degree + 1, dim_v)
 
static const unsigned int static_dofs_per_cell = static_dofs_per_cell_x * static_dofs_per_cell_v
 
static const unsigned int n_q_points_x = dealii::Utilities::pow(n_points, dim_x)
 
static const unsigned int n_q_points_v = dealii::Utilities::pow(n_points, dim_v)
 
static const unsigned int n_q_points = n_q_points_x * n_q_points_v
 
static const unsigned int n_vectors = PARENT::n_vectors
 
static const unsigned int n_vectors_v = PARENT::n_vectors_v
 
- Static Public Attributes inherited from hyperdeal::FEEvaluationBase< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >
static const int dim = dim_x + dim_v
 
static const unsigned int n_vectors = MF::VectorizedArrayTypeX::size()
 
static const unsigned int n_vectors_v = MF::VectorizedArrayTypeV::size()
 
static const int static_dofs
 

Protected Attributes

dealii::FEEvaluation< dim_x, degree, n_points, 1, Number, typename PARENT::VectorizedArrayTypeX > phi_x
 
dealii::FEEvaluation< dim_v, degree, n_points, 1, Number, typename PARENT::VectorizedArrayTypeV > phi_v
 
- Protected Attributes inherited from hyperdeal::FEEvaluationBase< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >
dealii::AlignedVector< VectorizedArrayType > data
 
const MFmatrix_free
 
const dealii::MatrixFree< dim_x, Number, VectorizedArrayTypeX > & matrix_free_x
 
const dealii::MatrixFree< dim_v, Number, VectorizedArrayTypeV > & matrix_free_v
 
unsigned int macro_cell_x
 
unsigned int macro_cell_v
 
unsigned int lane_y
 
unsigned int macro
 
const dealii::AlignedVector< Number > * shape_values
 
const dealii::AlignedVector< Number > * shape_gradients
 

Detailed Description

template<int dim_x, int dim_v, int degree, int n_points, typename Number, typename VectorizedArrayType>
class hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >

The class that provides all functions necessary to evaluate functions at quadrature points and cell integrations in phase space. It delegates the the actual evaluation task to two dealii::FEEvaluation objects and combines the result on-the-fly.

Member Typedef Documentation

◆ PARENT

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
using hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::PARENT
Initial value:
dim_v,
degree,
n_points,
Number,
VectorizedArrayType>
FEEvaluationBase(const MatrixFree< dim_x, dim_v, Number, VectorizedArrayType > &matrix_free)
Definition fe_evaluation_base.h:144

Constructor & Destructor Documentation

◆ FEEvaluation()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::FEEvaluation ( const MatrixFree< dim_x, dim_v, Number, VectorizedArrayType > & matrix_free,
const unsigned int dof_no_x,
const unsigned int dof_no_v,
const unsigned int quad_no_x,
const unsigned int quad_no_v )
inline

Constructor.

Parameters
matrix_freeData object that contains all data.
dof_no_xIf x-space matrix_free of matrix_free was set up with multiple DoFHandler objects, this parameter selects to which DoFHandler/AffineConstraints pair the given evaluator should be attached to.
dof_no_vSame as above but for v-space.
quad_no_xIf x-space matrix_free of matrix_free was set up with multiple Quadrature objects, this parameter selects the appropriate number of the quadrature formula.
quad_no_vSame as above but for v-space.

Member Function Documentation

◆ get_gradient_v()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, dim_v, VectorizedArrayType > hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::get_gradient_v ( const VectorizedArrayType *__restrict grad_in,
const unsigned int q,
const unsigned int qx,
const unsigned int qv ) const
inline

Get gradient in v-space

◆ get_gradient_x()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, dim_x, VectorizedArrayType > hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::get_gradient_x ( const VectorizedArrayType *__restrict grad_in,
const unsigned int q,
const unsigned int qx,
const unsigned int qv ) const
inline

Get gradient in x-space

◆ get_quadrature_point()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
dealii::Point< dim, VectorizedArrayType > hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::get_quadrature_point ( const unsigned int q) const
inline

Return location of quadrature point.

TODO: split up in x- and v-space.

◆ n_vectorization_lanes_filled()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
unsigned int hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::n_vectorization_lanes_filled ( ) const
inline

Return number of filled lanes.

◆ read_dof_values() [1/2]

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
void hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::read_dof_values ( const dealii::LinearAlgebra::distributed::Vector< Number > & src)
inline

Read values from src vector and write them into the internal buffer.

◆ read_dof_values() [2/2]

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
void hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::read_dof_values ( const dealii::LinearAlgebra::distributed::Vector< Number > & src,
VectorizedArrayType * data ) const
inline

Read values from src vector and write it into the buffer data.

◆ reinit()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
void hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::reinit ( typename MatrixFree< dim_x, dim_v, Number, VectorizedArrayType >::ID cell_index)
inlineoverridevirtual

◆ set_dof_values()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
void hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::set_dof_values ( dealii::LinearAlgebra::distributed::Vector< Number > & dst) const
inline

Sum values from the internal buffer into the vector dst.

◆ submit_gradient_v()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
DEAL_II_ALWAYS_INLINE void hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::submit_gradient_v ( VectorizedArrayType *__restrict data_ptr_out,
const VectorizedArrayType *__restrict grad_in,
const unsigned int q,
const unsigned int qx,
const unsigned int qv ) const
inline

Submit gradient in v-space

◆ submit_gradient_x()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
DEAL_II_ALWAYS_INLINE void hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::submit_gradient_x ( VectorizedArrayType *__restrict data_ptr_out,
const VectorizedArrayType *__restrict grad_in,
const unsigned int q,
const unsigned int qx,
const unsigned int qv ) const
inline

Submit gradient in x-space

◆ submit_inplace() [1/2]

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
VectorizedArrayType hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::submit_inplace ( const VectorizedArrayType data,
const unsigned int q ) const
inline

Return product of data and |J|xw.

TODO: remove.

◆ submit_inplace() [2/2]

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
VectorizedArrayType hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::submit_inplace ( const VectorizedArrayType data,
unsigned int ,
unsigned int qx,
unsigned int qv ) const
inline

Return product of data and |J|xw.

TODO: remove.

◆ submit_value()

template<int dim_x, int dim_v, int degree, int n_points, typename Number , typename VectorizedArrayType >
template<bool do_add>
DEAL_II_ALWAYS_INLINE void hyperdeal::FEEvaluation< dim_x, dim_v, degree, n_points, Number, VectorizedArrayType >::submit_value ( VectorizedArrayType *__restrict data_ptr_out,
const VectorizedArrayType values_in,
const unsigned int q,
const unsigned int qx,
const unsigned int qv ) const
inline

Submit value for integration along x-v coordinates


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