Extract coordinates of the nodes of a subboundary

That is not an error, that is simply from the fact that a CG-2 space has degrees of freedom on facets as well as vertices. Therefore you cannot use vertex_to_dof_map for this operation. You need a more complex operation for this:

from dolfin import *
import numpy as np

class boundary(SubDomain):
    def inside(self, x, on_boundary):
        return x[0] <= 0.5 and on_boundary


mesh = UnitSquareMesh(5,5)
V = FunctionSpace(mesh, "CG", 2)
v = project(Expression("x[0]+x[1]", degree=2),V)
mf = MeshFunction("size_t", mesh, mesh.topology().dim()-1, 0)
boundary().mark(mf, 1)

mesh.init(2,1)
dofs = []
cell_to_facets = mesh.topology()(2,1)
for cell in cells(mesh):
    facets = cell_to_facets(cell.index())
    for facet in facets:
        if mf[facet] == 1:
            dofs_ = V.dofmap().entity_closure_dofs(mesh, 1, [facet])
            for dof in dofs_:
                dofs.append(dof)
1 Like