Finish breaking the decompositions: deprecate chol and chol! #27249

mbauman · 2018-05-24T18:29:53Z

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) and chol(::UniformScaling) methods since we can't construct Choleskys for those — I believe that using sqrt is just fine here and cannot construct a case where this simplification would change any existing support.

KristofferC · 2018-05-24T18:32:03Z

stdlib/LinearAlgebra/src/cholesky.jl

@@ -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


thecholesky

fredrikekre · 2018-05-24T19:01:54Z

stdlib/LinearAlgebra/src/cholesky.jl

 if d == :U
- return UpperTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors'))
+ return @assertposdef UpperTriangular(Symbol(Cuplo) == d ? Cfactors : copy(Cfactors')) info


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?

I might err to following the existing *fact methods' behavior and forgoing these checks. Thoughts @andreasnoack?

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.

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?

Sacha0 · 2018-05-24T19:24:31Z

stdlib/LinearAlgebra/src/deprecated.jl

+# 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))


x -> A or vice versa?

Sacha0 · 2018-05-24T19:28:57Z

stdlib/LinearAlgebra/src/cholesky.jl

@@ -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


Errant comma following "cholesky"?

mschauer · 2018-05-24T19:37:00Z

Why should chol/chol! disappear? The factorization represents the original matrix, the root is of interests of its own. For example in Gaussian computations the cholesky root of the covariance Σ is applied to Gaussian white noise, with typical code:

chol(Σ)'*randn(d,n)

For a matrix A I'd say that chol(A) is as useful and justified as sqrt(A) (both are possible generalizations of sqrt(::Number). You would most likely not want to write
rootfact(A).root for sqrt(A).

Sacha0

Modulo whether to add the checks in the last commit, looks great! Much thanks Fredrik and Matt! :)

Sacha0 · 2018-05-25T20:17:00Z

(CI shows fourteen failures in the LinearAlgebra/cholesky test suite.)

StefanKarpinski · 2018-05-25T20:22:41Z

@mschauer has a point here. Perhaps we should keep chol and just make cholesky the name for the function that returns the factorization object.

Sacha0 · 2018-05-25T20:28:35Z

On the other hand, writing cholesky[!](A).U in place of chol[!](A) where necessary does not seem onerous, and there is a meaningful cost to retaining four functions instead of two, in both cognitive overhead / confusion for users and maintenance costs for developers. The tradeoffs strike me as favoring the two-function choice. ~~Another option is to generalize sqrt~~. Best!

Sacha0 · 2018-05-25T20:55:58Z

Observations from a chat with @mbauman: One troubling consideration is that multiple notions of matrix root exist for positive definite operators. Moreover, even fixing the definition of matrix root to traditional Cholesky factor, there remains a choice of factor: upper or lower by default? In other words, A == C'C or CC'?

StefanKarpinski · 2018-05-25T21:04:39Z

We should maybe have a chat with @alanedelman about the best approach.

andreasnoack · 2018-05-25T21:26:57Z

We settled upper a while ago.

StefanKarpinski · 2018-05-25T21:34:04Z

No one is proposing changing it. The question is really whether chol(A) is a useful and meaningful to want to compute on its own or if it really only makes sense as part of a factorization.

Sacha0 · 2018-05-26T03:04:42Z

stdlib/LinearAlgebra/test/cholesky.jl

@@ -87,7 +86,6 @@ end

 apos = apd[1,1]
 @test all(x -> x ≈ √apos, cholesky(apos).factors)
- @test_throws PosDefException cholesky(-one(eltya))


(cholesky returns a failed decomposition rather than throwing a PosDefException, hence the test failure following translation of the call from chol to cholesky.)

Sacha0 · 2018-05-26T03:05:25Z

stdlib/LinearAlgebra/test/cholesky.jl

@@ -64,7 +64,6 @@ end

 @testset "throw for non-square input" begin
 A = rand(eltya, 2, 3)
- @test_throws DimensionMismatch LinearAlgebra.chol!(A)


(This straggler generated a depwarn, and is now redundant with the test two lines below.)

Sacha0 · 2018-05-26T03:05:50Z

stdlib/LinearAlgebra/test/cholesky.jl

@@ -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))


(Likewise here, this straggler generated a depwarn, and is now redundant with the test a line above.)

Sacha0 · 2018-05-28T19:12:48Z

Now passing CI, this pull request needs only a decision :).

ViralBShah · 2018-05-28T19:57:57Z

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.λ)))

…], qrfact[!]

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.

Sacha0 · 2018-05-28T20:20:22Z

We should get travis to pass, and that will need a rebasing to master and resolving some merge conflicts.

Done! Time for another spin of the CI wheel of fortune.

andyferris · 2018-05-29T00:38:13Z

(Edit: I misunderstood one of the comments above).

I don't have a strong opinion on whether to deprecate chol or not... for eigen I can still use destructuring to get the components, like before, which is cute. (Would (c,) = cholesky(matrix) work?)

andyferris · 2018-05-29T00:54:22Z

Thinking about this some more, getproperty overloads for cholesky(matrix).upper and cholesky(matrix).lower would be pretty neat. Though longer to type than chol(matrix) it would be quite unambiguous, explicit and readable code (i.e. no confusion of whether chol returns an upper- or lower-triangular matrix, someone less mathematical can Google "cholesky", etc)

Sacha0 · 2018-05-29T01:00:10Z

(AV i686 timed out at three hours, though after passing the main test suite.)

JeffBezanson · 2018-05-29T04:23:50Z

Ready to merge?

fredrikekre · 2018-05-29T06:37:38Z

Did we reach a conclusion to #27249 (comment) ?

fredrikekre · 2018-05-29T06:38:48Z

Thinking about this some more, getproperty overloads for cholesky(matrix).upper and cholesky(matrix).lower would be pretty neat.

This already exists and is what we deprecate chol to in this PR.

mbauman · 2018-05-29T18:20:25Z

So there are two decisions remaining here:

Do we want to preserve chol because, unlike the other decompositions, it can be viewed as a matrix transformation unto itself? I know Sacha is against this, and I'm leaning against, too, on the basis of the difference between chol(A) \ X and cholesky(A) \ X.
Should we check for success in the property accessors? Looks like folks in the hidden thread were leaning against doing so. Finish breaking the decompositions: deprecate chol and chol! #27249 (comment). I don't have as good of a sense on this point.

mschauer · 2018-05-29T18:45:45Z

I always had the impression that Julia favours functional forms length(x) over element access x.len, and that U = cholesky(matrix).U is a second class interface. (Thats a style question and my impression could be off.) But I think U = chol(A) or maybe U = cholesky(A).upper would be nicer in this regard, I am happy as long as U has an interface.

fredrikekre · 2018-05-29T20:33:54Z

Why do you prefer U = cholesky(A).upper over U = cholesky(matrix).U?

Sacha0 · 2018-05-29T20:34:27Z

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!

simonbyrne · 2018-05-29T20:48:30Z

In either case, the changes vs 0.6 should be noted in NEWS.md.

mschauer · 2018-05-29T21:04:04Z

Why do you prefer U = cholesky(A).upper over U = cholesky(matrix).U?

I think this makes it more visible that .upper is an interface and not an ad hoc access. Compare upper(cholesky(A)) with U(cholesky(A)).

andyferris · 2018-05-29T23:43:03Z

Why do you prefer U = cholesky(A).upper over U = cholesky(matrix).U?

Plus, it's way clearer to read (if you're not familiar with all of our idioms).

ViralBShah · 2018-05-29T23:51:50Z

Should be good to merge, CI-wise.

mbauman · 2018-05-30T15:35:56Z

Alright, sounds like we're in agreement here. Let's merge this and then tackle the new error check flag in a followup PR.

andreasnoack · 2018-05-30T15:46:55Z

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 U vs upper and throwing behavior in separate PRs.

mschauer · 2018-08-30T10:37:03Z

Is this now locked in at .U?

mbauman added domain:linear algebra Linear algebra kind:deprecation This change introduces or involves a deprecation labels May 24, 2018

KristofferC reviewed May 24, 2018

View reviewed changes

mbauman force-pushed the fe/breakmore branch 2 times, most recently from 13442cd to 329e605 Compare May 24, 2018 18:50

fredrikekre reviewed May 24, 2018

View reviewed changes

Sacha0 reviewed May 24, 2018

View reviewed changes

Sacha0 reviewed May 26, 2018

View reviewed changes

martinholters mentioned this pull request May 28, 2018

Make Cholesky factorization compatible with 0.7 JuliaArrays/StaticArrays.jl#398

Merged

fredrikekre and others added 6 commits May 28, 2018 13:18

deprecate chol(A) to cholesky(A).U

64b9fe3

deprecate chol:

71ab63a

- chol(A::AbstractMatrix) -> cholesky(A).U - chol(A::Number) -> sqrt(A) - chol(A::UniformScaling) -> UniformScaling(sqrt(J.λ)))

deprecate chol!(A) to cholesky!(A).U

7112e64

update Val and uplo deprecations for cholfact[!], bkfact[!], lufact[!…

2e46f08

…], qrfact[!]

RFC: check the info code on getproperty

bed991d

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.

Remove latent chol[!] tests to address CI failures.

d6bc505

Sacha0 force-pushed the fe/breakmore branch from 0053809 to d6bc505 Compare May 28, 2018 20:19

Sacha0 mentioned this pull request May 29, 2018

AppVeyor CI failures in starting workers #27274

Closed

JeffBezanson added this to the 0.7 milestone May 29, 2018

andreasnoack merged commit 36bcbb4 into master May 30, 2018

Sacha0 deleted the fe/breakmore branch May 30, 2018 15:53

Finish breaking the decompositions: deprecate chol and chol! #27249

Finish breaking the decompositions: deprecate chol and chol! #27249

Conversation

mbauman commented May 24, 2018

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Sacha0 May 25, 2018 • edited Loading

Choose a reason for hiding this comment

Sacha0 May 24, 2018 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mschauer commented May 24, 2018 • edited Loading

Sacha0 left a comment

Choose a reason for hiding this comment

Sacha0 commented May 25, 2018

StefanKarpinski commented May 25, 2018

Sacha0 commented May 25, 2018 • edited Loading

Sacha0 commented May 25, 2018

StefanKarpinski commented May 25, 2018

andreasnoack commented May 25, 2018 via email • edited by StefanKarpinski Loading

StefanKarpinski commented May 25, 2018

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Sacha0 commented May 28, 2018

ViralBShah commented May 28, 2018 • edited Loading

Sacha0 commented May 28, 2018

andyferris commented May 29, 2018 • edited Loading

andyferris commented May 29, 2018

Sacha0 commented May 29, 2018

JeffBezanson commented May 29, 2018

fredrikekre commented May 29, 2018

fredrikekre commented May 29, 2018

mbauman commented May 29, 2018

mschauer commented May 29, 2018

fredrikekre commented May 29, 2018

Sacha0 commented May 29, 2018

simonbyrne commented May 29, 2018

mschauer commented May 29, 2018

andyferris commented May 29, 2018

ViralBShah commented May 29, 2018

mbauman commented May 30, 2018

andreasnoack commented May 30, 2018

mschauer commented Aug 30, 2018

Sacha0 May 25, 2018 •

edited

Loading

Sacha0 May 24, 2018 •

edited

Loading

mschauer commented May 24, 2018 •

edited

Loading

Sacha0 commented May 25, 2018 •

edited

Loading

andreasnoack commented May 25, 2018 via email •

edited by StefanKarpinski

Loading

ViralBShah commented May 28, 2018 •

edited

Loading

andyferris commented May 29, 2018 •

edited

Loading