Error with LagrangeInterpolator in parallel

Hi,

I construct a mesh with two subdomains. I create two meshviews from theses subdomains and I want to interpolate a function from one meshview to the other. Using LagrangeInterpolator the code runs fine in serial but I got the following error when running the code in parallel (with mpirun -np 5):

Traceback (most recent call last):
  File "interpolation_meshview.py", line 42, in <module>
    LagrangeInterpolator.interpolate(u2, u1)
RuntimeError: 

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
***     fenics-support@googlegroups.com
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to initialize mapping of degrees of freedom.
*** Reason:  Missing entities of dimension 0. Try calling mesh.init(0).
*** Where:   This error was encountered inside DofMap.cpp.
*** Process: 4
*** 
*** DOLFIN version: 2019.2.0.dev0
*** Git changeset:  f0a90b8fffc958ed612484cd3465e59994b695f9
*** -------------------------------------------------------------------------


Even more strange, when I run the code with 3 processors, I get no error but the resulting interpolated function is wrong…

Here is a minimal example:

from dolfin import *

# Parameters
N = 30
xc, yc = 0.5, 0.5
rc = 0.12

# Construct the mesh and mark subdomains (a disk in a square)
mesh = UnitSquareMesh(N,N)
subdomains = MeshFunction("size_t", mesh, 2, 1)
for cell in cells(mesh):
    if (cell.midpoint()[0]-xc)**2 + (cell.midpoint()[1]-yc)**2 < rc**2:
        subdomains.set_value(cell.index(),2)

# create meshviews
submesh1 = MeshView.create(subdomains, 1)
submesh2 = MeshView.create(subdomains, 2)

# Define functon spaces, functions and expression to interpolate
V1 = FunctionSpace(submesh1, "CG", 1)
V2 = FunctionSpace(submesh2, "CG", 1)
u1 = Function(V1)
u2 = Function(V2)
expr = Expression("x[0]*x[1] + 2*x[1]", degree=1)
LagrangeInterpolator.interpolate(u1, expr)

# Interpolation with LagrangeInterpolator
LagrangeInterpolator.interpolate(u2, u1)

# Save results
with XDMFFile(MPI.comm_world, "u1.xdmf") as file:
    file.write(u1)
with XDMFFile(MPI.comm_world, "u2.xdmf") as file:
    file.write(u2)

Any help would be appreciated.

Cheers,
Fabien

Hi Fabien,

I don’t have access to 2019.2.0.dev0, and (apparently) therefore don’t have access to MeshView.

However, it seems to me that you should not be able to interpolate u1 into u2, since these two functions are defined on non-intersecting domains. If you replace

LagrangeInterpolator.interpolate(u2, u1)

with

LagrangeInterpolator.interpolate(u2, expr)

does this remove the error?

Hi @conpierce8 ,

Thank you for your answer.

According to the documentation LagrangeInterpolator is for interpolation between non-matching meshes. And the code runs perfectly with only one processor, so I guess it should also work in parallel…

If I do what you suggest, It removes the error, but it does not solve the interpolation issue between u1 and u2.