DOLFINx 0.9.0
DOLFINx C++ interface
Loading...
Searching...
No Matches
Expression< T, U > Class Template Reference

Represents a mathematical expression evaluated at a pre-defined set of points on the reference cell. More...

#include <Expression.h>

Public Types

using scalar_type = T
 Scalar type.
 
using geometry_type = U
 Geometry type of the points.
 

Public Member Functions

 Expression (const std::vector< std::shared_ptr< const Function< scalar_type, geometry_type > > > &coefficients, const std::vector< std::shared_ptr< const Constant< scalar_type > > > &constants, std::span< const geometry_type > X, std::array< std::size_t, 2 > Xshape, std::function< void(scalar_type *, const scalar_type *, const scalar_type *, const geometry_type *, const int *, const uint8_t *)> fn, const std::vector< int > &value_shape, std::shared_ptr< const FunctionSpace< geometry_type > > argument_function_space=nullptr)
 Create an Expression.
 
 Expression (Expression &&e)=default
 Move constructor.
 
virtual ~Expression ()=default
 Destructor.
 
std::shared_ptr< const FunctionSpace< geometry_type > > argument_function_space () const
 Get argument function space.
 
const std::vector< std::shared_ptr< const Function< scalar_type, geometry_type > > > & coefficients () const
 Get coefficients.
 
const std::vector< std::shared_ptr< const Constant< scalar_type > > > & constants () const
 Get constants.
 
std::vector< int > coefficient_offsets () const
 Offset for each coefficient expansion array on a cell.
 
void eval (const mesh::Mesh< geometry_type > &mesh, std::span< const std::int32_t > entities, std::span< scalar_type > values, std::array< std::size_t, 2 > vshape) const
 Evaluate Expression on cells or facets.
 
const std::function< void(scalar_type *, const scalar_type *, const scalar_type *, const geometry_type *, const int *, const uint8_t *)> & get_tabulate_expression () const
 Get function for tabulate_expression.
 
int value_size () const
 Get value size.
 
const std::vector< int > & value_shape () const
 Get value shape.
 
std::pair< std::vector< geometry_type >, std::array< std::size_t, 2 > > X () const
 Evaluation points on the reference cell.
 

Detailed Description

template<dolfinx::scalar T, std::floating_point U>
class dolfinx::fem::Expression< T, U >

Represents a mathematical expression evaluated at a pre-defined set of points on the reference cell.

This class closely follows the concept of a UFC Expression.

The functionality can be used to evaluate a gradient of a Function at quadrature points in all cells. This evaluated gradient can then be used as input in to a non-FEniCS function that calculates a material constitutive model.

Template Parameters
TThe scalar type
UThe mesh geometry scalar type

Member Typedef Documentation

◆ scalar_type

template<dolfinx::scalar T, std::floating_point U>
using scalar_type = T

Scalar type.

Field type for the Expression, e.g. double, std::complex<float>, etc.

Constructor & Destructor Documentation

◆ Expression()

template<dolfinx::scalar T, std::floating_point U>
Expression ( const std::vector< std::shared_ptr< const Function< scalar_type, geometry_type > > > & coefficients,
const std::vector< std::shared_ptr< const Constant< scalar_type > > > & constants,
std::span< const geometry_type > X,
std::array< std::size_t, 2 > Xshape,
std::function< void(scalar_type *, const scalar_type *, const scalar_type *, const geometry_type *, const int *, const uint8_t *)> fn,
const std::vector< int > & value_shape,
std::shared_ptr< const FunctionSpace< geometry_type > > argument_function_space = nullptr )
inline

Create an Expression.

Note
Users should prefer the create_expression factory functions.
Parameters
[in]coefficientsCoefficients in the Expression.
[in]constantsConstants in the Expression
[in]XPoints on the reference cell, shape=(number of / points, tdim) and storage is row-major.
[in]XshapeShape of X.
[in]fnFunction for tabulating the Expression.
[in]value_shapeShape of Expression evaluated at single point.
[in]argument_function_spaceFunction space for Argument.

Member Function Documentation

◆ argument_function_space()

template<dolfinx::scalar T, std::floating_point U>
std::shared_ptr< const FunctionSpace< geometry_type > > argument_function_space ( ) const
inline

Get argument function space.

Returns
The argument function space, nullptr if there is no argument.

◆ coefficient_offsets()

template<dolfinx::scalar T, std::floating_point U>
std::vector< int > coefficient_offsets ( ) const
inline

Offset for each coefficient expansion array on a cell.

Used to pack data for multiple coefficients in a flat array. The last entry is the size required to store all coefficients.

Returns
The offsets.

◆ coefficients()

template<dolfinx::scalar T, std::floating_point U>
const std::vector< std::shared_ptr< const Function< scalar_type, geometry_type > > > & coefficients ( ) const
inline

Get coefficients.

Returns
Vector of attached coefficients.

◆ constants()

template<dolfinx::scalar T, std::floating_point U>
const std::vector< std::shared_ptr< const Constant< scalar_type > > > & constants ( ) const
inline

Get constants.

Returns
Vector of attached constants with their names. Names are used to set constants in user's c++ code. Index in the vector is the position of the constant in the original (nonsimplified) form.

◆ eval()

template<dolfinx::scalar T, std::floating_point U>
void eval ( const mesh::Mesh< geometry_type > & mesh,
std::span< const std::int32_t > entities,
std::span< scalar_type > values,
std::array< std::size_t, 2 > vshape ) const
inline

Evaluate Expression on cells or facets.

Parameters
[in]meshCells on which to evaluate the Expression.
[in]entitiesList of entities to evaluate the expression on. This could be either a list of cells or a list of (cell, local
[out]valuesA 2D array to store the result. Caller is responsible for correct sizing which should be (num_cells, / num_points * value_size * num_all_argument_dofs columns). facet index) tuples. Array is flattened per entity.
[in]vshapeThe shape of values (row-major storage).

◆ get_tabulate_expression()

template<dolfinx::scalar T, std::floating_point U>
const std::function< void(scalar_type *, const scalar_type *, const scalar_type *, const geometry_type *, const int *, const uint8_t *)> & get_tabulate_expression ( ) const
inline

Get function for tabulate_expression.

Returns
fn Function to tabulate expression.

◆ value_shape()

template<dolfinx::scalar T, std::floating_point U>
const std::vector< int > & value_shape ( ) const
inline

Get value shape.

Returns
The value shape.

◆ value_size()

template<dolfinx::scalar T, std::floating_point U>
int value_size ( ) const
inline

Get value size.

Returns
The value size.

◆ X()

template<dolfinx::scalar T, std::floating_point U>
std::pair< std::vector< geometry_type >, std::array< std::size_t, 2 > > X ( ) const
inline

Evaluation points on the reference cell.

Returns
Evaluation points.

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