evaluate() silently returns wrong values in parallel for points the rank does not own — including mesh nodes

[lmoresi]

In an MPI run, uw.function.evaluate silently returns wrong values for any
query point the calling rank does not own — including points that ARE mesh
points. Serial is exact, global_evaluate is exact, and nothing in evaluate's
signature, docstring or return value says the answer is rank-local.

Measured

UnstructuredSimplexBox, cellSize=0.2, a P2 MeshVariable holding
x^2 + 2 y^2 — exact in P2, so any error is location, not approximation. The
query set is every rank's own P2 DOF coordinates, allgathered, so every
query point is a node of the mesh. True values span [0, 3].

np	sampler	query set	extrapolated	max abs error
1	evaluate	all coords	0 / 153	8.9e-16
2	evaluate	all coords	69 / 166, 71 / 166	1.48, 0.54
4	evaluate	all coords	117-126 / 180	up to 2.59
2	evaluate	own coords	0 / 84	6.7e-16
4	evaluate	own coords	0 / 46	8.9e-16
2	global_evaluate	all coords	0 / 166	8.9e-16
4	global_evaluate	all coords	0 / 180	8.9e-16

Two thirds of the points are wrong at np=4, by up to 2.59 on a field whose whole
range is 3. The ranks also disagree with each other — 1.48 on rank 0 against
0.54 on rank 1 at np=2, for the identical query set — so the answer depends on
who asked.

Why

evaluate classifies points against the rank-local mesh. A point owned by
another rank is not in any local cell, so it is extrapolated from the nearest
local cell, which may be most of a domain away. The extrapolation flag is set
correctly — check_extrapolated=True reports every one of them — but the flag is
off by default and evaluate returns bare values.

This is not the curved-boundary case in #605, where the point is genuinely
outside the meshed domain and extrapolation is the only option. Here the point is
inside the domain, is a node of the mesh, and another rank holds the exact value.

Why it matters more than it looks

The natural thing to write is the thing that breaks:

line = np.column_stack([np.linspace(0, 1, 200), np.full(200, 0.5)])
profile = uw.function.evaluate(T.sym[0], line)      # every rank asks for all 200

In serial that is correct. On four ranks each rank gets a different, mostly
wrong profile, with no error, no warning and no NaN. Anything that samples a
fixed set of points — a profile, a probe location, a set of observation stations,
a comparison against an analytic solution — has this shape.

Options

1. Delegate. When uw.mpi.size > 1 and any point falls outside the local
   domain, route the whole call through global_evaluate. Correct by default,
   at the cost of making a fast local call sometimes collective — and a
   collective call from a conditional is its own hazard, so the condition would
   have to be reduced across ranks.
2. Refuse. Raise when a point is foreign in an MPI run, naming
   global_evaluate. Loud, cheap, and it breaks working rank-local code only if
   that code was already asking for points it does not own.
3. Warn. Emit once per call with a count of foreign points. Least
   disruptive, keeps the current contract, and makes the failure visible.
4. Document it. Necessary regardless: the docstring currently does not
   mention parallelism at all, and has no See Also pointing at
   global_evaluate.

Rank-locality is what makes evaluate fast for the common case where each rank
asks about its own points, and that is worth keeping. The defect is answering
confidently for points it cannot see.

Related: #551 measures this same path and records the rank-local classification
in one clause, but frames it as anti-scaling; the correctness consequence is not
stated there. #605 is the curved-boundary case. #604 is a located point
returning NaN.

Underworld development team with AI support from Claude Code

[lmoresi]

Folding this together with #606 (and what #605 left behind)

Bundling rather than tracking these separately, since they want one change.

Three findings, three mechanisms, one reason they are all invisible:

• A located point can interpolate to NaN, and #602 makes that invisible #604 — a point the mask reports as LOCATED returns NaN, and Recover non-finite values in parallel global evaluation #602 replaces
  it with an RBF extrapolation.
• evaluate() silently returns wrong values in parallel for points the rank does not own — including mesh nodes #606 — evaluate in parallel answers from the rank-local mesh, so a point
  another rank owns comes back extrapolated from up to a domain away.
• Points on a convex curved boundary are always extrapolated, silently, with O(h) error #605 (closed as not-a-defect) — a point on a convex curved boundary is
  outside the meshed domain and is extrapolated, correctly and unavoidably.

The mechanisms are unrelated. What they share is that evaluate and
global_evaluate both COMPUTE the extrapolation mask and then discard it unless
the caller passes check_extrapolated=True, which almost nobody does because the
default signature returns bare values.

One change addresses the reporting half of all three: count what was
extrapolated or recovered, and say so once per call. Then #605's curved-boundary
case explains itself instead of being a silent 6%, #604's substitution becomes
observable, and #606's rank-local answer announces itself rather than looking
like data.

That does not fix #604's underlying defect (a located point should not
interpolate to NaN) or #606's (the answer should be right, or refused). Those
stay as they are. But the reporting is the part that is cheap, shared, and worth
doing first — and it is what turns the other two from silent wrong answers into
things somebody notices.

Suggested split, so this does not sprawl:

1. Report — the shared change, here or in whichever of the two is picked up
   first.
2. evaluate() silently returns wrong values in parallel for points the rank does not own — including mesh nodes #606's correctness — delegate to global_evaluate, refuse, or document,
   as that issue sets out.
3. A located point can interpolate to NaN, and #602 makes that invisible #604's root cause — still undiagnosed; the RBF path succeeding where the
   located path returned NaN remains the strongest clue.

Underworld development team with AI support from Claude Code