sarah-quinones/faer-rs

`selfadjoint_eigendecomposition` sometimes produces incorrect results

jbncode opened this issue · 2 comments

With certain matrices, Mat::selfadjoint_eigendecomposition appears to produce incorrect results (faer 0.18.2). In the cases I've looked at it seems like the order of the eigenvalues is not consistent with the order of the eigenvectors.

Below is code that reproduces the bug with the simplest matrix I could find that demonstrates the effect. Changing the value of n makes the issue appear and disappear (e.g. n = 32 and n = 34 seem fine, but not n = 33).

use faer::{Mat, Side, solvers::SolverCore};

fn main() {
    let n = 33;

    let mut a = Mat::zeros(n, n);
    for i in 0 .. n {
        a[(i, i)] = (n - i) as f64;
    }

    let eig = a.selfadjoint_eigendecomposition(Side::Upper);
    let a_reconstructed = eig.reconstruct();

    for i in 0 .. n {
        println!("{:20} {:20} {:20}",
            a[(i, i)],
            a_reconstructed[(i, i)],
            a_reconstructed[(i, i)] - a[(i, i)]);
    }
}

thanks a lot for finding this! i pushed a commit to main that should fix the issue

Thank you!