[shinaoka]

Closes #112. Closes #144. Closes #145.

This PR combines #148 (#112), #151 (#144) and #152 (#145), as the maintainer directed on 2026-10-29 ("関連する修正は単一のPRにまとめてもいいよ"). All three change the nonsymmetric EVD path and conflicted with each other. The original commits are kept with their authorship and messages; only the upstream report numbers are renumbered (faer-029, -030, -031), because all three had taken faer-029. Each fix's regression tests, PATCHSET row and upstream report are kept.

1. Input scaling as ?GEEV (#112, from #148)

Stage-by-stage comparison with cgeev (c32 kinds 6/7, n = 63–65):

stage	as is	scaled to unit size
Hessenberg backward error / eps	16–20	16–20
complex Schur backward error / eps	226–258	67–93

?GEEV scales its input into [smlnum, bignum] (smlnum = sqrt(safmin)/eps) before the reduction; evd_imp did not.

Where the gap grows. A 372-matrix sweep against LAPACK ?geev (scipy/OpenBLAS) covered f32/f64/c32/c64, n = 16–257, the random, Grcar, near-Jordan and graded families, and scales up to the range limits. On random, Grcar and near-Jordan matrices:

• The GET22 ratio grew from about 2^-90 (32-bit).
• At 2^-100 (32-bit) and 2^-960 (64-bit) and below, it was 1e4 to 2e15, against LAPACK's 2–70.
• f32 also gave NoConvergence, and two inputs did not return.
• faer 0.24.4 returns NaN/inf eigenvectors with Ok there.

Fix. evd_imp scales by a power of two to just inside the bound, and scales the eigenvalues back.

• The scaling itself does not round.
• The bound is tested on the largest |re| or |im|, because faer's complex modulus is inf for 2^1023 + 0i and zero for 2^-1074 + 0i. This was a post-review BLOCKER on evd: scale the input into range as ?GEEV does (#112) #148.

Tests.

• geev_scaling_tests::evd_scales_its_input_into_range fails without the scaling (f32 2^-120: 1.7e6).
• evd_scaling_at_the_complex_range_ends fails with a modulus-based scan.

The known-failures entry for #112 is removed.

Upstream report: faer-029-evd-input-not-scaled.

2. The real Schur deflation-window sort always ends (#144, from #151)

aggressive_early_deflation's bubble sort, ported from dlaqr3, had two faults:

• It sorted the deflated part ns..jw instead of the undeflated part infqr..ns.
• It set the next window from the last swap.

Since schur_swap reported success for a rejected 2×2 swap, the sort compared the same pair forever. Mat::eigen of f32 J_n + 1e-6 G at 2^-100 hung at n = 128 and 257.

Fix. The sort follows dlaqr3:

• It sorts the undeflated part.
• Each pass ends one before where the previous pass stopped, so the sort ends after at most ns passes.
• It uses schur_swap_checked, and a rejected swap leaves the pair in place.

Test: aed_sort_terminates_when_a_swap_fails, run through hessenberg_schur_real (which, like ?HSEQR, does not scale its input) under a 120 s timeout. On main it times out at n = 128.

Upstream report: faer-030-real-schur-aed-sort-hang. Neither size returns within 60 s on 0.24.4.

3. The 2×2 block solves of the real eigenvector back substitution follow ?LALN2 (#145, from #152)

Stage-by-stage on f64 graded n = 64:

• The Schur backward error is 8 eps|A|.
• The eigenvectors of faer's own Schur form T have a residual of 3.4e7, where dtrevc3 on the same T gives 2.2e-4.
• All 24 failing blocks are complex pairs.

solve_cplx_shifted_upper_quasi_triangular_system solved each shifted 2×2 block through its inverse, and compared the determinant (a product of two entries) with the entry-sized clamp. So every block with entries below sqrt(clamp) got a floored determinant with unchanged numerators.

Fix. Both 2×2 solves use ?LALN2's method:

• an all-small block is taken as clamp * I;
• otherwise, Gaussian elimination with complete pivoting, with the second pivot raised to the clamp.

Balancing (?GEBAL) is not the cause. Unbalanced dgeev gives 0.58 on the same matrix, and the fixed faer gives 0.45. Adopting balancing is a separate decision.

Tests:

• lapack_nonsymmetric_evd_graded_pairs_keep_their_residual (f64 before: 5.9e8 at n = 16; 1.5e12 on a coupled 4×4);
• ..._pair_solve_with_a_small_determinant;
• ..._real_shift_block_below_the_square_root_of_min_positive (NoConvergence before).

Upstream report: faer-031-real-evd-pair-clamp. 0.24.4 gives 8.97e4 on a 4×4 where dgeev gives 9.6e-14.

Verification of the combined branch (primerose, 4 cores)

• scripts/run-lapack-validation.py: every family exits 0, including nonsymmetric_evd 40/40, schur 12/12, hessenberg 8/8, svd 88/88, self_adjoint_evd 20/20, lu, qr, cholesky, indefinite, triangular and lstsq. The known-failures list is empty.
• cargo test -p t4a-faer --lib linalg::evd: 119 passed.
• The LAPACK validation: the nonsymmetric EVD fixtures exceed the GET22-style residual threshold (4 dtypes) #112 sweep inputs, rerun on this branch:

inputs	before	this branch	LAPACK
real graded, n = 16–257 (D A D^-1, D = 1e±6)	up to 2.6e9	0.20–0.87	below 5e-4 (balanced)
complex graded	up to 1.2	0.19–1.19	below 5e-4
tiny (2^-100 32-bit, 2^-900 64-bit), all types and families	up to 4e6, plus NoConvergence and hangs	2.4–78; at most 1.68× LAPACK	1.9–102
f32 near-Jordan n = 128 / 257 at 2^-100 (hung before)	no return	42 / 60	44 / 69

• Each fix had its own independent post-review (gpt-6.1-sol) on its original PR, and the findings were fixed there. The conflicts were only in the ledgers (PATCHSET rows, the upstream-reports table and report numbers); the code merged without a conflict.
• Performance: not timed. The changes are one pass over the input, a sort bound, and a 2×2 solve.

LAPACK sources were read only; no code was copied. Nothing was posted upstream.

🤖 Generated with Claude Code