$$ \newcommand{\C}{{\mathbb{{C}}}} \newcommand{\R}{{\mathbb{{R}}}} \newcommand{\Q}{{\mathbb{{Q}}}} \newcommand{\Z}{{\mathbb{{Z}}}} \newcommand{\N}{{\mathbb{{N}}}} \newcommand{\uu}[1]{{\boldsymbol{{#1}}}} \newcommand{\uuuu}[1]{{\symbb{{#1}}}} \newcommand{\uv}[1]{{\underline{{#1}}}} \newcommand{\ve}[1]{{\uv{{e}}_{{#1}}}} \newcommand{\x}{{\uv{{x}}}} \newcommand{\n}{{\uv{{n}}}} \newcommand{\eps}{{\uu{{\varepsilon}}}} \newcommand{\E}{{\uu{{E}}}} \newcommand{\sig}{{\uu{{\sigma}}}} \newcommand{\Sig}{{\uu{{\Sigma}}}} \newcommand{\cod}{{\uv{{\symscr{b}}}}} \newcommand{\trans}[1]{{{}^{t}{#1}}} \newcommand{\sotimes}{{\stackrel{s}{\otimes}}} \newcommand{\sboxtimes}{\stackrel{s}{\boxtimes}} \newcommand{\norm}[1]{{\lVert{{#1}}\rVert}} \newcommand{\ud}{{\,\mathrm{d}}} \newcommand{\mat}{\mathsf} \DeclareMathOperator{\arcosh}{arcosh} \DeclareMathOperator{\divz}{div} \DeclareMathOperator{\divu}{\uv{div}} \DeclareMathOperator{\hess}{hess} \DeclareMathOperator{\gradu}{\uv{grad}} \DeclareMathOperator{\graduu}{\uu{grad}} \DeclareMathOperator{\Mat}{Mat} \DeclareMathOperator{\tr}{tr} \DeclareMathOperator{\ISO}{ISO} \newcommand{\jump}[1]{[\hspace*{-.15em}[\hspace*{.1em}{#1}% \hspace*{.1em}]\hspace*{-.15em}]} $$

4  The tensor object

Objectives

This tutorial presents the object tensor which is the main structure of echoes allowing to represent symmetric second-order or fourth-order tensors both in matrix and synthetic forms and containing information about anisotropy.

Imports

import numpy as np
from echoes import *
import math

np.set_printoptions(precision=6, suppress=True)
# to display only 6 significant digits of array components

Definition of the tensor object

The tensor object is a structure designed to represent symmetric second-order or fourth-order tensors. It gathers four members

• the material symmetry: ISO (isotropic), TI (transversely isotropic), ORTHO (orthotropic) or ANISO (anisotropic) [other material symmetries exist but have not been implemented]
• a condensed vector of parameters [as shown below the size of this vector unambiguously defines the order of the tensor]
• a vector of angles (in radians) [if required]
• a matrix [\(3×3\) for a second-order or \(6×6\) for a fourth-order tensor]

Warning

It is important to note that the tensor object is designed only to contain or retrieve information about a symmetric second-order or fourth-order tensor (for instance the material symmetry) and not to perform calculation. Although simple calculations such as addition, subtraction, multiplication by a scalar, inversion are available, most of other operations (products, contractions…) have not been implemented. To do so, it is necessary to extract the parameters or matrix and perform operations using usual numpy or scipy operations (dot, einsum…) before eventually building a new tensor.

If C is a python variable of tensor type, its members can be accessed by

• material symmetry: C.sym
• parameters: C.param or C.p
• angles: C.angles
• matrix: C.array or C.a

Second-order tensors

A symmetric second-order tensor object is designed by means of its three eigenvalues and the orientation of the eigenbasis through three Euler angles \(\theta\), \(\phi\) and \(\psi\) (see Fig. A.1). It is built using one of the following constructors

• T = tensor(T₁, T₂, T₃, θ=0, ϕ=0, ψ=0) # default values of the angles are 0
• T = tensor([T₁, T₂, T₃], angles=[0, 0, 0])

Warning

• The eigenvalues and eigenvectors are systematically reordered in decreasing order of eigenvalues and angles are recomputed accordingly.
• The eigenvalues are analyzed in order to characterize the symmetry type between ISO, TI and ORTHO.
• The Euler angles are not unique since unit eigenvectors are determined up to a multiplicative factor of ±1.
• Note that the six parameters of the constructor tensor should not be confused with the components of the Kelvin-Mandel representation of symmetric second-order tensors (see Section A.2).

T = tensor(2., 1., 3.) ; print(T)

Order 2 ORTHO tensor | Param(size=3)=[ 3 2 1 ] | Angles(size=3)=[ 1.5708 4.71239 3.14159 ]
[ 2 0 0 
  0 1 0 
  0 0 3 ]

π = math.pi
T = tensor(2., 1., 2., π/3, π/4, π/5) ; print(T)

Order 2 TI tensor | Param(size=3)=[ 2 2 1 ] | Angles(size=3)=[ 2.10486 5.84624 0.619825 ]
[ 1.3918 0.284068 0.396985 
  0.284068 1.86732 -0.185416 
  0.396985 -0.185416 1.74088 ]

T = tensor([2., 1., 3.], angles=[π/3, π/4, π/5]) ; print(T)

Order 2 ORTHO tensor | Param(size=3)=[ 3 2 1 ] | Angles(size=3)=[ 2.10486 5.84624 2.19062 ]
[ 1.7668 0.659068 0.703171 
  0.659068 2.24232 0.120771 
  0.703171 0.120771 1.99088 ]

A symmetric second-order tensor can also be built from a symmetric \(3×3\) matrix. In this case a diagonalization is performed so as to find the param vector (the eigenvalues) and the angles determining the orientation of the eigenvectors (see Fig. A.1)

M = np.array([ [1.7668, 0.659068, 0.703171], 
               [0.659068, 2.24232, 0.120771], 
               [0.703171, 0.120771, 1.99088] ])
T = tensor(M) ; print(T)

Order 2 ORTHO tensor | Param(size=3)=[ 3 2 1 ] | Angles(size=3)=[ 1.03673 2.70464 0.950971 ]
[ 1.7668 0.659068 0.703171 
  0.659068 2.24232 0.120771 
  0.703171 0.120771 1.99088 ]

Isotropic fourth-order tensors (ISO)

The special isotropic tensors \(\uuuu{I}\), \(\uuuu{J}\) and \(\uuuu{K}\) are available under the global variables tId4, tJ4 and tK4.

for T in [tId4, tJ4, tK4]:
    print(T)

Order 4 ISO tensor | Param(size=2)=[ 1 1 ] | Angles(size=0)=[ ]
[ 1 0 0 0 0 0 
  0 1 0 0 0 0 
  0 0 1 0 0 0 
  0 0 0 1 0 0 
  0 0 0 0 1 0 
  0 0 0 0 0 1 ]

Order 4 ISO tensor | Param(size=2)=[ 1 0 ] | Angles(size=0)=[ ]
[ 0.333333 0.333333 0.333333 0 0 0 
  0.333333 0.333333 0.333333 0 0 0 
  0.333333 0.333333 0.333333 0 0 0 
  0 0 0 0 0 0 
  0 0 0 0 0 0 
  0 0 0 0 0 0 ]

Order 4 ISO tensor | Param(size=2)=[ 1.11022e-16 1 ] | Angles(size=0)=[ ]
[ 0.666667 -0.333333 -0.333333 0 0 0 
  -0.333333 0.666667 -0.333333 0 0 0 
  -0.333333 -0.333333 0.666667 0 0 0 
  0 0 0 1 0 0 
  0 0 0 0 1 0 
  0 0 0 0 0 1 ]

Any fourth-order tensor depends on two parameters \(\alpha\) and \(\beta\) such that \(\uuuu{T}=\alpha \uuuu{J} + \beta \uuuu{K}\) and can be built by one of these constructors

• T = tensor(α, β)
• T = tensor([α, β])
• T = tensor(np.array([α, β]))

α, β = 7.2, 6.1
𝕋 = tensor(α, β)
print("𝕋 =\n", 𝕋)
print("𝕋.sym =", 𝕋.sym)
print("𝕋.a =\n", 𝕋.a)
print("𝕋.angles =", 𝕋.angles)
print("𝕋.p =", 𝕋.p)

𝕋 =
 Order 4 ISO tensor | Param(size=2)=[ 7.2 6.1 ] | Angles(size=0)=[ ]
[ 6.46667 0.366667 0.366667 0 0 0 
  0.366667 6.46667 0.366667 0 0 0 
  0.366667 0.366667 6.46667 0 0 0 
  0 0 0 6.1 0 0 
  0 0 0 0 6.1 0 
  0 0 0 0 0 6.1 ]

𝕋.sym = ISO
𝕋.a =
 [[6.466667 0.366667 0.366667 0.       0.       0.      ]
 [0.366667 6.466667 0.366667 0.       0.       0.      ]
 [0.366667 0.366667 6.466667 0.       0.       0.      ]
 [0.       0.       0.       6.1      0.       0.      ]
 [0.       0.       0.       0.       6.1      0.      ]
 [0.       0.       0.       0.       0.       6.1     ]]
𝕋.angles = []
𝕋.p = [7.2 6.1]

𝕀, 𝕁, 𝕂 = tId4, tJ4, tK4 # just to make it look nice!

𝕋2 = α*𝕁 + β*𝕂
print("𝕋2 =\n", 𝕋2)
assert np.allclose(𝕋.a, 𝕋2.a), "error 𝕋.a and 𝕋2.a should be equal"
assert np.allclose(𝕋.p, 𝕋2.p), "error 𝕋.p and 𝕋2.p should be equal"

𝕋2 =
 Order 4 ISO tensor | Param(size=2)=[ 7.2 6.1 ] | Angles(size=0)=[ ]
[ 6.46667 0.366667 0.366667 0 0 0 
  0.366667 6.46667 0.366667 0 0 0 
  0.366667 0.366667 6.46667 0 0 0 
  0 0 0 6.1 0 0 
  0 0 0 0 6.1 0 
  0 0 0 0 0 6.1 ]

A symmetric fourth-order tensor representing a stiffness tensor can also be built by means of

• bulk and shear moduli stiff_kmu(k, µ) → \(3\,k \uuuu{J} + 2\,\mu \uuuu{K}\)
• Young modulus and Poisson ratio stiff_Enu(E, ν) → \(\frac{E}{1-2\,\nu} \uuuu{J} + \frac{E}{1+\nu} \uuuu{K}\)
• Lamé moduli stiff_lambdamu(λ, μ) → \(3\,\lambda \uuuu{J} + 2\,\mu \uuuu{I}\)

Some converters are also available (see conversion table)

• E = E_from_kmu(k, µ)
• ν = nu_from_kmu(k, µ)
• k = k_from_Enu(E, ν)
• µ = mu_from_Enu(E, ν)
• E, ν = Enu_from_kmu(k, µ)
• k, µ = kmu_from_Enu(E, ν)

k, µ = 72.5, 32.7
ℂ = stiff_kmu(k, µ)
print("ℂ =\n", ℂ)
print(f"E = {round(ℂ.E,2)} ; ν = {round(ℂ.nu,2)} ; λ = {round(ℂ.lamelambda,2)} ; μ = {round(ℂ.mu,2)}")
E, ν = Enu_from_kmu(k, µ)
assert np.allclose([k, μ], [*kmu_from_Enu(E, ν)]), "error kmu_from_Enu(E, ν) should return k, μ"

ℂ =
 Order 4 ISO tensor | Param(size=2)=[ 217.5 65.4 ] | Angles(size=0)=[ ]
[ 116.1 50.7 50.7 0 0 0 
  50.7 116.1 50.7 0 0 0 
  50.7 50.7 116.1 0 0 0 
  0 0 0 65.4 0 0 
  0 0 0 0 65.4 0 
  0 0 0 0 0 65.4 ]

E = 85.28 ; ν = 0.3 ; λ = 50.7 ; μ = 32.7

Symmetric transversely isotropic fourth-order tensors (TI)

Such a tensor depends on 5 parameters and 2 angles \(\theta\) and \(\phi\) defining the orientation of the normal to the isotropy plane (see Fig. A.1).

The decomposition on the symmetric Walpole tensors writes by means of five parameters \((t_i)_{i=1,\ldots,5}\). In the frame where \(\uv{n}=\ve{3}\) the matrix is given by

\[ \uuuu{T}=\sum_{i=1}^5 t_i \uuuu{W}^s_i \mapsto \left( \begin{array}{cccccc} \frac{t_2 + t_4}{2} & \frac{t_2 - t_4}{2} & \frac{\sqrt{2} \, t_3}{2} & 0 & 0 & 0 \\ \frac{t_2 - t_4}{2} & \frac{t_2 + t_4}{2} & \frac{\sqrt{2} \, t_3}{2} & 0 & 0 & 0 \\ \frac{\sqrt{2} \, t_3}{2} & \frac{\sqrt{2} \, t_3}{2} & t_1 & 0 & 0 & 0 \\ 0 & 0 & 0 & t_5 & 0 & 0 \\ 0 & 0 & 0 & 0 & t_5 & 0 \\ 0 & 0 & 0 & 0 & 0 & t_4 \\ \end{array} \right) \]

A fourth-order transversely isotropic tensor can be built by one of these constructors

• T = tensor(t1, t2, t3, t4, t5, θ=0, ϕ=0)
• T = tensor(param, angles=[0, 0])

where param is a list or a numpy.ndarray of 5 items and angles is a list or a numpy.ndarray of 2 items θ, ϕ (if angles is omitted the angles are null and the axis is oriented along \(\ve{3}\))

param = [3., 4., math.sqrt(2), 2., 1.]
angles = [π/3, π/4]
print(tensor(param))
print(tensor(param, angles))

Order 4 TI tensor | Param(size=5)=[ 3 4 1.41421 2 1 ] | Angles(size=2)=[ 0 0 ]
[ 3 1 1 0 0 0 
  1 3 1 0 0 0 
  1 1 3 0 0 0 
  0 0 0 1 0 0 
  0 0 0 0 1 0 
  0 0 0 0 0 2 ]

Order 4 TI tensor | Param(size=5)=[ 3 4 1.41421 2 1 ] | Angles(size=2)=[ 1.0472 0.785398 ]
[ 2.53125 1.28125 1.1875 0.32476 -0.108253 -0.132583 
  1.28125 2.53125 1.1875 -0.108253 0.32476 -0.132583 
  1.1875 1.1875 2.625 -0.216506 -0.216506 0.265165 
  0.32476 -0.108253 -0.216506 1.75 3.33067e-16 0.153093 
  -0.108253 0.32476 -0.216506 3.33067e-16 1.75 0.153093 
  -0.132583 -0.132583 0.265165 0.153093 0.153093 1.8125 ]

It is also possible to build a TI tensor by means of 5 \(C_{ijkl}\) moduli:

\[ (T_{1 1 1 1}, T_{1 1 2 2}, T_{1 1 3 3}, T_{3 3 3 3}, T_{2 3 2 3}) \mapsto \left( \begin{array}{cccccc} T_{1 1 1 1} & T_{1 1 2 2} & T_{1 1 3 3} & 0 & 0 & 0 \\ T_{1 1 2 2} & T_{1 1 1 1} & T_{1 1 3 3} & 0 & 0 & 0 \\ T_{1 1 3 3} & T_{1 1 3 3} & T_{3 3 3 3} & 0 & 0 & 0 \\ 0 & 0 & 0 & 2 \, T_{2 3 2 3} & 0 & 0 \\ 0 & 0 & 0 & 0 & 2 \, T_{2 3 2 3} & 0 \\ 0 & 0 & 0 & 0 & 0 & T_{1 1 1 1} - T_{1 1 2 2} \\ \end{array} \right) \]

by

T = stiff_TI(T₁₁₁₁, T₁₁₂₂, T₁₁₃₃, T₃₃₃₃, T₂₃₂₃, θ=0, ϕ=0)

T = stiff_TI(6.2, 5.6, 1.3, 2.1, 2.3)
assert all(v > 0 for v in np.linalg.eigvals(T.a)), "one eigenvalue is negative"
# check that T is positive definite
print(T)
print(stiff_TI(6.2, 5.6, 1.3, 2.1, 2.3, π/3, π/4))

Order 4 TI tensor | Param(size=5)=[ 2.1 11.8 1.83848 0.6 4.6 ] | Angles(size=2)=[ 0 0 ]
[ 6.2 5.6 1.3 0 0 0 
  5.6 6.2 1.3 0 0 0 
  1.3 1.3 2.1 0 0 0 
  0 0 0 4.6 0 0 
  0 0 0 0 4.6 0 
  0 0 0 0 0 0.6 ]

Order 4 TI tensor | Param(size=5)=[ 2.1 11.8 1.83848 0.6 4.6 ] | Angles(size=2)=[ 1.0472 0.785398 ]
[ 5.48281 1.88281 2.58437 -2.43028 -0.698233 -0.855157 
  1.88281 5.48281 2.58438 -0.698233 -2.43028 -0.855157 
  2.58437 2.58438 5.83125 -0.50879 -0.50879 -2.74446 
  -2.43028 -0.698233 -0.50879 2.44375 0.84375 0.421006 
  -0.698233 -2.43028 -0.50879 0.84375 2.44375 0.421006 
  -0.855157 -0.855157 -2.74446 0.421006 0.421006 2.61563 ]

A engineering construction based on directional moduli and Poisson ratio is also available:

\[ (E_1, E_3, \nu_{1 2}, \nu_{3 1}, G_{3 1}) \mapsto \left( \begin{array}{cccccc} \frac{1}{E_1} & \frac{ - \nu_{1 2}}{E_1} & \frac{ - \nu_{3 1}}{E_3} & 0 & 0 & 0 \\ \frac{ - \nu_{1 2}}{E_1} & \frac{1}{E_1} & \frac{ - \nu_{3 1}}{E_3} & 0 & 0 & 0 \\ \frac{ - \nu_{3 1}}{E_3} & \frac{ - \nu_{3 1}}{E_3} & \frac{1}{E_3} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{2 \, G_{3 1}} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{2 \, G_{3 1}} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1 + \nu_{1 2}}{E_1} \\ \end{array} \right) \]

by

S = comp_TI(E₁, E₃, ν₁₂, ν₃₁, G₃₁, θ=0, ϕ=0)

S = comp_TI(1., 4., 0.2, 0.1, 2.)
print(S)
print(inv(S))
print(comp_TI(1., 4., 0.2, 0.1, 2., π/3, π/4))

Order 4 TI tensor | Param(size=5)=[ 0.25 0.8 -0.0353553 1.2 0.25 ] | Angles(size=2)=[ 0 0 ]
[ 1 -0.2 -0.025 0 0 0 
  -0.2 1 -0.025 0 0 0 
  -0.025 -0.025 0.25 0 0 0 
  0 0 0 0.25 0 0 
  0 0 0 0 0.25 0 
  0 0 0 0 0 1.2 ]

Order 4 TI tensor | Param(size=5)=[ 4.02516 1.25786 0.177888 0.833333 4 ] | Angles(size=2)=[ 0 0 ]
[ 1.0456 0.212264 0.125786 0 0 0 
  0.212264 1.0456 0.125786 0 0 0 
  0.125786 0.125786 4.02516 0 0 0 
  0 0 0 4 0 0 
  0 0 0 0 4 0 
  0 0 0 0 0 0.833333 ]

Order 4 TI tensor | Param(size=5)=[ 0.25 0.8 -0.0353553 1.2 0.25 ] | Angles(size=2)=[ 1.0472 0.785398 ]
[ 0.53125 0.04375 -0.015625 0.205681 -0.205681 -0.251907 
  0.04375 0.53125 -0.015625 -0.205681 0.205681 -0.251907 
  -0.015625 -0.015625 0.6625 -0.248982 -0.248982 0.198874 
  0.205681 -0.205681 -0.248982 0.75625 -0.20625 -0.107165 
  -0.205681 0.205681 -0.248982 -0.20625 0.75625 -0.107165 
  -0.251907 -0.251907 0.198874 -0.107165 -0.107165 0.7125 ]

Symmetric orthotropic fourth-order tensors (ORTHO)

Such a tensor depends on 9 parameters and 3 angles \(\theta\), \(\phi\) and \(\psi\) defining the orientation of an orthonormal frame (see Fig. A.1).

The convention here relies on the components \(T_{1 1 1 1}, T_{1 1 2 2}, T_{1 1 3 3}, T_{2 2 2 2}, T_{2 2 3 3}, T_{3 3 3 3}, T_{2 3 2 3}, T_{3 1 3 1}, T_{1 2 1 2}\) as inputs. In the canonical frame the matrix is given by

\[ \small (T_{1 1 1 1}, T_{1 1 2 2}, T_{1 1 3 3}, T_{2 2 2 2}, T_{2 2 3 3}, T_{3 3 3 3}, T_{2 3 2 3}, T_{3 1 3 1}, T_{1 2 1 2}) \mapsto \left( \begin{array}{cccccc} T_{1 1 1 1} & T_{1 1 2 2} & T_{1 1 3 3} & 0 & 0 & 0 \\ T_{1 1 2 2} & T_{2 2 2 2} & T_{2 2 3 3} & 0 & 0 & 0 \\ T_{1 1 3 3} & T_{2 2 3 3} & T_{3 3 3 3} & 0 & 0 & 0 \\ 0 & 0 & 0 & 2 \, T_{2 3 2 3} & 0 & 0 \\ 0 & 0 & 0 & 0 & 2 \, T_{3 1 3 1} & 0 \\ 0 & 0 & 0 & 0 & 0 & 2 \, T_{1 2 1 2} \\ \end{array} \right) \]

A fourth-order isotropic tensor can be built by one of these constructors

• T = tensor(T₁₁₁₁, T₁₁₂₂, T₁₁₃₃, T₂₂₂₂, T₂₂₃₃, T₃₃₃₃, T₂₃₂₃, T₃₁₃₁, T₁₂₁₂, θ=0, ϕ=0, ψ=0)
• T = tensor(param, angles=[0, 0, 0])

where param is a list or a numpy.ndarray of 9 items and angles is a list or a numpy.ndarray of 3 items θ, ϕ, ψ (if angles is omitted the frame is the canonical one).

T = tensor([9.,3.,5.,6.,3.,7.,1.,4.,3.])
assert all(v > 0 for v in np.linalg.eigvals(T.a)), "one eigenvalue is negative" # check that T is positive definite
print(T)
print(tensor([9.,3.,5.,6.,3.,7.,1.,4.,3.], [π/3, π/4, π/5]))

Order 4 ORTHO tensor | Param(size=9)=[ 9 3 5 6 3 7 1 4 3 ] | Angles(size=3)=[ 0 0 0 ]
[ 9 3 5 0 0 0 
  3 6 3 0 0 0 
  5 3 7 0 0 0 
  0 0 0 2 0 0 
  0 0 0 0 8 0 
  0 0 0 0 0 6 ]

Order 4 ORTHO tensor | Param(size=9)=[ 9 3 5 6 3 7 1 4 3 ] | Angles(size=3)=[ 1.0472 0.785398 0.628319 ]
[ 5.83478 3.57431 3.79988 0.0416272 0.710013 -0.509607 
  3.57431 9.74488 2.26743 -1.24237 -0.0801039 0.934783 
  3.79988 2.26743 9.1371 -0.976365 1.31117 0.522827 
  0.0416272 -1.24237 -0.976365 3.25126 -0.192475 -0.544059 
  0.710013 -0.0801039 1.31117 -0.192475 5.36512 -1.53671 
  -0.509607 0.934783 0.522827 -0.544059 -1.53671 4.66686 ]

Symmetric anisotropic fourth-order tensors (ANISO)

Such a tensor depends on 21 parameters, by convention here the 21 components of the Kelvin-Mandel notation. Components are considered line after line of the upper triangular part of the Kelvin-Mandel matrix

\[ \begin{array}{cccccc} (T_{11},&T_{12},&T_{13},&T_{14},&T_{15},&T_{16}\\ &T_{22},&T_{23},&T_{24},&T_{25},&T_{26}\\ & &T_{33},&T_{34},&T_{35},&T_{36}\\ & & &T_{44},&T_{45},&T_{46}\\ & & & &T_{55},&T_{56}\\ & & & & &T_{66}) \end{array} \mapsto \left( \begin{array}{cccccc} T_{11}&T_{12}&T_{13}&T_{14}&T_{15}&T_{16}\\ T_{12}&T_{22}&T_{23}&T_{24}&T_{25}&T_{26}\\ T_{13}&T_{23}&T_{33}&T_{34}&T_{35}&T_{36}\\ T_{14}&T_{24}&T_{34}&T_{44}&T_{45}&T_{46}\\ T_{15}&T_{25}&T_{35}&T_{45}&T_{55}&T_{56}\\ T_{16}&T_{26}&T_{36}&T_{46}&T_{56}&T_{66} \end{array} \right) \]

Warning

It is important to recall here that the parameters \(T_{IJ}\) with \(1\leq I,J\leq 6\) do not all coincide with the components \(T_{ijkl}\). Indeed the convention presented in Section A.2 implies that \(T_{11}=T_{1111}\) but \(T_{14}=\sqrt{2}\,T_{1123}\) and \(T_{44}=2\,T_{2323}\).

Such a tensor is built by

T = tensor(param)

where param is a list or a numpy.ndarray of 21 items.

param = [ 2.66011, 1.26432, 0.662772, 1.9402, 1.54905, 1.10384, 3.6072, 1.78964, 2.0247, 1.2701, 1.1089, 2.743, 1.3367, 1.2962, 0.897632, 4.42684, 2.05632, 1.52686, 3.54431, 1.3445, 1.99356 ]
T = tensor(param)
assert all(v > 0 for v in np.linalg.eigvals(T.a)), "one eigenvalue is negative" # check that T is positive definite
print(T)

Order 4 ANISO tensor | Param(size=21)=[ 2.66011 1.26432 0.662772 1.9402 1.54905 1.10384 3.6072 1.78964 2.0247 1.2701 1.1089 2.743 1.3367 1.2962 0.897632 4.42684 2.05632 1.52686 3.54431 1.3445 1.99356 ] | Angles(size=0)=[ ]
[ 2.66011 1.26432 0.662772 1.9402 1.54905 1.10384 
  1.26432 3.6072 1.78964 2.0247 1.2701 1.1089 
  0.662772 1.78964 2.743 1.3367 1.2962 0.897632 
  1.9402 2.0247 1.3367 4.42684 2.05632 1.52686 
  1.54905 1.2701 1.2962 2.05632 3.54431 1.3445 
  1.10384 1.1089 0.897632 1.52686 1.3445 1.99356 ]

Construction of a fourth-order tensor from a \(6×6\) matrix and projections

Given a symmetric positive definite \(6×6\) matrix M (warning: these properties are not checked on construction), a tensor can be built by one of the function

• T = tensor(M)
• T = tensor(M, angles=[θ, ϕ])

By default, this function finds out whether M represents an ISO, TI or ORTHO tensor in the canonical frame if angles is not precised or in the frame defined by angles. If none of these symmetries is suitable, the tensor is considered ANISO.

angles = [π/3, π/4, π/5]
T = tensor([9.,3.,5.,6.,3.,7.,1.,4.,3.], angles)
M = T.a
print("M =\n", M) # the information about ORTHO symmetry has vanished in M
T2 = tensor(M, angles=angles)
print(T2)

M =
 [[ 5.834782  3.574311  3.799884  0.041627  0.710013 -0.509607]
 [ 3.574311  9.744881  2.267425 -1.242374 -0.080104  0.934783]
 [ 3.799884  2.267425  9.137098 -0.976365  1.311174  0.522827]
 [ 0.041627 -1.242374 -0.976365  3.25126  -0.192475 -0.544059]
 [ 0.710013 -0.080104  1.311174 -0.192475  5.36512  -1.536706]
 [-0.509607  0.934783  0.522827 -0.544059 -1.536706  4.666858]]
Order 4 ORTHO tensor | Param(size=9)=[ 9 3 5 6 3 7 1 4 3 ] | Angles(size=3)=[ 1.0472 0.785398 0.628319 ]
[ 5.83478 3.57431 3.79988 0.0416272 0.710013 -0.509607 
  3.57431 9.74488 2.26743 -1.24237 -0.0801039 0.934783 
  3.79988 2.26743 9.1371 -0.976365 1.31117 0.522827 
  0.0416272 -1.24237 -0.976365 3.25126 -0.192475 -0.544059 
  0.710013 -0.0801039 1.31117 -0.192475 5.36512 -1.53671 
  -0.509607 0.934783 0.522827 -0.544059 -1.53671 4.66686 ]

Given a \(6×6\) matrix supposed to represent a tensor in the canonical frame, it is possible to find the closest tensor (in the sense of least-squares relatively to the euclidean matrix norm) among a family of given material symmetry. The syntax is simply T = tensor(M, SYM) where SYM denotes the chosen class of symmetry (ISO, TI, ORTHO).

For example the closest ISO tensor is calculated by

T = tensor(M, ISO)

In this case it boils down to the “isotropisation” of the tensor \(\uuuu{T}\) given by (Bornert et al., 2001)

\[ \ISO(\uuuu{T})=\left(\uuuu{T}::\uuuu{J}\right)\,\uuuu{J}+\left(\frac{\uuuu{T}::\uuuu{K}}{5}\right)\,\uuuu{K} \]

M = np.array([ [2.66011, 1.26432, 0.662772, 1.9402, 1.54905, 1.10384],
               [1.26432, 3.6072, 1.78964, 2.0247, 1.2701, 1.1089], 
               [0.662772, 1.78964, 2.743, 1.3367, 1.2962, 0.897632], 
               [1.9402, 2.0247 ,1.3367, 4.42684, 2.05632, 1.52686], 
               [1.54905, 1.2701, 1.2962, 2.05632, 3.54431, 1.3445], 
               [1.10384, 1.1089, 0.897632, 1.52686, 1.3445, 1.99356] ])

T = tensor(M,ISO)
print("T =\n", T)

cont4=lambda T1,T2:np.einsum('ij,ij',T1,T2) # T1 and T2 in KM notation so 2 indices
T2 = tensor(cont4(M,J4), cont4(M,K4)/5)
print("T2 =\n", T2)

T =
 Order 4 ISO tensor | Param(size=2)=[ 5.48126 2.69875 ] | Angles(size=0)=[ ]
[ 3.62625 0.927502 0.927502 0 0 0 
  0.927502 3.62625 0.927502 0 0 0 
  0.927502 0.927502 3.62625 0 0 0 
  0 0 0 2.69875 0 0 
  0 0 0 0 2.69875 0 
  0 0 0 0 0 2.69875 ]

T2 =
 Order 4 ISO tensor | Param(size=2)=[ 5.48126 2.69875 ] | Angles(size=0)=[ ]
[ 3.62625 0.927502 0.927502 0 0 0 
  0.927502 3.62625 0.927502 0 0 0 
  0.927502 0.927502 3.62625 0 0 0 
  0 0 0 2.69875 0 0 
  0 0 0 0 2.69875 0 
  0 0 0 0 0 2.69875 ]

Check a well-known result stating that the inverse of the closest isotropic tensor is in general not the closest isotropic tensor of the inverse.

M = np.array([ [2.66011, 1.26432, 0.662772, 1.9402, 1.54905, 1.10384],
               [1.26432, 3.6072, 1.78964, 2.0247, 1.2701, 1.1089], 
               [0.662772, 1.78964, 2.743, 1.3367, 1.2962, 0.897632], 
               [1.9402, 2.0247 ,1.3367, 4.42684, 2.05632, 1.52686], 
               [1.54905, 1.2701, 1.2962, 2.05632, 3.54431, 1.3445], 
               [1.10384, 1.1089, 0.897632, 1.52686, 1.3445, 1.99356] ])
T = tensor(M, ISO)
print("inv(ISO(M)) =\n", inv(tensor(M, ISO)))
print("ISO(inv(M)) =\n", tensor(np.linalg.inv(M), ISO))

inv(ISO(M)) =
 Order 4 ISO tensor | Param(size=2)=[ 0.18244 0.370542 ] | Angles(size=0)=[ ]
[ 0.307841 -0.0627006 -0.0627006 0 0 0 
  -0.0627006 0.307841 -0.0627006 0 0 0 
  -0.0627006 -0.0627006 0.307841 0 0 0 
  0 0 0 0.370542 0 0 
  0 0 0 0 0.370542 0 
  0 0 0 0 0 0.370542 ]

ISO(inv(M)) =
 Order 4 ISO tensor | Param(size=2)=[ 0.416399 0.615133 ] | Angles(size=0)=[ ]
[ 0.548889 -0.0662448 -0.0662448 0 0 0 
  -0.0662448 0.548889 -0.0662448 0 0 0 
  -0.0662448 -0.0662448 0.548889 0 0 0 
  0 0 0 0.615133 0 0 
  0 0 0 0 0.615133 0 
  0 0 0 0 0 0.615133 ]

The closest TI tensor in the frame defined by given angles is calculated by

T = tensor(M, TI, angles=[θ, ϕ], epsrel=1.e-3)

where epsrel is a relative stop criterion for the iterative least-square algorithm.

print(tensor(M, TI, angles=[π/3,π/4]))

Order 4 TI tensor | Param(size=5)=[ 9.90471 1.93007 0.947144 1.68894 1.88118 ] | Angles(size=2)=[ 1.0472 0.785398 ]
[ 3.29543 1.46231 1.08371 1.31151 1.39476 1.70822 
  1.46231 3.29543 1.08371 1.39476 1.31151 1.70822 
  1.08371 1.08371 2.59348 1.03685 1.03685 1.16792 
  1.31151 1.39476 1.03685 3.04891 1.31191 1.57733 
  1.39476 1.31151 1.03685 1.31191 3.04891 1.57733 
  1.70822 1.70822 1.16792 1.57733 1.57733 3.69285 ]

One can also optimize with respect to the angles so as to find the best frame for a given symmetry

T = tensor(M, TI, epsrel=1.e-3, optiangles=True)

where epsrel is a relative stop criterion for the iterative least-square algorithm.

print(tensor(M, TI, epsrel=1.e-3, optiangles=True))

Order 4 TI tensor | Param(size=5)=[ 10.2696 1.81638 0.89854 1.70466 1.73988 ] | Angles(size=2)=[ 0.92234 0.891541 ]
[ 2.52626 1.12576 1.07821 1.27354 1.04336 1.07127 
  1.12576 3.30972 1.51116 1.80824 1.44496 1.49908 
  1.07821 1.51116 3.17755 1.73199 1.39846 1.42041 
  1.27354 1.80824 1.73199 3.77321 1.65983 1.70366 
  1.04336 1.44496 1.39846 1.65983 3.05771 1.38017 
  1.07127 1.49908 1.42041 1.70366 1.38017 3.13058 ]

Exercise

1. Find the closest ORTHO tensor T from the previous matrix M (also optimizing the angles)

2. Build a new tensor from the matrix of T (in the canonical frame)

3. Build a new tensor from the matrix of T and the angles of T

4. Find the tensor of best symmetry from the matrix of T with angle optimization

M = np.array([ [2.66011, 1.26432, 0.662772, 1.9402, 1.54905, 1.10384],
               [1.26432, 3.6072, 1.78964, 2.0247, 1.2701, 1.1089], 
               [0.662772, 1.78964, 2.743, 1.3367, 1.2962, 0.897632], 
               [1.9402, 2.0247 ,1.3367, 4.42684, 2.05632, 1.52686], 
               [1.54905, 1.2701, 1.2962, 2.05632, 3.54431, 1.3445], 
               [1.10384, 1.1089, 0.897632, 1.52686, 1.3445, 1.99356] ])
T = tensor(M, ORTHO, optiangles=True) ; print(T)
print(tensor(T.a))
print(tensor(T.a, angles = T.angles))
print(tensor(T.a, optiangles = True))

Order 4 ORTHO tensor | Param(size=9)=[ 10.2695 -0.375452 1.64478 2.06277 -0.00326333 1.57933 0.792497 1.03801 0.701183 ] | Angles(size=3)=[ 1.522 2.49696 0.924649 ]
[ 2.80357 1.20599 1.14979 1.63413 1.50234 1.16022 
  1.20599 3.38542 1.441 2.18362 1.51678 1.17623 
  1.14979 1.441 2.6612 1.44697 1.1142 0.974349 
  1.63413 2.18362 1.44697 4.52485 2.11525 1.42722 
  1.50234 1.51678 1.1142 2.11525 3.37899 1.16562 
  1.16022 1.17623 0.974349 1.42722 1.16562 2.22099 ]

Order 4 ANISO tensor | Param(size=21)=[ 2.80357 1.20599 1.14979 1.63413 1.50234 1.16022 3.38542 1.441 2.18362 1.51678 1.17623 2.6612 1.44697 1.1142 0.974349 4.52485 2.11525 1.42722 3.37899 1.16562 2.22099 ] | Angles(size=0)=[ ]
[ 2.80357 1.20599 1.14979 1.63413 1.50234 1.16022 
  1.20599 3.38542 1.441 2.18362 1.51678 1.17623 
  1.14979 1.441 2.6612 1.44697 1.1142 0.974349 
  1.63413 2.18362 1.44697 4.52485 2.11525 1.42722 
  1.50234 1.51678 1.1142 2.11525 3.37899 1.16562 
  1.16022 1.17623 0.974349 1.42722 1.16562 2.22099 ]

Order 4 ORTHO tensor | Param(size=9)=[ 10.2695 -0.375452 1.64478 2.06277 -0.00326333 1.57933 0.792497 1.03801 0.701183 ] | Angles(size=3)=[ 1.522 2.49696 0.924649 ]
[ 2.80357 1.20599 1.14979 1.63413 1.50234 1.16022 
  1.20599 3.38542 1.441 2.18362 1.51678 1.17623 
  1.14979 1.441 2.6612 1.44697 1.1142 0.974349 
  1.63413 2.18362 1.44697 4.52485 2.11525 1.42722 
  1.50234 1.51678 1.1142 2.11525 3.37899 1.16562 
  1.16022 1.17623 0.974349 1.42722 1.16562 2.22099 ]

Order 4 ORTHO tensor | Param(size=9)=[ 10.2695 -0.375452 1.64478 2.06277 -0.00326333 1.57933 0.792497 1.03801 0.701183 ] | Angles(size=3)=[ 1.522 2.49696 0.924649 ]
[ 2.80357 1.20599 1.14979 1.63413 1.50234 1.16022 
  1.20599 3.38542 1.441 2.18362 1.51678 1.17623 
  1.14979 1.441 2.6612 1.44697 1.1142 0.974349 
  1.63413 2.18362 1.44697 4.52485 2.11525 1.42722 
  1.50234 1.51678 1.1142 2.11525 3.37899 1.16562 
  1.16022 1.17623 0.974349 1.42722 1.16562 2.22099 ]

\(\,\)

Bornert, M., Bretheau, T., Gilormini, P., 2001. Homogénéisation en mécanique des matériaux. Hermes science.