# Maths

## The Rotation Library: rotation.hpp

The Rotation Library provides various tools to rotate matrices that represents stiffness, compliance and localisation tensors, and vectors that represents strain and stress tensors.

# The Rotation Library: rotation.hpp

### Rot_strain(vec, mat)

Rotates a strain vector, providing a 3×3 orientation matrix. This function computes the 6×6 rotation matrix to rotate the strain vector in Voigt notation, and rotate it accordingly.

vec E = randu(6);
double theta = 30*pi/180.;
mat DR;
DR << cos(theta) << sin(theta) << 0 << endr
<< -1.*sin(theta) << cos(theta) << 0 << endr
<< 0 << 0 << 1;
Rot_strain(E, DR);


E is rotated according to DR.

### Rot_stress(vec, mat)

Rotates a stress vector, providing a 3×3 orientation matrix. This function computes the 6×6 rotation matrix to rotate the stress vector in Voigt notation, and rotate it accordingly.

vec S = randu(6);
double theta = 30*pi/180.;
mat DR;
DR << cos(theta) << sin(theta) << 0 << endr
<< -1.*sin(theta) << cos(theta) << 0 << endr
<< 0 << 0 << 1;
Rot_stress(S, DR);


S is rotated according to DR.

### fillQE(double, int)

Fills the 6×6 strain rotation matrix from a given angle and rotation axis (1, 2 or 3)

double theta = 30*pi/180.;
int dir = 1;
fillQE(theta, dir);


QE is filled for a angle of 30° around the axis 1.

### fillQS(double, int)

Fills the 6×6 stress rotation matrix from a given angle and rotation axis (1, 2 or 3)

double theta = 30*pi/180.;
int dir = 1;
fillQS(theta, dir);


QS is filled for a angle of 30° around the axis 1.

### rotateL(mat, double, int)

Rotates the 6×6 stiffness matrix from a given angle and rotation axis (1, 2 or 3)

mat L = randu(6,6);
double theta = 30*pi/180.;
int dir = 1;
rotateL(L, theta, dir);


L is rotated for a angle of 30° around the axis 1.

### rotateL(mat, double, int)

Rotates the 6×6 compliance matrix from a given angle and rotation axis (1, 2 or 3)

mat M = randu(6,6);
double theta = 30*pi/180.;
int dir = 1;
rotateL(L, theta, dir);


M is rotated for a angle of 30° around the axis 1.

### rotateA(mat, double, int)

Rotates the 6×6 strain concentration tensor from a given angle and rotation axis (1, 2 or 3)

mat A = randu(6,6);
double theta = 30*pi/180.;
int dir = 1;
rotateA(A, theta, dir);


A is rotated for a angle of 30° around the axis 1.

# The Statistics Library: stats.hpp

### normal_distrib(double, const double, const double)

provides an approximation (near $$10^{-7}$$) of the normal law distribution function of the variable x with $$mean$$, the mean of the distribution and $$dev$$, it’s standard deviation.

$$N=1-\frac{1}{ exp(-0.5*X_{norm}²)\sqrt{2pi}}(b_1k+b_2k²+b_3k³+b_4k⁴+b_5k⁵)$$

with $$X_{norm}=(x-mean)/dev$$,   $$k=\frac{1}{1+0.2316419X_{norm}}$$

and $$b_1=\small{0.319381530}$$, $$b_2=\small{0.356563782}$$, $$b_3=\small{1.781477937}$$, $$b_4=\small{-1.821255978}$$, $$b_5=\small{1.330274429})$$

double x = (double)rand();
const double mean = (double)rand();
const double dev = (double)rand();
double N= normal_distrib(x, mean, dev)


### ODF_sd(const double, const double, const double, const double, const double, const double)

provides a standard orientation distribution function of $$\theta$$ defined by :

$$F= a_1 * cos(\theta-m)^{(2*p_1)} + a_2 * cos(\theta-m)^{(2*p_2 + 1)} * sin(\theta)^{(2*p_2)}$$

example:

const double theta = (double)rand();
const double p1 = (double)rand();
const double p2 = (double)rand();
const double m =(double)rand();
const double a1 =(double)rand();
const double a2 =(double)rand();
double F= ODF_sd(theta, m, a1, a2, p1, p2);

### Gaussian(const double, const double, const double, const double)

provides the Gaussian distribution function of the variable X, with $$mean$$, the mean of the distribution and $$sd$$, it’s standard deviation.
the law is given by:

$$G=\frac{ampl}{(sd\sqrt{2 pi})} * exp(\small{-0.5} (\frac{X – mean}{sd})²)$$

example:

const double X = (double)rand();
const double mean = (double)rand();
const double sd = (double)rand();
const double ampl =1;
double G= Gaussian(x, mean, sd, ampl);

### Mult_Gaussian(const double, const int, const vec, const vec, const vec)

provides the cumulative distribution of the gaussian function of the variable X,with $$mean$$, the means vector of the distribution, $$sd$$, the standard deviations vector and Npeak the numberof peaks.

$$MG=\Sigma^{Npeak}_{i=1}\frac{ampl(i)}{(sd(i)\sqrt{2 pi})} * exp(\small{-0.5} (\frac{X – mean(i)}{sd(i)})²)$$

example:

const double X = (double)rand();
const int Npeak = (int)rand();
const vec mean = randu(Npeak);
const vec sd = randu(Npeak);
const vec ampl =randu(Npeak);
double Mg= Mult_Gaussian(x,Npeak, mean, sd, ampl);



### Lorentzian(const double, const double, const double, const double)

provides the Lorentzian distribution function of the variable X, with $$mean$$ is the location parameter wiche specifies the peak location of the distribution and $$Widht$$ is the full widht at the half maximum.the function is given by :

$$L=\frac{Widht .Ampl}{2pi(x-mean)²+(\frac{Widht}{2})²}$$
example:

const double X = (double)rand();
const double mean = (double)rand();
const double Widht = (double)rand();
const double ampl =1;
double L= Lorentzian(x, mean, Widht, ampl);

### Mult_Lorentzian(const double, const int, const vec, const vec, const vec)

provides the cumulative distribution of the Lorentzian function of the variable X,with $$mean$$, the location parameters vector of the distribution, $$Widht$$, the full widhts at the half maximum vector and Npeak the number of peaks. the function is given by:

$$ML=\Sigma^{Npeak}_{i=1}\frac{Widht(i) .Ampl(i)}{2pi(x-mean(i))²+(\frac{Widht(i)}{2})²}$$

example:

const double X = (double)rand();
const int Npeak = (int)rand();
const vec mean = randu(Npeak);
const vec Widht = randu(Npeak);
const vec ampl =randu(Npeak);
double Mg= Mult_Lorentzian(x,Npeak, mean, Widht, ampl);

### PseudoVoigt(const double, const double, const double, const double, const double, const double)

Provides the Pseudo-Voigt function $$P$$of the variable X. This function is the sum of a gaussian and a lorentzian function wich have the same position and the same aera. the expression is given by :
$$P= \eta L+(1-\eta)G$$

$$=\eta\frac{Widht_Lor.Ampl}{2pi(x-mean)²+(\frac{Widht_Lor}{2})²}+(1-\eta)\frac{ampl}{(sd_gau\sqrt{2 pi})} * exp(\small{-0.5}(\frac{X – mean}{sd_gau})²)$$

with $$\eta$$ the lorentzian factor.

example:

const double X = (double)rand();
const double eta = (double)rand();
const double mean = (double)rand();
const double Widht_Lor = (double)rand();
const double sd_gau = (double)rand();
const double ampl =1;
double Pv= PseudoVoigt(x, eta,mean, Widht_Lor,sd_gau, ampl);

### Mult_PseudoVoigt(const double, const int, const vec, const vec, const vec, const vec, const vec)

provides the cumulative distribution of the Pseudo-Voigt function of the variable X given by :
$$MP=\Sigma^{Npeak}_{i=1}\eta(i) L(i)+(1-\eta(i))G(i)$$

$$=\Sigma^{Npeak}_{i=1}\eta(i)\frac{Widht_Lor(i).Ampl(i)}{2pi(x-mean(i))²+(\frac{Widht_Lor(i)}{2})²}+(1-\eta(i))\frac{ampl(i)}{(sd_gau(i)\sqrt{2 pi})} * exp(\small{-0.5}(\frac{X – mean(i)}{sd_gau(i)})²)$$

example:

const double X = (double)rand();
const int Npeak = (int)rand();
const vec eta = randu(Npeak);
const vec mean = randu(Npeak);
const vec Widht_Lor = randu(Npeak);
const vec sd_gau = randu(Npeak);
const vec ampl =randu(Npeak);
double Pv= PseudoVoigt(x,Npeak, eta,mean, Widht_Lor,sd_gau, ampl);

### Pearson7(const double, const double, const double, const double, const double)

Provides the $$7^(th)$$ Pearson function $$P7$$of the variable X with $$invWidth$$ a normalization factor and $$shape$$ the shape parameter.The expression is given by:

$$P_7= max(1+\frac{(invWidth(X-mean))²}{shape})^{shape}$$

example:

const double X = (double)rand();
const double invWidth = (double)rand();
const double mean = (double)rand();
const double shape = (double)rand();
const double max = 1;
double P_7= Pearson7(x,mean, inv_width,shape, max);

### Mult_Pearson7(onst double, const int, const vec, const vec, const vec, const vec)

Provides the cumulative distribution of the $$7^(th)$$ Pearson function $$P7$$of the variable X with $$invWidth$$ a normalization factor, $$shape$$ the shape parameter and Npeak the number of peaks. The expression is given by:

$$MP_7=\Sigma^{Npeak}_{i=1}max(i)(1+\frac{(invWidth(i)(X-mean(i)))²}{shape(i)})^{shape(i)}$$

example:

const double X = (double)rand();
const int Npeak = (int)rand();
const vec invWidth = vec randu(Npeak);
const vec mean = vec randu(Npeak);
const vec shape = vec randu(Npeak);
const vec max = vec randu(Npeak);
double MP_7= Mult_Pearson7(x,Npeak, mean, inv_width,shape, max);