-
-
Notifications
You must be signed in to change notification settings - Fork 530
Implement absolute threshold for early deflation in Schur algorithm #1565
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Conversation
|
Hi everyone, could I get a feedback if this is implemented correctly or if I should have used another solution? I would love to tackle other issues in this repository :) |
|
@alexhroom I am tagging you since you handled the original issue report |
|
hey there, just approved the CI runs for you. I've got a pretty big stockpile of matrices that crash with the issue so next week I'll do some more robust checking that the issue's been sorted. Thanks for figuring it out by the way, I had a go myself ages ago but couldn't quite figure out the problem... |
Would love to help you on that :) I fixed the formatting issue preventing the pipeline to build also, first time contributor mistake |
|
I've had a play around and it's looking good. only thing is if you have a source/citation for where the threshold criterion comes from? |
|
I looked at LAPACK DLAHQR subroutine, and they use SMLNUM as their absolute criteria : SMLNUM is defined as the smallest possible value you can represent using the machine precision * number of rows of the matrix / epsilon. The issue is that depends on the machine precision and I didn't find an easy way to do that in rust. So I choose eps * eps, because the relative test for deflation was already eps * ( t[n,n] + t[m,m] ) |
|
great - could I ask you to add that citation as a comment to the routine, and then we should be all good to go? |
|
Nice suggestion, comment has been added :) |
Fix #1528.
The issue was that, during the deflation step, the diagonal element of the matrix were all zero so the deflation never triggered based on this line :
t[(n, m)].clone().norm1() <= eps.clone() * (t[(n, n)].clone().norm1() + t[(m, m)].clone().norm1())added a check to compare an off-diag element to an absolute very small value
(eps.clone() * eps.clone())the eigenvalues found are correctly :
Complex { re: -0.8611363115940531, im: 0.0 },
Complex { re: 0.8611363115940531, im: 0.0 },
Complex { re: -0.33998104358485637, im: 0.0 },
Complex { re: 0.33998104358485637, im: 0.0 }