If the grid that you’re evaluating basis functions on is (or can be) represented as a DOLFIN Mesh, then you can use the PETScDMCollection functionality to do this. I wrote up a short example of this in the following thread:
Each row of the matrix A gives the evaluations of all basis functions from the space V1 at a node from the space V2.