PhsycisSim.py

from dolfinx import *
import mshr

# Load STL file and create mesh
stl_file = "your_model.stl"
domain = mshr.Mesh3D(stl_file)
mesh = Mesh(domain)

# Define function space (Displacement)
V = VectorFunctionSpace(mesh, "Lagrange", 1)

# Define boundary condition
tol = 1E-14
# Apply zero displacement to a specific region (supports)
def boundary(x, on_boundary):
    return on_boundary and x[0] < tol  # Adjust this based on your supports

bc = DirichletBC(V, Constant((0.0, 0.0, 0.0)), boundary)

# Define force (load) applied to another region
# You can adjust the coordinates according to your problem
force = Constant((0, 0, -10))  # Example force in the z-direction

# Define the elasticity problem (displacement, stress, strain), we need to calculate these from the material properites of the plastic
mu = Constant(10.0)  # Lame's first parameter
lambda_ = Constant(100.0)  # Lame's second parameter

def epsilon(u):
    return sym(grad(u))

def sigma(u):
    return lambda_ * div(u) * Identity(3) + 2 * mu * epsilon(u)

# Define variational problem
u = TrialFunction(V)
d = u.geometric_dimension()  # number of space dimensions
v = TestFunction(V)
a = inner(sigma(u), epsilon(v)) * dx
L = dot(force, v) * ds

# Solve the problem
u = Function(V)
solve(a == L, u, bc)

# Post-processing: Compute stress and strain
# Strain
eps = epsilon(u)
# Stress
sig = sigma(u)

# Save the results
File("displacement.pvd") << u
File("strain.pvd") << project(eps, TensorFunctionSpace(mesh, "Lagrange", 1))
File("stress.pvd") << project(sig, TensorFunctionSpace(mesh, "Lagrange", 1))

print("FEM calculation complete. Results saved.")