Elasticity Tensor

The setup

In one dimension, Hooke's law reads σ = E ε — stress equals stiffness times strain. In three dimensions both stress and strain are tensors, and the linear relationship between them needs a tensor of higher rank to keep the indices balanced. That tensor is Cᵢⱼₖₗ, the elasticity tensor or stiffness tensor.

σᵢⱼ = Σ_{k,l} Cᵢⱼₖₗ · εₖₗ   (Einstein summation: σᵢⱼ = Cᵢⱼₖₗ εₖₗ)

Stress σᵢⱼ has units of pressure (Pa) and describes the i-component of force per unit area on a face whose outward normal is the j-direction. Strain εₖₗ is dimensionless and describes the symmetric part of the displacement gradient. C therefore has units of pressure: Pa, equivalently N/m².

From 81 to 21 to 2

Three universal symmetries reduce the count:

Reduction	Origin	Components left
None	Raw 4th-rank tensor	81
σᵢⱼ = σⱼᵢ	Angular-momentum balance	54
εₖₗ = εₗₖ	Definition of strain	36 (Voigt 6×6)
Cᵢⱼₖₗ = Cₖₗᵢⱼ	Existence of strain energy U(ε)	21 (triclinic max)

Crystal symmetries cut further:

Symmetry class	Independent constants	Example
Triclinic	21	Plagioclase feldspar
Monoclinic	13	Gypsum, augite
Orthorhombic	9	Olivine, fibre composites, wood
Tetragonal	6	Zircon, white tin
Trigonal	6 or 7	Quartz, calcite
Hexagonal	5	Graphite, ice Ih, magnesium
Cubic	3	Silicon, salt, copper, diamond
Isotropic	2	Glass, polycrystalline metal

Voigt 6×6 matrix form

Engineers prefer the 6×6 contracted form. Map index pairs (1,1)→1, (2,2)→2, (3,3)→3, (2,3)→4, (1,3)→5, (1,2)→6:

⎡σ₁⎤   ⎡C₁₁ C₁₂ C₁₃ C₁₄ C₁₅ C₁₆⎤ ⎡ε₁⎤
⎢σ₂⎥   ⎢C₁₂ C₂₂ C₂₃ C₂₄ C₂₅ C₂₆⎥ ⎢ε₂⎥
⎢σ₃⎥ = ⎢C₁₃ C₂₃ C₃₃ C₃₄ C₃₅ C₃₆⎥ ⎢ε₃⎥
⎢σ₄⎥   ⎢C₁₄ C₂₄ C₃₄ C₄₄ C₄₅ C₄₆⎥ ⎢2ε₂₃⎥
⎢σ₅⎥   ⎢C₁₅ C₂₅ C₃₅ C₄₅ C₅₅ C₅₆⎥ ⎢2ε₁₃⎥
⎣σ₆⎦   ⎣C₁₆ C₂₆ C₃₆ C₄₆ C₅₆ C₆₆⎦ ⎣2ε₁₂⎦

The factor of 2 on the strain shears is Voigt's quirk; alternatives (Mandel notation) use √2 instead to preserve the tensor norm. The matrix is always symmetric because of the strain-energy constraint.

Isotropic form — two numbers say it all

For an isotropic solid the only constants that survive are the Lamé parameters λ and μ:

Cᵢⱼₖₗ = λ δᵢⱼ δₖₗ + μ (δᵢₖ δⱼₗ + δᵢₗ δⱼₖ)
σᵢⱼ = λ tr(ε) δᵢⱼ + 2μ εᵢⱼ

Engineering forms:

Pair	Formula in (λ, μ)	Comment
Young's modulus E	μ (3λ + 2μ) / (λ + μ)	Tensile stiffness
Poisson's ratio ν	λ / (2(λ + μ))	Transverse contraction
Bulk modulus K	λ + 2μ / 3	Volumetric stiffness
Shear modulus G	μ	Pure shear stiffness
C₁₁	λ + 2μ	Diagonal Voigt
C₁₂	λ	Off-diagonal Voigt
C₄₄	μ	Shear Voigt

Worked example — silicon vs steel vs wood

Material	Symmetry	Key constants (GPa)	E (GPa)	ν
Silicon (Si)	Cubic	C₁₁=166, C₁₂=64, C₄₄=80	130 [100] · 169 [111]	0.27
Steel (mild)	Isotropic	λ≈115, μ≈80	200	0.30
Aluminium	Isotropic	λ≈56, μ≈26	69	0.33
Glass (fused silica)	Isotropic	λ≈16, μ≈31	73	0.17
Spruce wood (along grain)	Orthotropic	9 constants, E_L≈11, E_T≈0.7	11 along · 0.7 across	0.30 along · 0.43 across
Graphite (basal plane)	Hexagonal	C₁₁≈1060, C₃₃≈37	1020 in-plane · 35 out-of-plane	0.16

Silicon's stiffness along the body diagonal [111] is about 30 % higher than along the cube axis [100] — that anisotropy is exploited by MEMS designers who align resonators along the stiffest direction. Graphite is 30 times stiffer in the basal plane than perpendicular to it.

Variants

• Compliance tensor S = C⁻¹. Same physics, inverted view: ε = Sσ. Engineering moduli (E, ν, G) appear as 1/S₁₁, −S₁₂/S₁₁, 1/S₄₄ in isotropic form. Geophysicists usually report C; metallurgists usually report S.
• Anisotropic damping. Add an imaginary part Cᵢⱼₖₗ → C′ᵢⱼₖₗ + iC″ᵢⱼₖₗ for viscoelastic and seismic-attenuation modelling. The strain-energy symmetry still applies to both parts.
• Nonlinear elasticity. For large strains, σ becomes a nonlinear function of ε, and you Taylor-expand: σᵢⱼ = Cᵢⱼₖₗ εₖₗ + ½ Cᵢⱼₖₗₘₙ εₖₗ εₘₙ + … Third-order constants (six independent ones in isotropic case) describe acoustoelasticity and finite-strain ultrasound — used to detect stress in steel pipelines.
• Cosserat / micropolar elasticity. Materials with internal rotational degrees of freedom (foams, lattice metamaterials, bone) need an extended tensor that couples couple-stresses to micro-rotations — 18 independent constants in 3D, but rarely measured.
• Auxetic materials. Some carefully designed lattices have negative Poisson's ratio (ν < 0), so they get fatter when stretched. The full tensor remains positive definite — it's compatible with the same Hooke's-law framework, just with unusual off-diagonal entries.
• Effective tensors for composites. A composite's homogenized elasticity tensor C* depends on the constituents' moduli and microstructure (volume fraction, alignment, shape). Standard bounds: Reuss (lower) and Voigt (upper), with Hashin-Shtrikman bounds tighter for isotropic mixtures.

JavaScript — Voigt matrix utilities

// Build the isotropic 6×6 Voigt stiffness matrix from (E, nu)
function isotropicVoigt(E, nu) {
  const mu = E / (2 * (1 + nu));
  const lam = E * nu / ((1 + nu) * (1 - 2 * nu));
  const a = lam + 2 * mu, b = lam;
  return [
    [a, b, b, 0, 0, 0],
    [b, a, b, 0, 0, 0],
    [b, b, a, 0, 0, 0],
    [0, 0, 0, mu, 0, 0],
    [0, 0, 0, 0, mu, 0],
    [0, 0, 0, 0, 0, mu],
  ];
}

// Cubic Voigt matrix from 3 constants (C₁₁, C₁₂, C₄₄)
function cubicVoigt(C11, C12, C44) {
  return [
    [C11, C12, C12, 0,   0,   0],
    [C12, C11, C12, 0,   0,   0],
    [C12, C12, C11, 0,   0,   0],
    [0,   0,   0,   C44, 0,   0],
    [0,   0,   0,   0,   C44, 0],
    [0,   0,   0,   0,   0,   C44],
  ];
}

// Apply σ = C ε for Voigt 6-vectors
function applyVoigt(C, eps) {
  const sig = new Array(6).fill(0);
  for (let i = 0; i < 6; i++)
    for (let j = 0; j < 6; j++) sig[i] += C[i][j] * eps[j];
  return sig;
}

// Silicon at room temperature (GPa)
const Si = cubicVoigt(166, 64, 80);

// Uniaxial strain ε₁ = 0.001, rest zero
let eps = [0.001, 0, 0, 0, 0, 0];
console.log('Silicon stress for 0.1% uniaxial strain:', applyVoigt(Si, eps));
// → σ₁ = 0.166 GPa, σ₂ = σ₃ = 0.064 GPa

// Pure shear ε₆ = 0.001 (= 2 ε₁₂)
eps = [0, 0, 0, 0, 0, 0.001];
console.log('Silicon stress for 0.1% shear strain:', applyVoigt(Si, eps));
// → σ₆ = 0.080 GPa

Where the elasticity tensor matters

• Finite-element structural analysis. Every commercial FEA package stores C either implicitly (isotropic E, ν) or explicitly (orthotropic, anisotropic) and applies it billions of times per simulation to compute stress from strain.
• Seismic anisotropy. Earth's upper mantle is anisotropic — olivine crystals align under flow. P-wave and S-wave speeds vary with direction; inverting them maps mantle convection.
• Semiconductor MEMS. Silicon resonators are cut along low-coupling crystal directions to maximise Q. Designers consult the full Voigt matrix of Si.
• Carbon-fibre composites. An aerospace ply has very different stiffness along vs across fibres; designers stack plies at different angles to tune the effective elasticity tensor for the load case.
• Bone biomechanics. Cortical bone is approximately orthotropic; trabecular bone is more anisotropic still. Hip replacements and dental implants are designed with patient-specific elasticity tensors derived from CT scans.
• Game-engine soft bodies. Cloth and finite-element solvers in Unreal, Unity, Houdini, and Marvelous Designer parameterise stretch and shear stiffness — effectively four numbers of an orthotropic Voigt matrix per cloth.
• Ultrasonic non-destructive testing. Wave speeds in welds and rolled metals depend on the local C; deviations from baseline reveal residual stress or texture changes.

Common mistakes

• Confusing stiffness C with compliance S. C maps strain to stress; S = C⁻¹ maps stress to strain. Engineering moduli use S; finite-element solvers usually use C.
• Mixing Voigt and Mandel notation. Voigt has factors of 2 on shears, Mandel has √2. Mixing them silently scales shear moduli — a classic finite-element bug.
• Using isotropic E, ν for an anisotropic material. Composites, single crystals, and rolled metals all need at least an orthotropic tensor — using a single E underestimates or overestimates the off-axis response by factors of 2 to 30.
• Forgetting positive-definiteness. A valid Voigt matrix must have all six eigenvalues positive (strain energy ½ εᵀ C ε > 0 for all non-zero ε). Negative eigenvalues mean the material is unstable — useful as a sanity check.
• Using engineering shear strain instead of tensor shear strain. γᵢⱼ = 2εᵢⱼ. Voigt uses engineering shear in the strain vector but tensor shear in the components of ε. Sign of factor 2 is the difference between right and wrong.
• Assuming Cᵢⱼₖₗ is symmetric in i,j or k,l separately enough. The full Voigt symmetry depends on three pairwise symmetries — losing any one of them (in damage mechanics, for instance, where stress-strain symmetry breaks) drops you back to 81 components.

Performance — FEA cost scaling

A finite-element solver multiplies a 6×6 stiffness matrix by a 6-vector at every Gauss point of every element — typically 8 to 27 Gauss points per hexahedral element. A modern aerospace structural simulation with 10⁷ elements does ~10⁹ such multiplications per timestep, costing roughly 36 floating-point operations each (216 multiplications + additions). On a 1 TFLOPS workstation, that is about 40 ms per stiffness pass — comparable to the matrix solve itself for explicit dynamics. Caching the per-material Voigt matrix in L1 is essential; recomputing it from anisotropic E, ν tables on the fly costs 2–3× overall.