# function domain is not correct after changing mesh coordinates

I am trying to solve a free boundary problem. To do so I want to adjust the mesh coordinates to represent the domain as it changes. The following code illustrates a problem that I am having though. I define a dolfin Function (I call it "u") on a FunctionSpace with domain [0,1] then change the mesh coordinates so I think that the domain will become [0,2]. But when I try to evaluate my Function at 1, I get an error message saying "The point is not inside the domain."

------------------- python code -------------------

from dolfin import *

from numpy import linspace

mesh = IntervalMesh(10, 0.0, 1.0)

V = FunctionSpace(mesh, 'Lagrange', 1)

u = interpolate(

print 'u(0.5) = ', u(0.5)

print 'u(1.0) = ', u(1.0)

mesh.coordinate

print 'u(0.5) = ', u(0.5)

print 'u(1.0) = ', u(1.0)

------------------- output -------------------

u(0.5) = 0.5

u(1.0) = 1.0

u(0.5) = 0.25

u(1.0) = Evaluating at x = <Point x = 1 y = 0 z = 0>

Traceback (most recent call last):

File "trash.py", line 14, in <module>

print 'u(1.0) = ', u(1.0)

File "/usr/lib/

self.

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:/

***

*** 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_

*** Where: This error was encountered inside Function.cpp.

*** Process: 0

*** -------

------------------- end of output -------------------

If I include the statement

parameters[

as suggested by the error message the program works (with a warning message) but it seems like this should not be necsessary; the point (x=1) is inside the domain.

------------------- output after setting 'allow_

u(0.5) = 0.5

u(1.0) = 1.0

u(0.5) = 0.25

u(1.0) = Extrapolating function value at x = <Point x = 1 y = 0 z = 0> (not inside domain).

0.5

------------------- end of output -------------------

Is this a bug in FEniCS/DOLFIN or am I doing something wrong?

## Question information

- Language:
- English Edit question

- Status:
- Answered

- For:
- DOLFIN Edit question

- Assignee:
- No assignee Edit question

- Last query:
- 2013-03-08

- Last reply:
- 2013-03-23

Glen D. Granzow (ggranzow) said : | #1 |

The following code leads me to believe that my dolfin Function "u" doesn't work like I anticipated. I expected that the last two numbers numbers on each line of output would be the same (namely 0.09) but they are not. In addition I thought that all four plots would have the dashed black curve coinciding with the solid red curve but this is only true for the first two plots. I see that Mesh objects have a move function. Should I be using it instead of trying to assign new values to the mesh.coordinates? How would I do this?

Note that I have set the appropriate parameter to 'allow_

The real point of adding this additional information to my question is that what I originally thought was a minor irritant (that FEniCS thought it was extrapolating when it wasn't) is a more serious problem (the function "u" is not evaluating as expected.)

------------------- python code -------------------

from dolfin import *

from matplotlib import pyplot

parameters[

set_log_

mesh = IntervalMesh(10, 0.0, 1.0)

V = FunctionSpace(mesh, 'Lagrange', 1)

u = interpolate(

pyplot.clf()

for i in range(4):

x = mesh.coordinate

print x, u(x), u.vector(

pyplot.

pyplot.

pyplot.

mesh.

pyplot.show()

------------------- output -------------------

0.3 0.09 0.09

0.6 0.09 0.09

1.2 -0.03 0.09

2.4 -0.33 0.09

------------------- end of output -------------------

Marie Rognes (meg-simula) said : | #2 |

I don't have a good answer for you, but this doesn't seem to work as expected no. (Unless there is some magic with the mesh.coordinates that I'm missing) Could you please consider filing a bug?

One comment: you should not assume that the ordering of u.vector() directly corresponds to the mesh.coordinates() (see other questions regarding ordering of degrees of freedom).

Glen D. Granzow (ggranzow) said : | #3 |

Thank you Marie. As you suggested I created a bug report (#1155271).

Andre Massing (massing) said : | #4 |

Hi! I just gave a detailed answer in the bug report you submitted. Long story short,

use http://

to clear and rebuild the underlying search data structures for the mesh.

Andre Massing (massing) said : | #5 |

Hi! I just gave a detailed answer in the bug report you submitted. Long story short,

use http://

to clear and rebuild the underlying search data structures for the mesh.

## Can you help with this problem?

Provide an answer of your own, or ask Glen D. Granzow for more information if necessary.