I’m currently trying to extract a certain 2D portion of the boundary of a 3D mesh. As a minimum working example, consider trying to extract the bottom surface of a unit cube.
from dolfin import * dim = 3 N = 2**5 mesh = UnitCubeMesh(N,N,N) bdim = dim-1 bmesh = BoundaryMesh(mesh, "exterior") mapping = bmesh.entity_map(bdim) part_of_bot = MeshFunction("size_t", bmesh, bdim) for cell in cells(bmesh): curr_facet_normal = Facet(mesh, mapping[cell.index()]).normal() if near(curr_facet_normal.y(), -1.0): #On bot boundary part_of_bot[cell] = 1 bot_boundary = SubMesh(bmesh, part_of_bot, 1) File('bot_boundary.pvd') << bot_boundary
When directly examining the output mesh files, I can see that it preserved the y coordinate (which makes sense). However, I need the 2D mesh, with the zero y coordinate pruned. How might I go about getting rid of that coordinate from bot_boundary, or better yet, is there a way to build bot_boundary (say an option that I can pass to SubMesh) that doesn’t preserve the y coordinate in the first place? I can’t find a relevant function in the docs (Restriction seemed relevant, but that’s been removed since).
Thanks for any input.