[ParadaCarleton]

test_wave.jl reports a different list of special points on each run for the same travelling-wave branch, and the point that line 198 continues as a fold is typed :hopf on some runs and :nd on others. The branch itself is the same on every run, with the same parameters and finite eigenvalues at each step, but at random steps the spectrum has one more eigenvalue of about 1e17 to 1e19. That eigenvalue is the infinite eigenvalue of the pencil (J, B) with B = diag(1, …, 1, 0): LAPACK returns its β at rounding level instead of exactly zero, so α / β is finite, with a random sign. gev(::DefaultEig, A, B, nev) dropped only non-finite values, so a positive one was counted as unstable and created a spurious bifurcation point; the existing test at test_linear.jl:720 checks only all(isfinite), which such an eigenvalue passes. This PR drops an eigenpair (λ, v) when ‖B v‖₁ ≤ n eps ‖B‖₁ ‖v‖₁: an exact pair satisfies β A v = α B v, so |β| / |α| = ‖B v‖ / ‖A v‖, and the test reads β from the pair's own eigenvector. (schur(A, B) exposes α and β, but LAPACK's gges orders them differently from the ggev call that eigen makes, and matching the two would go through λ, the value that rounding corrupts.) The new test uses the singular Q diag(1, …, 1, 0) Q' over 20 random orthogonal Q, and checks that gev returns exactly n - 1 eigenvalues and that they match those of the Schur complement of Q' J Q; the test fails with the old filter.

https://claude.ai/code/session_015sE9Latuiqsp6zSHQ1cshz

[rveltz]

Thank you for this.
I was about to reject this because of the use of opnorm as it is problematic for matrix-free or large systems. But for DefautlEig, it does not matter