I cant see anything wrong with your code. However, as i do not have the meshes, i cant really test it.
First suggestion is to not use dolfin-convert, as it is no longer maintained, see: Transitioning from mesh.xml to mesh.xdmf, from dolfin-convert to meshio
for the appropriate syntax using the latest version of meshio. How to read it into dolfin is found in earlier posts of that topic.
Secondly, check that you are integrering over the correct boundaries using a mesh file.
One way to do that is to visualize the mesh function in paraview. The second is to Calculate the area of each neumann segment, printing assemble(1*ds(subdomain_id=i, domain=mesh, subdomain_data=boundaries))
1 Like