RFC: Extract Bunch-Kaufman factors and use them for printing #22601

andreasnoack · 2017-06-28T20:40:30Z

As part of adding tests, I've also cleaned up the signatures of bkfact similarly to the changes we have made for cholfact. That should actually go in a separate PR so I'll open that shortly.

I've also wrapped the LAPACK routine for extracting the factors when using rook pivoting but I realized that they were only added in LAPACK 3.7.0 so we'll have to wait until we bump OpenBLAS before we can enable that functionality. For now, the call to those routines are commented out and an error is thrown.

Fixes JuliaLang/LinearAlgebra.jl#441

Update: The printing looks like

julia> A = [1 2 3; 2 1 2; 3 2 1]
3×3 Array{Int64,2}:
 1  2  3
 2  1  2
 3  2  1

julia> F = bkfact(Symmetric(A, :L))
Base.LinAlg.BunchKaufman{Float64,Array{Float64,2}}
D factor:
3×3 Tridiagonal{Float64}:
 1.0  3.0    ⋅
 3.0  1.0   0.0
  ⋅   0.0  -1.0
L factor:
3×3 Base.LinAlg.UnitLowerTriangular{Float64,Array{Float64,2}}:
 1.0  0.0  0.0
 0.0  1.0  0.0
 0.5  0.5  1.0
permutation:
3-element Array{Int64,1}:
 1
 3
 2
successful: true

tkelman · 2017-06-28T20:45:48Z

test/linalg/bunchkaufman.jl


-                @testset "(Automatic) Bunch-Kaufman factor of indefinite matrix" begin
-                    bc1 = factorize(asym)
+                # check thactor factorize gives a Bunch-Kaufman


check that ?

tkelman · 2017-06-28T20:46:49Z

base/linalg/bunchkaufman.jl

+            LUD, od = LAPACK.syconv!(B.uplo, copy(B.LD), B.ipiv)
+        end
+        if d == :D
+            D = diagm(diag(LUD))


Would it be worth returning D as a SymTridiagonal to save space?

Or maybe a general Tridiagonal if there's no such thing as HermTridiagonal...

tkelman · 2017-06-28T20:47:40Z

base/linalg/bunchkaufman.jl

+    getindex(B::BunchKaufman, d::Symbol)
+
+Extract the factors of the Bunch-Kaufman factorization `B`. The factorization can take the
+two forms `L*D*L.'` or `U*D*U.'` where `L` is a `UnitLowerTriangular` matrix, `U` is a


there's also a hermitian version where it's L*D*L' right?

Yes, and it is also tested. I'll add it to the docstring.

tkelman · 2017-06-29T19:13:57Z

base/linalg/bunchkaufman.jl

+        return eye(T, n)[:,invperm(B[:p])]
+    elseif d == :L || d == :U || d == :D
+        if B.rook
+            # syconvf_rook just added to LAPACK 3.7.0. Uncomment and remove error when we distribute LAPACK 3.7.0


xref OpenMathLib/OpenBLAS#1118

Improve BunchKaufman printing

Adjust tests

…ivoting

tkelman · 2017-06-30T23:59:12Z

base/linalg/bunchkaufman.jl

+    p = T[1:maxi;]
+    uploL = uplo == 'L'
+    i = uploL ? 1 : maxi
+    # if uplo == 'U' we construct the permution backwards


permutation

tkelman · 2017-06-30T23:59:50Z

base/linalg/bunchkaufman.jl

+- `:p`: permutation vector
+- `:P`: permutation matrix
+
+```jldoctest


# Examples ```jldoctest

tkelman · 2017-07-01T00:01:15Z

base/linalg/bunchkaufman.jl

+symmetric or Hermitian matrix with 1x1 or 2x2 blocks. The argument `d` can be
+- `:D`: the block diagonal matrix
+- `:L`: the lower triangular factor (if factorization is `L*D*L.'`)
+- `:U`: the lower triangular factor (if factorization is `U*D*U.'`)


makes more sense to list ' form before .' form

tkelman · 2017-07-01T00:03:09Z

base/linalg/bunchkaufman.jl

 function Base.show(io::IO, mime::MIME{Symbol("text/plain")}, B::BunchKaufman)
    println(io, summary(B))
-    print(io, "successful: $(issuccess(B))")
+    println(io, "D factor:")


test the show method and/or fix the relevant doctest

I've adjusted the doctest

there's an existing one that was showing the raw fields of a BunchKaufman, came up in one of Kristoffer's recent pr's - this will make it fail in a different way

Okay. I've fixed that one as well.

tkelman · 2017-07-02T01:16:56Z

doc/src/manual/linear-algebra.md

  2.0  -1.0  -3.0
 -4.0  -3.0   5.0

+julia> B = [1.5 2 -4; 2 -1 -3; -4 -3 5]


unintentionally duplicated?

tkelman reviewed Jun 28, 2017

View reviewed changes

kshyatt added display and printing Aesthetics and correctness of printed representations of objects. linear algebra Linear algebra labels Jun 28, 2017

tkelman reviewed Jun 28, 2017

View reviewed changes

andreasnoack mentioned this pull request Jun 29, 2017

Adjust signatures for bkfact to use Symmetric and Hermitian #22605

Merged

tkelman reviewed Jun 29, 2017

View reviewed changes

andreasnoack added 3 commits June 30, 2017 16:26

Add getindex method for extraction of factors in BunchKaufman

e29cb2d

Improve BunchKaufman printing

Adjust bkfact signatures to use Symmetric and Hermitian

5878202

Adjust tests

Wrap LAPACK functions for reconstruction of Bunch-Kaufman with rook p…

8a20788

…ivoting

andreasnoack force-pushed the anj/bkfact branch from e255f1c to dbb2153 Compare June 30, 2017 23:56

tkelman reviewed Jun 30, 2017

View reviewed changes

base/linalg/bunchkaufman.jl

- `:p`: permutation vector

- `:P`: permutation matrix

```jldoctest

Copy link

Contributor

tkelman Jun 30, 2017

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

# Examples ```jldoctest

tkelman reviewed Jul 1, 2017

View reviewed changes

andreasnoack force-pushed the anj/bkfact branch 2 times, most recently from b3228a9 to fe131e4 Compare July 2, 2017 01:13

tkelman reviewed Jul 2, 2017

View reviewed changes

Update documentation and use Tridiagonal for storing D

19aa410

andreasnoack force-pushed the anj/bkfact branch from fe131e4 to 19aa410 Compare July 2, 2017 01:19

andreasnoack merged commit 709d65e into master Jul 3, 2017

andreasnoack deleted the anj/bkfact branch July 3, 2017 01:34

Uh oh!

RFC: Extract Bunch-Kaufman factors and use them for printing #22601

RFC: Extract Bunch-Kaufman factors and use them for printing #22601

Uh oh!

Conversation

andreasnoack commented Jun 28, 2017 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

4 participants

andreasnoack commented Jun 28, 2017 •

edited

Loading