Functions
Marmot::ContinuumMechanics::VoigtNotation::Invariants Namespace Reference

Functions

Eigen::Vector3d principalStrains (const Marmot::Vector6d &strain)
 
Eigen::Vector3d principalStresses (const Marmot::Vector6d &stress)
 
Eigen::Vector3d sortedPrincipalStrains (const Marmot::Vector6d &strain)
 
Eigen::Matrix3d principalStressesDirections (const Marmot::Vector6d &stress)
 
double vonMisesEquivalentStress (const Marmot::Vector6d &stress)
 
double vonMisesEquivalentStrain (const Marmot::Vector6d &strain)
 
template<typename T >
normStrain (const Eigen::Matrix< T, 6, 1 > &strain)
 
double normStress (const Marmot::Vector6d &stress)
 
double StrainVolumetricNegative (const Marmot::Vector6d &strain)
 
template<typename T >
I1 (const Eigen::Matrix< T, 6, 1 > &stress)
 
double I1Strain (const Marmot::Vector6d &strain)
 
template<typename T >
I2 (const Eigen::Matrix< T, 6, 1 > &stress)
 
double I2Strain (const Marmot::Vector6d &strain)
 
template<typename T >
I3 (const Eigen::Matrix< T, 6, 1 > &stress)
 
double I3Strain (const Marmot::Vector6d &strain)
 
template<typename T >
J2 (const Eigen::Matrix< T, 6, 1 > &stress)
 
double J2Strain (const Marmot::Vector6d &strain)
 
template<typename T >
J3 (const Eigen::Matrix< T, 6, 1 > &stress)
 
double J3Strain (const Marmot::Vector6d &strain)
 
std::pair< Eigen::Vector3d, Eigen::Matrix< double, 3, 6 > > principalValuesAndDerivatives (const Eigen::Matrix< double, 6, 1 > &S)
 

Function Documentation

◆ principalStrains()

Vector3d Marmot::ContinuumMechanics::VoigtNotation::Invariants::principalStrains ( const Marmot::Vector6d strain)

Computes the principal strains by solving the eigenvalue problem.

\[ \displaystyle |\varepsilon_{ij} - \lambda\, \delta_{ij}| = 0 \hspace{.5cm} \Rightarrow \hspace{.5cm} \lambda^{(1)},\lambda^{(2)},\lambda^{(3)}\hspace{0.3cm} \widehat{=}\hspace{0.3cm} \varepsilon_1,\, \varepsilon_2,\, \varepsilon_3 \]

The resulting principal strains are NOT sorted.

◆ principalStresses()

Vector3d Marmot::ContinuumMechanics::VoigtNotation::Invariants::principalStresses ( const Marmot::Vector6d stress)

Computes the principal stresses by solving the eigenvalue problem.

\[ \displaystyle |\sigma_{ij} - \lambda\, \delta_{ij}| = 0 \hspace{.5cm} \Rightarrow \hspace{.5cm} \lambda^{(1)},\lambda^{(2)},\lambda^{(3)}\hspace{0.3cm} \widehat{=}\hspace{0.3cm}\sigma_1,\, \sigma_2,\, \sigma_3 \]

The resulting principal stresses are NOT sorted.

◆ sortedPrincipalStrains()

Vector3d Marmot::ContinuumMechanics::VoigtNotation::Invariants::sortedPrincipalStrains ( const Marmot::Vector6d strain)

Calculates the principal strains from its corresponding haigh westergaard coordinates.

\[ \displaystyle \begin{bmatrix} \varepsilon_{1}\\ \varepsilon_{2}\\ \varepsilon_{3} \end{bmatrix} = \frac{1}{\sqrt{3}} \begin{bmatrix} \xi\\ \xi\\ \xi \end{bmatrix} + \sqrt{\frac{2}{3}}\,\rho\,\begin{bmatrix} \cos(\theta)\\ -\sin\left(\frac{\pi}{6} - \theta\right)\\ -\sin\left(\frac{\pi}{6} + \theta\right) \end{bmatrix} \]

The computation of \( \xi\) ,\ \( \rho\) and \(\theta\) can be found in haighWestergaardFromStrain()

◆ principalStressesDirections()

Eigen::Matrix3d Marmot::ContinuumMechanics::VoigtNotation::Invariants::principalStressesDirections ( const Marmot::Vector6d stress)

Computes the principal stress directions \(\boldsymbol{x}^{(k)}\) of the eigenvalues \( \sigma_k \) by solving

\[ \displaystyle \left(\boldsymbol{\sigma} - \sigma_k \cdot \boldsymbol{I}\right) \cdot \boldsymbol{x}^{(k)} = 0 \]

The resulting principal stress directions are NOT sorted.

◆ vonMisesEquivalentStress()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::vonMisesEquivalentStress ( const Marmot::Vector6d stress)

Computes the equivalent von Mises stress.

\[ \displaystyle \sigma^{(eq)} = \sqrt{3 \cdot J_2} \]

Wherein \( J_2 \) denotes the second invariant of the deviator stress tensor (see J2()).

◆ vonMisesEquivalentStrain()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::vonMisesEquivalentStrain ( const Marmot::Vector6d strain)

Computes the equivalent von Mises strain from deviatoric part of the strain tensor \( e_{ij} \)

\[ \displaystyle \varepsilon^{(eq)} = \sqrt{ \frac{2}{3} \cdot e_{ij}\,e_{ij}} \]

◆ normStrain()

template<typename T >
T Marmot::ContinuumMechanics::VoigtNotation::Invariants::normStrain ( const Eigen::Matrix< T, 6, 1 > &  strain)

Computes the euclidian norm of the strain tensor \( ||\boldsymbol{\varepsilon}|| \)

◆ normStress()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::normStress ( const Marmot::Vector6d stress)

Computes the euclidian norm of the stress tensor \( ||\boldsymbol{\sigma}|| \)

◆ StrainVolumetricNegative()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::StrainVolumetricNegative ( const Marmot::Vector6d strain)

Computes the volumetric plastic strains in compression

\[ \displaystyle \varepsilon^{vol}_{\ominus} = \sum^{3}_{i = 1} \left\langle -\varepsilon_i \right\rangle \]

using the Macaulay brackets \( \left\langle \bullet \right\rangle \) and the principal values of the strain tensor \(\varepsilon_i \)

◆ I1()

template<typename T >
T Marmot::ContinuumMechanics::VoigtNotation::Invariants::I1 ( const Eigen::Matrix< T, 6, 1 > &  stress)

Computes the first invariant \( I_1 \) of the stress tensor \( \boldsymbol{\sigma} \).

\[ \displaystyle I_1 = tr(\boldsymbol{\sigma}) \]

◆ I1Strain()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::I1Strain ( const Marmot::Vector6d strain)

Computes the first invariant \( I^{(\varepsilon)}_1 \) of the strain tensor \( \boldsymbol{\varepsilon} \).

\[ \displaystyle I^{(\varepsilon)}_1 = tr(\boldsymbol{\varepsilon}) \]

◆ I2()

template<typename T >
T Marmot::ContinuumMechanics::VoigtNotation::Invariants::I2 ( const Eigen::Matrix< T, 6, 1 > &  stress)

Computes the second invariant \( I_2 \) of the stress tensor \( \boldsymbol{\sigma} \).

\[ \displaystyle I_2 = \frac{1}{2} \left(tr(\boldsymbol{\sigma})^2 - tr(\boldsymbol{\sigma}^2)\right) \]

◆ I2Strain()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::I2Strain ( const Marmot::Vector6d strain)

Computes the second invariant \( I^{(\varepsilon)}_2 \) from a voigt notated strain vector \( \boldsymbol{\varepsilon} \).

\[ \displaystyle I^{(\varepsilon)}_2 = \varepsilon_{11}\,\varepsilon_{22} + \varepsilon_{22}\,\varepsilon_{33} + \varepsilon_{11}\,\varepsilon_{33} - \frac{1}{4}(\gamma^2_{12} - \gamma^2_{13} - \gamma^2_{23}) \]

◆ I3()

template<typename T >
T Marmot::ContinuumMechanics::VoigtNotation::Invariants::I3 ( const Eigen::Matrix< T, 6, 1 > &  stress)

Computes the third invariant \( I_3 \) of the stress tensor \( \boldsymbol{\sigma} \).

\[ \displaystyle I_3 = det(\boldsymbol{\sigma}) \]

◆ I3Strain()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::I3Strain ( const Marmot::Vector6d strain)

Computes the third invariant \( I^{(\varepsilon)}_3 \) from a voigt notated strain vector \( \boldsymbol{\varepsilon} \) by calling voigtToStrain() and calculating the determinant.

◆ J2()

template<typename T >
T Marmot::ContinuumMechanics::VoigtNotation::Invariants::J2 ( const Eigen::Matrix< T, 6, 1 > &  stress)

Computes the second invariant \( J_2 \) of the deviatoric part of the stress tensor \( \boldsymbol{s} \).

\[ \displaystyle J_2 = \frac{1}{3} I^2_1 - I_2 \]

◆ J2Strain()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::J2Strain ( const Marmot::Vector6d strain)

Computes the second invariant \( J^{(\varepsilon)}_2 \) of the deviatoric part of the strain tensor \( \boldsymbol{\varepsilon} \).

\[ \displaystyle J^{(\varepsilon)}_2 = \frac{1}{3} I^{(\varepsilon) 2}_1 - I^{(\varepsilon)}_2 \]

◆ J3()

template<typename T >
T Marmot::ContinuumMechanics::VoigtNotation::Invariants::J3 ( const Eigen::Matrix< T, 6, 1 > &  stress)

Computes the third invariant \( J_3 \) of the deviatoric part of the stress tensor \( \boldsymbol{s} \).

\[ \displaystyle J_3 = \frac{2}{27} I^3_1 - \frac{1}{3} I_1\cdot I_2\cdot I_3 \]

◆ J3Strain()

double Marmot::ContinuumMechanics::VoigtNotation::Invariants::J3Strain ( const Marmot::Vector6d strain)

Computes the third invariant \( J^{(\varepsilon)}_3 \) of a voigt notated deviatoric strain vector \( \boldsymbol{e} \) by calling voigtToStrain() and calculating the determinant.

◆ principalValuesAndDerivatives()

std::pair< Eigen::Vector3d, Eigen::Matrix< double, 3, 6 > > Marmot::ContinuumMechanics::VoigtNotation::Invariants::principalValuesAndDerivatives ( const Eigen::Matrix< double, 6, 1 > &  S)