On looking deeper, I see that LAPACK also has a Jacobi method for the SVD (as mentioned in Calling Lapack's Jacobi SVD - dgesvj) This approach is not mentioned in the Users’ Guide, but (often?) gives singular values with elementwise relative accuracy close to unit rounding error. It works for all of the s.v. for the companion matrices described above until N>64 where things are so bad that NaNs are generated.
1 Like