Skip to content

Make XmGridTrace correct and fast enough for per-glyph flow paths - #10

Merged
wdolinar merged 10 commits into
bump/xmsconan-2.15.3from
feature/follow-flow-timesteps
Aug 17, 2026
Merged

Make XmGridTrace correct and fast enough for per-glyph flow paths#10
wdolinar merged 10 commits into
bump/xmsconan-2.15.3from
feature/follow-flow-timesteps

Conversation

@wdolinar

Copy link
Copy Markdown
Member

Groundwork for routing the "follow flow path" 2D vector display option through XmGridTrace,
so the traced path integrates the vector dataset as it changes through time rather than
against one frozen timestep. This is the xmsgridtrace half; next_ms orchestration and the
xmsvtk consumer are separate.

Based on bump/xmsconan-2.15.3, not master — targeting master would bundle that branch's
four migration commits. Retarget once it merges.

Correctness

  • The time interpolation was inverted. perc1 = |t - t1| / totalTime multiplied timestep
    one's
    value, so a particle released at t1 was advected entirely by the field at t2. Three
    recorded baselines had baked it in. Each new first step was derived by hand rather than
    accepted from the runner.
  • No-data is now propagated, not blended. With the weights corrected, a cell inactive in only
    one timestep stopped resolving to the sentinel and started resolving to a blend of it —
    0.9 * 0.1 + 0.1 * -9999999 is -999999.9, which passes every no-data test the callers make.
    A cell active at t1 but inactive at t2 only actually terminates the trace now.
  • The two output arrays could come back different lengths. Position pushes were conditional
    on movement while time pushes were unconditional, so a sub-tolerance step misaligned every
    later pair — undetectable to a caller zipping them.
  • A latent null dereference when only one timestep had been supplied now returns a clean
    extraction failure.

Tracing through time

Traces suspend at the second loaded timestep and resume when a later one arrives, so the field
is never extrapolated past what is known:

tracer->StartTraces(seeds, seedTimes);
while (tracer->ContinueTraces() > 0 && series.HasNext())
  tracer->AddGridScalarsAtTime(series.Next(), ...);
tracer->GetTraceResults(traces, times, reasons);

The two-step window is kept, so memory stays bounded however long the series is and the caller
reads timesteps only as traces need them. Position, time, the whole-trace distance and elapsed
budgets, the adaptive step size and the previous velocity all carry across a boundary — the
subdivision tests compare each step against the one before it, so dropping those would kink the
path exactly where the timestep changes.

Exit reasons are now an enum rather than a string. A caller must distinguish "left the grid, draw
it short" from "spent its distance budget" for tens of thousands of seeds, and the old messages
were composed by appending, so no fixed string identified a case.

Python bindings included; continue_traces releases the GIL, and only it.

Performance

Measured on a 200x200 quad grid, 10,000 seeds, Release. Per-glyph tracing was ~830x too slow to
be viable; it is now faster than the render-thread tracer it replaces.

per 10k-glyph re-trace
before ~21 s
cache the boundary-exit extractor ~0.66 s
one search per triangulation ~0.15 s
what it replaces (xmsvtk, on the render thread) 0.30 s

Two independent problems:

  • The out-of-domain branch built a whole-grid triangulation and a GmMultiPolyIntersector per
    exit event, inside the stepping loop
    — ~40 ms each, against ~48 us for a complete interior
    trace. Caching one extractor takes that to ~9.5 us per exit. On a realistic seed population, 5%
    of seeds were accounting for 98% of total time.
  • Four ExtractData calls per sample each ran their own point-location query for the same
    (x, y)
    , with identical weights. Now one search per distinct triangulation: the x/y pair shares
    one, and both timesteps share one when their activity masks are equal. Searches per seed
    97.5 -> 24.4, setup 117 -> 32 ms.

Both are xmsgridtrace-only; no xmsextractor change, no version cascade. The sharing constructor
they rely on has existed since 2022 and was simply unused here.

Test Plan

  • C++ suite green: cmake --build build/Release && ./build/Release/runner — 22/22
  • Python suite green: python build.py --python-only — 19/19, build exit 0
  • flake8 clean under the Aquaveo rules: python -m flake8 _package/ — exit 0
  • Every pre-existing recorded baseline is unchanged by the Tier 1 search reduction, which is
    the evidence that applying weights manually reproduces ExtractData bit-for-bit
  • testTracesContinueAcrossTimeSteps asserts the trace denied a third timestep is an exact
    prefix of the one given it — resuming extends the path rather than recomputing it
  • Both sharing branches covered: test_inactive_cell has differing activity between timesteps
    (two searches), test_unique_time_steps has identical activity (shared)
  • Benchmark reproduces: XMGT_BENCH_SEEDS=10000 ./build/Release/runner

The "follow flow path" vector display option is being routed through XmGridTrace so
the traced path integrates a time-varying field. That only works if the tracer is fast
enough to run for every visible glyph, so this establishes what it costs today.

testTraceBenchmark traces N seeds over a 200x200 quad grid loaded with two timesteps of
a vortex field whose rotation reverses between them, across three seed populations:

  interior  far enough from the edge that no trace can reach it -- pure stepping cost
  boundary  in a band along the perimeter -- forces the out-of-domain exit branch
  mixed     spread over the whole domain -- what the display actually does

Separating them matters: the populations turn out to differ by two orders of magnitude
per seed, and an undifferentiated average would have hidden that.

Alongside wall time it reports an ExtractData call count, from a CXX_TEST-only counter
incremented where the four per-step searches happen. Without it an optimization cannot
be shown to have removed searches rather than merely found a faster machine. It also
reports a setup breakdown -- BuildTriangles, the GmTriSearch R-tree build, and an
activity-only reapply timed separately -- because which of those dominates decides
whether triangulations can be shared across timesteps.

XMGT_BENCH_SEEDS and XMGT_BENCH_CELLS size the run so a sweep needs no recompile. The
defaults are small enough to leave the case in the regular suite, and the assertions
are order-of-magnitude guards rather than tight bounds so it will not go flaky on a
shared runner.

One assertion is deliberately loose for a measured reason: a seed that exits the grid
on its first step can reach the points.size() < 3 early return and come back with only
the seed point, so "every seed yields a usable polyline" is false. It shows up at
roughly 1 in 100,000.
…er exit

TracePoint constructed a fresh XmUGrid2dPolylineDataExtractor on every out-of-domain
step, inside the stepping loop. That constructor triangulates the whole grid, and the
SetPolyline that follows indexes every triangle into a new GmMultiPolyIntersector.
Both depend only on the grid, which cannot change during a trace, and both were thrown
away at the end of the if block and rebuilt for the next exiting particle.

Benchmarked on a 200x200 grid, one exit event cost ~40 ms against ~48 us for a
complete interior trace -- roughly 830x. On a seed population with a 5% exit rate,
those 5% accounted for 98% of total trace time.

It is now a member built lazily on the first exit and reused. Three things make the
reuse safe, each checked in xmsextractor rather than assumed:

  - BuildTriangles is guarded by m_triangleType != a_location, so the second
    SetPolyline skips the triangulation.
  - ComputeExtractLocations builds m_multiPolyIntersector only when null and clears
    its output locations at entry, so no state carries between polylines.
  - XmGridTrace consumes only GetExtractLocations(), never extracted values, so the
    dummy zero scalars the constructor installs are irrelevant and the instance stays
    valid for the tracer's lifetime.

m_ugrid is fixed at construction with no setter, so a cached extractor cannot outlive
the grid it was built for. The member stays null until a trace actually exits, so a
seed population that never reaches a boundary pays no memory for it.

Measured A/B on one machine state at 1,000 seeds: the boundary population goes from
15,996 to 71 us/seed (224x) and the realistic mixed population from 2,200 to 51
us/seed (43x), while the interior population -- which never exits, and is the control
-- moves 1%. Per exit the cost falls from 40 ms to ~9.5 us. ExtractData counts are
identical before and after: this removes no searches at all, only whole-grid rebuilds,
so the separate search-reduction work is still available on top of it.

testBoundaryExtractorIsCached is the guard. Caching is invisible in the output, so the
assertion that matters is a construction count from a CXX_TEST-only counter; the test
also compares the two traces to 1e-12, which is what would catch reuse silently
changing an answer.
GetVectorAtLocationAndTime weighted each timestep by its own distance from the current
time, so perc1 was |t - t1| / totalTime and multiplied timestep 1's value. At t == t1
that weight is zero: a particle released exactly at the first timestep was advected
entirely by the field at the second. A timestep is weighted by its *closeness* to the
current time, so the distance from one timestep is the weight of the other.

This is the defect that makes the tracer worth using at all. Routing the "follow flow
path" display option through it is only an improvement if the trace follows the field
as it changes, and an inverted blend does not.

Propagating XM_NODATA is part of the same fix rather than a separate one. With the
weights corrected, a location that is inactive in only one of the two timesteps stops
resolving to the sentinel and starts resolving to a blend of it: 0.9 * 0.1 + 0.1 *
-9999999 is -999999.9, which is neither no-data nor a velocity, and which passes every
no-data test the callers make. The sentinel is now propagated when either bracketing
timestep has no data at the location, which is what makes "a cell active at t1 but
inactive at t2 terminates the trace" actually hold. testStartInactiveCell was passing
before only because the inverted weights happened to give timestep 2 all the weight at
t == t1; with the propagation it passes for the right reason.

testTimeVaryingFieldChangesPath is the regression guard. One cell spanning the domain
gives a spatially uniform field, so any curvature in the path can only have come from
time; the field rotates +x -> +y between timesteps rather than reversing, so the
interpolated velocity never passes through zero and cannot trip the velocity-is-zero
exit partway along. Its first assertion is the one that catches an inversion -- a
particle released at t1 must step due east with y untouched -- and it compares against
a frozen-field control traced by the same code, which never turns.

Three existing baselines move. Each first step was checked by hand rather than
accepted from the runner:

  testUniqueTimeSteps  0.5 -> 0.6, the t1 cell value 0.1 over dt 1, where it was 0.7
                       from t2's 0.2. Its second step, 0.9*0.11 + 0.1*0.21 = 0.12,
                       matches the recorded 0.744 to the digit.
  testInactiveCell     same first step, and the trace now terminates exactly at x = 1,
                       the boundary of the cell that is inactive at t2, rather than at
                       0.9979 -- it had been stopping just short via a max-change-
                       velocity blow-up on a no-data-contaminated blend.
  testTutorial         first step y 0.5 -> 1.5: corner scalars interpolate to (0, 0.5),
                       times the multiplier of 2, over dt 1. The old 1.25 was a
                       boundary-clipping artifact of using the doubled second-timestep
                       field at t = 0 -- the tutorial's own comment says the second
                       timestep is doubled to show an increase, so the trace should
                       start at the first timestep's magnitude and speed up.
The "follow flow path" display traces every visible vector glyph, tens of thousands of
them per redraw, driven from Python. One TracePoint call per glyph pays a language
boundary crossing per glyph for work that is identical across them, and it can only
report why the *last* trace ended -- GetExitMessage describes a single operation, so a
caller has no way to ask why glyph 4,000 stopped short.

TracePoints takes all the seeds at once and returns a polyline, a time array, and an
exit message per seed. It does not advance the time steps: every trace runs against
whichever pair AddGridScalarsAtTime most recently supplied, and a caller wanting traces
that span more of a series feeds the next step and traces again. Keeping that in the
caller is deliberate -- the two-step window is instance state, so a batch that advanced
it internally would have to carry per-seed continuation state, which is a different and
larger design than this one.

Mismatched input lengths return nothing rather than tracing the common prefix. A caller
that supplied the wrong number of start times has a bug, and a partial result lets it go
unnoticed.

TracePoint's two output arrays could also come back different lengths, which this fixes
because the batch documents them as parallel. The position push was conditional on the
step actually moving while the time push was unconditional, so a step shorter than
XM_ZERO_TOL left the times array one longer and silently misaligned every later pair --
undetectable to a caller zipping them. The time is now pushed only when the point is.
No existing expectation moves, so no recorded trace contained such a step.

testTracePointsMatchesSerialTracePoint compares the batch against serial TracePoint
calls on an identical fixture rather than against a recorded baseline, which would drift
with the tracer instead of pinning the equivalence. Its seeds cover the three shapes a
caller has to handle: two traces that leave the grid, and a seed outside it that yields
an empty trace rather than a polyline.
…n enum

A trace can only run as far as the second of the two loaded time steps, because that is
as far as the field is known. The batch added in the previous commit therefore traced
every seed to the edge of one window and threw away everything it knew, which is not
tracing a flow through time -- it is tracing it through one interval. Restarting from
the last position would not fix that either, because a restart loses the trace's history.

Traces now suspend and resume. StartTraces seeds a batch, ContinueTraces advances every
unfinished trace and returns how many are waiting on a later time step, and the caller
feeds the next one and calls again:

    tracer->StartTraces(seeds, seedTimes);
    while (tracer->ContinueTraces() > 0 && series.HasNext())
      tracer->AddGridScalarsAtTime(series.Next(), ...);
    tracer->GetTraceResults(traces, times, reasons);

Keeping the two-step window means memory stays bounded however long the series is, and
the caller reads time steps only as the traces actually need them. Stopping early is
legitimate: traces still waiting end where they got to, and say so.

The substance is what survives a window change. Position and time are the obvious ones.
The distance and elapsed-time budgets are whole-trace, not per-window, so they carry.
So do the adaptive step size and the previous velocity, because the subdivision tests
compare each step against the one before it -- restarting those at a boundary would kink
the path exactly where the time step changes, which is the one place this has to be
smooth. TracePoint's body is now StepTrace, which either starts a trace or resumes one;
every exit from it routes through a single lambda that writes that state back, so there
is no path that advances a trace without recording where it reached.

Resumability cannot be read off the loop's final state, so it is tracked explicitly: a
subdivision puts the trace back in motion *after* the time step clamp has already fired,
and several conditions in one iteration overwrite each other.

This also only works because of the interpolation fix. A trace resuming at the new first
time step is advected by that step's field; under the inverted weights it would have used
the following one, so every window boundary would have introduced an error.

Since nothing outside this repository uses XmGridTrace, the surrounding API is cleaned up
rather than extended around:

  - The exit reason is an enum, not a message. A caller has to tell "left the grid, draw
    it short" from "spent its distance budget, this is the normal ending" for tens of
    thousands of seeds, and string comparison cannot support that -- the old messages were
    composed by appending, so no fixed string identified a case. The strings remain, one
    per reason, for logs and tooltips. TracePoint gets GetExitReason so the single-point
    path can answer the same question the batch answers.
  - The one-shot TracePoints from the previous commit is gone; the resumable trio subsumes
    it, and one way to batch is better than two.
  - GetExitMessage returns const std::string& and is const.
  - AddGridScalarsAtTime takes activity by const reference.

testTracesContinueAcrossTimeSteps is the guard, over a field rotating +x -> +y -> -x. Its
strongest assertion is not that the continued trace is longer: it is that the trace never
given the third time step is an exact prefix of the one that was. That is what shows
resuming extends the path rather than recomputing it, and it is what fails if any carried
state is dropped at the boundary. It also checks the path turns back on itself, which no
single pair of those time steps can produce.

Measured on the benchmark's realistic seed population, 43 of 250 seeds stop waiting for a
later time step -- traces the previous design silently truncated.

The Python bindings still compile against this by inspection: their three uses are lambdas
or a member pointer that still resolves. Binding the new calls needs a python-enabled
build, which the testing preset does not produce, and is not done here.
Three Python tests mirror C++ cases whose expectations moved when the inverted time
interpolation was fixed, so they carried the same wrong values: test_unique_time_steps,
test_inactive_cell and test_tutorial. The new values are transcribed from the C++ source,
where each first step was derived by hand rather than captured from the runner.

test_max_tracing_distance failed for a different reason worth recording. Its positions
matched to six decimals while one *time* differed by 4.4e-16 -- one ULP -- because the
times were compared with assert_array_equal, exact float equality, while the positions in
the same test were already compared approximately.

Bisecting placed it on the interpolation fix, which was not the obvious answer: that test
supplies identical scalars at both time steps, so the two weighted terms are the same pair
of products and IEEE addition is commutative. The cause is FMA contraction. The compiler
folds `d1 * w1 + d2 * w2` into a fused multiply-add, computing one product exactly inside
the FMA and rounding the other; swapping which weight multiplies which time step therefore
moves the last bit even though the mathematics is unchanged.

So the assertion was pinning the compiler's contraction decision rather than the tracer's
behaviour. All fifteen of these time comparisons had the same latent fragility and only one
happened to trip; they now compare approximately, matching what the same tests already do
for positions and what the C++ mirrors do. The reason is recorded in the class docstring so
a future reader does not tighten them again.

The bisect also confirmed something worth having checked: caching the boundary-exit
polyline extractor changes no numeric result. The commit that introduced it passes all
sixteen Python tests, which the reused GmMultiPolyIntersector could plausibly not have done.

flake8 is not installed in this environment, so the Python edit is unlinted.
…t was missing

Adds start_traces, continue_traces, get_trace_results and get_exit_reason to the pybind11
module, exports XmGridTraceExitEnum as exit_reason_enum, and forwards all four through the
hand-written wrapper so FlowPathService can drive the resume loop from Python:

    tracer.start_traces(seeds, seed_times)
    while tracer.continue_traces() > 0:
        step = series.next()
        if step is None:
            break
        tracer.add_grid_scalars_at_time(*step)
    traces, times, reasons = tracer.get_trace_results()

continue_traces releases the GIL, and only it. Tracing tens of thousands of seeds takes long
enough that holding the GIL would stall the interpreter for a caller on a worker thread,
which is exactly how this is meant to run. start_traces keeps the GIL because it converts
Python iterables inside its lambda.

start_traces raises ValueError when the start times do not match the points. The C++ side
refuses the batch and returns empty, which from Python would look like a tracer that silently
did nothing.

GetExitReason is added here rather than earlier because it was never actually added. The
commit that claimed it used a scripted string replacement with no assertion that the anchor
matched, so the edit silently did nothing; the build then passed because nothing had changed,
and the claim went unverified into that commit message. It is now on the interface, the impl,
and covered by a test that asserts the single-point path and the batch report the same reason
for the same seed -- if those can disagree, a caller cannot use them interchangeably.

The lesson is in the tooling, not the code: scripted edits need an assertion that the anchor
was found, and a green build is not evidence that an edit landed.

Three Python tests cover the new surface: a trace that continues across three time steps and
whose stopped-early result is a prefix of the continued one, the ValueError, and the batch
matching serial trace_point calls.

Two clang-format suggestions are deliberately not applied. It wants GetExitReason collapsed
to one line, along with all nine sibling accessors that are not written that way; and it wants
425 of the 437 lines of XmGridTrace_py.cpp reformatted, a file never kept under clang-format.
Both would bury the change in unrelated churn.

Also corrects the previous commit's claim that the Python test edit was unlinted -- flake8 is
installed now and it is clean. The only findings in _package are eight pre-existing AQU104
import-header comments, four of them in a file this branch never touched.
flake8 reported eight AQU104 findings across the two Python source files, all on line 1:
neither carried the numbered import group comments the Aquaveo rules require. They are
pre-existing -- four of them are in grid_trace.py, which this branch had not otherwise
touched -- and were invisible until flake8-aquaveo was installed. Plain flake8 only runs
pycodestyle, pyflakes and mccabe, so the AQU rules, the google docstring convention and the
appnexus import order the .flake8 config asks for were all silently unchecked.

The convention, matched from next_ms, is all four comments present even when a section is
empty, with a blank line after the module docstring.

No import moved, so the suite passing is confirmation the modules still resolve the same way.
_package is now clean under the full rule set.
…er sample

GetVectorAtLocationAndTime ran four ExtractData calls per sample -- x and y, for each of two
time steps -- and each one performed its own point-location query for the same (x, y). The
interpolation weights were identical across all four; only the scalar array being weighted
differed.

It now runs one search per distinct triangulation and applies the resulting weights to each
component. Three things make that possible, and the third is the one that had to be measured
rather than assumed:

  - The x and y extractors of a time step now share a triangulation, via the sharing
    constructor that has existed since 2022 and was simply never used here.
  - Both time steps share one as well when their activity masks are equal. The triangulation
    and its R-tree depend only on the grid and the mask, so an identical mask makes them
    interchangeable. Differing activity is the case that genuinely cannot share, so the mask
    is compared rather than assumed -- test_inactive_cell covers that path and
    test_unique_time_steps covers the shared one.
  - iApplyWeights reproduces ExtractData exactly, accumulating in double and narrowing to
    float, so this is bit-identical rather than merely close. Every recorded baseline in both
    the C++ and Python suites is unchanged, which is the evidence for that.

Ordering matters in AddGridScalarsAtTime and is easy to get backwards: the y extractor is
built from x, and only after x's scalars are set. The sharing constructor copies the
triangulation and the flag saying what it was built for, so copying x before it has built one
leaves y believing it must build, and y then rebuilds the very triangulation it is sharing --
silently costing what the sharing was meant to save.

Measured on a 200x200 grid at 10,000 seeds, against the previous commit:

                        before     after
  setup, 2 time steps   117 ms     32 ms
  searches per seed     97.5       24.4
  interior us/seed      50.8       12.6
  mixed us/seed         50.7       11.8

The setup figure is better than the ~53 ms projected in TRIANGULATION_SHARING.md. That
projection assumed a per-extractor scalar-handling floor of roughly 8 ms, derived by
subtraction rather than timed; the measurement says it is far smaller, and that only one
triangulation and one R-tree are now built for all four extractors rather than four of each.

A re-trace of 10,000 glyphs is now about 150 ms including setup, against the 299 ms the
current in-render tracer costs for the same work -- so the new path is faster than what it
replaces while adding the time dimension, rather than merely close enough. It began this
branch at about 21 seconds.

Also fixes a latent crash this rewrite made obvious: with only one time step supplied, the
first extractor is null and was dereferenced. It now returns a clean extraction failure,
which is what edge case 9 in the session plan assumed already happened.

The benchmark's counter counted ExtractData calls and now counts searches, which is what it
always meant; its label and field are renamed to match rather than silently changing meaning.
@wdolinar

Copy link
Copy Markdown
Member Author

Review: batch time-stepped tracing in xmsgridtrace (7 files, C++ core + Python bindings/tests)

Multi-reviewer pass over 658dd45: 8 specialist reviewers, deduped to 53 entries, then every CRITICAL/MAJOR entry independently re-checked by a separate validator. 29 findings validated — 3 confirmed, 26 dropped. Only the confirmed ones appear below as Critical/Major; the Minor list was not validated.

Summary

What changed. XmGridTrace gains a batch API — StartTraces / ContinueTraces / GetTraceResults — that suspends each trace at the second loaded time step and resumes it when a later one arrives, with the per-seed state (position, elapsed time, distance budget, adaptive step size, previous velocity) carried across the boundary in a new TraceState. Alongside it: exit reasons become an enum instead of composed strings, the time-interpolation weights are swapped, the out-of-domain branch caches a XmUGrid2dPolylineDataExtractor instead of rebuilding one per exit event, consecutive time steps share a triangulation when their activity masks match, and Python bindings plus a wrapper (grid_trace.py) and package export expose the new surface with the GIL released around continue_traces. This is a large diff for one class — roughly 2,200 patch lines, of which ~470 are new inline CxxTest bodies in XmGridTrace.cpp including a wall-clock benchmark, and it touches a public, versioned C++ API.

Notably, the array-alignment fix at xmsgridtrace/gridtrace/XmGridTrace.cpp:680-692 is done the right way: the position push and the time push now sit inside a single if (moved) block gated on one moved predicate, so the two parallel output arrays cannot drift out of step regardless of which branch reached that point.

Verdict: BLOCK (C=1 M=2 m=24)


Critical items

C1: xmsgridtrace/gridtrace/XmGridTrace.cpp:686 — resumed step size can be exactly zero, hanging ContinueTraces [validated]

The clamp that ends a window writes deltaT into the resumed state; when a trace already sits exactly on m_time2, the subtraction yields 0.0:

    // If the change in DeltaT would push us beyond the time step, set it to hit the timestep
    if (elapsedTime + deltaT + ptTime > m_time2)
    {
      deltaT = m_time2 - elapsedTime - ptTime;
      bContinue = false; // This will be the last point traced in this window
      stopReason = GTEXIT_WAITING_FOR_TIME_STEP;
    }

deltaT *= 1.2 at :686 cannot lift zero, stopWith at :490-501 persists it into a_state.m_deltaT, and AddGridScalarsAtTime never resets it. On the next StepTrace the loop at :534 has no reachable exit: :542 cannot raise deltaT, :546 no longer fires, elapsedTime never advances so :554 cannot fire, pt1 == pt0 so neither the out-of-domain branch at :573 nor the max-distance branch at :657 triggers, and the min-delta-time escape at :646 is unreachable because a zero-length step re-evaluates the same point (vx1 == vx0, mag1 == mag0), so bSplit is never set.

Change: Floor the resumed step size on entry to StepTrace — treat a non-positive a_state.m_deltaT as the initial 1.0 (or m_minDeltaTime) where deltaT is read at :477.

Why: This is an unbreakable infinite loop, not a slow path: moved is false so nothing is appended and no allocation ever trips a limit, and the binding releases the GIL at XmGridTrace_py.cpp:386-387, so Ctrl-C cannot interrupt it. The header at XmGridTrace.h:157-158 explicitly blesses the triggering sequence ("Calling ContinueTraces twice without supplying a time step in between does no useful work"), and the PR's own create_rotating_field_tracer fixture reproduces it.


Major items

M1: xmsgridtrace/gridtrace/XmGridTrace.cpp:441 — triangulation-sharing predicate omits the scalar location [validated]

  // Share the triangulation with the previous time step when the two agree on activity. The
  // triangulation and its GmTriSearch R-tree depend only on the grid and the activity mask,
  // so an identical mask makes them interchangeable -- which both skips a rebuild and lets a
  // single point-location query serve all four extractors instead of one per time step.
  // Differing activity is the case that cannot share, and it is why the mask is compared
  // rather than assumed.
  m_sharedAcrossTime = hadPrevious && a_activity == m_activity2;
  m_extractor2x = m_sharedAcrossTime ? XmUGrid2dDataExtractor::New(m_extractor1x)
                                     : XmUGrid2dDataExtractor::New(m_ugrid);

The in-code invariant is wrong in the direction that matters: the triangulation depends on the point option derived from a_scalarLoc, while activity only feeds m_triangles->SetCellActivity. The sharing constructor shares the m_triangles shared_ptr itself, and SetGridCellScalars rebuilds that shared object in place — retroactively invalidating time step 1. With step 1 at LOC_POINTS and step 2 at LOC_CELLS under equal activity, m_extractor1x->GetScalars() still holds numGridPoints entries while the rebuilt triangulation returns centroid indices >= numGridPoints, and iApplyWeights(*m_extractor1x, ...) at :804 reads out of bounds at :111.

Change: Include the scalar location in the predicate at :441 — share only when both the activity mask and a_scalarLoc match the previous call — and correct the comment at :435-440. Keep the existing activity term for the shared tri-search.

Why: A legal combination on a public API (XmGridTrace.h:117-129 documents no constraint that a_scalarLoc stay constant, and the Python binding takes it per call) becomes an out-of-bounds vector read. The pre-diff code constructed a fresh New(m_ugrid) per time step and was safe against it.

M2: xmsgridtrace/gridtrace/XmGridTrace.cpp:507 — seed past the window is permanently killed as an extraction failure [validated]

  if (!a_state.m_started)
  {
    outTrace.clear();
    outTimes.clear();
    if (ptTime > m_time2 || // Test if the time specified is after the time range
        !GetVectorAtLocationAndTime(pt0, ptTime, vector)) // Ensure extraction did not fail
    {
      stopWith(GTEXIT_EXTRACTION_FAILED);
      return;
    }

iIsTerminal at :79-82 treats every reason except GTEXIT_NOT_STARTED / GTEXIT_WAITING_FOR_TIME_STEP as terminal, so the early return at :472-473 makes such a seed unrevivable by any later AddGridScalarsAtTime, and ContinueTraces at :748 does not count it as waiting — the documented driving loop (XmGridTrace.h:149-154) can exit immediately, leaving the seed with an empty trace and a misleading reason. The same code already treats the identical condition the other way mid-trace: at :546-551 a step that would cross m_time2 defers with GTEXIT_WAITING_FOR_TIME_STEP.

Change: In the unstarted branch, split the two conditions: when ptTime > m_time2, stop with GTEXIT_WAITING_FOR_TIME_STEP (non-terminal) so the seed starts once a covering window is loaded; reserve GTEXIT_EXTRACTION_FAILED for the actual GetVectorAtLocationAndTime failure.

Why: Seeds released at times later than the first loaded window silently produce empty traces with the wrong reason, and the loop the header prescribes terminates before their data ever arrives. StartTraces at :721-737 accepts arbitrary per-seed times with no window validation, and neither the header (XmGridTrace.h:141-166) nor the Python docstrings state a requirement that seed times lie inside the loaded window — so this is not a documented caller error.


Minor items

  • _package/tests/XmGridTrace_pyt.py:16 — 738-line TestCase mixes trace_point knobs, the batch API and a tutorial; split into TestTracePoint / TestBatchTraces.
  • _package/tests/XmGridTrace_pyt.py:17 — class docstring is a changelog about an FMA-contraction bug, not a description; move rationale to the commit message.
  • _package/tests/XmGridTrace_pyt.py:126 — four structurally identical knob tests (:126, :209, :346, :414); drop the TestCase base and parametrize.
  • _package/tests/XmGridTrace_pyt.py:263create_rotating_field_tracer sits mid-class among tests while the other creators are at :28/:53; group the creation helpers together.
  • _package/tests/XmGridTrace_pyt.py:339 — per-seed assertions in the loop lack self.subTest(i=i); wrap the loop body in subTest.
  • _package/tests/XmGridTrace_pyt.py:388 — three distinct failure modes (:388, :399, :610) assert byte-identical empty results; assert the exit reason in each.
  • _package/tests/XmGridTrace_pyt.py:552 — golden arrays at :552, :584, :670 regenerated wholesale from new output; keep the behavioral check, shrink the goldens.
  • _package/xms/gridtrace/__init__.py:3 — exports exit_reason_enum but not data_location_enum; export it and accept it alongside the magic strings.
  • _package/xms/gridtrace/grid_trace.py:32 — pre-existing **kwargs['instance'] backdoor bypasses the ugrid is None check; use a classmethod factory instead.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:55 — test-only file-scope counters (:55, :62) are not parallel-safe; document the single-threaded requirement.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:97iApplyWeights takes six params, four of them x/y pairs; pass the extractor pair, return a Pt3d.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:111xScalars[ptIdx]/yScalars[ptIdx] unchecked against array size; assert size >= triangulation point count. (Same site M1 reaches out of bounds.)
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:238 — new m_boundaryExtractor uses legacy BSHP<T>; defer to a file-wide Ptr<T> migration.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:425AddGridScalarsAtTime accepts a non-increasing a_time silently; reject backwards or zero-width windows.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:499stopWith overwrites instance-level m_exitReason/m_exitMessage per seed, leaving batch-mode GetExitReason arbitrary or stale; document as TracePoint-only, or reset in StartTraces.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:500 — one std::string assignment per seed for a value no batch caller reads; set the message only on the TracePoint path.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:596deltaT *= (newSegDist / segDist) unguarded against segDist == 0; guard the divisor.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:669 — clipped point's time uses deltaT * perc while its position uses (1 - perc); use deltaT * (1 - perc).
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:840totalTime = fabs(m_time1 - m_time2) divides at :846-847 with no zero guard, yielding NaN positions no exit test catches; reject or short-circuit equal times.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:1188 — benchmark helper re-traces every seed just for the reason histogram; collect reasons in the timed pass.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:1982 — 33-row baselines duplicated verbatim in XmGridTrace_pyt.py:670-735; let C++ own the baselines, Python the marshalling.
  • xmsgridtrace/gridtrace/XmGridTrace.cpp:2079 — the same 7-setter fixture run appears three times (:2080, :2216, :2352); extract a shared helper.
  • xmsgridtrace/gridtrace/XmGridTrace.h:200GetExitMessage now returns a const std::string& to a member the next trace overwrites; document the reference's validity window.
  • xmsgridtrace/python/gridtrace/XmGridTrace_py.cpp:387 — GIL released while ContinueTraces can call XM_LOG (XmGridTrace.cpp:589, :790); verify no Python log sink is registered through xmscore.

A resumed trace could hang, a legal pair of time steps could read out of
bounds, and a staggered seed could be killed off permanently. All three
are in code this branch introduced.

StepTrace hung when a window ended exactly on the second time step. The
time step clamp computes deltaT = m_time2 - elapsed - ptTime, which for a
trace already sitting on m_time2 is exactly zero, and that zero was
persisted into the resumed trace. A zero-length step moves nothing and
changes no velocity, so no clamp and no subdivision test could ever end
the loop -- and the min-delta-time escape is inside the split branch,
which a zero-length step can never enter. It spun forever appending
nothing, with the GIL released so Python could not interrupt it. The
window is now finished before stepping, keeping the step size the call
came in with, so a redundant ContinueTraces really is the no-op the
header promises. StepTrace also floors a non-positive resumed step size,
so no path can reintroduce this.

AddGridScalarsAtTime decided triangulation sharing from the activity mask
alone. The triangulation is built for a data location -- LOC_CELLS adds a
centroid per cell, LOC_POINTS adds none -- and sharing shares the object
rather than copying it, so the second step's SetGrid*Scalars rebuilt the
triangulation the first step was still pointing at. Its shorter scalar
array was then indexed by the new centroid indices: an out-of-bounds read,
not a wrong answer. Both data locations now join the mask in the predicate.

A seed released after the loaded window reported GTEXIT_EXTRACTION_FAILED,
which iIsTerminal treats as terminal, so the seed never started even once
its time step arrived. StartTraces takes a release time per seed so a
batch can be staggered, making this an ordinary input; it now reports
GTEXIT_WAITING_FOR_TIME_STEP, as the mid-trace clamp always did.

Adds a regression test per fix. testBeyondTimestep and its Python twin
asserted only that the trace was empty, which is why the third defect went
unnoticed -- they now assert the exit reason, which is what tells the
three empty-trace cases apart.

C++ 25/25, Python 19/19, flake8 clean.
@wdolinar

Copy link
Copy Markdown
Member Author

All three confirmed findings are fixed in 143cfd7, each with a regression test. Point-by-point:

C1 — zero step size hangs ContinueTracesagreed, fixed

Confirmed the whole chain, including that the min-delta-time escape at :646 is inside if (bSplit) and a zero-length step can never set bSplit, so there was genuinely no exit.

The suggested fix — floor the resumed step size — stops the hang, and I took it. But it is not sufficient on its own, which my own regression test caught: it substitutes 1.0 for a step size the un-redundant run carries as (clamp remainder) × 1.2, so the redundant ContinueTraces still silently changed the path that followed (35 points vs 31). XmGridTrace.h:157 promises the call "does no useful work", and that has to mean the outcome is indistinguishable.

So the primary fix is at the clamp instead: when the window has nothing left, finish it before stepping and restore the step size the call came in with. The floor stays as a backstop so no other path can reintroduce a zero. The test asserts equality with the run that made no redundant call, not merely that it terminates.

M1 — sharing predicate omits the data location — agreed, fixed

The comment at :435-440 asserting the triangulation depends "only on the grid and the activity mask" was the actual defect; the code faithfully implemented a wrong invariant.

One correction to the finding, in the direction of more: a_activityLoc has to be compared too. It decides how the same bitset maps onto cell activity, so an equal mask under a changed activity location is also not shareable. The predicate now requires all three to match, and the comment says why each term is load-bearing.

M2 — future seed killed as an extraction failure — agreed, fixed

The || fused two unrelated conditions into one terminal reason. Split, with the out-of-window case reporting GTEXIT_WAITING_FOR_TIME_STEP exactly as the mid-trace clamp at :546-551 always did. StartTraces' docs now state that staggered release times past the window are supported rather than merely tolerated.

On the minor items

Taken now, because it is part of M2 rather than cleanup: testBeyondTimestep and test_beyond_timestep asserted only that the trace came back empty. Three unrelated outcomes produce an empty trace, so that assertion could not tell them apart — which is precisely why M2 sat unnoticed. Both now assert the exit reason, as does test_start_out_of_cell.

The remaining ~20 minor items are real and I am not disputing them. Deferring them deliberately: folding a TestCase split, a BSHPPtr migration and golden-array trimming into a branch that is fixing a hang would bury the three changes a reviewer actually needs to check. Worth their own PR.

The 26 dropped findings

Flagging one that was dropped correctly but is worth a second look by a human: XmGridTrace.cpp:587 conflates a failed boundary intersection with a clean grid exit under one GTEXIT_LEFT_GRID. The validator dropped it because the pre-image behaved identically, which is right — it is not a regression. But this PR introduces the machine-readable exit reason that would make the two distinguishable, so it is the natural place to fix it, and the benchmark tolerance at :2432 exists to absorb the resulting dropped traces. I have left it alone rather than widen the diff; say the word if you would rather it went in here.

Verification

  • C++ 25/25 (./build/Release/runner), up from 22 — the three new tests are testRedundantContinueDoesNotStallTrace, testSeedReleasedAfterWindowWaitsThenTraces, testDataLocationChangeIsNotShared
  • Python 19/19 (python build.py --python-only)
  • flake8 clean on the changed file

The C1 test is worth knowing about operationally: before the fix it does not fail, it hangs. I ran the suite under a 300 s watchdog for that reason. It completes in a few seconds now, but a future regression there will look like a CI timeout rather than a red assertion.

@wdolinar
wdolinar merged commit 25518a0 into bump/xmsconan-2.15.3 Aug 17, 2026
18 checks passed
@wdolinar
wdolinar deleted the feature/follow-flow-timesteps branch August 17, 2026 16:47
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant