Interpolation nonmatching mesh in parallel

Hi,

I have some code I have been using for a while that is part of fenicstools and that has been requested by several. The code is basically a stand-alone function that interpolates a dolfin Function from one mesh to a Function with the same FunctionSpace defined on another mesh.

mesh0 = UnitSquareMesh(10, 10)
mesh1 = UnitSquareMesh(20, 20)
V0 = FunctionSpace(mesh0, "CG", 1)
V1 = FunctionSpace(mesh1, "CG", 1)
u0 = interpolate(Expression("x[0]*x[1]"), V0)
u0.update()
u1 = interpolate_nonmatching_mesh(u0, V1)

The new feature is that it works in parallel. I think this functionality should be part of dolfin. It could be added to the Function class.

I attach a cpp file that contains the implementation. It can be compiled using instant, or compile_extension_module. A complete test is

from dolfin import *

interpolation_code = open('interpolation_nonmatching.cpp', 'r').read()
compiled_fem_module = compile_extension_module(code=interpolation_code)

def interpolate_nonmatching_mesh(u0, V):
    """Interpolate from GenericFunction u0 to FunctionSpace V.

    The FunctionSpace V can have a different mesh than that of u0, if u0 
    has a mesh.

    """
    u = Function(V)
    compiled_fem_module.interpolate_nonmatching_mesh(u0, u)
    return u

mesh0 = UnitSquareMesh(10, 10)
mesh1 = UnitSquareMesh(20, 20)
V0 = FunctionSpace(mesh0, "CG", 1)
V1 = FunctionSpace(mesh1, "CG", 1)
u0 = interpolate(Expression("x[0]*x[1]"), V0)
u0.update()
u1 = interpolate_nonmatching_mesh(u0, V1)
plot(u0, title="mesh0")
plot(u1, title="mesh1")
Q0 = FunctionSpace(mesh0, "CR", 1)
Q1 = FunctionSpace(mesh1, "CR", 1)
Mixed0 = V0*Q0
Mixed1 = V1*Q1
m0 = interpolate(Expression(("x[0]", "x[1]")), Mixed0)
m0.update()
m1 = interpolate_nonmatching_mesh(m0, Mixed1)
plot(m0[1], title="mesh0")
plot(m1[1], title="mesh1")
interactive()

Comments (10)