-
-
Notifications
You must be signed in to change notification settings - Fork 5.5k
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
Finish breaking the decompositions: deprecate chol and chol! #27249
Conversation
stdlib/LinearAlgebra/src/cholesky.jl
Outdated
@@ -4,7 +4,7 @@ | |||
# Cholesky Factorization # | |||
########################## | |||
|
|||
# The dispatch structure in the chol!, chol, cholesky, and cholesky! methods is a bit | |||
# The dispatch structure in thecholesky, and cholesky! methods is a bit |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
thecholesky
13442cd
to
329e605
Compare
if d == :U | ||
return UpperTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors')) | ||
return @assertposdef UpperTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors')) info |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The reason I added this last commit is that chol(A)
used to error when the factorization failed. Now that chol(A)
is replaced with cholesky(A).U
we need this check here to keep the same behavior. If you really want to access the failed factor you can use getfield(A, :U)
.
On that note, we should probably do the same thing for lu(A).L
etc, although lu(A)
happily returned before without checking if the factorization failed. Thoughts?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I might err to following the existing *fact
methods' behavior and forgoing these checks. Thoughts @andreasnoack?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
We should probably check. It is a bit unfortunate that expert users would need to use getfield
but the alternative is that we'd have to teach all linear algebra users to always use issuccess
before doing anything with factorizations.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I assume that users call cholesky either: 1) knowing that their matrices are positive definite; or 2) not being certain that their matrices are positive definite, and in that case checking for success after the call. Otherwise I assume they would use another, more appropriate decomposition or, absent knowledge of which decomposition to use, call factorize
? And as such not checking in cholfact
seems to have been fine so far?
# deprecate chol to cholesky with getproperty | ||
@deprecate(chol(A::RealHermSymComplexHerm), cholesky(A).U) | ||
@deprecate(chol(A::AbstractMatrix), cholesky(A).U) | ||
@deprecate(chol(x::Number, args...), sqrt(A)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
x
-> A
or vice versa?
stdlib/LinearAlgebra/src/cholesky.jl
Outdated
@@ -4,7 +4,7 @@ | |||
# Cholesky Factorization # | |||
########################## | |||
|
|||
# The dispatch structure in the chol!, cholesky, and cholesky! methods is a bit | |||
# The dispatch structure in the cholesky, and cholesky! methods is a bit |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Errant comma following "cholesky"?
Why should chol(Σ)'*randn(d,n) For a matrix |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Modulo whether to add the checks in the last commit, looks great! Much thanks Fredrik and Matt! :)
(CI shows fourteen failures in the |
@mschauer has a point here. Perhaps we should keep |
On the other hand, writing |
Observations from a chat with |
We should maybe have a chat with @alanedelman about the best approach. |
We settled upper a while ago.
|
No one is proposing changing it. The question is really whether |
@@ -87,7 +86,6 @@ end | |||
|
|||
apos = apd[1,1] | |||
@test all(x -> x ≈ √apos, cholesky(apos).factors) | |||
@test_throws PosDefException cholesky(-one(eltya)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
(cholesky
returns a failed decomposition rather than throwing a PosDefException
, hence the test failure following translation of the call from chol
to cholesky
.)
@@ -64,7 +64,6 @@ end | |||
|
|||
@testset "throw for non-square input" begin | |||
A = rand(eltya, 2, 3) | |||
@test_throws DimensionMismatch LinearAlgebra.chol!(A) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
(This straggler generated a depwarn, and is now redundant with the test two lines below.)
@@ -263,7 +261,6 @@ end | |||
for A in (R, C) | |||
@test !LinearAlgebra.issuccess(cholesky(A)) | |||
@test !LinearAlgebra.issuccess(cholesky!(copy(A))) | |||
@test_throws PosDefException LinearAlgebra.chol!(copy(A)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
(Likewise here, this straggler generated a depwarn, and is now redundant with the test a line above.)
Now passing CI, this pull request needs only a decision :). |
We should get travis to pass, and that will need a rebasing to master and resolving some merge conflicts. |
- chol(A::AbstractMatrix) -> cholesky(A).U - chol(A::Number) -> sqrt(A) - chol(A::UniformScaling) -> UniformScaling(sqrt(J.λ)))
now that we have replaced e.g. chol(A) with cholesky(A).U we don't check for positive definiteness anymore. This adds a check for that.
Done! Time for another spin of the CI wheel of fortune. |
(Edit: I misunderstood one of the comments above). I don't have a strong opinion on whether to deprecate |
Thinking about this some more, |
(AV i686 timed out at three hours, though after passing the main test suite.) |
Ready to merge? |
Did we reach a conclusion to #27249 (comment) ? |
This already exists and is what we deprecate |
So there are two decisions remaining here:
|
I always had the impression that Julia favours functional forms |
Why do you prefer |
Regarding Matt's second bullet, a slack/#linalg discussion generated another potential approach: In decomposition methods that can fail, throw on failure by default and provide an escape hatch for expert users. In terms of user experience, this approach prevents non-expert users from contacting ill-formed decompositions, while allowing expert users to work also with ill-formed decompositions without too much hoop-jumping. Additionally, this approach has the software engineering upside of allowing failure checks to be consolidated in the decomposition methods and removed from the wider-spread and more-numerous use sites. If we went this route, I would advocate that we merge this pull request with the last commit in place, and change the above post-alpha; the impact of that change (on functioning, non-test code) should be minimal. Best! |
In either case, the changes vs 0.6 should be noted in |
I think this makes it more visible that |
Plus, it's way clearer to read (if you're not familiar with all of our idioms). |
Should be good to merge, CI-wise. |
Alright, sounds like we're in agreement here. Let's merge this and then tackle the new error check flag in a followup PR. |
I agree. I also support the proposal in #27249 (comment) and there so far no objections, so I'll merge and suggest that we sort out the |
Is this now locked in at |
This branch simply rebases @fredrikekre's work on top of the work done in #27212.
There was a slight question about what to do with the
chol(::Number)
andchol(::UniformScaling)
methods since we can't constructCholesky
s for those — I believe that usingsqrt
is just fine here and cannot construct a case where this simplification would change any existing support.