### See a field only on a complex boundary

thanks in advance for the help.

I want to simulate a flow around a complex 3D geometry, in particular a car mirror. I solved a the Navier-Stokes equation on a complex 3D geometry:

and I get the solution:

In my Fenics code I defined the boundary condition as follow:

#set boundaries

class Inlet(SubDomain):

def inside(self, x, on_boundary):

return on_boundary and abs(x[0]-x_min) < DOLFIN_EPS

class Outlet(SubDomain):

def inside(self, x, on_boundary):

return on_boundary and abs(x[0] - x_max) < DOLFIN_EPS

class Walls(SubDomain):

def inside(self, x, on_boundary):

return on_boundary and \

(abs(x[1] - y_max) < DOLFIN_EPS or abs(x[2] - z_min) < DOLFIN_EPS or abs(x[2] - z_max) < DOLFIN_EPS)

class AllBoundary(SubDomain):

def inside(self, x, on_boundary):

return on_boundary

boundaries = FacetFunction("size_t", mesh)

boundaries.set_all(0)

circle_ID = 4 # this will mark all the boundary, but it will be overwritten later

circle = AllBoundary()

circle.mark(boundaries, circle_ID)

walls_ID = 1

walls = Walls()

walls.mark(boundaries, walls_ID)

outlet_ID = 2

outlet = Outlet()

outlet.mark(boundaries, outlet_ID)

inlet_ID = 3

inlet = Inlet()

inlet.mark(boundaries, inlet_ID)

Now I want to see into paraview the pressure distribution only on the object (the mirror). The paraview commands don't seem to allow this operation. Is there a way to make this kind of plot?

Thank you very much for the help.

M Pesarin

### 2 Answers

```
from fenics import *
mesh = UnitCubeMesh(10, 10, 10)
mesh_markers = FacetFunction("size_t", mesh, 0)
bottom = CompiledSubDomain('on_boundary && near(x[2], 0.0)')
bottom.mark(mesh_markers, 4)
bmesh = BoundaryMesh(mesh, 'exterior')
cellmap = bmesh.entity_map(2)
bmesh_markers = CellFunction("size_t", bmesh, 0)
for c in cells(bmesh):
if mesh_markers[cellmap[c.index()]] == 4:
bmesh_markers[c] = 1
smesh = SubMesh(bmesh, bmesh_markers, 1)
plot(smesh, interactive=True)
```

If you have the degrees of freedom marked on the mirror, you can extract the boundary mesh using "fenics::BoundaryMesh()". With this in hand, create a P1 FunctionSpace defined over this mesh, then use "fenics::LagrangeInterpolator()" to interpolate the 3D data onto the 2D mesh. Then just save the 2D Function:

```
from fenics import *
mesh = UnitCubeMesh(20,20,20)
Q = FunctionSpace(mesh, 'CG', 1)
code = 'pow(x[0], 2) + pow(x[1], 2) + pow(x[2], 2)'
u_from = interpolate(Expression(code, degree=1), Q)
bmesh = BoundaryMesh(mesh, 'exterior')
cellmap = bmesh.entity_map(2)
pb = CellFunction("size_t", bmesh, 0)
for c in cells(bmesh):
if Facet(mesh, cellmap[c.index()]).normal().z() > 1e-3:
pb[c] = 1
submesh = SubMesh(bmesh, pb, 1)
Q_bmsh = FunctionSpace(bmesh, 'CG', 1)
u_to = Function(Q_bmsh)
lg = LagrangeInterpolator()
lg.interpolate(u_to, u_from)
File('u_to.pvd') << u_to
```

Here, just the top of the domain is extracted for example. You can use your own FacetFunction or whatever else to define "bmsh".
We solved successfully the problem.

M Pesarin

Did Adam's answer solve your problem? Cause i am facing the same issue... i am trying to extract the values of my solution on different subdomains in fenics. Unfortunately submesh is not working in parallel yet... and despite the fact that Adam's answer seems really helpful, it is not convenient when solving in parallel. Thank you all for your effort anyways.