uw.function.evaluate returns wrong values for P1 fields at points exactly on cell edges (3D)

[lmoresi]

Found while validating an MG prolongation against evaluate as the reference — the reference was the thing that was wrong.

Symptom

Evaluating a degree-1 scalar MeshVariable at points lying exactly on coarse-cell edges (3D tetrahedral mesh) returns values that are not the FE interpolant at those points.

Independent check: each failing point lies on exactly one coarse edge (a,b), so the P1 value there is unambiguously 0.5*(u_a + u_b). Comparing that against evaluate:

point	independent edge-truth	evaluate
504	+0.700344	+0.923709
565	+0.323892	−0.728195
2175	−0.287598	−2.759645
2177	−2.122083	−6.128067
2178	−0.119434	+0.782410
2205	+0.361674	+0.739306

Errors are O(1), not round-off. Affected points are a small fraction (2–12 per few thousand) but the wrong values are badly wrong.

Why it matters

These are not exotic points. Every vertex introduced by edge bisection lies exactly on a parent edge, so any workflow evaluating a field at refined-mesh vertices, or at points shared between cells, is exposed. It is silent — no warning, no NaN.

Reproducer sketch

base  = uw.meshing.UnstructuredSimplexBox(minCoords=(0.,0.,0.), maxCoords=(1.,1.,1.),
                                          cellSize=0.5, refinement=1, qdegree=2)
child = base.adapt(metric, max_levels=2)          # child vertices lie on parent edges
uc = uw.discretisation.MeshVariable("uc", coarse_mesh, 1, degree=1)
uc.data[:, 0] = np.random.standard_normal(uc.data.shape[0])
fx  = child.dm.getCoordinatesLocal().array.reshape(-1, 3)
val = uw.function.evaluate(uc.sym[0], fx)         # wrong at points on coarse edges

Compare against 0.5*(u_a + u_b) for the coarse edge each point sits on.

Full script: /tmp/who_is_wrong.py pattern — locate the coarse edge whose midpoint matches the query point, take the P1 value along it.

Suspected cause

Point location for a query lying exactly on a cell boundary. The point is in the closure of several cells; if location picks a cell and the basis is then evaluated with reference coordinates slightly outside that cell, the result extrapolates rather than interpolates. This repo has prior form here — see the quad TOP-face silent-zeros issue.

Notes

• 2D appears unaffected in the same test (0 disagreements across three refinement passes); 3D shows 2–12 per pass.
• A correct-by-construction MG prolongation agrees with the edge-truth in every case, so the discretisation and the field data are fine — only evaluate disagrees.

[lmoresi]

This is a RECURRENCE, not a new bug — maintainer: "Every time I think we have eliminated the element-edge bug in evaluate, it seems to come back." A session will be devoted to killing it properly.

That reframes the issue. If it keeps returning, the defect is not only in the location/evaluation code but in the absence of a permanent regression test that pins the on-boundary case. Whatever the fix is this time, it will come back again unless a test makes it impossible to reintroduce.

Prior incarnation for whoever takes it: #390, evaluate returning silent ZEROS on a quad TOP face, diagnosed then as a UW3 hint-policy regression, not PETSc — with the lesson recorded that byte-identical PETSc behaviour does not exonerate UW3. Worth starting from that diagnosis rather than re-deriving it.

What this instance adds:

• It is not confined to quads or to boundary faces of the domain. This is a 3D tetrahedral mesh, and the failing points are on interior coarse-cell edges.
• It is not silent zeros this time — the values are plausible but wrong (O(1) errors), which is worse: nothing looks anomalous downstream.
• 2D was unaffected in the same test (0 disagreements over three passes) while 3D showed 2-12 per pass. If the fix is dimension-general, 2D passing is not evidence.
• A ready-made oracle exists: any point on a coarse edge has unambiguous P1 value (1-t) u_a + t u_b along that edge. tests/test_0753_nested_mg_prolongation.py::test_reproduces_an_arbitrary_coarse_field already computes exactly this and would make a good basis for the regression test — it deliberately avoids evaluate as its reference.

Cost of not fixing it, from this session: I used evaluate as the reference to validate an MG prolongation, concluded my own (correct) code was broken, removed a working component and committed a detailed and entirely wrong explanation of why it was unsound. Recovering took several rounds of measurement. A wrong oracle is more expensive than a missing one.