Pluggable IOD layer, multi-root picker, residual prefilter, and IAS15 min-dt floor#324
Open
matthewholman wants to merge 6 commits into
Open
Pluggable IOD layer, multi-root picker, residual prefilter, and IAS15 min-dt floor#324matthewholman wants to merge 6 commits into
matthewholman wants to merge 6 commits into
Conversation
484a85b to
4bb0c43
Compare
… 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>
4bb0c43 to
34621d9
Compare
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Refactors
do_fit's root-handling and LM machinery into a coherent pipeline:This was developed iteratively with Claude code.
Pluggable IOD layer (
src/layup/iod.py)register_iod(name, callable),get_iod(name),iod_methods(). An IOD is any(observations, seq) -> list[FitResult] | None.motion-rate prior, etc.) drop in without modifying
do_fit.do_fit(..., iod=...)accepts either a registered name or acallable directly.
Multi-root picker — LM from every IOD candidate at
screen_iter_max=80, picks smallest-χ² converged root withr > 0.3 AU. Falls back to
full_iter_max=100if nothing convergesat the cheap tier. Closes the historical pitfall where
do_fitcommitted to
solns[0](largest r) even when that was a phantomouter-SS solution for an inner-SS target.
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).
IAS15 min-dt floor —
sim->ri_ias15.min_dtis floored at1e-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:classical_42AU_arc_002.00dcentaur_25AU_arc_002.00dTest suite went from ~9 min to ~3 s wallclock as a side effect.
Dependencies
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
tests/layup/test_iod_picker.py(8 tests)diagnostic/scan