How to set the pressure boundary condition？

hello,I was trying to do a simulation of a biweekly compression test (2D), and the part that added confining pressure to the left and right sides of the model stumped me. I don’t know how to properly add these pressures, I tried to write a part, but it did not work

``````from dolfin import *

L = 0.3
H = 0.6
Nx = 50
Ny = 25
mesh = RectangleMesh(Point(0., 0.), Point(L, H), Nx, Ny, "crossed")

def eps(v):

E = Constant(1e5)
nu = Constant(0.3)
model = "plane_stress"

mu = E/2/(1+nu)
lmbda = E*nu/(1+nu)/(1-2*nu)
if model == "plane_stress":
lmbda = 2*mu*lmbda/(lmbda+2*mu)

def sigma(v):
return lmbda*tr(eps(v))*Identity(2) + 2.0*mu*eps(v)

rho_g = 10.
f = Constant((0, -rho_g))

class Top(SubDomain):
def inside(self,x,on_boundary):
return near(x[1],H) and on_boundary

class Left(SubDomain):
def inside(self,x,on_boundary):
return near(x[0],0) and on_boundary
class Right(SubDomain):
def inside(self,x,on_boundary):
return near(x[0],0.3) and on_boundary

facets = MeshFunction('size_t',mesh,1)

Top().mark(facets,1)
Left().mark((facets,2))
Right().mark(facets,3)
ds = Measure('ds',subdomain_data=facets)

V = VectorFunctionSpace(mesh, 'Lagrange', degree=2)

#Bottom fixing
def bottom(x,on_boundary):
return  near(x[1],0) and on_boundary
bc = DirichletBC(V,Constant((0,0,0)),bottom)

T2 = Constant(6000.,0.)  #the confine pressure of the left boundary
T3 = Constant(-6000.,0.) #the confine pressure of the right boundary
du = TrialFunction(V)
u_ = TestFunction(V)
a = inner(sigma(du), eps(u_))*dx
l = inner(f, u_)*dx + dot(T1,u_)*ds(1)+dot(T2,u_)*ds(2)+dot()+dot(T3,u_)*ds(3)

u = Function(V, name="Displacement")
solve(a == l, u, bc)
``````

Could you be more specific as to what “did not work”? Are you receiving an error in the command line or are the results not what you expect. As posted you have a ufl.dot() with no arguments, are you receiving a TypeError? Please provide more information as per read-before-posting-how-do-i-get-my-question-answered.

1 Like