Rob Siegwart - Principal Stresses in 3D

When a stress element is characterized by its six stress components and these components are put in stress tensor form, the principal normal stresses are simply the eigenvalues of the stress tensor, and the principal planes are defined by the eigenvectors. This is easily computed using a mathematical program such as MathCAD, Matlab, etc. Here we will use Python and the library NumPy.

Principal Normal Stress

Put the stresses into the stress tensor form:

\begin{equation*} S=\begin{bmatrix}\sigma_x & \tau_{xy} & \tau_{zx} \\ \tau_{xy} & \sigma_y & \tau_{yz} \\ \tau_{zx} & \tau_{yz} & \sigma_z \end{bmatrix} \end{equation*}

import numpy as np
S = np.array([ [S_xx, S_xy, S_zx],
               [S_xy, S_yy, S_yz],
               [S_zx, S_yz, S_zz] ])

Then, compute the eigenvalues and eigenvectors. In Python:

e_val, e_vec = np.linalg.eig(S)

To store the eigenvalues in the variable e_val and the eigenvectors in the variable \(e_vec\).

The principal stresses can be assigned by sorting, based on the assignment relation of \(\sigma_1 > \sigma_2 > \sigma_3\).

p3, p2, p1 = np.sort(e_val)   # sort smallest to largest

The eigenvectors are stored as column vectors in the e_vec array. To access these, we first need to associate the principal stresses with their positions in the original result in order to assign them the correct vector. To do this we will convert the original eigenvalue array to a regular list and use standard Python indexing. Then, the vectors are extracted from the eigenvector array using these indices.

e_val_l = e_val.tolist()
p1_index, p2_index, p3_index = e_val_l.index(p1), e_val_l.index(p2), e_val_l.index(p3)
p1_vec, p2_vec, p3_vec = e_vec[:,p1_index], e_vec[:,p2_index], e_vec[:,p3_index]

Principal Shear Stress

The principal shearing stresses can be inferred from the usage of a Mohr circle. In constructing a Mohr circle for a 3D stress state, the maximum shearing stress can be shown to be one half of the difference between the maximum and minimum principal normal stresses. From the Mohr circle the principal shearing stresses can be determined as:

\begin{align*} \tau_1=(\sigma_1-\sigma_3)/2 \\ \tau_2=(\sigma_1-\sigma_2)/2 \\ \tau_3=(\sigma_2-\sigma_3)/2 \end{align*}

In Python,

tau1 = (p1-p3)/2
tau2 = (p1-p2)/2
tau3 = (p2-p3)/2

The principal shearing planes bisect the two principal normal planes used in determining the principal shearing stress. The \(\tau_1\) plane bisects the \(\sigma_1\) and \(\sigma_3\) planes, and since all principal normal planes are orthogonal to each other, the principal shearing planes are at 45 degrees to each principal normal plane.

To find the vectors for these planes, we compute the vector that bisects the vectors for the normal principal stresses. Add them vectorally to obtain a rhombus, the diagonal of which splits the rhombus into two equal segments.

In Python we add our normal principal vectors and divide by the magnitude to normalize:

tau1_vec = (p1_vec + p3_vec)/np.linalg.norm(p1_vec + p3_vec)
tau2_vec = (p1_vec + p2_vec)/np.linalg.norm(p1_vec + p2_vec)
tau3_vec = (p2_vec + p3_vec)/np.linalg.norm(p2_vec + p3_vec)