Pvar-dissip: NonInertialFrameForce rectangular Jacobians in C#935
Open
jobovy wants to merge 1 commit into
Open
Pvar-dissip: NonInertialFrameForce rectangular Jacobians in C#935jobovy wants to merge 1 commit into
jobovy wants to merge 1 commit into
Conversation
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## feat/pvar-dissip-struct #935 +/- ##
============================================================
- Coverage 99.92% 61.89% -38.03%
============================================================
Files 225 225
Lines 36140 36325 +185
Branches 778 810 +32
============================================================
- Hits 36113 22485 -13628
- Misses 27 13840 +13813 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
43451e7 to
85497c1
Compare
76fd6b8 to
640807e
Compare
85497c1 to
fa982b7
Compare
640807e to
093e2fb
Compare
fa982b7 to
862c9e5
Compare
093e2fb to
8b4bc77
Compare
862c9e5 to
2e23404
Compare
Wire the rectangular force Jacobian (dF/dx, dF/dv) of NonInertialFrameForce into the RectDissipativeForceJacobian machinery of the 3D variational equations (Orbit.integrate_dxdv): the frame force F = -2 Omega x (v+v0) - Omega x (Omega x [r+x0]) - Omegadot x [r+x0] - a0(t) is linear in (r, v), so the Jacobian blocks are EXACT closed forms dF/dv = -2 [Omega]_x dF/dx = |Omega|^2 I - Omega Omega^T - [Omegadot]_x with the translation terms (x0, v0, a0) contributing zero -- for EVERY configuration the C force supports: constant scalar Omega about z, constant vector Omega, constant Omegadot (Omega + Omegadot t), and Omega as function(s) of time through both the tfunc callbacks (pot_type 39) and the cinterp GSL splines (pot_type 45; Omegadot as the spline derivative, time queries clamped to [tmin,tmax] exactly as in the force; Omega/Omegadot and |Omega|^2 evaluated branch-for-branch as in the force). No configuration is unwireable, so hasC_dxdv3d=True unconditionally on the Python class (aggregated by CompositePotential as before). Because dF/dv = -2 [Omega]_x is antisymmetric, tr(dF/dv) = 0: rotating (even arbitrarily time-dependent, translating) frames are phase-volume PRESERVING, det M(t) = 1 -- exercising the velocity block of the variational Jacobian nontrivially, unlike the contracting friction case. Validation (orbit-level tests in tests/test_orbit.py, mirroring the Chandrasekhar dissipative dxdv battery; they validate the C Jacobian end-to-end against code paths independent of it): - FD-of-the-flow STM test (MWPotential2014 in a spinning-up frame, fixed- args path, and in a vector-Omega(t) + translating frame via the cinterp splines): all 6 columns < 2.7e-6 (tolerance 1e-4). - KEY physics test: det M(t) = 1 along the whole orbit to <= 9.2e-11 (tolerance 1e-8) for both flow configurations -- asserted directly on the STM determinant, no trace integral needed since tr(dF/dv) = 0 exactly. - Exact rotating-vs-inertial cross-check (Plummer sphere, constant vector Omega): r(t) = R(t)^T x(t) to 7.2e-12 and M_rot(t) = T(t) M_in(t) T(0)^-1 with T = [[R^T,0],[-[Omega]_x R^T,R^T]], R(t) = expm([Omega]_x t), to 2.4e-11 (tolerance 1e-9) -- validating all signs/factors of both Jacobian blocks at integrator precision, far beyond the 1e-4 FD-of-flow check. - Flag/gate test: hasC_dxdv3d True + CompositePotential aggregation; the pure-Python integrate_dxdv methods raise; a forced hasC_dxdv3d=False falls back to odeint with a warning and then raises (never silently wrong). Regression: pytest tests/test_orbit.py -k "dxdv or liouville" passes 130 (= 126 before + 4 new); tests/test_noninertial.py 53 passed unchanged. The numpy paths are untouched: the C force and the pure-Python force are unchanged (the new C Jacobian is only called by the dxdv path) and the only Python changes are the hasC_dxdv3d flag, docstrings, and error messages. Co-authored-by: Claude Fable 5 <noreply@anthropic.com>
8b4bc77 to
5849419
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.
Stacked on #933. Exact closed-form frame Jacobians: dF/dv = −2[Ω]ₓ, dF/dx = −[Ω]ₓ[Ω]ₓ − [Ω̇]ₓ (all Omega configurations the C force supports, incl. time-dependent/cinterp; translation terms zero).
Key physics test: tr(∂F/∂v)=0 ⇒ the rotating frame is phase-volume preserving — det M(t)=1 holds despite velocity dependence, exercising the velocity block in a way friction cannot. Plus Jacobian-vs-FD and FD-of-flow. Regression 126→131 (=+5 new); test_noninertial.py 53 unchanged. Verifier: pass.
numpy-path changes: none (dxdv path previously unavailable).
🤖 Generated with Claude Code