search

RalphAS / GenericSchur.jl

Julia package for Schur decomposition of matrices with generic element types
Other
26 stars 4 forks source link

better Hessenberg solvers #2

Closed stevengj closed 3 years ago

stevengj commented 5 years ago

I noticed that you implemented Hessenberg solvers in https://github.com/RalphAS/GenericSchur.jl/blob/a361308be3b52faee82a9db30ea87d40d82fb0c4/src/hessenberg.jl

In JuliaLang/julia#31853 we are adding such solvers to LinearAlgebra — the algorithm we are using seems to be significantly faster (5-10x) than what you are doing now because it combines the RQ factorization and the triangular backsolve, and is able to do them both without modifying H.

In the near future, if you already have an upper-Hessenberg matrix H, you will be able to do Hessenberg(H) \ b to use the LinearAlgebra solver. If you need any additional functionality, please comment.

RalphAS commented 5 years ago

I'm very pleased to see this development. I have also been wondering how Hessenberg types and methods should best be organized, and I think you are on the right track. Just yesterday I happily learnt that the the algorithm you implemented is stable so that I would gladly take advantage of it. I'd like to see negligible overhead for the simple case of zero shifts and trivial reflectors, but I gather you've already thought about that.

stevengj commented 5 years ago

Yes, the algorithm I implemented is backwards stable because it is equivalent to Givens RQ followed by backsubstitution. The code already has negligible overhead for zero shifts, and soon we will support trivial reflectors efficiently as well.

stevengj commented 5 years ago

In the latest PR, then if you already have an upper-Hessenberg matrix H (trivial reflectors), you can do UpperHessenberg(H) \ b to use the fast solver.

RalphAS commented 3 years ago

Closing since I dropped the inferior local method and have adopted the StdLib Hessenberg structs.