Give an algorithm to compute the $LDL^T$ decomposition.

4.1k Views Asked by At

Suppose that $A$ is a symmetric and positive definite tridiagonal matrix.

Question: how do I give an algorithm to compute the $LDL^T$ decomposition of $A$, where $l_{ii} = 1$ (all the diagonal elements of $L$ are equal to $1$)?

What I know:

  • A tridiagonal matrix is a band matrix that has non-zero elements only on the main diagonal.
  • Matrix $A$ is positive definite if for all $u\in\mathbb{R}^N\backslash\{0\}$ we have that $u'Au > 0$.

I have no clue how to handle/solve this, thanks in advance!

2

There are 2 best solutions below

1
On

Let

$$ L(i,j) = \begin{cases} 1 &,& i=j \\ l_{ij} &,& j<i \\ 0 \end{cases} $$

and $$ D(i,j) = \begin{cases} d_{i} &,& i=j \\ 0 \end{cases} $$

Thus $$ L^t(i,j) = L(j,i) = \begin{cases} 1 &,& i=j \\ l_{ij} &,& i<j \\ 0 \end{cases} $$

Hence $$ (LD)(i,j) = \begin{cases} d_j &,& i=j \\ d_jl_{ij} &,& j<i \\ 0 \end{cases} $$ So finnaly $$ (LDL^t)(i,j) = \begin{cases} d_i + \sum\limits_{k=1}^{j-1} d_k l_{ik}^2 &,& i=j \\ d_j l_{ij} + \sum\limits_{k=1}^{j-1} d_k l_{ik} l_{jk} &,& i < j \\ (LDL^t)(j,i) &,& j < i \end{cases} $$

As we want to make $A = LDL^t$ we equate elements and solve row by row or column by column starting by the diagonal element since both $A$ and $LDL^t$ are symmetric.

Illustratively, going row by row, in the first row we have $$ a_{1,1} = d_1 $$ $$ a_{1,2} = d_1 l_{2,1} \implies l_{2,1} = \frac{a_{1,2}}{d_1} $$ $$ 0 = a_{1,3} = d_1 l_{3,1} \implies l_{3,1} = 0 $$ and inductively $l_{i,1} = 0$ for $i > 2$.

Similarly for the second row we have $$ a_{2,2} = d_2 + d_1 l_{2,1}^2 \implies a_{2,2} - d_1 l_{2,1}^2 = d_2 $$ $$ a_{2,3} = d_1 l_{2,1} l_{3,1} + d_2 l_{3,2} \implies l_{3,2} = \frac{a_{2,3}}{d_2} $$ $$ 0 = a_{2,4} = d_1 l_{2,1} l_{4,1} + d_2 l_{4,2} \implies l_{4,2} = 0 $$ and inductively $l_{i,2} = 0$ for $i > 3$. And so on...

Because $A$ is positive-definite, it follows that $a_{ii} > 0$ and therefore we will not have branch decision problems involving zeros in the algorithm.

Now you need to work by yourself a few simple examples to grasp the algorithm. Let me know if you get stuck.

Furthermore, notice that $LDL^t$ is always symmetric positive-semi-definite and so we require $A$ to be symmetric positive-definite for if $LDL^t$ is symmetric positive-definite then both $L$ and $D$ are invertible, that is, the decomposition will be unique.

2
On

I will show how an $LDL^T$ factorization can be computed inductively. You will still need to explain why you never divide by zero, but this is the essence of my original hint.

Write $$A = \begin{bmatrix} A_{11} & A_{21}^T \\ A_{21} & A_{22} \end{bmatrix} = \begin{bmatrix} 1 & \\ v & L \end{bmatrix} \begin{bmatrix} \alpha & \\ & D \end{bmatrix} \begin{bmatrix} 1 & v^T \\ & L^T \end{bmatrix} = \begin{bmatrix} \alpha & \alpha v^T \\ \alpha v^T & LDL^T + \alpha vv^T \end{bmatrix}$$ Here $A_{11}$ has dimension $1$, $A_{21}$ and $v$ are column vectors of dimension $n-1$ and $A_{22}$, $D$ and $L$ have dimension $n-1$ and $\alpha$ is simply a scalar. The matrix $D$ is diagonal, $L$ is lower unit triangular. Evidently, $$\alpha = a_{11}, \quad v = \frac{1}{\alpha}A_{21}$$ and $L$ can be computed the matrix. $$A' = A_{22} - \alpha vv^T$$ which is symmetric and has dimension $n-1$.

The fact that your matrix has a particular sparsity pattern, i.e., is tridiagonal has not been used at all. It will allow you to greatly reduce the cost of computing $v$ and applying the symmetric update.