introduce PointEvaluator object

[leo-collins]

PointEvaluator is a convenience object for evaluating a Function at a set of points. The idea is for the user to avoid the boilerplate when point evaluating using a VertexOnlyMesh. Now that #4484 is merged we should be able to deprecate .at() in favour of this API.

The basic usage is:

evaluator = PointEvaluator(mesh, points)
function_at_points = evaluator.evaluate(function)

The kwargs you can pass to PointEvaluator are tolerance and missing_points_behaviour which are self-explanatory, and redundant. If redundant=True (the default) then only the points on rank 0 are evaluated, and we broadcast the result to all other ranks. If redundant=False then each rank evaluates the points it has been given and returns the results on that rank. This is useful for e.g. external data input.

[connorjward]

Awesome.

Could you amend the manual to recommend doing this instead of the full VoM ritual?

[pbrubeck]

This is just a suggestion that I think might align better with the rest of firedrake. How about something like this?

def VertexOnlyFunctionSpace(mesh, points):
    # TODO some caching
    return FunctionSpace(VertexOnlyMesh(mesh, points), "DG", 0)

DG_vom = VertexOnlyFunctionSpace(mesh, points)
f_at_points = assemble(interpolate(f, DG_vom))

[connorjward]

This is just a suggestion that I think might align better with the rest of firedrake. How about something like this?

def VertexOnlyFunctionSpace(mesh, points):
    # TODO some caching
    return FunctionSpace(VertexOnlyMesh(mesh, points), "DG", 0)

DG_vom = VertexOnlyFunctionSpace(mesh, points)
f_at_points = assemble(interpolate(f, DG_vom))

Personally I really dislike the functions masquerading as classes approach. But otherwise this seems elegant.

We could have something like

DG_vom = FunctionSpace.at_points(mesh, points)

[leo-collins]

This is just a suggestion that I think might align better with the rest of firedrake. How about something like this?

def VertexOnlyFunctionSpace(mesh, points):
    # TODO some caching
    return FunctionSpace(VertexOnlyMesh(mesh, points), "DG", 0)

DG_vom = VertexOnlyFunctionSpace(mesh, points)
f_at_points = assemble(interpolate(f, DG_vom))

Problem with this is if you want to points in order (which is what users expect) then you'd still need to create a DG_vom_input_ordering and interpolate f_at_points into that. In parallel you also end up with all the results on rank 0. If you want them on the ranks where the points are, you'd need to broadcast them, or set redundant=False when constructing the vom, but then you end up doing the same calculation on every rank.

[pbrubeck]

Problem with this is if you want to points in order (which is what users expect) then you'd still need to create a DG_vom_input_ordering and interpolate f_at_points into that. In parallel you also end up with all the results on rank 0. If you want them on the ranks where the points are, you'd need to broadcast them, or set redundant=False when constructing the vom, but then you end up doing the same calculation on every rank.

Seems like you want a VertexOnlyFunctionSpace(mesh, points, preserve_input_ordering=True) and let it compose the two interpolations

[dham]

need to record coordinate dat version and reimmerse the VOM when needed.

Also need a moving coordinate test.