Comment 2 for bug 1155271

Revision history for this message
Glen D. Granzow (ggranzow) wrote :

The problem is reproducible in 2D:

from dolfin import *

parameters['allow_extrapolation'] = True

mesh = RectangleMesh(0.0, 0.0, 1.0, 1.0, 10, 10)
V = FunctionSpace(mesh, 'Lagrange', 1)
u = interpolate(Expression('x[0]*x[0]'), V)

index_1 = [x and y for (x,y) in zip(mesh.coordinates()[:,0] == 0.3, mesh.coordinates()[:,1] == 0.5)].index(True)

x_coordinates = interpolate(Expression("x[0]"), V).vector().array()
y_coordinates = interpolate(Expression("x[1]"), V).vector().array()

index_2 = [x and y for (x,y) in zip(x_coordinates == 0.3, y_coordinates == 0.5)].index(True)

for i in range(4):
  x, y = mesh.coordinates()[index_1]
  print 'x =', x, 'and the x coordinates span', min(mesh.coordinates()[:,0]), max(mesh.coordinates()[:,0])
  print 'y =', y, 'and the y coordinates span', min(mesh.coordinates()[:,1]), max(mesh.coordinates()[:,1])
  print 'The point (%g, %g) is in cell %d' % (x, y, mesh.intersected_cell(Point(x,y)))
  u_of_x_and_y = u(x,y) # This statement generates an unexpected warning
  print 'u(x,y) =', u_of_x_and_y, 'when I expected', u.vector().array()[index_2],
  print '(the numbers'+('', ' DO NOT')[abs(u_of_x_and_y - u.vector().array()[index_2])>1.e-5], 'agree)'
  print '--------------------------------------------------------------'
  mesh.coordinates()[:,0] *= 2

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

x = 0.3 and the x coordinates span 0.0 1.0
y = 0.5 and the y coordinates span 0.0 1.0
The point (0.3, 0.5) is in cell 84
u(x,y) = 0.09 when I expected 0.09 (the numbers agree)
--------------------------------------------------------------
x = 0.6 and the x coordinates span 0.0 2.0
y = 0.5 and the y coordinates span 0.0 1.0
The point (0.6, 0.5) is in cell -1
Extrapolating function value at x = <Point x = 0.6 y = 0.5 z = 0> (not inside domain).
u(x,y) = 0.09 when I expected 0.09 (the numbers agree)
--------------------------------------------------------------
x = 1.2 and the x coordinates span 0.0 4.0
y = 0.5 and the y coordinates span 0.0 1.0
The point (1.2, 0.5) is in cell -1
Extrapolating function value at x = <Point x = 1.2 y = 0.5 z = 0> (not inside domain).
u(x,y) = -0.03 when I expected 0.09 (the numbers DO NOT agree)
--------------------------------------------------------------
x = 2.4 and the x coordinates span 0.0 8.0
y = 0.5 and the y coordinates span 0.0 1.0
The point (2.4, 0.5) is in cell -1
Extrapolating function value at x = <Point x = 2.4 y = 0.5 z = 0> (not inside domain).
u(x,y) = -0.33 when I expected 0.09 (the numbers DO NOT agree)
--------------------------------------------------------------

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