evaluate_basis_all() gives incorrect result unless passing xy.copy()

[Nico Schlömer]

I’m having a hard time evaluating basis functions here.

import numpy
from dolfin import (
    BoundingBoxTree, Cell, FunctionSpace, Mesh, MeshEditor, Point
)

import meshzoo


def create_dolfin_mesh(points, cells):
    # https://bitbucket.org/fenics-project/dolfin/issues/845/initialize-mesh-from-vertices
    editor = MeshEditor()
    mesh = Mesh()
    editor.open(mesh, "triangle", 2, 2)
    editor.init_vertices(points.shape[0])
    editor.init_cells(cells.shape[0])
    for k, point in enumerate(points):
        editor.add_vertex(k, point)
    for k, cell in enumerate(cells):
        editor.add_cell(k, cell)
    editor.close()
    return mesh


points, cells = meshzoo.triangle(
    n=8, corners=numpy.array([[0.0, 0.0], [1.0, 0.0], [0.0, 1.0]])
)
mesh = create_dolfin_mesh(points, cells)
V = FunctionSpace(mesh, "CG", 1)

bbt = BoundingBoxTree()
bbt.build(mesh)
el = V.element()

XY = numpy.array([[0.38764254, 0.29652301]])
print(XY.dtype)

for i, xy in enumerate(XY):
    cell_id = bbt.compute_first_entity_collision(Point(*xy))
    cell = Cell(mesh, cell_id)
    coordinate_dofs = cell.get_vertex_coordinates()

    print(el.signature())
    print(xy, coordinate_dofs, cell_id)
    print(type(xy))
    v = el.evaluate_basis_all(xy, coordinate_dofs, cell_id)

    print(v)
    exit(1)

This code correctly produces

float64
FiniteElement('Lagrange', triangle, 1)
[0.38764254 0.29652301] [0.375, 0.25, 0.5, 0.25, 0.375, 0.375] 31
<class 'numpy.ndarray'>
[0.5266756  0.10114032 0.37218408]

here. When I run it in larger production code, however, I’m getting

 [-0.72790868  0.10114035  1.62676833]

In the last line. This error can be avoided by passing xy.copy() to evaluate_basis_all!

Any idea what memory mix-up could be going on here?