# Shouldn't `bunchkaufman` be actually named `ldlt`?

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.

2 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`.

4 Likes

I like the first one (LBLᵀ) best as it explains what is being done and nicely explains how it relates to LDLᵀ.

1 Like

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.

• Matlab’s ldl:

`[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.

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