Note
Go to the end to download the full example code.
A square patch made of two triangular plate elements
- Basic implementation test with applied loads.
Testing the tangent stiffness computation for a
Triangle6()
(using quadratic shape functions).
Using
elements.linear.Triangle6
(see Triangle6 class)materials.PlaneStress
(see PlaneStress material class)
free free
^ ^
| |
3--6--2 -> free
|\ b | >
| \ | >
7 8 5 > (w = 1.0)
| \ | >
| a \| >
0--4--1 -> free
width: 10.
height: 10.
Material parameters: St. Venant-Kirchhoff, plane stress
E = 10.0
nu = 0.30
t = 1.0
Element loads:
node 0: [ 0.0, 0.0]
node 1: [ 10./6, 0.0]
node 2: [ 10./6, 0.0]
node 3: [ 0.0, 0.0]
node 4: [ 0.0, 0.0]
node 5: [ 10.*2/3, 0.0]
node 6: [ 0.0, 0.0]
node 7: [ 0.0, 0.0]
node 8: [ 0.0, 0.0]
2nd Piola-Kirchhoff stress:
S_XX = w = 1.000
S_YY = S_XY = S_YX = S_ZZ = 0.000
Green Lagrange strain:
eps_XX = (1/E) ((1.000) - (0.30)(0.000)) = 0.100
eps_YY = (1/E) ((0.000) - (0.30)(1.000)) = -0.030
eps_XY = eps_YX = 0.000
eps_ZZ = -nu * (eps_XX + eps_YY) = -0.021
Stretches:
lam_X = sqrt(1 + 2 eps_XX) = 1.095
lam_Y = sqrt(1 + 2 eps_YY) = 0.9695
Displacements:
ux = (lam_X - 1) * x, uy = (lam_Y - 1) * y
node 0: [ 0.000, 0.000 ]
node 1: [ 0.954, 0.000 ]
node 2: [ 0.954, -0.305 ]
node 3: [ 0.000, -0.305 ]
node 4: [ 0.477, 0.000 ]
node 5: [ 0.954, -0.1525 ]
node 6: [ 0.477, -0.305 ]
node 7: [ 0.954, -0.1525 ]
node 8: [ 0.477, -0.1525 ]
Author: Peter Mackenzie-Helnwein
from femedu.examples import Example
from femedu.domain import System, Node
from femedu.solver import NewtonRaphsonSolver
from femedu.elements.linear import Triangle6
from femedu.materials import PlaneStress
class ExamplePlate02b(Example):
def problem(self):
params = dict(
E = 10., # Young's modulus
nu = 0.3, # Poisson's ratio
t = 1.0, # thickness of the plate
fy = 1.e30 # yield stress
)
a = 10. # length of the plate in the x-direction
b = 10. # length of the plate in the y-direction
model = System()
model.setSolver(NewtonRaphsonSolver())
nd0 = Node( 0.0, 0.0)
nd1 = Node( a, 0.0)
nd2 = Node( a, b)
nd3 = Node( 0.0, b)
nd4 = Node( a/2, 0.0)
nd5 = Node( a, b/2)
nd6 = Node( a/2, b)
nd7 = Node( 0.0, b/2)
nd8 = Node( a/2, b/2)
nd0.fixDOF('ux', 'uy')
nd1.fixDOF('uy')
nd3.fixDOF('ux')
model.addNode(nd0, nd1, nd2, nd3, nd4, nd5, nd6, nd7, nd8)
elemA = Triangle6(nd0, nd1, nd3, nd4, nd8, nd7, PlaneStress(params))
elemB = Triangle6(nd2, nd3, nd1, nd6, nd8, nd5, PlaneStress(params))
model.addElement(elemA, elemB)
elemA.setSurfaceLoad(face=2, pn=1.0)
elemB.setSurfaceLoad(face=2, pn=1.0)
model.setLoadFactor(1.0)
nd0.setDisp([0.0, 0.0])
nd1.setDisp([5.0, 0.0])
nd2.setDisp([5.0, -5.0])
nd3.setDisp([0.0, -5.0])
nd4.setDisp([2.5, 0.0])
nd5.setDisp([5.0, -2.5])
nd6.setDisp([2.5, -5.0])
nd7.setDisp([0.0, -2.5])
nd8.setDisp([2.5, -2.5])
elemA.updateState()
elemB.updateState()
model.report()
model.plot(factor=1.0, filename="plate02b_deformed.png")
model.setLoadFactor(0.0)
# model.solver.assemble()
# model.solver.showKt()
#
model.solve()
model.report() # activate this line for lots of debug info
model.plot(factor=0.0, title="Undeformed system", filename="plate02b_undeformed.png", show_bc=1)
model.setLoadFactor(1.0)
model.solve(verbose=1)
model.plot(factor=1.0, filename="plate02b_deformed.png")
model.report()
# requires femedu-1.0.25 or newer
model.valuePlot('sxx', show_mesh=True)
model.valuePlot('syy', show_mesh=True)
model.valuePlot('sxy', show_mesh=True)
Run the example by creating an instance of the problem and executing it by calling Example.run()
if __name__ == "__main__":
ex = ExamplePlate02b()
ex.run()
System Analysis Report
=======================
Nodes:
---------------------
Node_601:
x: [0.000 0.000]
fix: ['ux', 'uy']
u: [0.000 0.000]
Node_602:
x: [10.000 0.000]
fix: ['uy']
u: [5.000 0.000]
Node_603:
x: [10.000 10.000]
u: [5.000 -5.000]
Node_604:
x: [0.000 10.000]
fix: ['ux']
u: [0.000 -5.000]
Node_605:
x: [5.000 0.000]
u: [2.500 0.000]
Node_606:
x: [10.000 5.000]
u: [5.000 -2.500]
Node_607:
x: [5.000 10.000]
u: [2.500 -5.000]
Node_608:
x: [0.000 5.000]
u: [0.000 -2.500]
Node_609:
x: [5.000 5.000]
u: [2.500 -2.500]
Elements:
---------------------
Triangle6_836: nodes ( Node_601 Node_602 Node_604 Node_605 Node_609 Node_608 )
material: PlaneStress
strain 0: xx=5.000e-01 yy=-5.000e-01 xy=0.000e+00 zz=-1.332e-16
stress 0: xx=3.846e+00 yy=-3.846e+00 xy=0.000e+00 zz=0.000e+00
strain 1: xx=5.000e-01 yy=-5.000e-01 xy=0.000e+00 zz=-0.000e+00
stress 1: xx=3.846e+00 yy=-3.846e+00 xy=0.000e+00 zz=0.000e+00
strain 2: xx=5.000e-01 yy=-5.000e-01 xy=0.000e+00 zz=-6.661e-17
stress 2: xx=3.846e+00 yy=-3.846e+00 xy=0.000e+00 zz=0.000e+00
Triangle6_837: nodes ( Node_603 Node_604 Node_602 Node_607 Node_609 Node_606 )
material: PlaneStress
strain 0: xx=5.000e-01 yy=-5.000e-01 xy=-1.243e-15 zz=1.665e-16
stress 0: xx=3.846e+00 yy=-3.846e+00 xy=-4.782e-15 zz=0.000e+00
strain 1: xx=5.000e-01 yy=-5.000e-01 xy=-1.776e-16 zz=-0.000e+00
stress 1: xx=3.846e+00 yy=-3.846e+00 xy=-6.832e-16 zz=0.000e+00
strain 2: xx=5.000e-01 yy=-5.000e-01 xy=-5.329e-16 zz=-0.000e+00
stress 2: xx=3.846e+00 yy=-3.846e+00 xy=-2.050e-15 zz=0.000e+00
+
System Analysis Report
=======================
Nodes:
---------------------
Node_601:
x: [0.000 0.000]
fix: ['ux', 'uy']
u: [0.000 0.000]
Node_602:
x: [10.000 0.000]
fix: ['uy']
u: [-0.000 0.000]
Node_603:
x: [10.000 10.000]
u: [-0.000 0.000]
Node_604:
x: [0.000 10.000]
fix: ['ux']
u: [0.000 0.000]
Node_605:
x: [5.000 0.000]
u: [-0.000 -0.000]
Node_606:
x: [10.000 5.000]
u: [-0.000 0.000]
Node_607:
x: [5.000 10.000]
u: [-0.000 0.000]
Node_608:
x: [0.000 5.000]
u: [-0.000 0.000]
Node_609:
x: [5.000 5.000]
u: [-0.000 0.000]
Elements:
---------------------
Triangle6_836: nodes ( Node_601 Node_602 Node_604 Node_605 Node_609 Node_608 )
material: PlaneStress
strain 0: xx=-9.897e-12 yy=3.290e-11 xy=-3.957e-17 zz=-6.900e-12
stress 0: xx=-3.070e-13 yy=3.289e-10 xy=-1.522e-16 zz=0.000e+00
strain 1: xx=-9.897e-12 yy=3.290e-11 xy=2.220e-16 zz=-6.900e-12
stress 1: xx=-3.070e-13 yy=3.289e-10 xy=8.540e-16 zz=0.000e+00
strain 2: xx=-9.897e-12 yy=3.290e-11 xy=8.398e-17 zz=-6.900e-12
stress 2: xx=-3.053e-13 yy=3.289e-10 xy=3.230e-16 zz=0.000e+00
Triangle6_837: nodes ( Node_603 Node_604 Node_602 Node_607 Node_609 Node_606 )
material: PlaneStress
strain 0: xx=-9.898e-12 yy=3.290e-11 xy=-7.105e-16 zz=-6.900e-12
stress 0: xx=-3.121e-13 yy=3.289e-10 xy=-2.733e-15 zz=0.000e+00
strain 1: xx=-9.897e-12 yy=3.290e-11 xy=-1.776e-16 zz=-6.900e-12
stress 1: xx=-3.077e-13 yy=3.289e-10 xy=-6.832e-16 zz=0.000e+00
strain 2: xx=-9.897e-12 yy=3.290e-11 xy=-3.553e-16 zz=-6.900e-12
stress 2: xx=-3.082e-13 yy=3.289e-10 xy=-1.366e-15 zz=0.000e+00
/Users/pmackenz/Development/Educational/FEM.edu/venv/lib/python3.13/site-packages/matplotlib/quiver.py:678: RuntimeWarning: Mean of empty slice.
amean = a.mean()
/Users/pmackenz/Development/Educational/FEM.edu/venv/lib/python3.13/site-packages/numpy/_core/_methods.py:145: RuntimeWarning: invalid value encountered in scalar divide
ret = ret.dtype.type(ret / rcount)
norm of the out-of-balance force: 9.7183e+00
norm of the out-of-balance force: 9.5873e-11
+
System Analysis Report
=======================
Nodes:
---------------------
Node_601:
x: [0.000 0.000]
fix: ['ux', 'uy']
u: [0.000 0.000]
Node_602:
x: [10.000 0.000]
fix: ['uy']
u: [1.000 0.000]
Node_603:
x: [10.000 10.000]
u: [1.000 -0.300]
Node_604:
x: [0.000 10.000]
fix: ['ux']
u: [0.000 -0.300]
Node_605:
x: [5.000 0.000]
u: [0.500 -0.000]
Node_606:
x: [10.000 5.000]
u: [1.000 -0.150]
Node_607:
x: [5.000 10.000]
u: [0.500 -0.300]
Node_608:
x: [0.000 5.000]
u: [0.000 -0.150]
Node_609:
x: [5.000 5.000]
u: [0.500 -0.150]
Elements:
---------------------
Triangle6_836: nodes ( Node_601 Node_602 Node_604 Node_605 Node_609 Node_608 )
material: PlaneStress
strain 0: xx=1.000e-01 yy=-3.000e-02 xy=6.716e-17 zz=-2.100e-02
stress 0: xx=1.000e+00 yy=9.865e-12 xy=2.583e-16 zz=0.000e+00
strain 1: xx=1.000e-01 yy=-3.000e-02 xy=1.332e-16 zz=-2.100e-02
stress 1: xx=1.000e+00 yy=9.866e-12 xy=5.124e-16 zz=0.000e+00
strain 2: xx=1.000e-01 yy=-3.000e-02 xy=6.607e-17 zz=-2.100e-02
stress 2: xx=1.000e+00 yy=9.869e-12 xy=2.541e-16 zz=0.000e+00
Triangle6_837: nodes ( Node_603 Node_604 Node_602 Node_607 Node_609 Node_606 )
material: PlaneStress
strain 0: xx=1.000e-01 yy=-3.000e-02 xy=7.745e-16 zz=-2.100e-02
stress 0: xx=1.000e+00 yy=9.865e-12 xy=2.979e-15 zz=0.000e+00
strain 1: xx=1.000e-01 yy=-3.000e-02 xy=5.329e-17 zz=-2.100e-02
stress 1: xx=1.000e+00 yy=9.864e-12 xy=2.050e-16 zz=0.000e+00
strain 2: xx=1.000e-01 yy=-3.000e-02 xy=3.197e-16 zz=-2.100e-02
stress 2: xx=1.000e+00 yy=9.867e-12 xy=1.230e-15 zz=0.000e+00
Total running time of the script: (0 minutes 0.751 seconds)