mirror of
https://github.com/storytold/effectcraft.git
synced 2026-11-07 02:51:24 +00:00
46f04d1a2d
* Subspace Warp: QL eigen solver and banded mesh solve (#2319) plan() spent 43-96 % of its time in the cyclic Jacobi sym_eigen, whose absolute stop test never triggers on Gram matrices of magnitude ~1e8, and then ran a dense Cholesky on a mesh normal matrix of bandwidth ~30. Use Householder tridiagonalisation plus implicit QL inside subspace.rs (only the leading subspace is used, so basis and sign do not matter), and accumulate and factor the mesh normal equations in lower-band storage, without a Vec per data row. camtrack/linalg.rs is untouched. Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com> * Subspace Warp: replace the QL eigen solver with a clean-room symmetric QR The previous sym_eigen_ql followed a well-known published tred2/tqli listing too closely for the clean-room rules. It is replaced by sym_eigen_qr, written from the textbook method (Golub & Van Loan, Matrix Computations, ch. 8) in its own structure: - householder_tridiagonal: explicit reflectors with the symmetric rank-2 update B -= v w^T + w v^T (p = beta B v, w = p - (beta p^T v / 2) v), Q accumulated as rows of Q^T. - wilkinson_qr_step: one implicit QR step with the Wilkinson shift, Givens rotations chasing the bulge from the top of the unreduced block down. - sym_eigen_qr: deflation of negligible off-diagonal entries and the step driver, split so each function stays under the complexity limit. It returns None for non-finite input, a size mismatch or no convergence within 30 n steps; smooth_trajectories then falls back to the Jacobi sym_eigen, so it can never loop for ever or be worse than before. Tests: random symmetric matrices (n = 1..61) against Jacobi eigenvalues, residual |A v - lambda v| and V^T V = I; identity, repeated 3x3 blocks, zero, 1x1, 0x0, rows scaled from 1e-8 to 1e8 and a 1-2-1 tridiagonal; NaN, infinity and size mismatch return None; the leading-subspace projector of Gram-like low-rank matrices matches Jacobi. An ignored benchmark times plan() on a 300-frame synthetic parallax shot. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> --------- Co-authored-by: Claude Sonnet 5.5 <noreply@anthropic.com> Co-authored-by: FredApps <153173854+FredApps@users.noreply.github.com>