How to interpolate between non-matching meshes in parallel?

I need to compute the harmonic extension u2 of a field u1 (defined in a domain Omega_1) on a domain Omega_2, sharing a common portion of the boundary. By calling u1.set_allow_extrapolation(True) and setting the BC as DirichletBC(W, u1, bdry()) it perfectly works (in serial).

However, if I run the code in parallel, even if the code end without errors, the results are bad: the Dirichlet condition is wrongly assigned on a portion of the boundary.

I’ve found old posts online saying that interpolation between different meshes doesn’t work in parallel. Has this been fixed in the latest version? Is there any workaround?

I guess that doing ALE for FSI the same problem shows up. How is it solved in that case?

Thank you!

A workaround is using FEniCSTools: https://github.com/mikaem/fenicstools/wiki/Interpolation-nonmatching-mesh-in-parallel

1 Like

Hi,

Sorry for bumping this thread. I have a similar problem with interpolation between two meshes. I am solving Poisson and Laplace equations defined on two meshes with matching boundaries (problem with two subdomains). Since subdomains do not work in parallel, I tried to find the alternative method to solve this problem. I found this topic, so I tried to use FEniCSTools. I installed it without the errors, but unfortunately, I cannot make it run. Namely, when I try to run a simple example code:

from dolfin import *
from mshr import *
from fenicstools import *
mesh1 = UnitSquareMesh(16, 16)
V1 = FunctionSpace(mesh1, 'CG', 2)
u1 = interpolate(Expression("sin(pi*x[0])*cos(pi*x[1])", degree = 1), V1)
# Create a new _different_ mesh and FunctionSpace
mesh2 = UnitSquareMesh(10, 10)
x = mesh2.coordinates()
x[:, :] = x[:, :] * 0.5 + 0.25
V2 = FunctionSpace(mesh2, 'CG', 1)
interpolate_nonmatching_mesh(u1, V2)

it reports the following error: name ‘interpolate_nonmatching_mesh’ is not defined. Code was ran on a machine with installed FEniCS 2019.1.0 version.

On FEniCSTools github page I found some similar issues reported. I tried suggested actions (installing cppimport, h5py and pyvtk), but nothing changed. Any suggestion how to fix this error would be helpful. Also, is there an alternative for interpolation between meshes in parallel?

You may be able to use PETScDMCollection from DOLFIN for this, without FEniCSTools. See the following modification of your example:

from dolfin import *

#from mshr import *
#from fenicstools import *
def transfer_function(fromFunc,toFunc):
    fromFunc.set_allow_extrapolation(True)
    A = PETScDMCollection.create_transfer_matrix\
        (fromFunc.ufl_function_space(),
         toFunc.ufl_function_space())
    toFunc.vector()[:] = A*fromFunc.vector()

mesh1 = UnitSquareMesh(16, 16)
V1 = FunctionSpace(mesh1, 'CG', 2)
u1 = interpolate(Expression("sin(pi*x[0])*cos(pi*x[1])", degree = 1), V1)
# Create a new _different_ mesh and FunctionSpace
mesh2 = UnitSquareMesh(10, 10)
x = mesh2.coordinates()
x[:, :] = x[:, :] * 0.5 + 0.25
V2 = FunctionSpace(mesh2, 'CG', 1)

#interpolate_nonmatching_mesh(u1, V2)
u2 = Function(V2)
transfer_function(u1,u2)

from matplotlib import pyplot as plt
plot(u1)
plt.show()
plot(u2)
plt.show()

(Depending on what you are doing, it may make sense to pre-compute and store the matrix A that is created in transfer_function instead of regenerating it with each call to transfer_function.)

3 Likes

Thank you for a quick reply and giving me an idea how to resolve the problem. It works for a given example. I’ll try to adapt it for my original problem.

Hi all.

I have tried to use the code proposed by Kamensky but I have a doubt. I have got an aproximated solution UP of Navier-Stokes equation in a fine mesh and I am trying to interpolate to a course mesh. So, I defined

up = Function(R)
UP_ref= transfer_function(UP, ur) 

In what follows, I need to do a split
(u,p)=UP_ref.split()
but I get AttributeError: ‘NoneType’ object has no attribute ‘split’. Is there any solution for that?

As you can see from Kamensky’s snippet, what transfer_function(u1, u2) does is that it interpolates the values of u1 into u2. therefore, you should therefore split ur, as it is you output (the transfer_function does not return any output)

Thanks. I will try!!