elasticity1.py 1.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748
  1. """
  2. Show mesh and displacement solution with arrows.
  3. Refer to original script for the detailed description:
  4. https://fenicsproject.org/docs/dolfin/2018.1.0/python/
  5. demos/hyperelasticity/demo_hyperelasticity.py.html
  6. """
  7. print(__doc__)
  8. ########################################################### dolfin
  9. from dolfin import *
  10. # Create mesh and define function space
  11. mesh = UnitCubeMesh(10, 10, 10)
  12. V = VectorFunctionSpace(mesh, "Lagrange", 1)
  13. # Mark boundary subdomains
  14. left = CompiledSubDomain("near(x[0], side) && on_boundary", side=0.0)
  15. right = CompiledSubDomain("near(x[0], side) && on_boundary", side=1.0)
  16. # Define Dirichlet boundary (x = 0 or x = 1)
  17. c = Constant((0.0, 0.0, 0.0))
  18. r = Expression((
  19. "scale*0.0",
  20. "scale*(y0 + (x[1]-y0)*cos(theta) - (x[2]-z0)*sin(theta) - x[1])",
  21. "scale*(z0 + (x[1]-y0)*sin(theta) + (x[2]-z0)*cos(theta) - x[2])",
  22. ),
  23. scale=0.5, y0=0.5, z0=0.5, theta=pi/3, degree=2 )
  24. bcl = DirichletBC(V, c, left)
  25. bcr = DirichletBC(V, r, right)
  26. w = TrialFunction(V) # Incremental displacement
  27. v = TestFunction(V) # Test function
  28. u = Function(V) # solution
  29. solve(inner(grad(w), grad(v)) * dx == inner(c, v) * dx, u, [bcl, bcr])
  30. bmesh = BoundaryMesh(mesh, "exterior")
  31. ########################################################### vedo
  32. from vedo.dolfin import *
  33. # ps = point size, only mesh vertices are shown
  34. plot(u, mode='mesh', ps=10, scalarbar='3d')
  35. # plot displacements as white arrows, lw controls the mesh visibility
  36. plot(u, mode='arrows', add=True,
  37. color='w', alpha=0.5, cmap='gist_earth', lw=1)