PhysicsBasedAnimationToolkit 0.0.10
Cross-platform C++20 library of algorithms and data structures commonly used in computer graphics research on physically-based simulation.
Loading...
Searching...
No Matches
pbat::fem::LoadVector< TMesh, QuadratureOrder > Struct Template Reference

A matrix-free representation of a finite element load vector \( \mathbf{f}_i = \int_\Omega \mathbf{f}(X) \phi_i(X) \) under Galerkin projection. More...

#include <LoadVector.h>

Public Types

using SelfType = LoadVector<TMesh, QuadratureOrder>
 Self type.
 
using MeshType = TMesh
 Mesh type.
 
using ElementType = typename TMesh::ElementType
 Element type.
 
using QuadratureRuleType
 Quadrature rule type.
 

Public Member Functions

template<class TDerived>
 LoadVector (MeshType const &mesh, Eigen::Ref< MatrixX const > const &detJe, Eigen::DenseBase< TDerived > const &fe, int dims=1)
 Construct a new LoadVector object.
 
SelfTypeoperator= (SelfType const &)=delete
 
VectorX ToVector () const
 Transforms this element-wise load vector representation into a global load vector.
 
template<class TDerived>
void SetLoad (Eigen::DenseBase< TDerived > const &fe)
 Set the external loading.
 
void CheckValidState () const
 Check if the state of this load vector is valid.
 

Public Attributes

MeshType const & mesh
 The finite element mesh.
 
MatrixX fe
 |dims|x|# elements| piecewise element constant load
 
MatrixX N
 
Eigen::Ref< MatrixX const > detJe
 
int dims
 Dimensionality of external loading. Should have dims >= 1.
 

Static Public Attributes

static int constexpr kOrder = ElementType::kOrder
 Polynomial order of the load vector.
 
static int constexpr kQuadratureOrder = QuadratureOrder
 Quadrature order.
 

Detailed Description

template<CMesh TMesh, int QuadratureOrder>
struct pbat::fem::LoadVector< TMesh, QuadratureOrder >

A matrix-free representation of a finite element load vector \( \mathbf{f}_i = \int_\Omega \mathbf{f}(X) \phi_i(X) \) under Galerkin projection.

Note
Assumes \( \mathbf{f}(X) \) is piecewise constant over the elements.
Todo
Link to my higher-level FEM crash course doc.
Template Parameters
TMeshType satisfying concept CMesh
QuadratureOrderQuadrature order for integrating the load vector

Member Typedef Documentation

◆ QuadratureRuleType

template<CMesh TMesh, int QuadratureOrder>
using pbat::fem::LoadVector< TMesh, QuadratureOrder >::QuadratureRuleType
Initial value:
typename ElementType::template QuadratureType<QuadratureOrder>

Quadrature rule type.

Constructor & Destructor Documentation

◆ LoadVector()

template<CMesh TMesh, int QuadratureOrder>
template<class TDerived>
pbat::fem::LoadVector< TMesh, QuadratureOrder >::LoadVector ( MeshType const & mesh,
Eigen::Ref< MatrixX const > const & detJe,
Eigen::DenseBase< TDerived > const & fe,
int dims = 1 )
inline

Construct a new LoadVector object.

Template Parameters
TDerivedEigen dense expression type
Parameters
meshFinite element mesh
detJe|# quad.pts.|x|# elements| affine element jacobian determinants at quadrature points
fe|dims|x|# elements| piecewise element constant load, or |dims|x1 constant load
dimsDimensionality of image of FEM function space
Precondition
dims >= 1 and fe.cols() == 1 || fe.cols() == mesh.E.cols() and fe.rows() == dims

Member Function Documentation

◆ SetLoad()

template<CMesh TMesh, int QuadratureOrder>
template<class TDerived>
void pbat::fem::LoadVector< TMesh, QuadratureOrder >::SetLoad ( Eigen::DenseBase< TDerived > const & fe)
inline

Set the external loading.

Template Parameters
TDerivedEigen dense expression type
Parameters
fe|dims|x|# elements| piecewise element constant load, or |dims|x1 constant load

◆ ToVector()

template<CMesh TMesh, int QuadratureOrder>
VectorX pbat::fem::LoadVector< TMesh, QuadratureOrder >::ToVector ( ) const
inline

Transforms this element-wise load vector representation into a global load vector.

Returns
Global load vector

Member Data Documentation

◆ detJe

template<CMesh TMesh, int QuadratureOrder>
Eigen::Ref<MatrixX const> pbat::fem::LoadVector< TMesh, QuadratureOrder >::detJe

|# element quadrature points|x|# elements| matrix of jacobian determinants at element quadrature points

◆ N

template<CMesh TMesh, int QuadratureOrder>
MatrixX pbat::fem::LoadVector< TMesh, QuadratureOrder >::N

|ElementType::kNodes|x|# elements| integrated element shape functions. See (IntegratedShapeFunctions()).


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