Hello dear FEniCS users,

I am solving a mixed formulation using DG0 and BDM1 mixed function spaces. Moreover, I have the boundary tags for different boundaries of my domain. Everything works fine and I did get the solution fields.

Now I want to extract/project solution fields only at certain boundaries of my domain in order to get the normal and tangential components of the solution fields at these boundaries. I tried to find a way of doing it on FEniCS forum, but could not find the exact post.

Any help in this regard will be highly appreciated! Below I am attaching a minimal working example of my code

```
from dolfin import *
mesh = Mesh("domain.xml")
domain = MeshFunction("size_t", mesh, "domain_physical_region.xml")
boundaries = MeshFunction("size_t", mesh, "domain_facet_region.xml")
bdm_elements = FiniteElement("BDM", mesh.ufl_cell(), 1)
dg_elements = FiniteElement("DG", mesh.ufl_cell(), 0)
W = FunctionSpace(mesh, bdm_elements * dg_elements)
(u, p) = TrialFunctions(W)
(v, q) = TestFunctions(W)
lhs = alpha*inner(u,v)*dx - p*div(v)*dx - q*div(u)*dx
rhs = sum([-coeff_1*inner(n,v)*ds(i+1) for i in range(6)]) - coeff_2*inner(n,v)*ds(7)
bc = DirichletBC(W.sub(0), G, boundaries, 8)
usol = Function(W)
solve(lhs == rhs, usol, bc)
(vsol,psol) = usol.split()
```

Now let’s say I want to extract the solution at boundaries 1-6 only. What should be the procedure?

Thanks and cheers!