Page 1 of 1

Implementing Low rank update of block-LDLᵀ decomposition

PostPosted: Thu Dec 29, 2016 3:30 pm
by pushpendre
Hi,

I was wondering if any stable algorithm for efficiently updating the block-LDLᵀ decomposition of symmetric indefinite matrices was already available in the LAPACK library? I know that the dsytf2* functions compute this decomposition but I was not able to find anything for updating the decomposition.

Professor Danny Sorensen's thesis "UPDATING THE SYMMETRIC INDEFINITE FACTORIZATION WITH APPLICATIONS IN A MODIFIED NEWTON'S METHOD"(https://www.osti.gov/scitech/servlets/purl/7220580) contains one such algorithm including fortran source code on page 140, but the thesis is available as OCR'd pdf therefore I thought of asking on the LAPACK forum before spending time manually typing something that might already exist.

Thanks.

Re: Implementing Low rank update of block-LDLᵀ decomposition

PostPosted: Mon Jan 02, 2017 7:32 am
by Julien Langou
I understand you want to do: Given symmetric indefinite n-by-n A. Given n-by-k B with k << n. Given P, L and D, the LDLᵀ decomposition of A. Compute (in a fast way) the LDLᵀ decomposition of A + BBᵀ (by updating somehow P, L and D).

There is no such low rank update of LDLᵀ decomposition in LAPACK.

Note: You mention DSYTF2 to obtain an LDLᵀ decomposition from LAPACK. This is Level 2 BLAS using Bunch-Kaufman diagonal pivoting. Note that you can use DSYTRF ( Level 3 BLAS, Bunch-Kaufman diagonal pivoting ) and DSYTRF_RK (since v3.7, Level 3 BLAS, bounded Bunch-Kaufman (rook) diagonal pivoting). There is also DSYTRF_AA which performs LTLᵀ where T is symmetric tridiagonal (since v3.7, Level 3 BLAS, Aasen's factorization).

Note: DSYTRF_ROOK (since v3.5) should be deprecated any time soon, so please do not use.

If you want to contribute some codes, please contact us. ;)

Cheers, Julien.