version 0.3.7
Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef > Class Template Reference

NonLinearElastic class represents a non-linear elastic finite element. More...

#include <ikarus/finiteelements/mechanics/nonlinearelastic.hh>

Inheritance diagram for Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >:
[legend]

Public Types

using Basis = Basis_
 
using Material = Material_
 
using FlatBasis = typename Basis::FlatBasis
 
using BasePowerFE = PowerBasisFE< FlatBasis >
 
using FERequirementType = FERequirements_
 
using ResultRequirementsType = ResultRequirements< FERequirementType >
 
using LocalView = typename FlatBasis::LocalView
 
using Element = typename LocalView::Element
 
using Geometry = typename Element::Geometry
 
using GridView = typename FlatBasis::GridView
 
using Traits = TraitsFromLocalView< LocalView, useEigenRef >
 
using LocalBasisType = decltype(std::declval< LocalView >().tree().child(0).finiteElement().localBasis())
 
using RootBasis = Basis_::FlatBasis
 Type of the root basis. More...
 
using GlobalIndex = typename LocalView::MultiIndex
 Type of the global index. More...
 
using GridElement = typename LocalView::Element
 Type of the grid element. More...
 

Public Member Functions

template<typename VolumeLoad = utils::LoadDefault, typename NeumannBoundaryLoad = utils::LoadDefault>
 NonLinearElastic (const Basis &globalBasis, const typename LocalView::Element &element, const Material &p_mat, VolumeLoad p_volumeLoad={}, const BoundaryPatch< GridView > *p_neumannBoundary=nullptr, NeumannBoundaryLoad p_neumannBoundaryLoad={})
 Constructor for the NonLinearElastic class. More...
 
template<typename ScalarType = double>
auto displacementFunction (const FERequirementType &par, const std::optional< const Eigen::VectorX< ScalarType > > &dx=std::nullopt) const
 Get the displacement function for the given FERequirementType. More...
 
template<typename ScalarType = double>
auto strainFunction (const FERequirementType &par, const std::optional< const Eigen::VectorX< ScalarType > > &dx=std::nullopt) const
 The strain function for the given FERequirementType. More...
 
template<typename ScalarType , int strainDim, bool voigt = true>
auto materialTangent (const Eigen::Vector< ScalarType, strainDim > &strain) const
 Get the material tangent for the given strain. More...
 
template<typename ScalarType , int strainDim>
auto getInternalEnergy (const Eigen::Vector< ScalarType, strainDim > &strain) const
 Get the internal energy for the given strain. More...
 
template<typename ScalarType , int strainDim, bool voigt = true>
auto getStress (const Eigen::Vector< ScalarType, strainDim > &strain) const
 Get the stress for the given strain. More...
 
double calculateScalar (const FERequirementType &par) const
 Calculate the scalar value associated with the given FERequirementType. More...
 
void calculateVector (const FERequirementType &par, typename Traits::template VectorType<> force) const
 Calculate the vector associated with the given FERequirementType. More...
 
void calculateMatrix (const FERequirementType &par, typename Traits::template MatrixType<> K) const
 Calculate the matrix associated with the given FERequirementType. More...
 
void calculateAt (const ResultRequirementsType &req, const Dune::FieldVector< double, Traits::mydim > &local, ResultTypeMap< double > &result) const
 Calculate specified results at a given local position. More...
 
constexpr size_t size () const
 Get the size of the local view. More...
 
void globalFlatIndices (std::vector< GlobalIndex > &globalIndices) const
 Get the global flat indices for the power basis. More...
 
const GridElementgridElement () const
 Get the grid element associated with the local view. More...
 
const LocalViewlocalView () const
 Get the const reference to the local view. More...
 
LocalViewlocalView ()
 Get the reference to the local view. More...
 

Public Attributes

std::shared_ptr< const Geometrygeo_
 
Dune::CachedLocalBasis< std::remove_cvref_t< LocalBasisType > > localBasis
 
std::function< Eigen::Vector< double, worldDim >(const Dune::FieldVector< double, worldDim > &, const double &)> volumeLoad
 
std::function< Eigen::Vector< double, worldDim >(const Dune::FieldVector< double, worldDim > &, const double &)> neumannBoundaryLoad
 
const BoundaryPatch< GridView > * neumannBoundary
 
Material mat
 
size_t numberOfNodes {0}
 
int order {}
 

Static Public Attributes

static constexpr int myDim = Traits::mydim
 
static constexpr int worldDim = Traits::worlddim
 
static constexpr auto strainType = StrainTags::greenLagrangian
 
static constexpr int num_children
 Number of children in the powerBasis. More...
 

Protected Member Functions

template<typename ScalarType >
auto calculateScalarImpl (const FERequirementType &par, const std::optional< const Eigen::VectorX< ScalarType > > &dx=std::nullopt) const -> ScalarType
 
template<typename ScalarType >
void calculateVectorImpl (const FERequirementType &par, typename Traits::template VectorType< ScalarType > force, const std::optional< const Eigen::VectorX< ScalarType > > &dx=std::nullopt) const
 

Detailed Description

template<typename Basis_, typename Material_, typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
class Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >
Template Parameters
Basis_The basis type for the finite element.
Material_The material type for the finite element.
FERequirements_The requirements for the finite element.
useEigenRefA boolean flag indicating whether to use Eigen references.

Member Typedef Documentation

◆ BasePowerFE

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::BasePowerFE = PowerBasisFE<FlatBasis>

◆ Basis

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::Basis = Basis_

◆ Element

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::Element = typename LocalView::Element

◆ FERequirementType

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::FERequirementType = FERequirements_

◆ FlatBasis

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::FlatBasis = typename Basis::FlatBasis

◆ Geometry

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::Geometry = typename Element::Geometry

◆ GlobalIndex

using Ikarus::PowerBasisFE< Basis_::FlatBasis >::GlobalIndex = typename LocalView::MultiIndex
inherited

◆ GridElement

using Ikarus::PowerBasisFE< Basis_::FlatBasis >::GridElement = typename LocalView::Element
inherited

◆ GridView

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::GridView = typename FlatBasis::GridView

◆ LocalBasisType

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::LocalBasisType = decltype(std::declval<LocalView>().tree().child(0).finiteElement().localBasis())

◆ LocalView

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::LocalView = typename FlatBasis::LocalView

◆ Material

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::Material = Material_

◆ ResultRequirementsType

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::ResultRequirementsType = ResultRequirements<FERequirementType>

◆ RootBasis

using Ikarus::PowerBasisFE< Basis_::FlatBasis >::RootBasis = Basis_::FlatBasis
inherited

◆ Traits

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
using Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::Traits = TraitsFromLocalView<LocalView, useEigenRef>

Constructor & Destructor Documentation

◆ NonLinearElastic()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename VolumeLoad = utils::LoadDefault, typename NeumannBoundaryLoad = utils::LoadDefault>
Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::NonLinearElastic ( const Basis globalBasis,
const typename LocalView::Element &  element,
const Material p_mat,
VolumeLoad  p_volumeLoad = {},
const BoundaryPatch< GridView > *  p_neumannBoundary = nullptr,
NeumannBoundaryLoad  p_neumannBoundaryLoad = {} 
)
inline
Template Parameters
VolumeLoadThe type for the volume load function.
NeumannBoundaryLoadThe type for the Neumann boundary load function.
Parameters
globalBasisThe global basis for the finite element.
elementThe element for which the finite element is constructed.
p_matThe material for the non-linear elastic element.
p_volumeLoadVolume load function (default is LoadDefault).
p_neumannBoundaryNeumann boundary patch (default is nullptr).
p_neumannBoundaryLoadNeumann boundary load function (default is LoadDefault).

Member Function Documentation

◆ calculateAt()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
void Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::calculateAt ( const ResultRequirementsType req,
const Dune::FieldVector< double, Traits::mydim > &  local,
ResultTypeMap< double > &  result 
) const
inline
Parameters
reqThe ResultRequirementsType object specifying the required results.
localThe local position for which results are to be calculated.
resultThe ResultTypeMap object to store the calculated results.

◆ calculateMatrix()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
void Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::calculateMatrix ( const FERequirementType par,
typename Traits::template MatrixType<>  K 
) const
inline
Template Parameters
ScalarTypeThe scalar type for the calculation.
Parameters
parThe FERequirementType object specifying the requirements for the calculation.
KThe matrix to store the calculated result.

◆ calculateScalar()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
double Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::calculateScalar ( const FERequirementType par) const
inline
Template Parameters
ScalarTypeThe scalar type for the calculation.
Parameters
parThe FERequirementType object specifying the requirements for the calculation.
Returns
The calculated scalar value.

◆ calculateScalarImpl()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename ScalarType >
auto Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::calculateScalarImpl ( const FERequirementType par,
const std::optional< const Eigen::VectorX< ScalarType > > &  dx = std::nullopt 
) const -> ScalarType
inlineprotected

◆ calculateVector()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
void Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::calculateVector ( const FERequirementType par,
typename Traits::template VectorType<>  force 
) const
inline
Template Parameters
ScalarTypeThe scalar type for the calculation.
Parameters
parThe FERequirementType object specifying the requirements for the calculation.
forceThe vector to store the calculated result.

◆ calculateVectorImpl()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename ScalarType >
void Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::calculateVectorImpl ( const FERequirementType par,
typename Traits::template VectorType< ScalarType >  force,
const std::optional< const Eigen::VectorX< ScalarType > > &  dx = std::nullopt 
) const
inlineprotected

◆ displacementFunction()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename ScalarType = double>
auto Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::displacementFunction ( const FERequirementType par,
const std::optional< const Eigen::VectorX< ScalarType > > &  dx = std::nullopt 
) const
inline
Template Parameters
ScalarTypeThe scalar type for the displacement function.
Parameters
parThe FERequirementType object.
dxOptional displacement vector.
Returns
A StandardLocalFunction representing the displacement function.
Here is the caller graph for this function:

◆ getInternalEnergy()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename ScalarType , int strainDim>
auto Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::getInternalEnergy ( const Eigen::Vector< ScalarType, strainDim > &  strain) const
inline
Template Parameters
ScalarTypeThe scalar type for the material and strain.
strainDimThe dimension of the strain vector.
Parameters
strainThe strain vector.
Returns
The internal energy calculated using the material's storedEnergy function.
Here is the caller graph for this function:

◆ getStress()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename ScalarType , int strainDim, bool voigt = true>
auto Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::getStress ( const Eigen::Vector< ScalarType, strainDim > &  strain) const
inline
Template Parameters
ScalarTypeThe scalar type for the material and strain.
strainDimThe dimension of the strain vector.
voigtA boolean indicating whether to use the Voigt notation for stress.
Parameters
strainThe strain vector.
Returns
The stress vector calculated using the material's stresses function.
Here is the caller graph for this function:

◆ globalFlatIndices()

void Ikarus::PowerBasisFE< Basis_::FlatBasis >::globalFlatIndices ( std::vector< GlobalIndex > &  globalIndices) const
inlineinherited

The global indices are collected in a FlatInterLeaved order. I.e. x_0, y_0, z_0, ..., x_n, y_n, z_n

Parameters
globalIndicesOutput vector to store global indices.

◆ gridElement()

const GridElement & Ikarus::PowerBasisFE< Basis_::FlatBasis >::gridElement ( ) const
inlineinherited
Returns
The grid element.

◆ localView() [1/2]

LocalView & Ikarus::PowerBasisFE< Basis_::FlatBasis >::localView ( )
inlineinherited
Returns
The reference to the local view.

◆ localView() [2/2]

const LocalView & Ikarus::PowerBasisFE< Basis_::FlatBasis >::localView ( ) const
inlineinherited
Returns
The const reference to the local view.

◆ materialTangent()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename ScalarType , int strainDim, bool voigt = true>
auto Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::materialTangent ( const Eigen::Vector< ScalarType, strainDim > &  strain) const
inline
Template Parameters
ScalarTypeThe scalar type for the material and strain.
strainDimThe dimension of the strain vector.
voigtFlag indicating whether to use Voigt notation.
Parameters
strainThe strain vector.
Returns
The material tangent calculated using the material's tangentModuli function.
Here is the caller graph for this function:

◆ size()

constexpr size_t Ikarus::PowerBasisFE< Basis_::FlatBasis >::size ( ) const
inlineconstexprinherited
Returns
The size of the local view.

◆ strainFunction()

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
template<typename ScalarType = double>
auto Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::strainFunction ( const FERequirementType par,
const std::optional< const Eigen::VectorX< ScalarType > > &  dx = std::nullopt 
) const
inline
Template Parameters
ScalarTypeThe scalar type for the strain function.
Parameters
parThe FERequirementType object.
dxOptional displacement vector.
Returns
The strain function calculated using greenLagrangeStrains.
Here is the caller graph for this function:

Member Data Documentation

◆ geo_

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
std::shared_ptr<const Geometry> Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::geo_

◆ localBasis

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
Dune::CachedLocalBasis<std::remove_cvref_t<LocalBasisType> > Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::localBasis

◆ mat

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
Material Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::mat

◆ myDim

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
constexpr int Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::myDim = Traits::mydim
staticconstexpr

◆ neumannBoundary

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
const BoundaryPatch<GridView>* Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::neumannBoundary

◆ neumannBoundaryLoad

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
std::function<Eigen::Vector<double, worldDim>(const Dune::FieldVector<double, worldDim>&, const double&)> Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::neumannBoundaryLoad

◆ num_children

constexpr int Ikarus::PowerBasisFE< Basis_::FlatBasis >::num_children
staticconstexprinherited

◆ numberOfNodes

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
size_t Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::numberOfNodes {0}

◆ order

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
int Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::order {}

◆ strainType

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
constexpr auto Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::strainType = StrainTags::greenLagrangian
staticconstexpr

◆ volumeLoad

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
std::function<Eigen::Vector<double, worldDim>(const Dune::FieldVector<double, worldDim>&, const double&)> Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::volumeLoad

◆ worldDim

template<typename Basis_ , typename Material_ , typename FERequirements_ = FERequirements<>, bool useEigenRef = false>
constexpr int Ikarus::NonLinearElastic< Basis_, Material_, FERequirements_, useEigenRef >::worldDim = Traits::worlddim
staticconstexpr

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