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]-0.5)**2 + (x[1]-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.