constitutive.hpp

Overview

Utility functions for constitutive modeling, including Voigt notation conversions, identity tensors, and stiffness/compliance tensor construction for various material symmetries (isotropic, cubic, orthotropic, transversely isotropic).

API Reference

arma::mat Ireal()

Returns the fourth order identity tensor Ireal written in Voigt notation.

Example:

mat Ir = Ireal();

Returns:

The following 6x6 mat (arma::mat)

\[\begin{split} I_{real} = \left( \begin{array}{cccccc} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0.5 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0.5 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0.5 \end{array} \right) \end{split}\]

arma::mat Ivol()

Returns the volumic of the identity tensor Ivol written in Voigt notation.

Example:

mat Iv = Ivol();

Returns:

The following 6x6 mat (arma::mat)

\[\begin{split} I_{vol} = \left( \begin{array}{cccccc} 1/3 & 1/3 & 1/3 & 0 & 0 & 0 \\ 1/3 & 1/3 & 1/3 & 0 & 0 & 0 \\ 1/3 & 1/3 & 1/3 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 \end{array} \right) \end{split}\]

arma::mat Idev()

Returns the deviatoric of the identity tensor Idev written in Voigt notation.

Example:

mat Id = Idev();

Returns:

The following 6x6 mat (arma::mat)

\[\begin{split} I_{dev} = I_{real} - I_{vol} = \left( \begin{array}{cccccc} 2/3 & -1/3 & -1/3 & 0 & 0 & 0 \\ -1/3 & 2/3 & -1/3 & 0 & 0 & 0 \\ -1/3 & -1/3 & 2/3 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0.5 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0.5 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0.5 \end{array} \right) \end{split}\]

arma::mat Ireal2()

Returns the fourth order identity tensor \( \widehat{I} \) written in Voigt notation.

For example, this tensor allows to obtain : \( L*\widehat{M}=I \) or \( \widehat{L}*M=I \), where a matrix \( \widehat{A} \) is set by \( \widehat{A}=\widehat{I}\,A\,\widehat{I} \)

mat Ir2 = Ireal2();

Returns:

The following 6x6 mat (arma::mat)

\[\begin{split} \widehat{I} = \left( \begin{array}{cccccc} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 2 & 0 & 0 \\ 0 & 0 & 0 & 0 & 2 & 0 \\ 0 & 0 & 0 & 0 & 0 & 2 \end{array} \right) \end{split}\]

arma::mat Idev2()

Returns the deviatoric part of the identity tensor, in the form of \( \widehat{I} \) considering the Voigt notation.

Example:

mat Id2 = Idev2();

Returns:

The following 6x6 mat (arma::mat)

\[\begin{split} I_{dev2} = \left( \begin{array}{cccccc} 2/3 & -1/3 & -1/3 & 0 & 0 & 0 \\ -1/3 & 2/3 & -1/3 & 0 & 0 & 0 \\ -1/3 & -1/3 & 2/3 & 0 & 0 & 0 \\ 0 & 0 & 0 & 2 & 0 & 0 \\ 0 & 0 & 0 & 0 & 2 & 0 \\ 0 & 0 & 0 & 0 & 0 & 2 \end{array} \right) \end{split}\]

arma::vec Ith()

Returns the expansion vector.

Example:

vec It = Ith();

Returns:

The following 6 vec (arma::vec)

\[\begin{split} I_{th} = \left( \begin{array}{c} 1 \\ 1 \\ 1 \\ 0 \\ 0 \\ 0 \end{array} \right) \end{split}\]

arma::vec Ir2()

Returns the operator from a stress to strain Voigt convention.

Example:

vec I2 = Ir2();

Returns:

The following 6 vec (arma::vec)

\[\begin{split} I_{r2} = \left( \begin{array}{ccc} 1 \\ 1 \\ 1 \\ 2 \\ 2 \\ 2 \end{array} \right) \end{split}\]

arma::vec Ir05()

Returns the operator from a strain to stress Voigt convention.

Example:

vec I05 = Ir05();

Returns:

The following 6 vec (arma::vec)

\[\begin{split} I_{r05} = \left( \begin{array}{ccc} 1 \\ 1 \\ 1 \\ 0.5 \\ 0.5 \\ 0.5 \end{array} \right) \end{split}\]

arma::mat L_iso(const double &C1, const double &C2, const std::string &conv = "Enu")

Provides the elastic stiffness tensor for an isotropic material.

The two first arguments are a couple of Lamé coefficients. The third argument specify which couple has been provided and the order of coefficients. Exhaustive list of possible third argument : ‘Enu’,’nuE,’Kmu’,’muK’, ‘KG’, ‘GK’, ‘lambdamu’, ‘mulambda’, ‘lambdaG’, ‘Glambda’.

Example:

double E = 210000;
double nu = 0.3;
mat Liso = L_iso(E, nu, "Enu");

Parameters:

C1, C2, conv

Returns:

The 6x6 stiffness matrix considering a Voigt notation (arma::mat)

arma::mat M_iso(const double &C1, const double &C2, const std::string &conv = "Enu")

Provides the elastic compliance tensor for an isotropic material.

The two first arguments are a couple of Lamé coefficients. The third argument specify which couple has been provided and the order of coefficients. Exhaustive list of possible third argument : ‘Enu’,’nuE,’Kmu’,’muK’, ‘KG’, ‘GK’, ‘lambdamu’, ‘mulambda’, ‘lambdaG’, ‘Glambda’.

double E = 210000;
double nu = 0.3;
mat Miso = M_iso(E, nu, "Enu");

Parameters:

C1, C2, conv

Returns:

The 6x6 compliance matrix considering a Voigt notation (arma::mat)

arma::mat L_cubic(const double &C1, const double &C2, const double &C3, const std::string &conv = "EnuG")

Provides the elastic stiffness tensor for a cubic material.

Arguments are the stiffness coefficients C11, C12 and C44 or E, nu and G. The fourth argument specify which couple has been provided and the order of coefficients: ‘EnuG’ or ’Cii’.

double C11 = alead(10000., 100000.);
double C12 = alead(10000., 100000.);
double C44 = alead(10000., 100000.);
mat Lcubic = L_cubic(C11, C12, C44, "Cii");

Parameters:

C1, C2, C3, conv

Returns:

The 6x6 stiffness matrix considering a Voigt notation (arma::mat)

arma::mat M_cubic(const double &C1, const double &C2, const double &C3, const std::string &conv = "EnuG")

Provides the elastic compliance tensor for a cubic material.

Arguments are the stiffness coefficients C11, C12 and C44 or E, nu and G. The fourth argument specify which couple has been provided and the order of coefficients: ‘EnuG’ or ’Cii’.

double C11 = alead(10000., 100000.);
double C12 = alead(10000., 100000.);
double C44 = alead(10000., 100000.);
mat Mcubic = M_cubic(C11,C12,C44);

Parameters:

C1, C2, C3, conv

Returns:

The 6x6 compliance matrix considering a Voigt notation (arma::mat)

arma::mat L_ortho(const double &C11, const double &C12, const double &C13, const double &C22, const double &C23, const double &C33, const double &C44, const double &C55, const double &C66, const std::string &conv = "EnuG")

Provides the elastic stiffness tensor for an orthotropic material.

Arguments could be all the stiffness coefficients ( \( C_{ij} \)) or the material parameter. In this last case for an orthotropic material the material parameters should be : \( E_x, E_y, E_z, \nu_{xy}, \nu_{yz}, \nu_{xz}, G_{xy}, G_{yz}, G_{xz} \). The last argument must be set to “Cii” if the inputs are the stiffness coefficients or to “EnuG” if the inputs are the material parameters.

double C11 = alead(10000., 100000.);
double C12 = alead(10000., 100000.);
double C13 = alead(10000., 100000.);
double C22 = alead(10000., 100000.);
double C23 = alead(10000., 100000.);
double C33 = alead(10000., 100000.);
double C44 = alead(10000., 100000.);
double C55 = alead(10000., 100000.);
double C66 = alead(10000., 100000.);
mat Lortho = L_ortho(C11, C12, C13, C22, C23, C33, C44, C55, C66,"Cii");

Parameters:
  • C11 – first stiffness coefficient (or Ex)

  • C12 – second stiffness coefficient (or Ey)

  • C13 – third stiffness coefficient (or Ez)

  • C22 – fourth stiffness coefficient (or nu_xy)

  • C23 – fifth stiffness coefficient (or nu_yz)

  • C33 – sixth stiffness coefficient (or nu_xz)

  • C44 – seventh stiffness coefficient (or G_xy)

  • C55 – eighth stiffness coefficient (or G_yz)

  • C66 – ninth stiffness coefficient (or G_xz)

  • conv – convention string: “Cii” for stiffness coefficients or “EnuG” for material parameters

Returns:

The 6x6 stiffness matrix considering a Voigt notation (arma::mat)

arma::mat M_ortho(const double &C11, const double &C12, const double &C13, const double &C22, const double &C23, const double &C33, const double &C44, const double &C55, const double &C66, const std::string &conv = "EnuG")

Provides the elastic compliance tensor for an orthotropic material.

Arguments could be all the stiffness coefficients ( \( C_{ij} \)) or the material parameter. In this last case for an orthotropic material the material parameters should be : \( E_x, E_y, E_z, \nu_{xy}, \nu_{yz}, \nu_{xz}, G_{xy}, G_{yz}, G_{xz} \). The last argument must be set to “Cii” if the inputs are the stiffness coefficients or to “EnuG” if the inputs are the material parameters.

double C11 = alead(10000., 100000.);
double C12 = alead(10000., 100000.);
double C13 = alead(10000., 100000.);
double C22 = alead(10000., 100000.);
double C23 = alead(10000., 100000.);
double C33 = alead(10000., 100000.);
double C44 = alead(10000., 100000.);
double C55 = alead(10000., 100000.);
double C66 = alead(10000., 100000.);
mat Mortho = M_ortho(C11, C12, C13, C22, C23, C33, C44, C55, C66,"Cii");

Parameters:
  • C11 – first stiffness coefficient (or Ex)

  • C12 – second stiffness coefficient (or Ey)

  • C13 – third stiffness coefficient (or Ez)

  • C22 – fourth stiffness coefficient (or nu_xy)

  • C23 – fifth stiffness coefficient (or nu_yz)

  • C33 – sixth stiffness coefficient (or nu_xz)

  • C44 – seventh stiffness coefficient (or G_xy)

  • C55 – eighth stiffness coefficient (or G_yz)

  • C66 – ninth stiffness coefficient (or G_xz)

  • conv – convention string: “Cii” for stiffness coefficients or “EnuG” for material parameters

Returns:

The 6x6 compliance matrix considering a Voigt notation (arma::mat)

arma::mat L_isotrans(const double &EL, const double &ET, const double &nuTL, const double &nuTT, const double &GLT, const int &axis)

Provides the elastic stiffness tensor for an isotropic transverse material.

Arguments are longitudinal Young modulus \( E_L \) (EL), transverse Young modulus \( E_T \) (ET), Poisson’s ratio for loading along the longitudinal axis \( \nu_{TL} \) (nuTL), Poisson’s ratio for loading along the transverse axis \( \nu_{TT} \) (nuTT), shear modulus \( G_{LT} \) (GLT) and the axis of symmetry.

double EL = alead(10000., 100000.);
double ET = alead(10000., 100000.);
double nuTL = alead(0., 0.5);
double nuTT = alead(0.5, 0.5);
double GLT = alead(10000., 100000.);
double axis = 1;
mat Lisotrans = L_isotrans(EL, ET, nuTL, nuTT, GLT, axis);

Parameters:

EL, ET, nuTL, nuTT, GLT, axis

Returns:

The 6x6 stiffness matrix considering a Voigt notation (arma::mat)

arma::mat M_isotrans(const double &EL, const double &ET, const double &nuTL, const double &nuTT, const double &GLT, const int &axis)

Provides the elastic compliance tensor for an isotropic transverse material.

Arguments are longitudinal Young modulus \( E_L \) (EL), transverse Young modulus \( E_T \) (ET), Poisson’s ratio for loading along the longitudinal axis \( \nu_{TL} \) (nuTL), Poisson’s ratio for loading along the transverse axis \( \nu_{TT} \) (nuTT), shear modulus \( G_{LT} \) (GLT) and the axis of symmetry.

double EL = alead(10000., 100000.);
double ET = alead(10000., 100000.);
double nuTL = alead(0., 0.5);
double nuTT = alead(0.5, 0.5);
double GLT = alead(10000., 100000.);
double axis = 1;
mat Misotrans = M_isotrans(EL, ET, nuTL, nuTT, GLT, axis);

Parameters:

EL, ET, nuTL, nuTT, GLT, axis

Returns:

The 6x6 compliance matrix considering a Voigt notation (arma::mat)

arma::mat H_iso(const double &etaB, const double &etaS)

Provides the viscous tensor H for an isotropic material.

Arguments are the Bulk viscosity \( \eta_B \) (etaB) and shear viscosity \( \eta_S \) (etaS).

double etaB = alead(0., 0.1);
double etaS = alead(0., 0.1);
mat Hiso = H_iso(etaB, etaS);

Parameters:
  • etaB – Bulk viscosity \( \eta_B \)

  • etaS – Shear viscosity \( \eta_S \)

Returns:

The 6x6 viscous matrix considering a Voigt notation (arma::mat)

\[\begin{split} H_{iso} = \left( \begin{array}{cccccc} \eta_B + \frac{4}{3} \eta_S & \eta_B - \frac{2}{3} \eta_S & \eta_B - \frac{2}{3} \eta_S & 0 & 0 & 0 \\ \eta_B - \frac{2}{3} \eta_S & \eta_B + \frac{4}{3} \eta_S & \eta_B - \frac{2}{3} \eta_S & 0 & 0 & 0 \\ \eta_B - \frac{2}{3} \eta_S & \eta_B - \frac{2}{3} \eta_S & \eta_B + \frac{4}{3} \eta_S & 0 & 0 & 0 \\ 0 & 0 & 0 & \eta_S & 0 & 0 \\ 0 & 0 & 0 & 0 & \eta_S & 0 \\ 0 & 0 & 0 & 0 & 0 & \eta_S \end{array} \right) \end{split}\]

arma::vec el_pred(const arma::vec &sigma_start, const arma::mat &L, const arma::vec &Eel, const int &ndi = 3)

Provides the stress tensor (elastic prediction), from the previous stress increment, providing the elastic stiffness tensor and the trial elastic strain increment:

The input parameters are : sigma_start: The previous stress, L : Stiffness matrix; Eel : elastic strain vector, ndi (optional, default = 3): number of dimensions

vec sigma_start = zeros(6);
sigma_start.randu(6);
mat L = L_iso(70000, 0.3,"Enu");
vec Eel;
Eel.randu(6);
int ndi = 3;
vec sigma =  el_pred(sigma_start,L, Eel, ndi);

Parameters:

sigma_start, L, Eel, ndi

Returns:

The predicted 6 vec stress vector (arma::vec)

arma::vec el_pred(const arma::mat &L, const arma::vec &Eel, const int &ndi = 3)

Provides the stress tensor (elastic prediction), from the elastic stiffness tensor and the trial elastic strain:

The input parameters are : L : Stiffness matrix; Eel ; elastic strain vector, ndi (optional, default = 3): number of dimensions

mat L = L_iso(70000, 0.3,"Enu");
vec Eel;
Eel.randu(6);
int ndi = 3;
vec sigma =  el_pred(L, Eel, ndi);

Parameters:

L, Eel, ndi

Returns:

The predicted 6 vec stress vector (arma::vec)

arma::mat Isotropize(const arma::mat &Lt)

Provides the isotropized tangent modulus from the spectral decomposition (Bornet et al., 2001)

The returned isotropic tensor is called consistent since for any given strain it return the same stress as the anisotropic version.

double EL = (double)rand();
double ET = (double)rand();
double nuTL = (double)rand();
double nuTT = (double)rand();
double GLT = (double)rand();
double axis = 1;
mat L_isotrans = L_isotrans(EL, ET, nuTL, nuTT, GLT, axis);
mat L_iso = Isotropize(Lisotrans);

Parameters:

Lt

Returns:

The 6x6 isotropic matrix considering a Voigt notation (arma::mat)