C_C(lambda1,lambda3,rho,vL,C1,C2,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Carroll model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Carroll strain energy density function is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3\right) + \sqrt{3} C_2 \left(\sqrt{\frac{I_2}{J^{4/3}}} - \sqrt{3}\right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_CSH(lambda1,lambda3,rho,vL,C1,C2,C3,N,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Carroll model with a$I_1$ strain hardening term subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Gent-Carroll strain energy density function is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3 \right) + \sqrt{3} C_2 \left(\sqrt{\frac{I_2}{J^{4/3}}} - \sqrt{3}\right) + \frac{3^{1-N}}{2N} C_3 \left(\left(\frac{I_1}{J^{2/3}}\right)^N - 3^N \right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3,N: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_DCGT(lambda1,lambda3,rho,vL,C1,C2,C3,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Dobrynin & Carrillo-Gent-Thomas model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Dobrynin & Carrillo-Mooney-Rivlin strain energy density function is writen as $W = \frac{C_1}{6} \left(\frac{I_1}{J^{2/3}} - 3 + \frac{2C_3}{1 - \left(I_1/J^{2/3} - 3\right)/C_3} - 2C_3\right) + \frac{3}{2} C_2 ln \left(\frac{I_2}{3J^{4/3}}\right) + \frac{\kappa}{2}(J-1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_DCMR(lambda1,lambda3,rho,vL,C1,C2,C3,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Dobrynin & Carrillo-Mooney-Rivlin model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Dobrynin & Carrillo-Mooney-Rivlin strain energy density function is writen as $W = \frac{1 - C_3^2}{2(3 - 2C_3 + C_3^2)} C_1 \left(\frac{I_1}{J^{2/3}} - 3 + \frac{6}{C_3}\left(1-C_3\frac{I_1}{3J^{2/3}}\right)^{-1} - \frac{6}{C_3(1-C_3)}\right) + \frac{1}{2} C_2 \left(\frac{I_2}{J^{4/3}} - 3\right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

NB: This form of W does not give correct results. Use DCMR2 instead.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_DCMR2(lambda1,lambda3,rho,vL,C1,C2,C3,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Dobrynin & Carrillo-Mooney-Rivlin model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Dobrynin & Carrillo-Mooney-Rivlin strain energy density function is writen as $W = \frac{C_1}{6} \left( \frac{I_1}{J^{2/3}} - 3 + \frac{2C_3}{1 - \left(I_1/J^{2/3} - 3\right)/C_3} - 2C_3\right) + \frac{C_2}{2} \left(\frac{I_2}{J^{4/3}} - 3\right) + \frac{\kappa}{2}(J-1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_GC(lambda1,lambda3,rho,vL,C1,C2,C3,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Gent-Carroll model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Gent-Carroll strain energy density function is writen as $W = - \frac{1}{2} C_1 C_3 ln \left(1 - \frac{I_1 / J^{2/3} - 3}{C_3}\right) + \sqrt{3} C_2 \left(\sqrt{\frac{I_2}{J^{4/3}}} - \sqrt{3}\right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_GG(lambda1,lambda3,rho,vL,C1,C2,C3,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Gent-Gent model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Gent-Gent strain energy density function is writen as $W = - \frac{1}{2} C_1 C_3 ln \left(1 - \frac{I_1 / J^{2/3} - 3}{C_3}\right) + \frac{3}{2} C_2 ln\left(\frac{I_2}{3 J^{4/3}}\right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_GMR(lambda1,lambda3,rho,vL,C1,C2,C3,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Gent-Mooney-Rivlin model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Gent-Mooney-Rivlin strain energy density function is writen as $W = - \frac{1}{2} C_1 C_3 ln \left(1 - \frac{I_1 / J^{2/3} - 3}{C_3}\right) + \frac{1}{2} C_2 \left(\frac{I_2}{J^{4/3}} - 3\right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_GT(lambda1,lambda3,rho,vL,C1,C2,f,tau,n,betaprime)

Returns the 4-th order elastic tensor for a nearly incompressible hyperelastic medium described by the Gent-Thomas model for a large deformation specified through the elongations in two principal directions lambda1 and lambda3. The stretch ratio in the third direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Gent-Thomas strain energy density function is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3 \right) + \frac{3}{2} C_2 ln \left(\frac{I_2}{3J^{4/3}}\right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_GTSH(lambda1,lambda3,rho,vL,C1,C2,C3,N,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Gent-Thomas model with a $I_1$ strain hardening term subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Gent-Carroll strain energy density function is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3 \right) + \frac{3}{2} C_2 ln \left(\frac{I_2}{3J^{4/3}}\right) + \frac{3^{1-N}}{2N} C_3 \left(\left(\frac{I_1}{J^{2/3}}\right)^N - 3^N \right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3,N: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_GTSH2(lambda1,lambda3,rho,vL,C1,C2,C3,N,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Gent-Thomas model with a $I_2$ strain hardening term subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Mooney-Rivlin strain energy density function with a $I_2$ strain hardening term is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3\right) + \frac{3}{2} C_2 ln \left(\frac{I_2}{3J^{4/3}}\right) + \frac{3^{1-N}}{2N} C_3 \left(\left(\frac{I_2}{J^{4/3}}\right)^N - 3^N \right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3, N: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_MR(lambda1,lambda3,rho,vL,C1,C2,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Mooney-Rivlin model subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Mooney-Rivlin strain energy density function is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3\right) + \frac{1}{2} C_2 \left(\frac{I_2}{J^{4/3}} - 3\right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_MRSH(lambda1,lambda3,rho,vL,C1,C2,C3,N,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Mooney-Rivlin model with a $I_1$strain hardening term subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Mooney-Rivlin strain energy density function with a $I_1$ strain hardening term is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3\right) + \frac{1}{2} C_2 \left(\frac{I_2}{J^{4/3}} - 3\right) + \frac{3^{1-N}}{2N} C_3 \left(\left(\frac{I_1}{J^{2/3}}\right)^N - 3^N \right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3, N: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_MRSH2(lambda1,lambda3,rho,vL,C1,C2,C3,N,f,tau,n,betaprime)

Returns the 4-th order incremental elastic tensor for a nearly incompressible hyperelastic medium described by the Mooney-Rivlin model with a $I_2$ strain hardening term subjected to a large deformation specified by the elongations in the principal directions 1, lambda1, and 3, lambda3. The stretch ratio in the remaining direction is given by imposing the stress free condition in direction 2. We also consider viscoelastic material properties through the fractional Kelvin-Voigt model.

The Mooney-Rivlin strain energy density function with a $I_2$ strain hardening term is writen as $W = \frac{1}{2} C_1 \left(\frac{I_1}{J^{2/3}} - 3\right) + \frac{1}{2} C_2 \left(\frac{I_2}{J^{4/3}} - 3\right) + \frac{3^{1-N}}{2N} C_3 \left(\left(\frac{I_2}{J^{4/3}}\right)^N - 3^N \right) + \frac{K}{2} (J - 1)^2$, where $I_1$, $I_2$ and $J^2$ are the first three invariants of the left Cauchy-Green tensor.

Input arguments:

• lambda1: stretch ratio in direction 1
• lambda3: stretch ratio in direction 3
• rho: material density
• vL: longitudinal wave velocity
• C1, C2, C3, N: hyperelastic parameters
• f: frequency
• tau: relaxation time
• n: fractional derivative order
• betaprime: coupling between viscoelasticity and hyperelasticity

C_elastic(λ,μ)

Returns the 4-th order elastic tensor for a material with Lamé parameters λ and μ.

Celastodynamics2D(C,x2,D2,D3,D22,D33,A,udof)

Discretize the elastodynamics equation in 2D given the elasticity tensor C, , the collocation points x2, the derivation matrices D2, D3, D22, and D33, and the displacement degrees of freedom udof. In polar coordinates, we additionally need the matrix A for the derivation of the basis vectors. Returns the four matrices to build the polynomial eigenvalue problem (ik)^2 Lkk + (ik) Lk + L0 + ω^2 M = 0.

assemblelayer(Lkk,Lk,L0,Md,Lfkk,Lfk,Lf0,Mfd,bc,bcf,D0b,Bfkt,Bf0t)

Returns the matrices Fkk, Fk, F0 and Fd which define the polynomial eigenproblem for a bilayer with PEP Lkk, Lk, L0 and Md (Lfkk, Lfk, Lf0 and Mfd, respectively). The Matrices D0b, Bfkt, and Bf0t couple the normal displacement and stress at the interface.

cheb1extrema(N)

Computes the extrema of the Chebyshev polynomial of the first kind of degree N. AKA, Chebyshev-Gauss-Lobatto (CGL) Points, or Chebyshev points of the second kind.

cheb1roots(N)

Computes the roots of the Chebyshev polynomial of the first kind of degree N. AKA, Chebyshev-Gauss (CG) Points, or Chebyshev points of the first kind.

chebdif(N, M)

Computes the differentiation matrices D1, D2, ..., DM on Chebyshev nodes.

Input arguments:

• N: size of differentiation matrix.
• M: number of derivatives required [integer]. Note: 0 < M <= N-1.

Outputs

• x: the Chebyshev nodes.
• DM: DM[1:N,1:N,ell] contains ell-th derivative matrix, ell=1,...,M.

Details

The code implements two strategies for enhanced accuracy suggested by W. Don & S. Solomonoff in SIAM J. Sci. Comp. Vol. 6, pp. 1253–1268 (1994). The two strategies are (a) the use of trigonometric identities to avoid the computation of differences x[k]-x[j] & (b) the use of the "flipping trick" which is necessary since $sin(t)$ can be computed to high relative precision when t is small whereas $sin(pi-t)$ cannot.

elastodynamics1D(C,D2,D22,udof)

Discretize the elastodynamics equation in 1D given the elasticity tensor C, the derivation matrices D2 and D22, and the displacement degrees of freedom udof. Returns the four matrices to build the polynomial eigenvalue problem (ik)^2 Lkk + (ik) Lk + L0 + ω^2 M = 0.

elastodynamics1D_tube(C,x2,D2,D22,AA,udof)

Discretize the elastodynamics equation in 1D for a tube on given the elasticity tensor C, the collocation points x2, derivation matrices D2 and D22, and the displacement degrees of freedom udof. In polar coordinates, we additionally need the matrix AA for the derivation of the basis vectors. Returns the four matrices to build the polynomial eigenvalue problem (ik)^2 Lkk + (ik) Lk + L0 + ω^2 M = 0.

elastodynamics2D(C,D2d,D3d,D22d,D33d,D23d,udof)

Discretize the elastodynamics equation in 2D given the elasticity tensor C, the derivation matrices D2d, D3d, D22d, D33d, and D23d and the displacement degrees of freedom udof. Returns the four matrices to build the polynomial eigenvalue problem (ik)^2 Lkk + (ik) Lk + L0 + ω^2 M = 0.

lgrdiffleft(n,x)

Compute the first order differentiation matrix of size n × n associated to the Legendre Gauss Radau points vector x with x(1) = -1. Refer to L. N. Trefethen: Spectral Methods in Matlab(2000). Based on lgrdiff.m (SPECTRAL METHODS: Algorithms, Analysis and Applications, Jie Shen et al.).

lgrpointsleft(n)

Compute the Legendre n Legendre Gauss Radau points with x(1) = -1. Refer to L. N. Trefethen: Spectral Methods in Matlab(2000). Based on legsrd.m (SPECTRAL METHODS: Algorithms, Analysis and Applications, Jie Shen et al.).

lpoly(n,x)

Compute the Legendre polynomial of degree n at the value given in the vector x. Refer to L. N. Trefethen: Spectral Methods in Matlab(2000). Based on lepoly.m (SPECTRAL METHODS: Algorithms, Analysis and Applications, Jie Shen et al.).

modedirection(p::result1D)

Compute an array of integers idx associating each eigenvalue of a SCM 1D computation p to its principal displacement direction (i.e. 1, 2, or 3).

modedirection(s::result2D)

Compute an array of integers idx associating each eigenvalue of a SCM 2D computation s to its principal displacement direction (i.e. 1, 2, or 3).

modepolarisation(p::result1D)

Compute the average displacement polarisation pol ∈ [0, 1] where 0 is a purely in-plane polarisation and 1 a purely out-of-plane polarisation for the displacement field associated to each eigenvalue of a SCM 1D computation p.

modepolarisation(s::result2D)

Compute the average displacement polarisation pol ∈ [0, 1] where 0 is a purely in-plane polarisation and 1 a purely out-of-plane polarisation for the displacement field associated to each eigenvalue of a SCM 2D computation s.

modesurfacepolarisation(p::result1D)

Compute the surface displacement polarisation pol ∈ [0, 1] where 0 is a purely in-plane polarisation and 1 a purely out-of-plane polarisation for the displacement field associated to each eigenvalue of a SCM 1D computation p.

modesymmetry(p::result1D)

Determine the symmetry sym (1 if symmetric, 0 if antisymmetric) of the displacement field associated to each eigenvalue of a SCM 1D computation p.

modesymmetry(s::result2D; dir)

Determine the symmetry sym (1 if symmetric, 0 if antisymmetric) of the displacement field associated to each eigenvalue of a SCM 2D computation s. Specify the Keyword dir="ip" (respectively dir="oop") when looking to determine the symmetry of in-plane (respectively out-of-plane) modes.