I feel like I am going a bit crazy here. Pretty simple set-up, I have a unit square mesh, and I would like to define a circle subdomain. The code I’ve used is below.
from dolfin import * import numpy as np msh = UnitSquareMesh(500, 500) def circle(x): return (x-0.5)**2 + (x-0.5)**2 class Circle(SubDomain): def inside(self, x, on_boundary): return near(circle(x) , 0.2**2) domains = MeshFunction("size_t", msh, msh.topology().dim()-1, 0) circ = Circle() circ.mark(domains,1)
However, the circle is still marked 0 when I try and plot domains. Refining the mesh does not work and neither does increasing or decreasing the tolerance in the
near function. I cannot think of any reason this does not work. If I try and set Circle to return
circle(x) > 0.2**2 or
circle(x) < 0.2**2, it works perfects fine, so I really do not understand why this doesn’t work.
Any help would be fantastic.