Linear Elasticity Equations

Summary of the Equations of Linear Elasticity

The stiffness matrix of linear elastic isotropic material contains two parameters:

E, Young's modulus (elastic modulus)
ν, Poisson’s ratio

Define the following quantities.

$\begin{array}{l} σ = stress \\ f = body force \\ ε = strain \\ u = displacement \end{array}$

The equilibrium equation is

$- \nabla \cdot σ = f$

The linearized, small-displacement strain-displacement relationship is

$ε = \frac{1}{2} (\nabla u + \nabla u^{T})$

The balance of angular momentum states that stress is symmetric:

$σ_{i j} = σ_{j i}$

The Voigt notation for the constitutive equation of the linear isotropic model is

$[\begin{matrix} σ_{11} \\ σ_{22} \\ σ_{33} \\ σ_{23} \\ σ_{13} \\ σ_{12} \end{matrix}] = \frac{E}{(1 + ν) (1 - 2 ν)} [\begin{matrix} 1 - ν & ν & ν & 0 & 0 & 0 \\ ν & 1 - ν & ν & 0 & 0 & 0 \\ ν & ν & 1 - ν & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 - 2 ν & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 - 2 ν & 0 \\ 0 & 0 & 0 & 0 & 0 & 1 - 2 ν \end{matrix}] [\begin{matrix} ε_{11} \\ ε_{22} \\ ε_{33} \\ ε_{23} \\ ε_{13} \\ ε_{12} \end{matrix}]$

The expanded form uses all the entries in σ and ε takes symmetry into account.

[\begin{matrix} σ_{11} \\ σ_{12} \\ σ_{13} \\ σ_{21} \\ σ_{22} \\ σ_{23} \\ σ_{31} \\ σ_{32} \\ σ_{33} \end{matrix}] = \frac{E}{(1 + ν) (1 - 2 ν)} [\begin{matrix} 1 - ν & 0 & 0 & 0 & ν & 0 & 0 & 0 & ν \\ • & 1 - 2 ν & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ • & • & 1 - 2 ν & 0 & 0 & 0 & 0 & 0 & 0 \\ • & • & • & 1 - 2 ν & 0 & 0 & 0 & 0 & 0 \\ • & • & • & • & 1 - ν & 0 & 0 & 0 & ν \\ • & • & • & • & • & 1 - 2 ν & 0 & 0 & 0 \\ • & • & • & • & • & • & 1 - 2 ν & 0 & 0 \\ • & • & • & • & • & • & • & 1 - 2 ν & 0 \\ • & • & • & • & • & • & • & • & 1 - ν \end{matrix}] [\begin{matrix} ε_{11} \\ ε_{12} \\ ε_{13} \\ ε_{21} \\ ε_{22} \\ ε_{23} \\ ε_{31} \\ ε_{32} \\ ε_{33} \end{matrix}]

(1)

In the preceding diagram, • means the entry is symmetric.

3D Linear Elasticity Problem

The toolbox form for the equation is

$- \nabla \cdot (c \otimes \nabla u) = f$

But the equations in the summary do not have ∇u alone, it appears together with its transpose:

$ε = \frac{1}{2} (\nabla u + \nabla u^{T})$

It is a straightforward exercise to convert this equation for strain ε to ∇u. In column vector form,

$\nabla u = [\begin{matrix} \partial u_{x} / \partial x \\ \partial u_{x} / \partial y \\ \partial u_{x} / \partial z \\ \partial u_{y} / \partial x \\ \partial u_{y} / \partial y \\ \partial u_{y} / \partial z \\ \partial u_{z} / \partial x \\ \partial u_{z} / \partial y \\ \partial u_{z} / \partial z \end{matrix}]$

Therefore, you can write the strain-displacement equation as

$ε = [\begin{matrix} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & \frac{1}{2} & 0 & \frac{1}{2} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & \frac{1}{2} & 0 & 0 & 0 & \frac{1}{2} & 0 & 0 \\ 0 & \frac{1}{2} & 0 & \frac{1}{2} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{2} & 0 & \frac{1}{2} & 0 \\ 0 & 0 & \frac{1}{2} & 0 & 0 & 0 & \frac{1}{2} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{2} & 0 & \frac{1}{2} & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \end{matrix}] \nabla u \equiv A \nabla u$

where A stands for the displayed matrix. So rewriting Equation 1, and recalling that • means an entry is symmetric, you can write the stiffness tensor as

$\begin{matrix} σ = \frac{E}{(1 + ν) (1 - 2 ν)} [\begin{matrix} 1 - ν & 0 & 0 & 0 & ν & 0 & 0 & 0 & ν \\ • & 1 - 2 ν & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ • & • & 1 - 2 ν & 0 & 0 & 0 & 0 & 0 & 0 \\ • & • & • & 1 - 2 ν & 0 & 0 & 0 & 0 & 0 \\ • & • & • & • & 1 - ν & 0 & 0 & 0 & ν \\ • & • & • & • & • & 1 - 2 ν & 0 & 0 & 0 \\ • & • & • & • & • & • & 1 - 2 ν & 0 & 0 \\ • & • & • & • & • & • & • & 1 - 2 ν & 0 \\ • & • & • & • & • & • & • & • & 1 - ν \end{matrix}] A \nabla u \\ = \frac{E}{(1 + ν) (1 - 2 ν)} [\begin{matrix} 1 - ν & 0 & 0 & 0 & ν & 0 & 0 & 0 & ν \\ 0 & 1 / 2 - ν & 0 & 1 / 2 - ν & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 / 2 - ν & 0 & 0 & 0 & 1 / 2 - ν & 0 & 0 \\ 0 & 1 / 2 - ν & 0 & 1 / 2 - ν & 0 & 0 & 0 & 0 & 0 \\ ν & 0 & 0 & 0 & 1 - ν & 0 & 0 & 0 & ν \\ 0 & 0 & 0 & 0 & 0 & 1 / 2 - ν & 0 & 1 / 2 - ν & 0 \\ 0 & 0 & 1 / 2 - ν & 0 & 0 & 0 & 1 / 2 - ν & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 1 / 2 - ν & 0 & 1 / 2 - ν & 0 \\ ν & 0 & 0 & 0 & ν & 0 & 0 & 0 & 1 - ν \end{matrix}] \nabla u \end{matrix}$

Make the definitions

$\begin{array}{l} μ = \frac{E}{2 (1 + ν)} \\ λ = \frac{E ν}{(1 + ν) (1 - 2 ν)} \\ \frac{E (1 - ν)}{(1 + ν) (1 - 2 ν)} = 2 μ + λ \end{array}$

and the equation becomes

$σ = [\begin{matrix} 2 μ + λ & 0 & 0 & 0 & λ & 0 & 0 & 0 & λ \\ 0 & μ & 0 & μ & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & μ & 0 & 0 & 0 & μ & 0 & 0 \\ 0 & μ & 0 & μ & 0 & 0 & 0 & 0 & 0 \\ λ & 0 & 0 & 0 & 2 μ + λ & 0 & 0 & 0 & λ \\ 0 & 0 & 0 & 0 & 0 & μ & 0 & μ & 0 \\ 0 & 0 & μ & 0 & 0 & 0 & μ & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & μ & 0 & μ & 0 \\ λ & 0 & 0 & 0 & λ & 0 & 0 & 0 & 2 μ + λ \end{matrix}] \nabla u \equiv c \nabla u$

If you are solving a 3-D linear elasticity problem by using PDEModel instead of StructuralModel, use the elasticityC3D(E,nu) function (included in your software) to obtain the c coefficient. This function uses the linearized, small-displacement assumption for an isotropic material. For examples that use this function, see StationaryResults.

Plane Stress

Plane stress is a condition that prevails in a flat plate in the x-y plane, loaded only in its own plane and without z-direction restraint. For plane stress, σ₁₃ = σ₂₃ = σ₃₁ = σ₃₂ = σ₃₃ = 0. Assuming isotropic conditions, the Hooke's law for plane stress gives the following strain-stress relation:

$[\begin{matrix} ε_{11} \\ ε_{22} \\ 2 ε_{12} \end{matrix}] = \frac{1}{E} [\begin{matrix} 1 & - ν & 0 \\ - ν & 1 & 0 \\ 0 & 0 & 2 + 2 ν \end{matrix}] [\begin{matrix} σ_{11} \\ σ_{22} \\ σ_{12} \end{matrix}]$

Inverting this equation, obtain the stress-strain relation:

$(\begin{matrix} σ_{11} \\ σ_{22} \\ σ_{12} \end{matrix}) = \frac{E}{1 - ν^{2}} (\begin{matrix} 1 & ν & 0 \\ ν & 1 & 0 \\ 0 & 0 & \frac{1 - ν}{2} \end{matrix}) (\begin{matrix} ε_{11} \\ ε_{22} \\ 2 ε_{12} \end{matrix})$

Convert the equation for strain ε to ∇u.

$ε = [\begin{matrix} 1 & 0 & 0 & 0 \\ 0 & \frac{1}{2} & \frac{1}{2} & 0 \\ 0 & \frac{1}{2} & \frac{1}{2} & 0 \\ 0 & 0 & 0 & 1 \end{matrix}] \nabla u \equiv A \nabla u$

Now you can rewrite the stiffness matrix as

$[\begin{matrix} σ_{11} \\ σ_{12} \\ σ_{21} \\ σ_{22} \end{matrix}] = [\begin{matrix} \frac{E}{1 - ν^{2}} & 0 & 0 & \frac{E ν}{1 - ν^{2}} \\ 0 & \frac{E}{2 (1 + ν)} & \frac{E}{2 (1 + ν)} & 0 \\ 0 & \frac{E}{2 (1 + ν)} & \frac{E}{2 (1 + ν)} & 0 \\ \frac{E ν}{1 - ν^{2}} & 0 & 0 & \frac{E}{1 - ν^{2}} \end{matrix}] \nabla u = [\begin{matrix} \frac{2 μ (μ + λ)}{2 μ + λ} & 0 & 0 & \frac{2 λ μ}{2 μ + λ} \\ 0 & μ & μ & 0 \\ 0 & μ & μ & 0 \\ \frac{2 λ μ}{2 μ + λ} & 0 & 0 & \frac{2 μ (μ + λ)}{2 μ + λ} \end{matrix}] \nabla u$

Plane Strain

Plane strain is a deformation state where there are no displacements in the z-direction, and the displacements in the x- and y-directions are functions of x and y but not z. The stress-strain relation is only slightly different from the plane stress case, and the same set of material parameters is used.

For plane strain, ε₁₃ = ε₂₃ = ε₃₁ = ε₃₂ = ε₃₃ = 0. Assuming isotropic conditions, the stress-strain relation can be written as follows:

$(\begin{matrix} σ_{11} \\ σ_{22} \\ σ_{12} \end{matrix}) = \frac{E}{(1 + ν) (1 - 2 ν)} (\begin{matrix} 1 - ν & ν & 0 \\ ν & 1 - ν & 0 \\ 0 & 0 & \frac{1 - 2 ν}{2} \end{matrix}) (\begin{matrix} ε_{11} \\ ε_{22} \\ 2 ε_{12} \end{matrix})$

Convert the equation for strain ε to ∇u.

$ε = [\begin{matrix} 1 & 0 & 0 & 0 \\ 0 & \frac{1}{2} & \frac{1}{2} & 0 \\ 0 & \frac{1}{2} & \frac{1}{2} & 0 \\ 0 & 0 & 0 & 1 \end{matrix}] \nabla u \equiv A \nabla u$

Now you can rewrite the stiffness matrix as

$[\begin{matrix} σ_{11} \\ σ_{12} \\ σ_{21} \\ σ_{22} \end{matrix}] = [\begin{matrix} \frac{E (1 - ν)}{(1 + ν) (1 - 2 ν)} & 0 & 0 & \frac{E ν}{(1 + ν) (1 - 2 ν)} \\ 0 & \frac{E}{2 (1 + ν)} & \frac{E}{2 (1 + ν)} & 0 \\ 0 & \frac{E}{2 (1 + ν)} & \frac{E}{2 (1 + ν)} & 0 \\ \frac{E ν}{(1 + ν) (1 - 2 ν)} & 0 & 0 & \frac{E (1 - ν)}{(1 + ν) (1 - 2 ν)} \end{matrix}] \nabla u = [\begin{matrix} 2 μ + λ & 0 & 0 & λ \\ 0 & μ & μ & 0 \\ 0 & μ & μ & 0 \\ λ & 0 & 0 & 2 μ + λ \end{matrix}] \nabla u$

Axisymmetric Analysis

Axisymmetric analysis speeds up simulations by simplifying 3-D solids using their symmetry around the axis of rotation and analyzing only the 2-D axisymmetric section. Use polar coordinates r,θ,z for radial, circumferential, and axial directions. If z is the axis of rotation, the symmetry around the z-axis means that the stress components are independent of the θ coordinate. The stress equilibrium equations for an axisymmetric structural mechanics are as follows:

$\begin{array}{l} \frac{\partial σ_{r r}}{\partial r} + \frac{\partial σ_{r z}}{\partial z} + \frac{1}{r} (σ_{r r} - σ_{θ θ}) + f_{r} = 0 \\ \frac{\partial τ_{r z}}{\partial r} + \frac{\partial σ_{z z}}{\partial z} + \frac{τ_{r z}}{r} + f_{z} = 0 \end{array}$

τ is the shear stress, and γ is the shear strain. Assuming isotropic conditions, the stress-strain relation can be written as follows:

$(\begin{matrix} σ_{r} \\ σ_{z} \\ \begin{array}{l} σ_{θ} \\ τ_{r z} \end{array} \end{matrix}) = \frac{E}{(1 + ν) (1 - 2 ν)} (\begin{matrix} 1 - ν & ν & 0 & 0 \\ ν & 1 - ν & 0 & 0 \\ 0 & 0 & 1 - ν & 0 \\ 0 & 0 & 0 & \frac{1 - 2 ν}{2} \end{matrix}) (\begin{matrix} ε_{r} \\ ε_{z} \\ \begin{array}{l} ε_{θ} \\ γ_{r z} \end{array} \end{matrix})$