Hi, I am new to FE and trying to model the following problem.
I have an elastic square with the following topology:
where blue and orange domains have different 'pressures', such that the force on the boundary between two domains is
Currently I define pressure as a discontinuous Lagrange element of degree 0, as in The FEniCS book Chapter 1.5.
My uff code looks like this:
V = VectorElement("Lagrange", triangle, 2)
W = FiniteElement("DG", triangle, 0)
u = TrialFunction(V)
v = TestFunction(V)
p = Coefficient(W)
E = Coefficient(W)
nu = Coefficient(W)
n = FacetNormal(triangle)
mu = E / (2.0*(1.0 + nu))
lmbda = E*nu / ((1.0+nu) * (1.0 - 2.0*nu))
def sigma(v):
return 2.0 * mu * sym(grad(v)) + lmbda*tr(sym(grad(v))) * Identity(2)
a = inner( sigma(u), sym(grad(v)) ) * dx
L = dot(p*n,v)*ds(0)
Since I want boundary integral to be calculated on every element face I define the Neumann BC as
class NeumannBoundary: public SubDomain{
bool inside(const Array<double>& x, bool on_boundary) const{
return true;
}
};
I also use a simple Dirichlet BC:
class DirichletBoundary: public SubDomain{
bool inside(const Array<double>& x, bool on_boundary) const{
return (abs(x[0]) <= 10 && on_boundary);
}
};
Where I prescribe u = (0,0)
When I run the solver, with zero 'pressure' in the outside (blue) material and 100 on the inside (orange) material, I get u=(0,0) everywhere.
If I set pressure of both regions to 100 I get the following solution:
Which is the same if I set pressure of the outer region to 100 and inner to 0.
So, I guess my boundary condition definition is wrong. Hence is the question: How to define internal surface conditions? (Sorry for the long intro)