I would like to numerically find a mutual capacitance of two stripes of metal on the opposites sides of a cylinder. The problem is obviously a 2D Laplace equation. I would like to find the potential outside the cylinder as well. Therefore I have something like this:

mesh = UnitCircleMesh(50)
V = FunctionSpace(mesh, 'Lagrange', 1)

u_L = Constant(-1)
def left_boundary(x, on_boundary):
    r = math.sqrt(x[0] * x[0] + x[1] * x[1])
    return near(r, 0.7) and x[0] < 0 and between(x[1], (-0.5, 0.5))

u_R = Constant(1)
def right_boundary(x, on_boundary):
    r = math.sqrt(x[0] * x[0] + x[1] * x[1])
    return near(r, 0.7) and x[0] > 0 and between(x[1], (-0.5, 0.5))

leftPlate = DirichletBC(V, u_L, left_boundary)
rightPlate = DirichletBC(V, u_R, right_boundary)
bcs = [leftPlate, rightPlate]

u = TrialFunction(V)
v = TestFunction(V)
a = dot(grad(u), grad(v)) * dx

u = Function(V)
solve(a == Constant(0) * v * dx, u, bcs)

When I run it, I receive *** Warning: Found no facets matching domain for boundary condition.. And the found solution is 0. What is wrong?

  • $\begingroup$ Are you trying to hit boundary facets or facets on radius 0.7 (in unit circle)? $\endgroup$ – Jan Blechta May 24 '13 at 21:35
  • $\begingroup$ Are you trying to hit exterior or interior (on radius 0.7) boundary? If the latter is your intent consider wheter there are facets with their vertices and midpoint on radius cca 0.7 forming connected interior boundary. If yes, you probably need to increase tolerance by near(r, 0.7, tol) as target facets would be rough approximation to circle. $\endgroup$ – Jan Blechta May 24 '13 at 21:43
  • $\begingroup$ Also fit your mesh = UnitCircleMesh(n) where $n = \frac{2}{0.7} m$ both $n$, $m$ being natural. $\endgroup$ – Jan Blechta May 24 '13 at 21:54
  • $\begingroup$ I want to calculate the potential both in the exterior and the interior of the cylinder. $\endgroup$ – facetus May 24 '13 at 22:38
  • $\begingroup$ Then increase tolerance of near function very much and tweak n so that mesh hits better radius 0.7 as suggested above. $\endgroup$ – Jan Blechta May 24 '13 at 22:44

As the error message suggest the DirichletBC does not hit any mesh entities with corresponding dofs. You need to examine your subdomains. One way to debug these are to apply the DirichletBC to a Function and plot the result:

def left_boundary(x, on_boundary):
    r = math.sqrt(x[0] * x[0] + x[1] * x[1])
    return r < 0.5

leftPlate = DirichletBC(V, u_L, left_boundary)
u = Function(V)
u.vector()[:] = 0
plot(u, interactive=True)

This will plot a nice hat. Now you can start tweak your left_boundary function until you reach the result you want.

| cite | improve this answer | |
  • $\begingroup$ Thanks! I'll try this. What is a subdomain (excuse me for a silly question)? $\endgroup$ – facetus May 24 '13 at 22:41
  • 1
    $\begingroup$ @noxmetus: SubDomain is DOLFIN class for definition of some subdomain of your spatial domain, see fenicsproject.org/documentation/dolfin/1.2.0/python/…. Function left_boundary(x, on_boundary) is convenience way of defining subdomains without actually subclassing SubDomain. Notice that left_boundary has same interface as SubDomain.inside(). $\endgroup$ – Jan Blechta May 25 '13 at 0:05
  • 1
    $\begingroup$ Note that you can plot boundary conditions directly without creating a Function first: plot(leftPlate) $\endgroup$ – Anders Logg May 26 '13 at 21:55

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.