Skip to content
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

\ for square singular matrices #939

Open
mikmoore opened this issue Jul 25, 2022 · 6 comments
Open

\ for square singular matrices #939

mikmoore opened this issue Jul 25, 2022 · 6 comments
Labels
speculative Whether the change will be implemented is speculative

Comments

@mikmoore
Copy link
Contributor

From this Discourse discussion.

using LinearAlgebra

julia> ones(3,2) \ ones(3)
2-element Vector{Float64}:
 0.5
 0.4999999999999999

julia> ones(3,3) \ ones(3)
ERROR: SingularException(2)

julia> ones(3,4) \ ones(3)
4-element Vector{Float64}:
 0.2500000000000001
 0.25
 0.25
 0.25

The issue here is that \ applies the LU decomposition to a square input but a pivoted QR to a rectangular input. The LU solver cannot handle singular inputs while QR can.

julia> qr(ones(3,3),ColumnNorm()) \ ones(3)
3-element Vector{Float64}:
 0.33333333333333337
 0.33333333333333337
 0.33333333333333337

A naive fix might just be to check the lu(A) call within \ and re-try with qr(A,ColumnNorm()) if it fails via singularity.

@fredrikekre
Copy link
Member

I would much rather have an error thrown at me when the matrix is singular. If you want a least-squares solution you should, IMO, use qr explicitly.

@fredrikekre fredrikekre added the speculative Whether the change will be implemented is speculative label Jul 25, 2022
@mikmoore
Copy link
Contributor Author

I agree that a SingularException is often useful - and that sidestepping it could cause plenty of chaos. However, I don't think it's terribly consistent that this works for all systems EXCEPT those that happen to be square.

At the very least, this might be a good target for an error hint indicating that a manual qr call can provide a minimum-norm least-squares solution.

@fredrikekre
Copy link
Member

However, I don't think it's terribly consistent that this works for all systems EXCEPT those that happen to be square.

Agreed, but I would prefer \ to only apply to square matrices instead (but obviously that ship has sailed).

@stevengj
Copy link
Member

stevengj commented Aug 2, 2022

I don't think it is feasible to change the behavior of \ at this point. A more helpful error message is very possible, however.

However, we might need to add an additional field to SingularException to add more context, as this error can be thrown in other contexts than ldiv!/rdiv! (e.g. in inv or cond) where a suggestion to use qr would not be appropriate, similar to what we currently do with throw_complex_domainerror.

@martinholters
Copy link
Member

Ref. #554.

@stevengj
Copy link
Member

Should we close as a duplicate of 28827?

@KristofferC KristofferC transferred this issue from JuliaLang/julia Nov 26, 2024
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
speculative Whether the change will be implemented is speculative
Projects
None yet
Development

No branches or pull requests

4 participants