Hello,
I was wondering if there are routines in NGSolve to evaluate a FE-Function based on a higher order polynomial space (like the P2 in the Taylor-Hood Element) on the nodes of a mesh? I know, there is a VTK output routine, but this writes Files and I would like to evaluate the function in python and use it then as numpy-array for something else.

Hi, yes you can create a numpy array of evaluation - points by giving an array of coordinates to the mesh and then put this array into the coefficientfunction to evaluate it:

Yes, a fe solution GridFunction is derived from CoefficientFunction, so allows to evaluate itself like a CoefficientFunction. CoefficientFunction is any function on the mesh, like for example the coordinate-CFs x, y and z.

Also grad(u) is a coefficientfunction, also x * sin(u), …

All of these function can be evaluated on points in the mesh.