I am implementing a boundary condition that is based on precomputing quantities
based on geometric information on the boundary (2D domain, 1D boundary).
I want to implement this by extracting a submesh of the boundary, which I can then do computations on. I have tried this post: Extract parts of the boundary
which is written for 3D. Translated to 2D, my code looks like this:
import dolfin as dl import mshr as ms domain = ms.Rectangle(dl.Point(0,0),dl.Point(10,10)) sphere = ms.Circle(dl.Point(5, 5), 3) domain = domain - sphere mesh = ms.generate_mesh(domain,10) exterior = dl.BoundaryMesh(mesh, 'exterior') from numpy import array from numpy.linalg import norm def inSphere(x): v = x - array([5,5]) return norm(v)<3+1e2*dl.DOLFIN_EPS class SphereDomain(dl.SubDomain): def inside(self,x,on_boundary): return inSphere(x) class BoxDomain(dl.SubDomain): def inside(self,x,on_boundary): return not inSphere(x) print(dl.__version__) dl.plot(dl.SubMesh(exterior, SphereDomain()))
With the output:
(This is dolfin v 2019.1.0)
But as you see, the resulting mesh is not connected in the way you’d expect. Any help on this would be appreciated. How do I make the resulting submesh agree with the boundary?
Thanks to user dokken, I now know that the mesh is connected, it’s the plotting that misrepresents the mesh structure. The following manual plot:
bmesh = dl.SubMesh(exterior, SphereDomain()) for f in dl.cells(bmesh): pts = [v.point() for v in dl.vertices(f)] plt.plot([pts.x(), pts.x()], [pts.y(), pts.y()], 'r.-')
Shows a correct visualisation of the boundary mesh.