QuOc

 

Public Member Functions | Static Public Member Functions | Protected Types | Protected Attributes

qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType > Class Template Reference

This class computes via "apply(...)" or "applyAdd(...)" parts of a hyperelastic energy Hessian with respect to the displacement, i.e. $ (<\frac{\partial^2 \int_\Omega W(I_1,I_2,I_3,x) dx}{\partial \phi^2},\psi_i,\psi_j>)_{i,j} $, where the $ \psi_i $ and $ \psi_i $ are those vector-valued FE basis functions that are only nonzero in the kth and lth component, the hyperelastic invariants are given by $ I_1=||\nabla\phi||_F^2 $, $ I_2=||cof(\nabla\phi)||_F^2 $, $ I_3=det(\nabla\phi)$, and the function $ W(I_1,I_2,I_3,x) $ is passed to the constructor (e.g. HyerelasticEnergyDensityDefault). In its first argument, "apply" and "applyAdd" expect a Multivector which represents the displacements $ d_i $ in the different Cartesian directions; the deformation $ \phi $ is then automatically computed via $ \phi=d+$ identity. More...

#include <hyperelastic.h>

Inheritance diagram for qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >:
aol::FELinAsymMatrixWeightedStiffInterface< ConfiguratorType, HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType > > aol::FELinOpInterface< ConfiguratorType::RealType, ConfiguratorType, FELinAsymMatrixWeightedStiffInterface< ConfiguratorType, HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >, QUOC_GRID_INDEX_MODE >, QUOC_GRID_INDEX_MODE > aol::FEOpInterface< ConfiguratorType, aol::Vector< ConfiguratorType::RealType > > aol::Op< aol::Vector< ConfiguratorType::RealType >, aol::Vector< ConfiguratorType::RealType > >

List of all members.

Public Member Functions

 HyperelasticSubHessian (const typename ConfiguratorType::InitType &Grid, const HyperelasticEnergyDensityType &HyperelasticEnergyDensity, const aol::MultiVector< RealType > &Displacement, const int K, const int L)
void getCoeffMatrix (const typename ConfiguratorType::ElementType &El, int QuadPoint, const typename ConfiguratorType::DomVecType &, typename ConfiguratorType::MatType &Matrix) const

Static Public Member Functions

static void elasticityTensor (const typename ConfiguratorType::MatType &DPhi, const HyperelasticEnergyDensityType &EnergyDensity, const typename ConfiguratorType::ElementType &El, int QuadPoint, const int K, const int L, aol::Mat< ConfiguratorType::Dim, ConfiguratorType::Dim, RealType > &ElasticityTensor)
 Computes ElasticityTensor= $ \frac{\partial^2W}{\partial\nabla\phi_K\partial\nabla\phi_L}(I_1,I_2,I_3) $ for the invariants $ I_1=||\nabla\phi||_F^2 $, $ I_2=||cof(\nabla\phi)||_F^2$, $ I_3=det(\nabla\phi) $, where $ \nabla\phi $, $ W $, and the spatial position are passed as arguments. The result matrix $ C= $ ElasticityTensor is such that $ <\frac{\partial^2W}{\partial\nabla\phi_K\partial\nabla\phi_L},\theta_K,\psi_L>=\theta_K\cdot(C\psi_L) $ for a variation $ \theta_K $ of $ \nabla\phi_K $ and $ \psi_L $ of $ \nabla\phi_L $.

Protected Types

typedef ConfiguratorType::RealType RealType

Protected Attributes

const
HyperelasticEnergyDensityType & 
_hyperelasticEnergyDensity
const
aol::DiscreteVectorFunctionDefault
< ConfiguratorType,
ConfiguratorType::Dim > 
_displacement
const int _k
const int _l

Detailed Description

template<typename ConfiguratorType, typename HyperelasticEnergyDensityType>
class qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >

This class computes via "apply(...)" or "applyAdd(...)" parts of a hyperelastic energy Hessian with respect to the displacement, i.e. $ (<\frac{\partial^2 \int_\Omega W(I_1,I_2,I_3,x) dx}{\partial \phi^2},\psi_i,\psi_j>)_{i,j} $, where the $ \psi_i $ and $ \psi_i $ are those vector-valued FE basis functions that are only nonzero in the kth and lth component, the hyperelastic invariants are given by $ I_1=||\nabla\phi||_F^2 $, $ I_2=||cof(\nabla\phi)||_F^2 $, $ I_3=det(\nabla\phi)$, and the function $ W(I_1,I_2,I_3,x) $ is passed to the constructor (e.g. HyerelasticEnergyDensityDefault). In its first argument, "apply" and "applyAdd" expect a Multivector which represents the displacements $ d_i $ in the different Cartesian directions; the deformation $ \phi $ is then automatically computed via $ \phi=d+$ identity.

Author:
Wirth

Definition at line 375 of file hyperelastic.h.


Member Typedef Documentation

template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
typedef ConfiguratorType::RealType qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::RealType [protected]

Constructor & Destructor Documentation

template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::HyperelasticSubHessian ( const typename ConfiguratorType::InitType &  Grid,
const HyperelasticEnergyDensityType &  HyperelasticEnergyDensity,
const aol::MultiVector< RealType > &  Displacement,
const int  K,
const int  L 
) [inline]

Definition at line 387 of file hyperelastic.h.

      : aol::FELinAsymMatrixWeightedStiffInterface<ConfiguratorType, HyperelasticSubHessian<ConfiguratorType,HyperelasticEnergyDensityType> > ( Grid ),
        _hyperelasticEnergyDensity( HyperelasticEnergyDensity ),
        _displacement( Grid, Displacement ),
        _k( L ), // l and k are exchanged here, since FELinAsymMatrixStiffOp assembles the matrix in the way that \int \nabla\phi^T A \nabla\psi dx becomes \Psi^T L \Phi
        _l( K ) {}


Member Function Documentation

template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
static void qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::elasticityTensor ( const typename ConfiguratorType::MatType &  DPhi,
const HyperelasticEnergyDensityType &  EnergyDensity,
const typename ConfiguratorType::ElementType &  El,
int  QuadPoint,
const int  K,
const int  L,
aol::Mat< ConfiguratorType::Dim, ConfiguratorType::Dim, RealType > &  ElasticityTensor 
) [inline, static]

Computes ElasticityTensor= $ \frac{\partial^2W}{\partial\nabla\phi_K\partial\nabla\phi_L}(I_1,I_2,I_3) $ for the invariants $ I_1=||\nabla\phi||_F^2 $, $ I_2=||cof(\nabla\phi)||_F^2$, $ I_3=det(\nabla\phi) $, where $ \nabla\phi $, $ W $, and the spatial position are passed as arguments. The result matrix $ C= $ ElasticityTensor is such that $ <\frac{\partial^2W}{\partial\nabla\phi_K\partial\nabla\phi_L},\theta_K,\psi_L>=\theta_K\cdot(C\psi_L) $ for a variation $ \theta_K $ of $ \nabla\phi_K $ and $ \psi_L $ of $ \nabla\phi_L $.

Definition at line 418 of file hyperelastic.h.

References aol::Mat< numRows, numCols, _DataType >::addMultiple(), aol::Mat< numRows, numCols, _DataType >::makeTensorProduct(), aol::Mat< numRows, numCols, _DataType >::setZero(), aol::Sqr(), and aol::Mat< numRows, numCols, _DataType >::transposeFrom().

Referenced by qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::getCoeffMatrix().

                                                                                                                        {
    // the cofactor matrix of the deformation gradient
    typename ConfiguratorType::MatType cof;
    // auxiliary vectors and matrices
    typename ConfiguratorType::VecType auxVec1;
    typename ConfiguratorType::MatType auxMat1, auxMat2;

    ElasticityTensor.setZero();

    // compute the cofactor matrix of the deformation gradient, $cof(\nabla\phi)$
    cof.makeCofactorMatrix( DPhi );
    // compute the hyperelastic invariants $||\nabla\phi||_F^2,||cof(\nabla\phi)||_F^2,det(\nabla\phi)$
    RealType
      I1 = DPhi.normSqr(),
      I2 = cof.normSqr(),
      I3 = DPhi.det();

    // compute the outer first and second derivative of the hyperelastic energy
    aol::Vec3<RealType> outerDeriv ( EnergyDensity.evaluateDerivativeAtQuadPoint ( I1, I2, I3, El, QuadPoint ) );
    aol::Matrix33<RealType> outerSecondDeriv ( EnergyDensity.evaluateSecondDerivativeAtQuadPoint ( I1, I2, I3, El, QuadPoint ) );

    // compute auxiliary matrices h, H
    typename ConfiguratorType::MatType H( cof.transposed() ), h, hCof;
    h.makeProduct( cof, H );
    RealType auxScal1 = h.tr();
    auxMat1 = h;
    for ( int i = 0; i < ConfiguratorType::Dim; i++ )
      auxMat1[i][i] -= auxScal1;
    H.makeProduct( auxMat1, cof );
    H /= - I3;
    hCof.makeProduct( h, cof );

    // compute the Hessian matrix for the kth component multiplied from the left and the lth from the right
    if ( K == L )
      for ( int i = 0; i < ConfiguratorType::Dim; i++ )
        ElasticityTensor[i][i] = 2 * outerDeriv[0]; // note: \epsilon:\varespilon = \sum_{k,l=1..d}\delta_{kl}\epsilon_{k:}I\varepsilon_{l:}^T
    auxVec1.setMultiple( DPhi[K], 4 * outerSecondDeriv[0][0] );
    auxVec1.addMultiple( H[K], 4 * outerSecondDeriv[0][1] );
    auxVec1.addMultiple( cof[K], 2 * outerSecondDeriv[0][2] );
    auxMat1.makeTensorProduct( auxVec1, DPhi[L] );
    ElasticityTensor += auxMat1; // note: (A:\epsilon)(B:\varepsilon) = \sum_{k,l=1..d}\epsilon_{k:}A_{k:}^T B_{l:}\varepsilon_{l:}^T

    auxVec1.setMultiple( cof[L], - outerDeriv[2] / I3 );
    auxMat1.makeTensorProduct( auxVec1, cof[K] );
    ElasticityTensor += auxMat1; // note: tr(A\epsilon^T B\varepsilon^T) = \sum_{k,l=1..d}\epsilon_{k:}A_{l:}^T B_{k:}\varepsilon_{l:}^T
    auxVec1.setMultiple( cof[K], outerDeriv[2] / I3 );
    auxVec1.addMultiple( DPhi[K], 2 * outerSecondDeriv[2][0] );
    auxVec1.addMultiple( H[K], 2 * outerSecondDeriv[2][1] );
    auxVec1.addMultiple( cof[K], outerSecondDeriv[2][2] );
    auxMat1.makeTensorProduct( auxVec1, cof[L] );
    ElasticityTensor += auxMat1;

    auxVec1.setMultiple( cof[K], - 2 * outerDeriv[1] / aol::Sqr( I3 ) );
    auxMat1.makeTensorProduct( auxVec1, hCof[L] );
    ElasticityTensor += auxMat1;
    auxVec1.setMultiple( cof[L], 2 * outerDeriv[1] / aol::Sqr( I3 ) );
    auxMat1.makeTensorProduct( auxVec1, hCof[K] );
    ElasticityTensor += auxMat1;
    auxMat2.transposeFrom( cof );
    auxMat1.makeProduct( auxMat2, cof );
    auxMat1 *= 2 * outerDeriv[1] / aol::Sqr( I3 ) * h[L][K];
    ElasticityTensor += auxMat1; // note: tr(A\epsilon B\varepsilon^T) = \sum_{k,l=1..d}\epsilon_{k:}A_{lk}B\varepsilon_{l:}^T
    auxVec1.setMultiple( cof[K],  2 * outerDeriv[1] / aol::Sqr( I3 ) * h.tr() );
    auxVec1.addMultiple( hCof[K], - 4 * outerDeriv[1] / aol::Sqr( I3 ) );
    auxMat1.makeTensorProduct( auxVec1, cof[L] );
    ElasticityTensor += auxMat1;
    auxVec1.setMultiple( H[L], - 2 * outerDeriv[1] / I3 );
    auxMat1.makeTensorProduct( auxVec1, cof[K] );
    ElasticityTensor += auxMat1;
    auxVec1.setMultiple( cof[K], 2 * outerDeriv[1] / I3 );
    auxVec1.addMultiple( DPhi[K], 4 * outerSecondDeriv[1][0] );
    auxVec1.addMultiple( H[K], 4 * outerSecondDeriv[1][1] );
    auxVec1.addMultiple( cof[K], 2 * outerSecondDeriv[1][2] );
    auxMat1.makeTensorProduct( auxVec1, H[L] );
    ElasticityTensor += auxMat1;
  }

template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
void qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::getCoeffMatrix ( const typename ConfiguratorType::ElementType &  El,
int  QuadPoint,
const typename ConfiguratorType::DomVecType &  ,
typename ConfiguratorType::MatType &  Matrix 
) const [inline]

Member Data Documentation

template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
const aol::DiscreteVectorFunctionDefault<ConfiguratorType,ConfiguratorType::Dim> qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::_displacement [protected]
template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
const HyperelasticEnergyDensityType& qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::_hyperelasticEnergyDensity [protected]
template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
const int qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::_k [protected]
template<typename ConfiguratorType , typename HyperelasticEnergyDensityType >
const int qc::HyperelasticSubHessian< ConfiguratorType, HyperelasticEnergyDensityType >::_l [protected]

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

Generated on Fri Sep 9 2011 21:09:48 for QuocMesh by doxygen 1.7.1