TMM4175 Polymer Composites

Home About Python Links Next Previous Table of Contents

Layer results

Obtaining the strains and stresses in the individual layers

Once the laminate load-deformation relations are solved, the deformation is known.

Then, the strains in the laminate coordinate system at any position $z$ through the thickness of the laminate can be computed by:

\begin{equation} \begin{bmatrix} \varepsilon_x \\ \varepsilon_y \\ \gamma_{xy} \end{bmatrix}= \begin{bmatrix} \varepsilon_x^0 \\ \varepsilon_y^0 \\ \gamma_{xy}^0 \end{bmatrix}+ z\begin{bmatrix} \kappa_x \\ \kappa_y \\ \kappa_{xy} \end{bmatrix} \end{equation}

Example:

In [1]:
from laminatelib import laminateStiffnessMatrix, solveLaminateLoadCase, laminateThickness

m1={'E1':140000, 'E2':10000, 'v12':0.3, 'G12':5000}

layup1  =[ {'mat':m1 , 'ori':  0  , 'thi':1} , 
           {'mat':m1 , 'ori': 90  , 'thi':1} ,
           {'mat':m1 , 'ori':  0  , 'thi':1}  ]

ABD1 = laminateStiffnessMatrix(layup1)

result = solveLaminateLoadCase(ABD1,Nx=500,Mx=100)
d=result[1]

print(d[0:3])  # mid-plane strains
print(d[3:6])  # curvatures

[ 1.71604988e-03 -9.65278056e-05  2.96942846e-20]
[ 3.28128174e-04 -6.64459553e-05  3.41450199e-21]

The strains at the bottom of the laminate ( $z=-1.5$ ) and the top of the laminate ( $z=1.5$ ) are:

In [2]:
strn_xyz_botLam = d[0:3] -1.5*d[3:6]
strn_xyz_topLam = d[0:3] +1.5*d[3:6]
print(strn_xyz_botLam)
print(strn_xyz_topLam)

[1.22385762e-03 3.14112736e-06 2.45725316e-20]
[ 2.20824214e-03 -1.96196739e-04  3.48160376e-20]

The most relevant positions for asessing strains and stresses are the ones at bottom and top of the individual layers:

In [3]:
bot = -laminateThickness(layup1)/2    # The bottom coordinate of the laminate
counter=1
for layer in layup1:
    top = bot + layer['thi']
    strnXYZbot=d[0:3]+bot*d[3:6]
    strnXYZtop=d[0:3]+top*d[3:6]
    
    print('Layer',str(counter),'bot','z=',bot,strnXYZbot)
    print('Layer',str(counter),'top','z=',top,strnXYZtop)
    counter=counter+1
    bot=top  # for the next layer, top becomes the bot

Layer 1 bot z= -1.5 [1.22385762e-03 3.14112736e-06 2.45725316e-20]
Layer 1 top z= -0.5 [ 1.55198579e-03 -6.33048280e-05  2.79870336e-20]
Layer 2 bot z= -0.5 [ 1.55198579e-03 -6.33048280e-05  2.79870336e-20]
Layer 2 top z= 0.5 [ 1.88011396e-03 -1.29750783e-04  3.14015356e-20]
Layer 3 bot z= 0.5 [ 1.88011396e-03 -1.29750783e-04  3.14015356e-20]
Layer 3 top z= 1.5 [ 2.20824214e-03 -1.96196739e-04  3.48160376e-20]

The stresses at bottom and top of all layers can be obtained by a few more lines including the transformed stiffness matrix and Hooke's law: \begin{equation} \mathbf{Q{'}}_k = \mathbf{T}_{\sigma}^{-1}\mathbf{Q}_k \mathbf{T}_{\epsilon} \end{equation}

\begin{equation} \boldsymbol{\sigma}{'} = \mathbf{Q{'}}_k\boldsymbol{\epsilon}{'} \end{equation}
In [4]:
from laminatelib import Q2D, Q2Dtransform
import numpy as np

bot = -laminateThickness(layup1)/2    # The bottom coordinate of the laminate
counter=1
for layer in layup1:
    top = bot + layer['thi']
    strnXYZbot=d[0:3]+bot*d[3:6]
    strnXYZtop=d[0:3]+top*d[3:6]
    
    Q=Q2D(layer['mat'])         # the stiffness matrix of current layer
    Qt=Q2Dtransform(Q,layer['ori']) # the transformed stiffness matrix of current layer
    strsXYZbot = np.dot( Qt, strnXYZbot)   # stresses by Hooke's law
    strsXYZtop = np.dot( Qt, strnXYZtop)   # stresses by Hooke's law
    
    print('Layer',str(counter),'bot','z=',bot,strsXYZbot)
    print('Layer',str(counter),'top','z=',top,strsXYZtop)
    counter=counter+1
    bot=top  # for the next layer, top becomes the bot

Layer 1 bot z= -1.5 [1.72458149e+02 3.72694304e+00 1.22862658e-16]
Layer 1 top z= -0.5 [2.18492692e+02 4.04893798e+00 1.39935168e-16]
Layer 2 bot z= -0.5 [ 1.54291307e+01 -4.23393671e+00 -7.50002167e-17]
Layer 2 top z= 0.5 [ 1.85310153e+01 -1.26058051e+01 -5.18885475e-16]
Layer 3 bot z= 0.5 [2.64527235e+02 4.37093292e+00 1.57007678e-16]
Layer 3 top z= 1.5 [3.10561778e+02 4.69292785e+00 1.74080188e-16]

In order to find the stresses in the material coordinate system we must add a transformation:

\begin{equation} \boldsymbol{\sigma} = \mathbf{T}_{\sigma}\boldsymbol{\sigma}{'} \end{equation}
In [5]:
from laminatelib import T2Ds, T2De

bot = -laminateThickness(layup1)/2    # The bottom coordinate of the laminate
counter=1
for layer in layup1:
    top = bot + layer['thi']
    strnXYZbot=d[0:3]+bot*d[3:6]
    strnXYZtop=d[0:3]+top*d[3:6]
    
    Q=Q2D(layer['mat'])         # the stiffness matrix of current layer
    Qt=Q2Dtransform(Q,layer['ori']) # the transformed stiffness matrix of current layer
    strsXYZbot = np.dot( Qt, strnXYZbot)   # stresses by Hooke's law
    strsXYZtop = np.dot( Qt, strnXYZtop)   # stresses by Hooke's law
    
    T=T2Ds(layer['ori'])  # transformation matrix for stresses
    strs123bot = np.dot( T, strsXYZbot)   # strains in material coordinate system
    strs123top = np.dot( T, strsXYZtop)   # strains in materila coordinate system
    
    print('Layer',str(counter),'bot','z=',bot,strs123bot)
    print('Layer',str(counter),'top','z=',top,strs123top)
    counter=counter+1
    bot=top  # for the next layer, top becomes the bot

Layer 1 bot z= -1.5 [1.72458149e+02 3.72694304e+00 1.22862658e-16]
Layer 1 top z= -0.5 [2.18492692e+02 4.04893798e+00 1.39935168e-16]
Layer 2 bot z= -0.5 [-4.23393671e+00  1.54291307e+01 -1.12901541e-15]
Layer 2 top z= 0.5 [-1.26058051e+01  1.85310153e+01 -1.38769489e-15]
Layer 3 bot z= 0.5 [2.64527235e+02 4.37093292e+00 1.57007678e-16]
Layer 3 top z= 1.5 [3.10561778e+02 4.69292785e+00 1.74080188e-16]

Obtaining exposure factors of individual layers

When the stresses in the material coordinate system are found, the exposure factors can be computed using the previously defined plane stress failure criteria functions

In [6]:
from laminatelib import fE2DMS, fE2DME, fE2DTW

A structure for the layer results

The derived results for one layer includes stresses and strains in two coordinate systems, as well as stress exposure factors and other relevant data such as the z-coordinates. All of these results should generally be computed at two locations; at the bottom and at the top of the layer. A suggestion for a useful datastructure is illustrated:

                                    Layer k
                                       |
           ------------------------------------------------------------------
           |                   |                      |                     |
        stress              strain                   fail                   h
     -----------         -----------         -------------------          -----
     |         |         |         |         |        |        |          |   |    
    xyz       123       xyz       123        MS       ME       TW        bot top
   -----     -----     -----     -----      -----    -----    -----
   |   |     |   |     |   |     |   |      |   |    |   |    |   |
  bot top   bot top   bot top   bot top    bot top  bot top  bot top

The data structure can easily be created using dictionaries and then be organized into a list containing the results from all layers:

In [7]:
def layerResults(layup,deformations):
    # An empty list that shall be populated with results from all layers:
    results=[]
    # The bottom coordinate of the laminate:
    bot = -laminateThickness(layup)/2
    for layer in layup:
        # the top of current layer is just the bottom + the thickness:
        top=bot+layer['thi']
        # creating a dictionary with the two coordinates:
        h={'bot':bot, 'top':top}
        # computing the strains in laminate coordinate system:
        strnXYZbot=deformations[0:3]+bot*deformations[3:6]
        strnXYZtop=deformations[0:3]+top*deformations[3:6]
        # put the strains into a dictionary:
        strnXYZ={'bot':strnXYZbot, 'top':strnXYZtop}
        # the transformation matrix for strains:
        Te=T2De(layer['ori'])
        # transforming the strains into the material coordinate system:
        strn123bot=np.dot(Te, strnXYZbot)
        strn123top=np.dot(Te, strnXYZtop)
        # the strains at top and bottom as a dictionary:
        strn123={'bot':strn123bot, 'top':strn123top}
        # all strains as a new dictionary:
        strn={'xyz':strnXYZ, '123':strn123}
        # stiffness matrix of the material:
        Q=Q2D(layer['mat'])
        # Hooke's law: finding the stresses in material system:
        strs123bot=np.dot(Q,strn123bot)
        strs123top=np.dot(Q,strn123top)
        # the material stresses as a dictionary having bottom and topp results:
        strs123={'bot':strs123bot, 'top':strs123top}
        # transformation matrix for stress, negativ rotation applies now,
        # since this is from material system to laminat system (reversed..)
        Ts=T2Ds(-layer['ori'])
        # transforming stresses into laminate system:
        strsXYZbot=np.dot(Ts,strs123bot)
        strsXYZtop=np.dot(Ts,strs123top)
        # organizing the stresses for the two location in a dictionary
        strsXYZ={'bot':strsXYZbot, 'top':strsXYZtop}
        # all stresses as a new nested dictionary:
        strs={'xyz':strsXYZ, '123':strs123}
        # Failure criteria...follows the same logic as above
        MSbot=fE2DMS(strs123bot,layer['mat'])
        MStop=fE2DMS(strs123top,layer['mat'])
        MS={'bot':MSbot, 'top':MStop}
        MEbot=fE2DME(strs123bot,layer['mat'])
        MEtop=fE2DME(strs123top,layer['mat'])
        ME={'bot':MEbot, 'top':MEtop}
        TWbot=fE2DTW(strs123bot,layer['mat'])
        TWtop=fE2DTW(strs123top,layer['mat'])
        TW={'bot':TWbot, 'top':TWtop}
        fail={'MS':MS, 'ME':ME, 'TW':TW}
        # and finally put everything into a top level dictionary,
        # and add that to the list
        results.append( {'h':h , 'strain':strn, 'stress':strs, 'fail':fail }  )
        bot=top
    return results

Example:

In [8]:
m1={'E1':140000, 'E2':10000, 'v12':0.3, 'G12':5000, 'XT':1200, 'XC':800, 'YT':50, 'YC':120, 'S12':75, 'f12':-0.5}

layup1  =[ {'mat':m1 , 'ori':   0  , 'thi':1} , 
           {'mat':m1 , 'ori': -45  , 'thi':1} ,
           {'mat':m1 , 'ori':  90  , 'thi':1} ,
           {'mat':m1 , 'ori':  90  , 'thi':1} ,
           {'mat':m1 , 'ori': +45  , 'thi':1} ,
           {'mat':m1 , 'ori':   0  , 'thi':1}  ]

ABD1=laminateStiffnessMatrix(layup1)

load,defs = solveLaminateLoadCase(ABD=ABD1, Nx=500)

res = layerResults(layup1,defs)

z-coordinate at the bottom of the first layer:

In [9]:
print(  res[0]['h']['bot']        )

-3.0

Stress in the global coordinate system at the bottom of the first layer:

In [10]:
print(  res[0]['stress']['xyz']['bot']   )

[203.66497762   2.71737626   8.01833861]

Stress in the material coordinate system at the top of the second layer:

In [11]:
print(  res[1]['stress']['123']['top']   )

[55.56324148 10.28461825  8.06807109]

Exposure factors for all layers and all positions based on the Tsai-Wu failure criterion:

In [12]:
for r in res:
    print (r['fail']['TW']['bot'])
    print (r['fail']['TW']['top'])

0.19541827216084476
0.18062327743397455
0.2623044214314118
0.21343870531417133
0.2996638330079345
0.29666963435116495
0.29666963435116495
0.2996638330079345
0.21343870531417133
0.2623044214314118
0.18062327743397455
0.19541827216084476

Visualization

The function plotLayerStresses(results) found in Plot Gallery makes graphs of layer results through the thickness of a laminate where results are obtained by the function layerResults(layup,deformations).

In [13]:
from plotlib import plotLayerStresses
%matplotlib inline

Examples:

In [14]:
m1={'E1':140000, 'E2':10000, 'v12':0.3, 'G12':5000, 'XT':1200, 'XC':800, 'YT':50, 'YC':120, 'S12':75, 'f12':-0.5}

layup1  =[ {'mat':m1 , 'ori':   0  , 'thi':1} , 
           {'mat':m1 , 'ori':  90  , 'thi':1} ,
           {'mat':m1 , 'ori':   0  , 'thi':1}   ]

ABD1=laminateStiffnessMatrix(layup1)
load,defs = solveLaminateLoadCase(ABD=ABD1, Nx=500, Nxy=200)
res = layerResults(layup1,defs)    
plotLayerStresses(res)
In [15]:
m1={'E1':140000, 'E2':10000, 'v12':0.3, 'G12':5000, 'XT':1200, 'XC':800, 'YT':50, 'YC':120, 'S12':75, 'f12':-0.5}

layup1  =[ {'mat':m1 , 'ori':   0  , 'thi':1} , 
           {'mat':m1 , 'ori': -45  , 'thi':1} ,
           {'mat':m1 , 'ori': +45  , 'thi':1} ,
           {'mat':m1 , 'ori': +45  , 'thi':1} ,
           {'mat':m1 , 'ori': -45  , 'thi':1} ,
           {'mat':m1 , 'ori':   0  , 'thi':1}  ]

ABD1=laminateStiffnessMatrix(layup1)
load,defs = solveLaminateLoadCase(ABD=ABD1, Nx=500)
res = layerResults(layup1,defs)    
plotLayerStresses(res)
In [16]:
m1={'E1':140000, 'E2':10000, 'v12':0.3, 'G12':5000, 'XT':1200, 'XC':800, 'YT':50, 'YC':120, 'S12':75, 'f12':-0.5}

layup1  =[ {'mat':m1 , 'ori':   0  , 'thi':1} , 
           {'mat':m1 , 'ori':  90  , 'thi':1} ,
           {'mat':m1 , 'ori':   0  , 'thi':1} ,
           {'mat':m1 , 'ori':  90  , 'thi':1} ,
           {'mat':m1 , 'ori':   0  , 'thi':1}  ]

ABD1=laminateStiffnessMatrix(layup1)
load,defs = solveLaminateLoadCase(ABD=ABD1, Ny=500, Mx=200)
res = layerResults(layup1,defs)    
plotLayerStresses(res)

The function plotLayerFailure(results) found in Plot Gallery makes graphs of layer results through the thickness of a laminate where results are obtained by the function layerResults(layup,deformations).

In [17]:
from plotlib import plotLayerFailure
%matplotlib inline

Examples:

In [18]:
m1={'E1':140000, 'E2':10000, 'v12':0.3, 'G12':5000, 'XT':1200, 'XC':800, 'YT':50, 'YC':120, 'S12':75, 'f12':-0.5}

layup1  =[ {'mat':m1 , 'ori':   0  , 'thi':1} , 
           {'mat':m1 , 'ori':  90  , 'thi':1} ,
           {'mat':m1 , 'ori':   0  , 'thi':1} ,
           {'mat':m1 , 'ori':  90  , 'thi':1} ,
           {'mat':m1 , 'ori':   0  , 'thi':1}  ]

ABD1=laminateStiffnessMatrix(layup1)
load,defs = solveLaminateLoadCase(ABD=ABD1, Ny=500, Mx=200)
res = layerResults(layup1,defs)    
plotLayerFailure(res)

Next Previous

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

Copyright 2021, All right reserved, I guess.