TMM4175 Polymer Composites

Home About Python Links Table of Contents

CASE STUDY: Principal stresses and transformation¶

Problem statement: Find the pricipal stresses and directions for the plane stress state

\begin{equation} \boldsymbol{\sigma}= \begin{bmatrix} \sigma_{x} & \tau_{xy} & \tau_{xz} \\ \tau_{xy} & \sigma_{y} & \tau_{yz} \\ \tau_{xz} & \tau_{yz} & \sigma_{z} \end{bmatrix}= \begin{bmatrix} 200 & -75 & 0 \\ -75 & 50 & 0 \\ 0 & 0 & 0 \end{bmatrix} \end{equation}

by using the relations for stress transformation.

Solution¶

Write the stress state as a 1-dimensional array:

import numpy as np

sxyz=np.array([200,50,0,0,0,-75])

The transformation matrix for stress by rotation about the z-axis as implemented in Transformations:

from compositelib import T3Dsz

Then, a straight forward sweep through rotation angles from -90 to 90 in order to find where the shear stress is zero, which also implies that the normal stresses show maximum or minimum:

s_xyz=(200,50,0,0,0,-75)
a=np.linspace(-90,90)
s11,s22,s12=[],[],[]
for angle in a:
    s_123=np.dot(T3Dsz(angle),s_xyz)
    s11.append(s_123[0])
    s22.append(s_123[1])
    s12.append(s_123[5])
import matplotlib.pyplot as plt
%matplotlib inline
fig,ax=plt.subplots(figsize=(14,5))
ax.plot(a,s11,color='red', label='$\sigma_1$')
ax.plot(a,s22,color='blue', label='$\sigma_2$')
ax.plot(a,s12,color='green', label=r'$\tau_{12}$')
ax.set_xlim(-90,90)
ax.legend(loc='best')
ax.grid(True)
plt.show()

According to the result, the shear stress is zero for rotation angles of -22.5 and 67.5 (a couple of iterations were required to find those numbers):

T=T3Dsz(-22.5)
s=np.dot(T,sxyz)
print(np.round(s,1))

[231.1  18.9   0.    0.    0.    0. ]

T=T3Dsz(67.5)
s=np.dot(T,sxyz)
print(np.round(s,1))

[ 18.9 231.1   0.    0.    0.   -0. ]

Note that the difference between the angles implies a orthogonal system:

67.5-(-22.5)

90.0

Hence, the first principal stress is 231.1, the second is 18.9 and the third is zero. The direction of the first principal stress is inclined by -22.5 degrees in the x-y plane relatively to the x-axis.

Verify using numpy.linalg.eig() (see Stress)

ST=np.array ( [ [200.0, -75.0,   0.0 ],
              [  -75.0,  50.0,   0.0 ],
              [    0.0,   0.0,   0.0 ]  ])

L,v = np.linalg.eig(ST)

print('Eigenvalues:')
print(L)
print()
print('Eigenvectors:')
print(v)

Eigenvalues:
[231.06601718  18.93398282   0.        ]

Eigenvectors:
[[ 0.92387953  0.38268343  0.        ]
 [-0.38268343  0.92387953  0.        ]
 [ 0.          0.          1.        ]]

Remember that the eigenvectors are given by the columns:

v1 = v[:,0]
v2 = v[:,1]
v3 = v[:,2]

import matplotlib.pyplot as plt
%matplotlib inline
fig,ax = plt.subplots(figsize=(4,4))
ax.set_axis_off()
ax.set_xlim(-1.1,1.1)
ax.set_ylim(-1.1,1.1)

ax.arrow(0, 0, 1, 0.0, head_width=0.1, head_length=0.1, fc='black', ec='black')
ax.arrow(0, 0, 0.0, 1, head_width=0.1, head_length=0.1, fc='black', ec='black')
ax.text(0.9, 0.1, 'x')
ax.text(0.1, 0.9, 'y')

ax.arrow(0, 0, v1[0], v1[1], head_width=0.1, head_length=0.1, fc='red', ec='red')
ax.arrow(0, 0, v2[0], v2[1], head_width=0.1, head_length=0.1, fc='blue', ec='blue')

plt.show()

Disclaimer:This site is about polymer composites, designed for educational purposes. Consumption and use of any sort & kind is solely at your own risk.
Fair use: I spent some time making all the pages, and even the figures and illustrations are my own creations. Obviously, you may steal whatever you find useful here, but please show decency and give some acknowledgment if or when copying. Thanks! Contact me: nils.p.vedvik@ntnu.no www.ntnu.edu/employees/nils.p.vedvik