This page was generated from wta/elasticity3D.ipynb.

3D Solid Mechanics¶

[1]:

from netgen.occ import *
from netgen.webgui import Draw as DrawGeo
import ngsolve

[2]:

box = Box((0,0,0), (3,0.6,1))
box.faces.name="outer"
cyl = sum( [Cylinder((0.5+i,0,0.5), Y, 0.25,0.8) for i in range(3)] )
cyl.faces.name="cyl"
geo = box-cyl

DrawGeo(geo);

find edges between box and cylinder, and build chamfers (requires OCC 7.4 or newer):

[3]:

cylboxedges = geo.faces["outer"].edges * geo.faces["cyl"].edges
cylboxedges.name = "cylbox"
geo = geo.MakeChamfer(cylboxedges, 0.03)

name faces for boundary conditions:

[4]:

geo.faces.Min(X).name = "fix"
geo.faces.Max(X).name = "force"

DrawGeo(geo);

[5]:

from ngsolve import *
from ngsolve.webgui import Draw
mesh = Mesh(OCCGeometry(geo).GenerateMesh(maxh=0.1)).Curve(3)
Draw (mesh);

Linear elasticity¶

Displacement: \(u : \Omega \rightarrow R^3\)

Linear strain:

\[\varepsilon(u) := \tfrac{1}{2} ( \nabla u + (\nabla u)^T )\]

Stress by Hooke’s law:

\[\sigma = 2 \mu \varepsilon + \lambda \operatorname{tr} \varepsilon I\]

Equilibrium of forces:

\[\operatorname{div} \sigma = f\]

Displacement boundary conditions:

\[u = u_D \qquad \text{on} \, \Gamma_D\]

Traction boundary conditions:

\[\sigma n = g \qquad \text{on} \, \Gamma_N\]

Variational formulation:¶

Find: \(u \in H^1(\Omega)^3\) such that \(u = u_D\) on \(\Gamma_D\)

\[\int_\Omega \sigma(\varepsilon(u)) : \varepsilon(v) \, dx = \int_\Omega f v dx + \int_{\Gamma_N} g v ds\]

holds for all \(v = 0\) on \(\Gamma_D\).

[6]:

E, nu = 210, 0.2
mu  = E / 2 / (1+nu)
lam = E * nu / ((1+nu)*(1-2*nu))

def Stress(strain):
    return 2*mu*strain + lam*Trace(strain)*Id(3)

[7]:

fes = VectorH1(mesh, order=3, dirichlet="fix")
u,v = fes.TnT()
gfu = GridFunction(fes)

with TaskManager():
    a = BilinearForm(InnerProduct(Stress(Sym(Grad(u))), Sym(Grad(v))).Compile()*dx)
    pre = Preconditioner(a, "bddc")
    a.Assemble()

[8]:

force = CF( (1e-3,0,0) )
f = LinearForm(force*v*ds("force")).Assemble()

[9]:

from ngsolve.krylovspace import CGSolver
inv = CGSolver(a.mat, pre, printrates='\r', tol=1e-8)
gfu.vec.data = inv * f.vec

CG iteration 1, residual = 0.00018026680835987196

CG iteration 2, residual = 7.56651448702291e-05: </pre>
CG iteration 2, residual = 7.56651448702291e-05: end{sphinxVerbatim}

[2KCG iteration 2, residual = 7.56651448702291e-05

CG iteration 3, residual = 8.455269170458239e-05

CG iteration 4, residual = 6.401519640134442e-05: </pre>
CG iteration 4, residual = 6.401519640134442e-05: end{sphinxVerbatim}

[2KCG iteration 4, residual = 6.401519640134442e-05

CG iteration 5, residual = 5.364767130768229e-05

CG iteration 6, residual = 4.8415655080243797e-05: </pre>
CG iteration 6, residual = 4.8415655080243797e-05: end{sphinxVerbatim}

[2KCG iteration 6, residual = 4.8415655080243797e-05

CG iteration 7, residual = 3.18085266995706e-05

CG iteration 8, residual = 2.500044384709958e-05: </pre>
CG iteration 8, residual = 2.500044384709958e-05: end{sphinxVerbatim}

[2KCG iteration 8, residual = 2.500044384709958e-05

CG iteration 9, residual = 1.8919147949149747e-05

CG iteration 10, residual = 1.3320560343250076e-05: </pre>
CG iteration 10, residual = 1.3320560343250076e-05: end{sphinxVerbatim}

[2KCG iteration 10, residual = 1.3320560343250076e-05

CG iteration 11, residual = 1.0719430516097942e-05

CG iteration 12, residual = 7.345111519344642e-06: </pre>
CG iteration 12, residual = 7.345111519344642e-06: end{sphinxVerbatim}

[2KCG iteration 12, residual = 7.345111519344642e-06

CG iteration 13, residual = 5.620103679748098e-06

CG iteration 14, residual = 3.8723573951237415e-06: </pre>
CG iteration 14, residual = 3.8723573951237415e-06: end{sphinxVerbatim}

[2KCG iteration 14, residual = 3.8723573951237415e-06

CG iteration 15, residual = 3.0676127439615533e-06

CG iteration 16, residual = 2.132427470070547e-06: </pre>
CG iteration 16, residual = 2.132427470070547e-06: end{sphinxVerbatim}

[2KCG iteration 16, residual = 2.132427470070547e-06

CG iteration 17, residual = 1.6603637101490366e-06

CG iteration 18, residual = 1.1314359574074014e-06: </pre>
CG iteration 18, residual = 1.1314359574074014e-06: end{sphinxVerbatim}

[2KCG iteration 18, residual = 1.1314359574074014e-06

CG iteration 19, residual = 8.685484773102635e-07

CG iteration 20, residual = 6.400136038782427e-07: </pre>
CG iteration 20, residual = 6.400136038782427e-07: end{sphinxVerbatim}

[2KCG iteration 20, residual = 6.400136038782427e-07

CG iteration 21, residual = 4.830602416846575e-07

CG iteration 22, residual = 3.7299454675368264e-07: </pre>
CG iteration 22, residual = 3.7299454675368264e-07: end{sphinxVerbatim}

[2KCG iteration 22, residual = 3.7299454675368264e-07

CG iteration 23, residual = 2.5067668469398104e-07

CG iteration 24, residual = 1.9245101311214237e-07: </pre>
CG iteration 24, residual = 1.9245101311214237e-07: end{sphinxVerbatim}

[2KCG iteration 24, residual = 1.9245101311214237e-07

CG iteration 25, residual = 1.4045037293564497e-07

CG iteration 26, residual = 1.0298833212850939e-07: </pre>
CG iteration 26, residual = 1.0298833212850939e-07: end{sphinxVerbatim}

[2KCG iteration 26, residual = 1.0298833212850939e-07

CG iteration 27, residual = 7.603234147346014e-08

CG iteration 28, residual = 5.773480819122819e-08: </pre>
CG iteration 28, residual = 5.773480819122819e-08: end{sphinxVerbatim}

[2KCG iteration 28, residual = 5.773480819122819e-08

CG iteration 29, residual = 4.261973216124953e-08

CG iteration 30, residual = 3.0223706215661516e-08: </pre>
CG iteration 30, residual = 3.0223706215661516e-08: end{sphinxVerbatim}

[2KCG iteration 30, residual = 3.0223706215661516e-08

CG iteration 31, residual = 2.1880790005333037e-08

CG iteration 32, residual = 1.4270526359066452e-08: </pre>
CG iteration 32, residual = 1.4270526359066452e-08: end{sphinxVerbatim}

[2KCG iteration 32, residual = 1.4270526359066452e-08

CG iteration 33, residual = 1.011043670092474e-08

CG iteration 34, residual = 8.005974956583896e-09: </pre>
CG iteration 34, residual = 8.005974956583896e-09: end{sphinxVerbatim}

[2KCG iteration 34, residual = 8.005974956583896e-09

CG iteration 35, residual = 5.153363589259495e-09

CG iteration 36, residual = 4.161271410520076e-09: </pre>
CG iteration 36, residual = 4.161271410520076e-09: end{sphinxVerbatim}

[2KCG iteration 36, residual = 4.161271410520076e-09

CG iteration 37, residual = 3.0463852839213758e-09

CG iteration 38, residual = 1.9992452271656236e-09: </pre>
CG iteration 38, residual = 1.9992452271656236e-09: end{sphinxVerbatim}

[2KCG iteration 38, residual = 1.9992452271656236e-09

CG iteration 39, residual = 1.5715765197986897e-09

CG iteration 40, residual = 1.4792843210847952e-09: </pre>
CG iteration 40, residual = 1.4792843210847952e-09: end{sphinxVerbatim}

[2KCG iteration 40, residual = 1.4792843210847952e-09

CG iteration 41, residual = 8.785156006686644e-10

CG iteration 42, residual = 6.664434413268961e-10: </pre>
CG iteration 42, residual = 6.664434413268961e-10: end{sphinxVerbatim}

[2KCG iteration 42, residual = 6.664434413268961e-10

CG iteration 43, residual = 4.422222154437819e-10

CG iteration 44, residual = 3.1062963423086934e-10: </pre>
CG iteration 44, residual = 3.1062963423086934e-10: end{sphinxVerbatim}

[2KCG iteration 44, residual = 3.1062963423086934e-10

CG iteration 45, residual = 2.1443535849919867e-10

CG iteration 46, residual = 1.4780169642981532e-10: </pre>
CG iteration 46, residual = 1.4780169642981532e-10: end{sphinxVerbatim}

[2KCG iteration 46, residual = 1.4780169642981532e-10

CG iteration 47, residual = 1.047017840572877e-10

CG iteration 48, residual = 7.136892563595456e-11: </pre>
CG iteration 48, residual = 7.136892563595456e-11: end{sphinxVerbatim}

[2KCG iteration 48, residual = 7.136892563595456e-11

CG iteration 49, residual = 4.809978048866156e-11

CG iteration 50, residual = 3.5477143940252223e-11: </pre>
CG iteration 50, residual = 3.5477143940252223e-11: end{sphinxVerbatim}

[2KCG iteration 50, residual = 3.5477143940252223e-11

CG iteration 51, residual = 2.36468357705152e-11

CG iteration 52, residual = 1.6333236744045032e-11: </pre>
CG iteration 52, residual = 1.6333236744045032e-11: end{sphinxVerbatim}

[2KCG iteration 52, residual = 1.6333236744045032e-11

CG iteration 53, residual = 1.074067278811655e-11

CG iteration 54, residual = 7.637737766003765e-12: </pre>
CG iteration 54, residual = 7.637737766003765e-12: end{sphinxVerbatim}

[2KCG iteration 54, residual = 7.637737766003765e-12

CG iteration 55, residual = 7.374278223678638e-12

CG iteration 56, residual = 4.397200485927479e-12: </pre>
CG iteration 56, residual = 4.397200485927479e-12: end{sphinxVerbatim}

[2KCG iteration 56, residual = 4.397200485927479e-12

CG iteration 57, residual = 2.8993813773832056e-12

CG iteration 58, residual = 2.050711041964762e-12: </pre>
CG iteration 58, residual = 2.050711041964762e-12: end{sphinxVerbatim}

[2KCG iteration 58, residual = 2.050711041964762e-12

CG iteration 59, residual = 1.2333447206342998e-12

CG converged in 59 iterations to residual 1.2333447206342998e-12: </pre>
CG converged in 59 iterations to residual 1.2333447206342998e-12: end{sphinxVerbatim}

[2KCG converged in 59 iterations to residual 1.2333447206342998e-12

[10]:

with TaskManager():
    fesstress = MatrixValued(H1(mesh,order=3), symmetric=True)
    gfstress = GridFunction(fesstress)
    gfstress.Interpolate (Stress(Sym(Grad(gfu))))

[11]:

Draw (5e4*gfu, mesh);

[12]:

Draw (Norm(gfstress), mesh, deformation=1e4*gfu, draw_vol=False, order=3);

[ ]: