Interpolation between none-matching meshes in parallel gives zero value at some vertex

Question

Hi,

I am interpolating a function from a refine mesh to its parent mesh in 3D and parallel. (MPI)

I am using the LagrangeInterpolator :

 LagrangeInterpolator  lip;
 lip.interpolate(func_active_strain_mech_proj, func_active_strain_elec);

However, I notice that many vertices in the parent mesh have values of zero after the interpolation, which is far from good...

I am also considering the solutions here
http://fenicsproject.org/qa/3020/projection-interpolation-between-meshes
But the reading and writing could be taking too much time in solving a large problem.

So any ideas? Thanks for your help!

Comments

Do you remember to set parameters["allow_extrapolation"] to False? Otherwise, if you find an error with LagrangeInterpolator you should create a minimal example that reproduces the error and report it as an issue on bitbucket.

Interpolation in both directions from one mesh to another and back is possible with LagrangeInterpolator. I think perhaps that when you are interpolating, some points on the boundary of one mesh could fall outside the domain of the other mesh, just by a very small number due to roundoff. Then the interpolation point is not found inside any cell on the other mesh and the value will be zero. Would be easier to tell if you could provide complete code.

---

Thanks.

1. I did set the parameters["allow_extrapolation"] to False.

2. I am modelling the electrical activity in the atria, of which the structure is fairly irregular. It is highly possible that some of the nodes at the boundary of the parent mesh is outside the child mesh. I moved these nodes of the parent mesh giving zeros after the child mesh was obtained through refining, yet haven't got a better result.

3. I'll try to see if I can do a minimal example. Thanks again!
   I moved the nodes by:

   /*for (int i = 0; i < size_vertices_parent_mesh; ++i)
   {
   Point p(coord_parent_mesh[3 * i], coord_parent_mesh[3 * i + 1], coord_parentmesh[3 * i + 2] );

   bool inside = false;
   CellIterator cell_elec(mesh_child);
   for (; !cell_elec.end(); ++cell_elec)
   {
       if (cell_elec->contains(p))
       {
           inside = true;
           break;
       }
   }
   if (!inside)
   {
       // move the point until it is contained by a cell of the child mesh
   }
   

   }
   */

Answer 1

In theory, all the information should be there, if you use refine() to refine your mesh in parallel, (set redistribute=False). There is some mesh data,

mesh.data().array("parent_cell", mesh.topology().dim()

which will give you the parent cell index for each cell. It should be possible to use that to
get the values you need.

Note that the other way, interpolate from parent to child, will now work in master branch.

Comments

Thanks.

Set redistribute=False does help in removing the zero values interpolating from child to parent, though the working load among processes become unbalanced.

The interpolation of both directions are quite useful in multi-physics simulations. Glad to know this will be potentially working soon.