Skip to content

Pluggable IOD layer, multi-root picker, residual prefilter, and IAS15 min-dt floor#324

Open
matthewholman wants to merge 6 commits into
mainfrom
feat/iod-picker-pipeline
Open

Pluggable IOD layer, multi-root picker, residual prefilter, and IAS15 min-dt floor#324
matthewholman wants to merge 6 commits into
mainfrom
feat/iod-picker-pipeline

Conversation

@matthewholman
Copy link
Copy Markdown
Collaborator

@matthewholman matthewholman commented May 12, 2026

Refactors do_fit's root-handling and LM machinery into a coherent pipeline:

This was developed iteratively with Claude code.

  1. Pluggable IOD layer (src/layup/iod.py)

    • Registry-based: register_iod(name, callable), get_iod(name),
      iod_methods(). An IOD is any
      (observations, seq) -> list[FitResult] | None.
    • Gauss ships as the default. Future IOD methods (Lambert,
      motion-rate prior, etc.) drop in without modifying do_fit.
    • do_fit(..., iod=...) accepts either a registered name or a
      callable directly.
  2. Multi-root picker — LM from every IOD candidate at
    screen_iter_max=80, picks smallest-χ² converged root with
    r > 0.3 AU. Falls back to full_iter_max=100 if nothing converges
    at the cheap tier. Closes the historical pitfall where do_fit
    committed to solns[0] (largest r) even when that was a phantom
    outer-SS solution for an inner-SS target.

  3. Residual prefilter — drops candidates whose predicted
    positions miss the observations by more than threshold_sigma=1000σ
    before running LM. Candidates whose inertial trajectory passes
    within 0.1 AU of the observer are passed through unfiltered (real
    NEO close approaches need their own integrator-budget solution,
    addressed by Update README.md #4 below).

  4. IAS15 min-dt floorsim->ri_ias15.min_dt is floored at
    1e-3 days during the picker LM loop. Caps integrator wallclock on
    phantom roots whose trajectories pass close to Earth (1000× speedup
    on diagnostic/scan, bit-identical fits).

Validation (diagnostic/scan, 98 cases)

Cartesian-LM via picker is 400-1000× faster on close-encounter cases
and strictly more accurate on cases where Gauss-IOD's solns[0] was wrong:

case pre-picker dr picker dr change
classical_42AU_arc_002.00d 1.57e10 km 7.22e8 km −95%
centaur_25AU_arc_002.00d 2.38e9 km 3.68e8 km −85%

Test suite went from ~9 min to ~3 s wallclock as a side effect.

Dependencies

  • Stacks on feat/pybind11-eigen-includes (PR #N). Once that lands,
    this PR's diff collapses to its 4 new commits.

Review Checklist for Source Code Changes

  • Does pip install still work?
  • Have you written a unit test for any new functions? — tests/layup/test_iod_picker.py (8 tests)
  • Do all the units tests run successfully?
  • Does Layup run successfully on a test set of input files/databases? — full diagnostic/scan
  • Have you used black on the files you have updated?

@matthewholman matthewholman force-pushed the feat/iod-picker-pipeline branch from 484a85b to 4bb0c43 Compare May 12, 2026 01:29
matthewholman and others added 6 commits May 14, 2026 17:30
… Python

Without these, accessing rho_hat / a_vec / d_vec from Python raises
"Unable to convert function return value to a Python type" because the
binding can't see the Eigen and std::array conversions.

This is a minimal enabler (no behavior change); the A/D-vector
correctness fix on fix/tangent-basis-vectors is independent and lives
on its own branch.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
The previous binding hardcoded iter_max = 100 inside the C++ entry
point. Promote it to a parameter (default 100, so callers see no
behavior change) and rebind with explicit pybind11 keyword args so
Python callers can reduce the LM budget for cheap screening passes.

This is the enabler the multi-root IOD picker needs: it runs LM from
each Gauss root at a small iter_max first, falls back to the full
budget only if nothing converges.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Replaces the hardcoded "gauss" dispatch in do_fit with a registry-based
lookup. layup/iod.py exposes register_iod / get_iod / iod_methods, with
gauss_iod registered at import time. do_fit's `iod` parameter now
accepts either a registered name or any callable matching the
signature (observations, seq) -> list[FitResult] | None, so callers can
drop in Lambert, motion-rate-prior, or other future IOD methods
without subclassing.

Cartesian-LM path now runs LM from every IOD candidate at a cheap
screening budget (screen_iter_max=80) and picks the smallest-χ²
converged root with heliocentric distance > min_r_helio_AU. Falls
back to the full budget (full_iter_max=100) only if no candidate
converges at the cheap tier. Memory's NEO-scan finding (3 cases moved
from ~1e9 km to ~1e5–1e6 km, 0 regressions on 45 comparable cases)
reproduces here.

BK engine path is unchanged: BK does its own internal prelim, so the
IOD candidate list is not used when do_fit routes to BK.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
The multi-root picker runs LM from every Gauss root, but most cases
have one geometrically correct root and 2-3 phantom roots (often
inner-SS solutions that LM grinds on uselessly). The pre-filter drops
those before LM ever runs.

For each candidate:
  1. Cheap algebraic bounds — r in [0.05, 1000] AU. (No bound-energy
     check: Gauss's velocity can be wildly wrong for a geometrically-
     correct root, and LM walks those to convergence.)
  2. Close-Earth-approach detection: candidates whose inertial
     trajectory passes within 0.1 AU of the observer at any obs time
     are passed through unfiltered. Full ASSIST integration would
     get stuck on the close encounter (tens of seconds per
     propagation), and a 2-body shortcut would silently mishandle
     real NEO close passes. Deferred problem; LM may also be slow
     on these, which is the same deferred issue.
  3. Residual check on every observation via full ASSIST propagation
     (Sun + planets + asteroids, all the real physics). Reject
     candidates whose predicted RA/Dec misses any observation by
     more than threshold_sigma (default 1000σ — loose enough to
     never drop the right root, tight enough that phantoms miss by
     10⁵-10⁶ σ are eliminated).

In do_fit, the filter runs once between the IOD call and the picker
loop. Falls back to the unfiltered list if no candidate survives,
because LM should still try something.

Observed effects on the smoke set:
  - apollo NEO (3 roots → 2 survivors after physical bounds; 4ms filter)
  - classical 42 AU (3 roots → all 3 pass through, the 2 near-Earth
    phantoms via the close-approach skip; 5ms filter; LM converges to
    r=41.247 AU as before)
  - mainbelt 2.5 AU (3 → 3, all near-Earth pass-through; 4ms filter)
  - centaur 25 AU previously crashed: filter eliminates the bad roots
    in a few ms instead of 60s+ ASSIST integration on phantoms

The filter is a 2-4ms operation per case. The remaining slow cases are
the ones with close-Earth phantom roots that we pass through to LM —
LM grinds on them but eventually picks the right root. That's the
deferred problem mentioned in the close-approach pass-through.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Without a floor on the IAS15 adaptive step, LM grinds for minutes on
phantom Gauss roots whose trajectories pass close to Earth: the
integrator chases ever-smaller steps to resolve the close encounter,
and on the diagnostic/scan set we saw 110s (classical_42AU) and 1617s
(centaur_25AU) per do_fit call dominated by this. Setting
sim->ri_ias15.min_dt to a non-zero floor caps that adaptive shrinking.

C++ exposes set_ias15_min_dt / get_ias15_min_dt (the floor is a
process-wide static applied at every reb_simulation_create in
compute_residuals_sequence and predict::predict). Default 0 means
no behavior change for direct callers of
run_from_vector_with_initial_guess.

do_fit's picker loop now push/pops a default of 1e-3 days (~86 s):

  case               baseline → with min_dt   speedup   chi²Δ   dr Δ
  classical_42AU     110.6s   → 0.28s         400×      0       0
  centaur_25AU      1617s    → 1.57s         1030×     0       0
  mainbelt           0.19s   → 0.19s         (no encounter; floor never engages)
  apollo NEO         0.11s   → 0.10s         (same)

Real-physics close-Earth approaches lose a few arcsec of integrator
accuracy *during the encounter itself*, but astrometric observations
are essentially never taken at closest approach, so the LM fit is
unaffected. Users who need full precision can pass
picker_ias15_min_dt=0 to do_fit, or set_ias15_min_dt(0) directly.

The IOD-picker test suite went from ~9 min to ~4 s wallclock as a
side-effect.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Pure formatting (def line wrapping, trailing commas, long-string
splits).  No behavior change; all 8 tests in
tests/layup/test_iod_picker.py still pass.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
@matthewholman matthewholman force-pushed the feat/iod-picker-pipeline branch from 4bb0c43 to 34621d9 Compare May 14, 2026 21:31
@matthewholman matthewholman requested a review from maxwest-uw May 14, 2026 21:38
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