Fenics create a rectangle mesh with markers

fenics-create-rectangle-mesh-with-markers.py

Lx, Ly = 1., .1
ny = 5
nx = int(ny*Lx/Ly)
mesh = dolfin.RectangleMesh(dolfin.Point(0, 0), dolfin.Point(Lx, Ly), nx, ny)

# Mark boundary subdomains
left = dolfin.CompiledSubDomain("near(x[0],0) && on_boundary")
right = dolfin.CompiledSubDomain("near(x[0],Lx) && on_boundary",Lx=Lx)
bottom = dolfin.CompiledSubDomain("near(x[1],0) && on_boundary")
top = dolfin.CompiledSubDomain("near(x[1],Ly) && on_boundary",Ly=Ly)

# Mark facets for Neumann BCs
boundary_indices = {'left': 1,
                    'right': 2,
                    'top': 3,
                    'bottom': 4}
boundary_markers = dolfin.MeshFunction("size_t", mesh, dim=1, value=0)
left.mark(boundary_markers, boundary_indices["left"])
right.mark(boundary_markers, boundary_indices["right"])
top.mark(boundary_markers, boundary_indices["top"])
bottom.mark(boundary_markers, boundary_indices["bottom"])

# Redefine element of area to include information about the markers
ds = dolfin.ds(domain=mesh, subdomain_data=boundary_markers)
dx = dolfin.dx(domain=mesh)
dolfin.plot(mesh)
with dolfin.XDMFFile("mesh_functions.xdmf") as f:
    f.write(boundary_markers)

# A simple test
for i in range(1,5):
    length = dolfin.assemble(1.*ds(i))
    print("side {:d} is of length {:4.2f}".format(i,length))

A simple test

for i in range(4): -> this should be range(1,5):

As boundary markers are 1 to 4.