Setting values on specified interval in 2D domain as boundary condition

Question

I solve problem in circle and I want to set one condition on boundary and one on specified interval inside the domain. So I want something like circle with incision to set boundary conditions.
I tried next approaches:

1. To cut off the thick rectangle and set a DirichletBC on it. If thickness is less than 0.001 I don't get incision, but can see an interval on mesh. So this approach works only with large enough thickness. It leads to a large error relative to the exact solution.

   mesh = generate_mesh(Circle(Point(0., 0.), R, 150) - Rectangle(Point(-l, -w), Point(l, w)), count)
   

2. To cut off rectangle with thickness less than 0.001 to get an interval on mesh and set DirichletBC:

   class InnerBoundary(SubDomain):
       def inside(self, x, on_boundary):
           return abs(x[0]) <= l and near(abs(x[1]), 0, eps=10E-5)   
   
   bcs = [DirichletBC(V, Constant(1.), InnerBoundary()),
      DirichletBC(V, Constant(0.), OuterBoundary())]
   

   But when I check values in points of interval, I get them with error about 10^-2. As a result it leads to undesirable error in solution.

Is there an approach to set condition inside the domain, not only on boudary?

Answer 1

Consider the following example for generating a domain with a slit.

from dolfin import *
from mshr import *
from matplotlib import pyplot
parameters["plotting_backend"] = "matplotlib"

N = 40
r = 1.0
l = 1.0
w = 1e-2

D = Circle(Point(0.,0.), 1)

# define polygon to cut out, 
c = Polygon([Point(-l/2,0), Point(0,-w), Point(l/2,0), Point(0, w)])

domain = D - c
mesh = generate_mesh(domain, N)

# mark the different parts of the boundary
boundary = CompiledSubDomain("on_boundary")
inner_boundary = CompiledSubDomain("on_boundary&&(fabs(x[0])<l)&&(fabs(x[1])<w+t)",
                                   l = l, w = w, t = DOLFIN_EPS)
facet_domains = FacetFunction("size_t", mesh)
boundary.mark(facet_domains, 1) 
inner_boundary.mark(facet_domains, 2)

# move mesh to close the gap in the domain
W = VectorFunctionSpace(mesh, "CG", 1)
t = Function(W)
DirichletBC(W, Expression(("0", "-x[1]"), degree = 1), facet_domains, 2).apply(t.vector())
ALE.move(mesh, t)

pyplot.figure()
plot(mesh)
pyplot.show()

Solving Laplace equation with Dirichlet boundaries:

V = FunctionSpace(mesh, "CG", 1)
u, v = TrialFunction(V), TestFunction(V)
a = inner(grad(u), grad(v)) * dx
L = Constant(0.) * v * dx

bcs = [DirichletBC(V, Constant(0.0), facet_domains, 1),
       DirichletBC(V, Constant(1.0), facet_domains, 2)]
uh = Function(V)
solve(a == L, uh, bcs)

pyplot.figure()
plot(uh)
pyplot.show()

You can also have discontinuity across the slit, for example with a Neumann boundary condition on the slit:

bcs = DirichletBC(V, Expression("x[1]", degree = 1), facet_domains, 1)
solve(a == L, uh, bcs)

pyplot.figure()
plot(uh)
pyplot.show()

Comments

I can't see any significant difference between this example and my approach.