Demo: Principal Stress Directions in 2D
Colab note: This notebook is designed to run on Google Colab. The first code cell installs dependencies.
This notebook helps visualize the principal stress directions and principal stresses for a 2D stress tensor. It follows the worked example in which you compute the eigenvalues and eigenvectors of a 2D stress tensor from Shearer’s Chapter 2.
Goals¶
By working through this notebook, you should be able to:
see how a 2D stress tensor acts on planes of different orientation
identify the principal stress directions
understand why the principal directions are the eigenvectors
see that the principal stresses are the corresponding eigenvalues
verify that shear traction vanishes on principal planes
We use the 2D stress tensor
throughout.
1. Imports¶
Source
# Install dependencies (for Google Colab or missing packages)
import sys
# Check if running in Colab
try:
import google.colab
IN_COLAB = True
print("Running in Google Colab")
except:
IN_COLAB = False
print("Running in local environment")
# Install required packages if needed
required_packages = {
'numpy': 'numpy',
'matplotlib': 'matplotlib',
'ipywidgets': 'ipywidgets'
}
missing_packages = []
for package, pip_name in required_packages.items():
try:
__import__(package)
print(f"✓ {package} is already installed")
except ImportError:
missing_packages.append(pip_name)
print(f"✗ {package} not found")
if missing_packages:
print(f"\nInstalling missing packages: {', '.join(missing_packages)}")
import subprocess
subprocess.check_call([sys.executable, "-m", "pip", "install", "-q"] + missing_packages)
print("✓ Installation complete!")
else:
print("\n✓ All required packages are installed!")Running in local environment
✓ numpy is already installed
✓ matplotlib is already installed
✓ ipywidgets is already installed
✓ All required packages are installed!
Source
import numpy as np
import matplotlib.pyplot as plt
from ipywidgets import interact, FloatSlider, Checkbox
plt.style.use("default")
plt.rcParams.update({
"figure.facecolor": "white",
"axes.facecolor": "white",
"savefig.facecolor": "white",
"font.size": 12,
"lines.linewidth": 2,
})2. Stress Tensor Utilities¶
def stress_tensor(txx, txy, tyy):
return np.array([[txx, txy],
[txy, tyy]], dtype=float)
def principal_stresses_and_directions(T):
vals, vecs = np.linalg.eigh(T)
order = np.argsort(vals)[::-1]
vals = vals[order]
vecs = vecs[:, order]
return vals, vecs
def traction_on_plane(T, theta_deg):
theta = np.deg2rad(theta_deg)
n = np.array([np.cos(theta), np.sin(theta)])
f = np.array([-np.sin(theta), np.cos(theta)])
t = T @ n
t_normal = np.dot(t, n)
t_shear = np.dot(t, f)
return n, f, t, t_normal, t_shear
def angle_of_vector(v):
return np.rad2deg(np.arctan2(v[1], v[0]))3. Worked Example¶
This first cell reproduces the kind of calculation you might do by hand for a tensor such as
You can replace these numbers with your own example.
txx, txy, tyy = -40.0, -10.0, -60.0
T = stress_tensor(txx, txy, tyy)
vals, vecs = principal_stresses_and_directions(T)
print("Stress tensor T:")
print(T)
print("\nPrincipal stresses (eigenvalues):")
print(vals)
print("\nPrincipal directions (eigenvectors):")
print(vecs)
print("\nAngles of principal directions measured from +x (degrees):")
for i in range(2):
print(f"Direction {i+1}: {angle_of_vector(vecs[:, i]):.2f}°")Stress tensor T:
[[-40. -10.]
[-10. -60.]]
Principal stresses (eigenvalues):
[-35.85786438 -64.14213562]
Principal directions (eigenvectors):
[[-0.92387953 0.38268343]
[ 0.38268343 0.92387953]]
Angles of principal directions measured from +x (degrees):
Direction 1: 157.50°
Direction 2: 67.50°
4. Visualizing the Tensor and Its Principal Directions¶
The plot below shows:
the x and y axes
the two principal stress directions
an optional plane normal ( \hat{\mathbf{n}} )
the corresponding traction vector ( \mathbf{t} = \boldsymbol{\tau}\hat{\mathbf{n}} )
On a principal plane, the traction vector is parallel to the normal, so the shear traction is zero.
def plot_stress_state(txx=-40.0, txy=-10.0, tyy=-60.0, plane_angle=45.0,
show_plane=True, show_traction=True):
T = stress_tensor(txx, txy, tyy)
vals, vecs = principal_stresses_and_directions(T)
n, f, t, t_normal, t_shear = traction_on_plane(T, plane_angle)
fig = plt.figure(figsize=(7, 7))
ax = fig.add_subplot(111)
ax.axhline(0)
ax.axvline(0)
colors = ["C1", "C2"]
for i in range(2):
v = vecs[:, i]
ax.arrow(0, 0, v[0], v[1], head_width=0.05, length_includes_head=True, color=colors[i])
ax.arrow(0, 0, -v[0], -v[1], head_width=0.05, length_includes_head=True, color=colors[i])
ax.text(1.12*v[0], 1.12*v[1], rf"$\sigma_{i+1}={vals[i]:.2f}$", color=colors[i], ha="center")
if show_plane:
ax.arrow(0, 0, n[0], n[1], head_width=0.05, length_includes_head=True, color="black")
ax.arrow(0, 0, f[0], f[1], head_width=0.05, length_includes_head=True, color="gray")
ax.text(1.1*n[0], 1.1*n[1], r"$\hat{\mathbf{n}}$", ha="center")
ax.text(1.1*f[0], 1.1*f[1], r"$\hat{\mathbf{f}}$", ha="center")
if show_traction:
scale = max(np.linalg.norm(T), 1.0)
t_plot = t / scale
ax.arrow(0, 0, t_plot[0], t_plot[1], head_width=0.05, length_includes_head=True, color="C3")
ax.text(1.12*t_plot[0], 1.12*t_plot[1], r"$\mathbf{t}$", color="C3", ha="center")
ax.set_aspect("equal")
ax.set_xlim(-1.4, 1.4)
ax.set_ylim(-1.4, 1.4)
ax.set_title("2D Stress Tensor: Principal Directions and Traction")
ax.set_xlabel("x")
ax.set_ylabel("y")
ax.grid(True, alpha=0.3)
plt.show()
print("Stress tensor:")
print(T)
print("\nPrincipal stresses:")
print(vals)
print("\nPrincipal direction angles from +x:")
for i in range(2):
print(f" direction {i+1}: {angle_of_vector(vecs[:, i]):.2f}°")
print(f"\nSelected plane normal angle: {plane_angle:.2f}°")
print(f"Normal traction: {t_normal:.3f}")
print(f"Shear traction: {t_shear:.3f}")interact(
plot_stress_state,
txx=FloatSlider(value=-40.0, min=-100.0, max=100.0, step=1.0, description="τxx"),
txy=FloatSlider(value=-10.0, min=-100.0, max=100.0, step=1.0, description="τxy"),
tyy=FloatSlider(value=-60.0, min=-100.0, max=100.0, step=1.0, description="τyy"),
plane_angle=FloatSlider(value=45.0, min=0.0, max=180.0, step=1.0, description="θ (deg)"),
show_plane=Checkbox(value=True, description="show plane"),
show_traction=Checkbox(value=True, description="show traction"),
);5. Shear Traction as a Function of Plane Orientation¶
A key result is that the shear traction becomes zero on the principal planes.
The next plot shows how normal traction and shear traction vary with plane orientation.
def plot_traction_vs_angle(txx=-40.0, txy=-10.0, tyy=-60.0):
T = stress_tensor(txx, txy, tyy)
vals, vecs = principal_stresses_and_directions(T)
angles = np.linspace(0, 180, 721)
normal = []
shear = []
for ang in angles:
_, _, _, tn, ts = traction_on_plane(T, ang)
normal.append(tn)
shear.append(ts)
normal = np.array(normal)
shear = np.array(shear)
principal_angles = [angle_of_vector(vecs[:, i]) % 180 for i in range(2)]
fig = plt.figure(figsize=(9, 5))
ax = fig.add_subplot(111)
ax.plot(angles, normal, label="normal traction")
ax.plot(angles, shear, label="shear traction")
for pa in principal_angles:
ax.axvline(pa, linestyle="--", alpha=0.8)
ax.set_xlabel("Plane normal angle (degrees)")
ax.set_ylabel("Traction (same units as stress)")
ax.set_title("Normal and Shear Traction vs Plane Orientation")
ax.grid(True, alpha=0.3)
ax.legend()
plt.show()
print("Principal direction angles (mod 180°):")
for i, pa in enumerate(principal_angles, start=1):
print(f" direction {i}: {pa:.2f}°")
print("\nAt these angles, the shear traction should be zero.")interact(
plot_traction_vs_angle,
txx=FloatSlider(value=-40.0, min=-100.0, max=100.0, step=1.0, description="τxx"),
txy=FloatSlider(value=-10.0, min=-100.0, max=100.0, step=1.0, description="τxy"),
tyy=FloatSlider(value=-60.0, min=-100.0, max=100.0, step=1.0, description="τyy"),
);6. Rotating Into the Principal Coordinate System¶
If we build a matrix ( \mathbf{N} ) from the eigenvectors, we can rotate the tensor into the principal coordinate system:
[ \boldsymbol{\tau}_R = \mathbf{N}^T \boldsymbol{\tau} \mathbf{N} ]
In that coordinate system, the stress tensor is diagonal.
def show_rotation_to_principal(txx=-40.0, txy=-10.0, tyy=-60.0):
T = stress_tensor(txx, txy, tyy)
vals, vecs = principal_stresses_and_directions(T)
N = vecs
T_rot = N.T @ T @ N
print("Original tensor T:")
print(T)
print("\nMatrix of principal directions N:")
print(N)
print("\nRotated tensor N^T T N:")
print(T_rot)
print("\nNotice that the off-diagonal terms are approximately zero.")interact(
show_rotation_to_principal,
txx=FloatSlider(value=-40.0, min=-100.0, max=100.0, step=1.0, description="τxx"),
txy=FloatSlider(value=-10.0, min=-100.0, max=100.0, step=1.0, description="τxy"),
tyy=FloatSlider(value=-60.0, min=-100.0, max=100.0, step=1.0, description="τyy"),
);7. Things to Try¶
Change and see how the principal directions rotate.
Set . What happens to the principal directions?
Try with nonzero . What do the principal directions look like?
Try a hydrostatic case: and . What happens to the traction?
Find the plane angles at which the shear traction is zero and compare them to the eigenvector directions.
Key Takeaway¶
The principal stress directions are the directions in which traction is purely normal.
Mathematically, these are the eigenvectors of the stress tensor, and the corresponding principal stresses are the eigenvalues.
8. Visual Connection to Mohr’s Circle¶
Some of you may already have seen Mohr’s circle in other classes.
Mohr’s circle provides another way to visualize the same 2D stress state:
the center gives the mean normal stress
the radius gives the maximum shear stress
the points where the circle crosses the horizontal axis are the principal stresses
The plot below is meant as a quick visual connection rather than a full derivation.
def mohr_circle_2d(txx=-40.0, txy=-10.0, tyy=-60.0, theta_deg=45.0):
T = stress_tensor(txx, txy, tyy)
vals, vecs = principal_stresses_and_directions(T)
# Circle parameters
center = 0.5 * (txx + tyy)
radius = np.sqrt(((txx - tyy) / 2.0)**2 + txy**2)
# Stress on plane with normal angle theta
_, _, _, tn, ts = traction_on_plane(T, theta_deg)
# Parametric circle
ang = np.linspace(0, 2*np.pi, 400)
sigma = center + radius * np.cos(ang)
tau = radius * np.sin(ang)
fig = plt.figure(figsize=(7, 6))
ax = fig.add_subplot(111)
ax.plot(sigma, tau, label="Mohr's circle")
ax.axhline(0)
ax.axvline(0)
# Original coordinate planes
ax.plot([txx], [-txy], 'o', label=r'$(\tau_{xx},-\tau_{xy})$')
ax.plot([tyy], [ txy], 'o', label=r'$(\tau_{yy},\tau_{xy})$')
# Principal stresses
ax.plot([vals[0], vals[1]], [0, 0], 'o', label="principal stresses")
# Selected plane
ax.plot([tn], [ts], 'o', label=rf"plane at $\theta={theta_deg:.0f}^\circ$")
ax.set_xlabel("Normal stress")
ax.set_ylabel("Shear stress")
ax.set_title("Mohr's Circle for a 2D Stress Tensor")
ax.set_aspect("equal", adjustable="box")
ax.grid(True, alpha=0.3)
ax.legend(loc="best")
plt.show()
print(f"Mean normal stress (circle center): {center:.3f}")
print(f"Maximum shear stress (circle radius): {radius:.3f}")
print(f"Principal stresses: {vals[0]:.3f}, {vals[1]:.3f}")
print(f"Stress on selected plane: normal = {tn:.3f}, shear = {ts:.3f}")interact(
mohr_circle_2d,
txx=FloatSlider(value=-40.0, min=-100.0, max=100.0, step=1.0, description="τxx"),
txy=FloatSlider(value=-10.0, min=-100.0, max=100.0, step=1.0, description="τxy"),
tyy=FloatSlider(value=-60.0, min=-100.0, max=100.0, step=1.0, description="τyy"),
theta_deg=FloatSlider(value=45.0, min=0.0, max=180.0, step=1.0, description="θ (deg)"),
);Quick takeaway¶
Mohr’s circle and the eigenvector picture are two views of the same stress tensor.
In the physical x-y plane, the principal directions are the eigenvectors.
In Mohr’s circle space, the principal stresses are where the circle crosses the horizontal axis.
The circle radius is the maximum possible shear stress.
If you have seen Mohr’s circle before, this section is meant to help connect that picture to the tensor and eigenvector viewpoint.