Hi,
Consider the following code:
from dolfin import *
mesh = UnitIntervalMesh(100)
inside = (0.9,)
V = FunctionSpace(mesh, "CG", 1)
f = interpolate(Constant(1), V)
print "Value of f at %s: %s" % (inside, f(inside))
On one process, this works fine:
[pef@aislinn:/tmp]$ python eval.py
Computed bounding box tree with 199 nodes for 100 entities.
Value of f at (0.9,): 1.0
But on two processes, one process succeeds and the other fails:
[pef@aislinn:/tmp]$ mpiexec -n 2 python eval.py
Process 0: Number of global vertices: 101
Process 0: Number of global cells: 100
Process 1: Computed bounding box tree with 99 nodes for 50 entities.
Process 0: Computed bounding box tree with 99 nodes for 50 entities.
Value of f at (0.9,): 1.0
Traceback (most recent call last):
File "eval.py", line 9, in <module>
print "Value of f at %s: %s" % (inside, f(inside))
File "/usr/lib/python2.7/dist-packages/dolfin/functions/function.py", line 587, in __call__
self.eval(values, x)
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
***
*** fenics@fenicsproject.org
***
*** 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: 1
*** -------------------------------------------------------------------------
My question: what's the recommended way of handling this? Is there any parallel-safe way of evaluating functions at points? (I know I could check for this case manually with MPI, but I'm hoping there's a cleaner way than making the user take care of it.) If there's nothing currently in the code, what do you think of a parallel=True flag to function evaluation that would check if it belonged to any process, not just this one?
