# Continuum_Mechanics

For details about the tensor notation, see “tensors in SMART+.pdf” in the Continuum Mechanics library.

# The Continuum Mechanics Library: contimech.hpp

### tr(vec)

Provides the trace of a second order tensor written as a vector v in the SMART+ formalism.
Return a double.
Example:

vec v = randu(6);
double trace = tr(v);


### dev(vec)

Provides the deviatoric part of a second order tensor written as a vector v in the SMART+ formalism.
Return a vec.
Example:

vec v = randu(6);
vec deviatoric = dev(v);


### Mises_stress(vec)

Provides the Von Mises stress $$\sigma^{Mises}$$ of a second order stress tensor written as a vector v in the SMART+ formalism.
Return a double.

vec v = randu(6);
double Mises_sig = Mises_stress(v);


### Mises_strain(vec)

Provides the Von Mises strain $$\varepsilon^{Mises}$$ of a second order stress tensor written as a vector v in the SMART+ formalism.
Return a double.

vec v = randu(6);
double Mises_eps = Mises_strain(v);


### eta_stress(vec)

Provides the stress flow $$\eta_{stress}=\frac{3/2\sigma_{dev}}{\sigma_{Mises}}$$ from a second order stress tensor written as a vector v in the SMART+ formalism (i.e. the shear terms are multiplied by 2, providing shear angles).
Return a vec.

vec v = randu(6);
vec sigma_f = eta_stress(v);


### eta_strain(vec)

Provides the strain flow $$\eta_{strain}=\frac{2/3\varepsilon_{dev}}{\varepsilon_{Mises}}$$ from a second order strain tensor written as a vector v in the SMART+ formalism (i.e. the shear terms are multiplied by 2, providing shear angles).
Return a vec.

vec v = randu(6);
vec eps_f = eta_strain(v);


### v2t_strain(vec)

Converts a second order strain tensor written as a vector v in the SMART+ formalism into a second order strain tensor written as a matrix m.
Return a mat.

vec v = randu(6);
mat m = v2t_strain(v);


### t2v_strain(vec)

Converts a second order strain tensor written as a matrix m in the SMART+ formalism into a second order strain tensor written as a vector v.
Return a vec.

mat m = randu(6,6);
vec v = t2v_strain(m);


### v2t_stress(vec)

Converts a second order stress tensor written as a vector v in the SMART+ formalism into a second order stress tensor written as a matrix m.
Return a mat.

vec v = randu(6);
mat m = v2t_stress(v);


### t2v_stress(vec)

Converts a second order stress tensor written as a matrix m in the SMART+ formalism into a second order stress tensor written as a vector v.
Return a vec.

mat m = randu(6,6);
vec v = t2v_stress(m);


### J2_strain(vec)

Provides the second invariant of a second order strain tensor written as a vector v in the SMART+ formalism.
Return a vec.

vec v = randu(6);
double J2 = J2_strain(v);


### J2_stress(vec)

Provides the second invariant of a second order stress tensor written as a vector v in the SMART+ formalism.
Return a double.

vec v = randu(6);
double J2 = J2_stress(v);


### J3_strain(vec)

Provides the third invariant of a second order strain tensor written as a vector v in the SMART+ formalism.
Return a vec.

vec v = randu(6);
double J3 = J3_strain(v);


### J3_stress(vec)

Provides the third invariant of a second order stress tensor written as a vector v in the SMART+ formalism.
Return a double.

vec v = randu(6);
double J3 = J3_stress(v);


### normal_ellipsoid(double, double, double, double, double)

Provides the normalized vector to an ellipsoid with semi-principal axes of length $$a_{1}$$, $$a_{2}$$, $$a_{3}$$. The direction of the normalized vector is set by angles u and v. These 2 angles correspond to the rotations in the plan defined by the center of the ellipsoid, a1 and a2 directions for u, a1 and a3 ones for v. u = 0 corresponds to a1 direction and v = 0 correspond to a3 one. So the normal vector is set at the parametrized position :
\begin{align} x & = a_{1} cos(u) sin(v) \\ y & = a_{2} sin(u) sin(v) \\ z & = a_{3} cos(v) \end{align}

Return a vec.

const double Pi = 3.14159265358979323846

double u = (double)rand()/(double)(RAND_MAX) % 2*Pi - 2*Pi;
double v = (double)rand()/(double)(RAND_MAX) % Pi - Pi;
double a1 = (double)rand();
double a2 = (double)rand();
double a3 = (double)rand();
vec v = normal_ellipsoid(u, v, a1, a2, a3);


### sigma_int(vec , double, double, double, double, double)

Provides the normal and tangent components of a stress vector $$\sigma_{in}$$ in accordance with the normal direction n to an ellipsoid with axes $$a_{1}$$, $$a_{2}$$, $$a_{3}$$. The normal vector is set at the parametrized position :
\begin{align} x & = a_{1} cos(u) sin(v) \\ y & = a_{2} sin(u) sin(v) \\ z & = a_{3} cos(v) \end{align}
Return a vec.

vec sigma_in = randu(6);
double u = (double)rand()/(double)(RAND_MAX) % Pi - Pi/2;
double v = (double)rand()/(double)(RAND_MAX) % 2*Pi - Pi;
double a1 = (double)rand();
double a2 = (double)rand();
double a3 = (double)rand();
vec sigma_i = sigma_int(sigma_in, a1, a2, a3, u, v));


### p_ikjl(vec)

Provides the Hill interfacial operator according to a normal a (see papers of Siredey and Entemeyer Ph.D. dissertation).
Return a mat.

vec v = randu(6);
mat H = p_ikjl(v);


# The Constitutive Library: constitutive.hpp

### Ireal()

Provides the fourth order identity tensor written in Voigt notation $$I_{real}$$, where :
$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)$
Return a mat.
Example:

mat Ir = Ireal();


### Ireal2()

Provides the fourth order identity tensor $$\widehat{I}$$ written in the form. So :
$\widehat{I} = \left( \begin{array}{ccc} 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)$

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}$$
Return a mat.
Example:

mat Ir2 = Ireal2();


### Ivol()

Provides the volumic of the identity tensor $$I_{real}$$ written in the SMART+ formalism. So :
$I_{vol} = \left( \begin{array}{ccc} 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)$
Return a mat.
Example:

mat Iv = Ivol();


### Idev()

Provides the deviatoric of the identity tensor $$I_{real}$$ written in the SMART+ formalism. So :
$I_{dev} = \left( \begin{array}{ccc} 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)$
Return a mat.
Example:

mat Id = Idev();


### Idev2()

Provides the deviatoric of the identity tensor $$\widehat{I}$$ written in the SMART+ formalism. So :
$I_{dev2} = \left( \begin{array}{ccc} 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)$
Return a mat.
Example:

mat Id2 = Idev2();


### Ith()

Provides the vector $$I_{th} = \left( \begin{array}{ccc} 1 \\ 1 \\ 1 \\ 0 \\ 0 \\ 0 \end{array} \right)$$.
Return a vec.
Example:

vec It = Ith();


### Ir2()

Provides the vector $$I_{r2} = \left( \begin{array}{ccc} 1 \\ 1 \\ 1 \\ 2 \\ 2 \\ 2 \end{array} \right)$$.
Return a vec.
Example:

vec I2 = Ir2();


### Ir05()

Provides the vector $$I_{r05} = \left( \begin{array}{ccc} 1 \\ 1 \\ 1 \\ 0.5 \\ 0.5 \\ 0.5 \end{array} \right)$$.
Return a vec.
Example:

vec I05 = Ir05();


### L_iso(double, double, string)

Provides the elastic stiffness tensor for an isotropic material.
The two first arguments are a couple of elastic properties. The third argument specifies which couple has been provided and the nature and order of coefficients.
Exhaustive list of possible third argument :
‘Enu’,’nuE,’Kmu’,’muK’, ‘KG’, ‘GK’, ‘lambdamu’, ‘mulambda’, ‘lambdaG’, ‘Glambda’.
Return a mat.
Example :

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


### M_iso(double, double, string)

Provides the elastic compliance tensor for an isotropic material.
The two first arguments are a couple of elastic properties. The third argument specify which couple has been provided and the nature and order of coefficients.
Exhaustive list of possible third argument :
‘Enu’,’nuE,’Kmu’,’muK’, ‘KG’, ‘GK’, ‘lambdamu’, ‘mulambda’, ‘lambdaG’, ‘Glambda’.
Return a mat.
Example :

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


### L_cubic(double, double, double)

Provides the elastic stiffness tensor for a cubic material.
Arguments are the stiffness coefficients C11, C12 and C44.

Return a mat.
Example :

double C11 = (double)rand();
double C12 = (double)rand();
doubel C44 = (double)rand();
mat Liso = L_cubic(C11,C12,C44);


### M_cubic(double, double, double)

Provides the elastic compliance tensor for a cubic material.
Arguments are the stiffness coefficients C11, C12 and C44.

Return a mat.
Example :

double C11 = (double)rand();
double C12 = (double)rand();
double C44 = (double)rand();
mat Miso = M_cubic(C11,C12,C44);


### L_ortho(double, double, double, double, double, double, double, double, double, string)

Provides the elastic stiffness tensor for an orthotropic material.
Arguments could be all the stiffness coefficients or the material parameter. For an orthotropic material the material parameters should be : Ex,Ey,Ez,nuxy,nuyz,nxz,Gxy,Gyz,Gxz.

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.

Return a mat.
Example :

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


### M_ortho(double, double, double, double, double, double, double, double, double,string)

Provides the elastic compliance tensor for an orthotropic material.
Arguments could be all the stiffness coefficients or the material parameter. For an orthotropic material the material parameters should be : Ex,Ey,Ez,nuxy,nuyz,nxz,Gxy,Gyz,Gxz.

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.

return a mat
Example :

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


### L_isotrans(double, double, double, double, double, int)

Provides the elastic stiffness tensor for an isotropic transverse material.
Arguments are longitudinal Young modulus EL, transverse young modulus, Poisson’s ratio for loading along the longitudinal axis nuTL, Poisson’s ratio for loading along the transverse axis nuTT, shear modulus GLT and the axis of symmetry.

Return a mat.
Example :

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


### M_isotrans(double, double, double, double, double, int)

Provides the elastic compliance tensor for an isotropic transverse material.
Arguments are longitudinal Young modulus EL, transverse young modulus, Poisson’s ratio for loading along the longitudinal axis nuTL, Poisson’s ratio for loading along the transverse axis nuTT, shear modulus GLT and the axis of symmetry.

Return a mat.
Example :

double EL = (double)rand();
double ET = (double)rand();
double nuTL = (double)rand();
double nuTT = (double)rand();
double GLT = (double)rand();
double axis = 1;
mat Misotrans = M_isotrans(EL, ET, nuTL, nuTT, GLT, axis);


### Christoffel(mat, vec)

Provides the Christoffel tensor in the anisotropic case.
Parameters are the elastic stiffness matrix L and the wavenumber vector k.
Return a mat.
Example :

vec L = randu(6,6);
vec k = randu(3);
mat Christ = Christoffel(L, k);


# The Damage Library: damage.hpp

### weibull(vec, double, double, double, double, string)

Provides the damage evolution $$\delta D$$ considering a Weibull damage law.
It is given by : $$\delta D = (1-D_{old})*\Big(1-exp\big(-1(\frac{crit}{\beta})^{\alpha}\big)\Big)$$
Parameters of this function are: the stress vector $$\sigma$$, the old damage $$D_{old}$$, the shape parameter $$\alpha$$, the scale parameter $$\beta$$, the time increment $$\Delta T$$ and the criterion (which is a string).

The criterion possibilities are :
“vonmises” : $$crit = \sigma_{Mises}$$
“hydro” : $$crit = tr(\sigma)$$
“J3” : $$crit = J3(\sigma)$$.
Default value of the criterion is “vonmises”.
Return a scalar.
Example:

double varD = damage_weibull(stress, damage, alpha, beta, DTime, criterion);


### damage_kachanov(vec, vec, double, double, double, string)

Provides the damage evolution $$\delta D$$ considering a Kachanov’s creep damage law.
It is given by : $$\delta D = \Big(\frac{crit}{A_0(1-D_{old})}\Big)^r$$
Parameters of this function are: the stress vector $$\sigma$$, the strain vector $$\varepsilon$$, the old damage $$D_{old}$$, the material properties characteristic of creep damage $$(A_0,r)$$ and the criterion (which is a string).

The criterion possibilities are :
“vonmises” : $$crit = (\sigma*(1+\varepsilon))_{Mises}$$
“hydro” : $$crit = tr(\sigma*(1+\varepsilon))$$
“J3” : $$crit = J3(\sigma*(1+\varepsilon))$$.
Here, the criterion has no default value.
Return a scalar.
Example:

double varD = damage_kachanov(stress, strain, damage, A0, r, criterion);


### damage_miner(double, double, double, double, double, double, double)

Provides the constant damage evolution $$\Delta D$$ considering a Woehler- Miner’s damage law.
It is given by : $$\Delta D = \big(\frac{S_{Max}-S_{Mean}+Sl_0*(1-b*S_{Mean})}{S_{ult}-S_{Max}}\big)*\big(\frac{S_{Max}-S_{Mean}}{B_0*(1-b*S_{Mean})}\big)^\beta$$
Parameters of this function are: the max stress value $$\sigma_{Max}$$, the mean stress value $$\sigma_{Mean}$$, the “ult” stress value $$\sigma_{ult}$$, the $$b$$, the $$B_0$$, the $$\beta$$ and the $$Sl_0$$.

Default value of $$Sl_0$$ is 0.0.
Return a scalar.
Example:

double varD = damage_minerl(S_max, S_mean, S_ult, b, B0, beta, Sl_0);


### damage_manson(double, double, double)

Provides the constant damage evolution $$\Delta D$$ considering a Coffin-Manson’s damage law.
It is given by : $$\Delta D = \big(\frac{\sigma_{Amp}}{C_{2}}\big)^{\gamma_2}$$
Parameters of this function are: the “amp” stress value $$\sigma_{Amp}$$, the $$C_{2}$$ and the $$\gamma_2$$.

Return a scalar.
Example:

double varD = damage_manson(S_amp, C2, gamma2);