Plane stress¶
A material is said to be under plane stress if the stress vector is zero across a particular surface. The assumption of plane stress simplifies the analysis considerably, and is an approximation that is typically made when analysing thin plates and shells.
Composites are often found as relatively thin plates and shells, which suggest that the plane stress assumption can be sufficiently accurate, at least on the larger scale. However, detailed 3-dimensional state of stress analysis must often be supplemented at joints, transitions and boundaries of the shell structure.
Hooke's law for an orthotropic material is
\begin{equation} \begin{bmatrix} \varepsilon_1 \\ \varepsilon_2 \\ \varepsilon_3 \\ \gamma_{23} \\ \gamma_{13} \\ \gamma_{12} \end{bmatrix} = \begin{bmatrix} S_{11} & S_{12} & S_{13} & 0 & 0 & 0 \\ S_{12} & S_{22} & S_{23} & 0 & 0 & 0 \\ S_{13} & S_{23} & S_{33} & 0 & 0 & 0 \\ 0 & 0 & 0 & S_{44} & 0 & 0 \\ 0 & 0 & 0 & 0 & S_{55} & 0 \\ 0 & 0 & 0 & 0 & 0 & S_{66} \end{bmatrix} \begin{bmatrix} \sigma_1 \\ \sigma_2 \\ \sigma_3 \\ \tau_{23} \\ \tau_{13} \\ \tau_{12} \end{bmatrix} \tag{1} \end{equation}Assume now that the stress vector is zero across surfaces with normals alligned with the 3-axis. This assumption implies that the non-zero stresses exist only in the 1-2 plane, i.e.,
\begin{equation} \sigma_3 = \tau_{23} = \tau_{13} = 0 \tag{2} \end{equation}The set of equations (1) can therfore be reduced to
\begin{equation} \begin{bmatrix} \varepsilon_1 \\ \varepsilon_2 \\ \gamma_{12} \end{bmatrix} = \begin{bmatrix} S_{11} & S_{12} & 0 \\ S_{12} & S_{22} & 0 \\ 0 & 0 & S_{66} \end{bmatrix} \begin{bmatrix} \sigma_1 \\ \sigma_2 \\ \tau_{12} \end{bmatrix} \tag{3} \end{equation}where the compliance matrix expressed by the enginering constants is
\begin{equation} \mathbf{S}=\begin{bmatrix} 1/E_1 & -\nu_{12}/E_1 & 0 \\ -\nu_{12}/E_1 & 1/E_{2} & 0 \\ 0 & 0 & 1/G_{12} \end{bmatrix} \tag{4} \end{equation}Implementation of a function for the 2D compliance matrix:
import numpy as np
def S2D(m):
return np.array([[ 1/m['E1'], -m['v12']/m['E1'], 0],
[-m['v12']/m['E1'], 1/m['E2'], 0],
[ 0, 0, 1/m['G12']]], float)
The 2D stiffness matrix of the material is simply the inverse of the compliance matrix:
\begin{equation} \mathbf{Q} = \mathbf{S}^{-1} \tag{5} \end{equation}In the following implementation we call the previous function for the compliance matrix S with a material (m) as argument. Then, the inverse of the compliance matrix S is returned:
def Q2D(m):
S=S2D(m)
return np.linalg.inv(S)
Example:
import matlib
m1=matlib.get('Carbon/Epoxy(a)')
Q=Q2D(m1)
print('Plane stress stiffness matrix:')
print()
print(np.round(Q))
Observe that the components of Q is not equal to the corresponding components of C
from compositelib import C3D
C=C3D(m1)
print(np.round(C))
print(Q[0,0]/C[0,0])
print(Q[0,1]/C[0,1])
print(Q[1,1]/C[1,1])
print(Q[2,2]/C[5,5])
The 2D transformation matrices are obtained by reducing the 3D transformation matrices for roation about the 3-axis:
\begin{equation} \mathbf{T}_{\sigma}=\begin{bmatrix} c^2 & s^2 & 2cs \\ s^2 & c^2 & -2cs \\ -cs & cs & c^2-s^2 \end{bmatrix} \tag{6} \end{equation}and \begin{equation} \mathbf{T}_{\epsilon}=\begin{bmatrix} c^2 & s^2 & cs \\ s^2 & c^2 & -cs \\ -2cs & 2cs & c^2-s^2 \end{bmatrix} \tag{7} \end{equation}
Both functions created below takes the orientation angle rot as the only argument:
from math import cos,sin,radians
def T2Ds(rot):
c,s = cos(radians(rot)), sin(radians(rot))
return np.array([[ c*c , s*s , 2*c*s],
[ s*s , c*c , -2*c*s],
[-c*s, c*s , c*c-s*s]], float)
def T2De(rot):
c,s = cos(radians(rot)), sin(radians(rot))
return np.array([[ c*c, s*s, c*s ],
[ s*s, c*c, -c*s ],
[-2*c*s, 2*c*s, c*c-s*s ]], float)
Transformation of the 2D stiffness matrix is given by: \begin{equation} \mathbf{Q{'}} = \mathbf{T}_{\sigma}^{-1}\mathbf{Q} \mathbf{T}_{\epsilon} \tag{8} \end{equation}
The transformed stiffnes matrix has the components
\begin{equation} \mathbf{Q}'= \begin{bmatrix} Q_{xx} & Q_{xy} & Q_{xs} \\ Q_{xy} & Q_{yy} & Q_{ys} \\ Q_{xs} & Q_{ys} & Q_{ss} \end{bmatrix} \tag{9} \end{equation}The function Q2Dtransform which takes the material stiffness matrix Q and the orientation angle rot as arguments:
def Q2Dtransform(Q,rot):
return np.dot(np.linalg.inv( T2Ds(rot) ) , np.dot(Q,T2De(rot)) )
Example:
Q=Q2D(m1)
print(np.round(Q))
print()
Qt=Q2Dtransform(Q,60)
print(np.round(Qt))
Hooke's law expressed in the global coordinate system is written as
\begin{equation} \begin{bmatrix} \sigma_x \\ \sigma_y \\ \tau_{xy} \end{bmatrix}= \begin{bmatrix} Q_{xx} & Q_{xy} & Q_{xs} \\ Q_{xy} & Q_{yy} & Q_{ys} \\ Q_{xs} & Q_{ys} & Q_{ss} \end{bmatrix} \begin{bmatrix} \varepsilon_x \\ \varepsilon_y \\ \gamma_{xy} \end{bmatrix} \tag{10} \end{equation}Plane stress thermomechanical relations
In the material coordinate system of an orthotropic material:
\begin{equation} \boldsymbol{\sigma} = \mathbf{Q}(\boldsymbol{\varepsilon} - \boldsymbol{\alpha}\Delta T ) \tag{11} \end{equation}\begin{equation} \begin{bmatrix} \sigma_1 \\ \sigma_2 \\ \tau_{12} \end{bmatrix} = \begin{bmatrix} Q_{11} & Q_{12} & 0 \\ Q_{12} & Q_{22} & 0 \\ 0 & 0 & Q_{66} \end{bmatrix} \begin{bmatrix} \begin{bmatrix} \varepsilon_1 \\ \varepsilon_2 \\ \gamma_{12} \end{bmatrix} -\Delta T \begin{bmatrix} \alpha_1 \\ \alpha_2 \\ 0 \end{bmatrix} \end{bmatrix} \tag{12} \end{equation}In a transformed coordinate system:
\begin{equation} \boldsymbol{\sigma}' = \mathbf{Q}'(\boldsymbol{\varepsilon}' - \boldsymbol{\alpha}'\Delta T ) \tag{13} \end{equation}\begin{equation} \begin{bmatrix} \sigma_x \\ \sigma_y \\ \tau_{xy} \end{bmatrix} = \begin{bmatrix} Q_{xx} & Q_{xy} & Q_{xs} \\ Q_{xy} & Q_{yy} & Q_{ys} \\ Q_{xs} & Q_{ys} & Q_{ss} \end{bmatrix} \begin{bmatrix} \begin{bmatrix} \varepsilon_x \\ \varepsilon_y \\ \gamma_{xy} \end{bmatrix} -\Delta T \begin{bmatrix} \alpha_x \\ \alpha_y \\ \alpha_{xy} \end{bmatrix} \end{bmatrix} \tag{14} \end{equation}where
\begin{equation} \boldsymbol{\alpha}'= \mathbf{T}_{\varepsilon}^{-1} \boldsymbol{\alpha} \tag{15} \end{equation}
