The Continuum Mechanics Library

double tr(const vec &v)

Provides the trace of a second order tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
double trace = tr(v);
vec dev(const vec &v)

Provides the deviatoric part of a second order tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
vec deviatoric = dev(v);
double Mises_stress(const vec &v)

Provides the Von Mises stress \(\sigma^{Mises}\) of a second order stress tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
double Mises_sig = Mises_stress(v);
vec eta_stress(const vec &v)

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 ‘simcoon’ formalism (i.e. the shear terms are multiplied by 2, providing shear angles).

vec v = randu(6);
vec sigma_f = eta_stress(v);
double Mises_strain(const vec &v)

Provides the Von Mises strain \(\varepsilon^{Mises}\) of a second order stress tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
double Mises_eps = Mises_strain(v);
vec eta_strain(const vec &v)

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 ‘simcoon’ formalism (i.e. the shear terms are multiplied by 2, providing shear angles).

vec v = randu(6);
vec eps_f = eta_strain(v);
double J2_stress(const vec &v)

Provides the second invariant of a second order stress tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
double J2 = J2_stress(v);
double J2_strain(const vec &v)

Provides the second invariant of a second order strain tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
double J2 = J2_strain(v);
double J3_stress(const vec &v)

Provides the third invariant of a second order stress tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
double J3 = J3_stress(v);
double J3_strain(const vec &v)

Provides the third invariant of a second order strain tensor written as a vector v in the ‘simcoon’ formalism.

vec v = randu(6);
double J3 = J3_strain(v);
double Macaulay_p(const double &d)

This function returns the value if it’s positive, zero if it’s negative (Macaulay brackets <>+)

double Macaulay_n(const double &d)

This function returns the value if it’s negative, zero if it’s positive (Macaulay brackets <>-)

double sign(const double &d)

This function returns the value if it’s negative, zero if it’s positive (Macaulay brackets <>-)

vec normal_ellipsoid(const double &u, const double &v, const double &a1, const double &a2, const double &a3)

Provides the normalized vector to an ellipsoid with semi-principal axes of length a1, a2, a3. 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{split}\begin{align} x & = a_{1} cos(u) sin(v) \\ y & = a_{2} sin(u) sin(v) \\ z & = a_{3} cos(v) \end{align}\end{split}\]
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);
vec sigma_int(const vec &sigma_in, const double &a1, const double &a2, const double &a3, const double &u, const double &v)

Provides the normal and tangent components of a stress vector σin in accordance with the normal direction n to an ellipsoid with axes a1, a2, a3. The normal vector is set at the parametrized position :

\[\begin{split}\begin{align} x & = a_{1} cos(u) sin(v) \\ y & = a_{2} sin(u) sin(v) \\ z & = a_{3} cos(v) \end{align}\end{split}\]
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));
mat p_ikjl(const vec &a)

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

vec v = randu(6);
mat H = p_ikjl(v);
mat sym_dyadic(const mat &A, const mat &B)

Provides the dyadic product (in Voigt Notation) of two 2nd order tensors converted. The function returns a 6x6 matrix that correspond to a 4th order tensor. Note that such conversion to 6x6 matrices product correspond to a conversion with the component of the 4th order tensor correspond to the component of the matrix (such as stiffness matrices)

mat A = randu(3,3);
mat B = randu(3,3);
mat C = sym_dyadic(A,B);