I am just musing over the choise of the name bunchkaufman for the bunkaufman function. If I understand it correctly, the function provides an implementation of an algorithm for LDLᵀ factorization of a symmetric and in general indefinite matrix. The thing is that there is a dedicated function for performing the LDLᵀ factorization – it is named ldlt. It seems to me that the only difference is that ldlt function works for either dense tridiagonal or sparse matrices, while bunchkaufman works for dense (symmetric) matrices. Am I right? If yes, wouldn’t it be more appropriate to keep just ldlt function and dispatch on the type of a matrix?
It seems to me that it is more common to reflect in the name of the function the type of factorization (WHAT) rather than the method (HOW). We have lu, qr, svd, … It seems to me that bunchkaufman is the only one among the matrix factorization functions in the standard LinearAlgebra library that tells HOW instead of WHAT. Am I right? If find it a bit confusing, the more so that the other ldlt function exists.
3 Likes
In Bunch-Kaufman the D is not necessarily diagonal—to ensure numeric stability it can be block-2x2. I’ve sometimes seen it written LBLᵀ, LDLᵀ, or just Bunch-Kaufman. I prefer the first and third over the second because I find it too confusing to refer to a non-diagonal matrix as D.
5 Likes
I like the first one (LBLᵀ) best as it explains what is being done and nicely explains how it relates to LDLᵀ.
2 Likes
Just to complement the discussion, let me document here how LDLᵀ is defined elsewhere:
Every nonsingular symmetric matrix A can be factored as A = PLDLᵀPᵀ, where where P is a permutation matrix, L is lower triangular with positive diagonal elements, and D is block diagonal, with nonsingular 1 × 1 and 2 × 2 diagonal blocks. This is called an LDLT factorization of A. (The Cholesky factorization can be considered a special case of LDLT factorization, with P = I and D = I.)
Computes the LDLt or Bunch-Kaufman factorization of a symmetric/ hermitian matrix. This function returns a block diagonal matrix D consisting blocks of size at most 2x2 and also a possibly permuted unit lower triangular matrix L such that the factorization A = L D L^H or A = L D L^T holds.
[L,D] = ldl(A) stores a block diagonal matrix D and a permuted lower triangular matrix in L such that A = L*D*L' . The block diagonal matrix D has 1-by-1 and 2-by-2 blocks on its diagonal. Note, this syntax is not valid for sparse A .
Diagonal pivoting methods.
We next describe the computation of the block LDLᵀ factorization…
This is not to say that Julia has to mimic all these and rename anything. I am just sharing my own experience that as a BFU it took me some time to find the appropriate name under which the functionality (the matrix factorization) I was looking for is hidden. Perhaps just a matter of explicitly stating these connections in the documentation for bunchkaufman? And possibly factorize too? I can certainly try to prepare some PR for the docs unless it is viewed as bikeshedding.
2 Likes
I think there’s some merit to renaming bunchkaufman as ldlt. Consider the Schur factorization. Based on the spectrum of the input, the triangular Schur factor may be genuinely upper-triangular (real spectrum) or it may have 2 x 2 blocks straddling the main diagonal (complex spectrum, psychologically upper-triangular). This is in a similar spirit as the LDLt factorization, where D is diagonal when the input matrix is symmetric positive-definite and may be only block diagonal when it is symmetric indefinite.
help?> schur
search: schur schur! Schur ordschur ordschur! GeneralizedSchur schedule isdispatchtuple
schur(A::StridedMatrix) -> F::Schur
Computes the Schur factorization of the matrix A. The (quasi) triangular Schur factor can be obtained from the Schur object F with either F.Schur or F.T and the orthogonal/unitary Schur vectors
can be obtained with F.vectors or F.Z such that A = F.vectors * F.Schur * F.vectors'. The eigenvalues of A can be obtained with F.values.
For real A, the Schur factorization is "quasitriangular", which means that it is upper-triangular except with 2×2 diagonal blocks for any conjugate pair of complex eigenvalues; this allows the
factorization to be purely real even when there are complex eigenvalues. To obtain the (complex) purely upper-triangular Schur factorization from a real quasitriangular factorization, you can use
Schur{Complex}(schur(A)).
julia> schur([1 2 3; 4 5 6; 7 8 3])
Schur{Float64, Matrix{Float64}}
T factor:
3×3 Matrix{Float64}:
13.1956 2.08647 4.38594
0.0 -0.35519 0.640425
0.0 0.0 -3.84045
Z factor:
3×3 Matrix{Float64}:
-0.280212 -0.769762 -0.57354
-0.65233 0.591018 -0.474514
-0.704235 -0.241173 0.667749
eigenvalues:
3-element Vector{Float64}:
13.195636108300677
-0.3551897919163341
-3.8404463163843325
julia> schur([1 2 3; 4 5 2; 3 -3 4])
Schur{Float64, Matrix{Float64}}
T factor:
3×3 Matrix{Float64}:
-1.86745 -1.31242 1.48681
0.0 5.93372 3.64964
0.0 -1.35705 5.93372
Z factor:
3×3 Matrix{Float64}:
-0.78126 0.214717 0.586114
0.294819 0.954572 0.043281
0.550195 -0.206611 0.809072
eigenvalues:
3-element Vector{ComplexF64}:
-1.86744569673359 + 0.0im
5.933722848366795 + 2.225476029157484im
5.933722848366795 - 2.225476029157484im
Depending on how the type LDLt currently exists, this may require some modifications, allocating enough memory for D to be symmetric tridiagonal. By the documentation, “The main use of an LDLt factorization F = ldlt(S) is to solve the linear system of equations Sx = b with F\b.” I think this can still be achieved with block diagonal D.
As with LU, the pivoting in the Bunch-Kaufman LDLt factorization would be elided in the name.
1 Like
I don’t understand how you give the Schur decomposition as an example why Bunch-Kaufman should be renamed to let, yet Schur is not named quqt