Proper coordinates for eval()-method


I’m fairly new to fenics and I need to evaluate a Function-object on specific points of the mesh. For example, take the non-linear poisson example (code see there). What I would now like to do is to specify the points, at which I want to evaluate the solution, and basically get the values at that point. After reading some other posts I came with the following:

dkiCoordX = np.linspace(0,1, 100)
dkiCoordY = np.linspace(0,1, 100)

# basically the cartesian product of the coordinate vectors, format 10000x2
dkiGrid = np.array(np.meshgrid(dkiCoordX, dkiCoordY)).T.reshape(-1,2) 

values1 = np.ones(10000)
values2 = np.ones(10000)
u.eval(values1, dkiGrid.flatten())
u.eval(values2, mesh.coordinates().flatten())

Apparently, this is not what eval() is expecting from me. Though I get no errors, the only thing that happens to the values vectors is that the first entry changes, everything else stays 1. I’ve also tried a lower number (5x5), but it still doesn’t work.

I would be thankful for any hints and, of course, solutions. It doens’t even have to be the eval()-method, I really just need to evaluate the solution at given points.


value = u(Point(x,y,z))


I’m sorry I should’ve been more specific; it does indeed work, just as

value = u(x,y)

but I need a fast way to evaluate u at many points. I really don’t want to loop over the coordinates in plain python and evaluate u pointwise, that would just take way too much time.

There seems to be no other way:

Ah yes, I’ve already stumbled upon that post. I hoped that maybe in Fenics 2018 there would be some way. And fenicstools seems to be only offering a Python2 interface… At least that’s what I got when I tried it out, I didn’t spend that much time on it, though.

Anyway, thanks for your answers. I’ve written myself small wrapper function via loops now… But it looks like the performance bottleneck actually is not even the loop, but the evaluation process in itself. It still works surprisingly fast.