2D boundary Submesh from 3D mesh

mesh
1

Hi everyone, I’m trying to create a 2D mesh from a 3D one. The goal is to solve 2DHelmholtz equation on a particular 2D boundary of a global 3D domain.
Here is my code which is not working

import dolfin as dl
class Conn_domain(dl.SubDomain):
def inside(self, x, on_boundary):
test_x = np.isclose(coords_vertex_conn[:,0], x[0]).any()
test_y = np.isclose(coords_vertex_conn[:,1], x[1]).any()
test_z = np.isclose(coords_vertex_conn[:,2], x[2]).any()
test = test_x and test_y and test_z
#if test:
# print(f’coord {x} is on the boundary’) # That is correct, this function returns True is x is on the 2D boundary
return test

    conn_domain = Conn_domain()
    bmesh = dl.BoundaryMesh(mesh, 'exterior', True)
    boundary_markers = dl.MeshFunction('size_t', bmesh, dim = 2)
    boundary_markers.set_all(0)
    conn_domain.mark(boundary_markers, 1)
    mesh_conn = dl.SubMesh(bmesh, boundary_markers,1)
    print('mesh_conn.num_vertices() : ',mesh_conn.num_vertices())
    print('mesh_conn.num_cells() : ',mesh_conn.num_cells())

I expect as output 1030 cells and 578 vertices on this 2D boundary subdomain but I get

mesh_conn.num_vertices() : 0
mesh_conn.num_cells() : 0

2

Please use 3x` encapsulation when wrapping code, such that it renders as

```
def test(x):
    return x[0]
```
3

Ok, actually I solve the issue. The problem was on the definition of inside function. The argument x of this function will be the coordinates of each vertex AND each cell midpoint. So if the surface coordinates are in coords_vertex_conn and cell midpoints coordinates in cells_conn_midpoint, the correct inside function is :
‘’’
class Conn_domain(SubDomain):

        def inside(self, x, on_boundary):
            tol = 1e-4 # Can be less than this value but that works
            testv_x = np.isclose(coords_vertex_conn[:,0], x[0], atol=tol) 
            testv_y = np.isclose(coords_vertex_conn[:,1], x[1], atol=tol)
            testv_z = np.isclose(coords_vertex_conn[:,2], x[2], atol=tol)
            testv_xyz = testv_x * testv_y * testv_z

            # Test on cells midpoint
            testc_x = np.isclose(cells_conn_midpoint[:,0], x[0], atol=tol)
            testc_y = np.isclose(cells_conn_midpoint[:,1], x[1], atol=tol)
            testc_z = np.isclose(cells_conn_midpoint[:,2], x[2], atol=tol)
            testc_xyz = testc_x * testc_y * testc_z
            return (testv_xyz.any()) or (testc_xyz.any())
    '''