Interpolation (problem with current development version)

Question

Hello,

since recently I have a problem interpolating functions in nested meshes. The following example code triggers this problem:

from dolfin import *

elem = "CG"
deg = 1

mesh1 = UnitSquareMesh(10, 10)
mesh2 = refine(refine(mesh1))

space1 = FunctionSpace(mesh1, elem, deg)
space2 = FunctionSpace(mesh2, elem, deg)

fun = Function(space1)

interpolate(fun, space2)

This results in:

Traceback (most recent call last):
  File "/home/neumann/software/fenics/FEniCS/test.py", line 17, in <module>
    interpolate(fun, space2)
  File "/usr/lib/python2.7/dist-packages/dolfin/fem/interpolation.py", line 64, in interpolate
    Pv.interpolate(v)
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
***
***     https://answers.launchpad.net/dolfin
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to evaluate function at point.
*** Reason:  The point is not inside the domain. Consider setting "allow_extrapolation" to allow extrapolation.
*** Where:   This error was encountered inside Function.cpp.
*** Process: 0
*** -------------------------------------------------------------------------

Does anybody else have this problem.

parameters['allow_extrapolation'] = True # TODO: Bug Workaround?!?

solves the problem and the results seem to be ok so it might just be a corner case which is triggered. Probibly it has to do with rounding errrors (I have a 64-Bit machine) and was not happening until now. I used this function in this manner for a long time now. So the error was introduced not long ago.

I tested with the Ubuntu Version as well as PPA and PPA-dev.

Thanks,
Jo

My System:
Ubuntu 13.10 64-Bit, 8GB RAM, Intel Core i7

Answer 1

This works for me with development version from bitbucket built from source.

There were quite a few bugs related to introducing a new mesh intersection code that should now, hopefully, have been taken care of. Try to install Fenics from source, otherwise I think you just have to wait for the fix to enter the binaries. For your case, where the points probably miss intersection by close to machine precision, I guess allow_extrapolation will do the trick.

Comments

Thank you very much for the information. My further testing seems to support your guess. allow_extrapolation does not give wrong results and does not impact performance to much. The downside of this workaround is the missing check for domain bounds on nonrelated meshes obviously.