This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

Error with interpolate after using ALE.move

+1 vote

Hello,

I have two contiguous meshes and I would like to interpolate a function from one mesh to another, after having moved them with ALE.move. However, I get the following error using interpolate :

Error: Unable to create mesh entity.
Reason: Mesh entity index -1 out of range [0, 211] for entity of dimension 2.
Where: This error was encountered inside MeshEntity.cpp.
Process: 0

Here is a minimal exemple :

from dolfin import *
from mshr import *

N = 10 
xmin,xmax = 0.,1.
ymin,ymax = 0.,1.
xc = 0.5
length = 0.6
radius = 0.1

domain1_vertices = [Point(xmin,ymin), Point(xc-radius,ymin), Point(xc-radius,length), Point(xc+radius,length), Point(xc+radius,ymin),Point(xmax,ymin), Point(xmax,ymax), Point(xmin,ymax)]
domain1 = Polygon(domain1_vertices)

domain2_vertices = [Point(xc-radius,ymin), Point(xc+radius,ymin), Point(xc+radius,length), Point(xc-radius,length)]
domain2 = Polygon(domain2_vertices)

mesh1 = generate_mesh(domain1,N)
mesh2 = generate_mesh(domain2,N)

V2 = VectorFunctionSpace(mesh2,"Lagrange",1)
V1 = VectorFunctionSpace(mesh1,"Lagrange",1)

f = Expression(("scale*x[0]","scale*x[1]"),scale=1.,degree=1)
f1 = project(f,V1)
f1.set_allow_extrapolation(True)

d = Expression(("scale*x[1]","0."),scale=0.1,degree=1)
d2 = project(d,V2)
d1 = project(d,V1)
ALE.move(mesh2,d2)
ALE.move(mesh1,d1)

f2 = interpolate(f1,V2) 

The strange thing is that, if I define the expression f and project it directly on the moved mesh, it works fine.

Thank you for your help.

Best regards,
Fabien

asked Mar 16, 2017 by Fvergnet FEniCS Novice (200 points)

It seems that adding

mesh1.bounding_box_tree().build(mesh1)
mesh2.bounding_box_tree().build(mesh2)

works.

Does anybody can confirm and/or explain why?

Regards,
Fabien

Exactly what I needed +1! Unfortunately, I do not know why.

...