Number of degrees of freedom of constrained functionspace dependent on mesh partitioning

When creating a periodic functionspace from a constrained domain in fenics (version 2018.1.0),
I observe that the number of degrees of freedom of this periodic/constrained functionspace is not always the same.
At least, I observe that when the periodic functionspace is created in parallel, its number of degrees of freedom changes seemingly random from one code execution to another.
However, this only happens when the mesh is sufficiently large and distributed over many processes (i.e. when running mpirun -np 72 or mpirun -np 108).
Anyway, I suppose that this kind of behaviour is not desired, and perhaps it might indicate a bug.
Therefore, I would like to ask whether someone can help me figuring out what the cause is.

The output of the following working example illustrates my observations:

#  Import dolfin
import dolfin as df

# Define mesh
n = 15
Nx = 30 
Ny = 30			
Nz = 30 
mesh = df.BoxMesh(df.Point(0.0, 0.0, 0.0), df.Point(n, 1., 1.), n*Nx, Ny , Nz)

near_tol = 1e-7
map_tol = 1e-5

# Create sub domain for periodic boundary condition
class PeriodicDomain(df.SubDomain):

	def __init__(self,map_tol=map_tol):
		df.SubDomain.__init__(self,map_tol=map_tol) # Call base class constructor!

	def inside(self, x, on_boundary):
		return  df.near(x[1], 0,eps=near_tol) and on_boundary 

	def map(self, x, p):
		p[0] = x[0]
		p[1] = x[1] - 1.0
		p[2] = x[2]

# Create periodic/constrained functionspace 	
constrained_domain = PeriodicDomain(map_tol=map_tol)
V = df.FunctionSpace(mesh, 'CG', 1, constrained_domain=constrained_domain)
u0_funct = df.Function(V)

# Print number of degrees of freedom
print("u0_funct.vector().size()", u0_funct.vector().size())

The output of the last print-statement (for 72 processes) is one of the following lines:
u0_funct.vector().size() 419432,
u0_funct.vector().size() 419430, or
u0_funct.vector().size() 419431,
which indicates that the number of degrees of freedom changes each execution time (and hence with the mesh partitioning?).

I remark again that the latter behaviour only is observed for sufficiently large meshes (so n>=15) and for a sufficiently large number of MPI processes (72 for example), when run on our supercomputer.
The fenics environment on our supercomputer has been build from the easybuild scripts of https://github.com/easybuilders/easybuild-easyblocks/blob/master/easybuild/easyblocks/d/dolfin.py .
Furthermore, I have checked that the problem does not dissapear if the geometric tolerances (near_tol and map_tol are changed), so I don’t see how the map or inside functions would fail to detect some degrees of freedom.

4 Likes