Is there Fortran subroutine in LAPACK/BLAS or somewhere else to calculate LDL decomposition?

Question

Like the title says, I need to form cholesky LDL decomposition for my positive definite matrix A (Like normal cholesky, but there's ones one diagonal of L, and D is diagonal matrix). I have found only one function in Lapack which does that, but it says the matrix A has to be tridiagonal. Is there somekind of function which does that in some free subroutine libraries like lapack?

Answer 1

Check out the SSPTRF function from LAPACK:

Purpose
=======

SSPTRF computes the factorization of a real symmetric matrix A stored in packed format using the Bunch-Kaufman diagonal pivoting method:

A = U*D*U**T or A = L*D*L**T

where U (or L) is a product of permutation and unit upper (lower) triangular matrices, and D is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.

Also, Golub and van Loan's "Matrix Computations" book gives an algorithm for doing the decomposition. In my third edition, it's on page 138, Section 4.1.2 "Symmetry and the LDL^T Factorization".