(Deleted) Elastic Membrane Problems
1st question) Does the xml file, converted with dolfin keep exact the same mesh ?
I use Jupyter inside a Docker container and i converted the file from inside . I say this because i wrote a code for
a simple elastic problem as in the tutorial but it seems it has a few problems .
from fenics import *
from dolfin import *
from mshr import *
import numpy as np
#Geometries and Meshing with GMSH
mesh = Mesh("me.xml");
V = VectorFunctionSpace(mesh, "P", 1)
d = 799
r = d/2
mu = 1
rho = 1
beta = 1.25
lambda_ = beta
g = 9.81
#Def boundary condition, with tolerance
tol = 1E-14
def clamped_boundary(x, on_boundary):
return on_boundary and np.abs(x**2 + x**2 > r**2)
bc = DirichletBC(V, (0,0,0), clamped_boundary)
# Define Strain
return lambda_*nabla_div(u)*Identity(d) + 2*mu*epsilon(u)
# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Constant((0, 0, -rho*g))
d = u.geometric_dimension() #SPACE DIMENSION
T = Constant((0,0,0))
a = inner(sigma(u), epsilon(v))*dx
L = dot(f, v)*dx + dot(T, v)*ds
# Compute solution
u = Function(V)
solve(a == L, u, bc)
vtkfile = File("displacement.pvd")
vtkfile << (u)
s = sigma(u) - (1./3)*tr(sigma(u))*Identity(d) # deviatoric stress
von_Mises = sqrt(3./2*inner(s, s))
V = FunctionSpace(mesh, "Lagrange", 2)
von_Mises = project(von_Mises, V)
vtkfile = File("stresses.pvd")
vtkfile << (von_Mises)
This is the code and these are the results:
Geometry and mesh imported from gmsh and resaved as vtk file, it looks completely fine but...
Stresses and here it starts to show problems, it could be still fine but ...
And this is the displecemnt once i show WarpByVector, and it's not that i m watching only one axis , if i rotate it still keep being a spike like if the elements fall in themself. When if i use WarpByVector in the Linear Elastic Beam problem in the tutorial everything works fine.
I'm thinking there are some problem with the mesh but i'm not sure about it.
Thanks to everyone
The diameter of your disc is 799 meters while the Lamé constants are 1 Pa (should be in the range of 1e7-1e9)
You can set the scale factor really small in warp by vector. But I'd start using reasonable values for the constants.
If you use piecewise linear function space for the displacements, why would you project the mises stress onto a piecewise quadratic space? Naturally, the Mises stress should be of the same order as your strains, in this case come from a piecewise constant space:
DGSpace = FunctionSpace(mesh, FiniteElement('DG', mesh.ufl_cell, 0))
Maybe this would lift the problem of the awkward oscillations in your Mises plot.