Linear Systems of Algebraic Equations

This page presents some topics from Linear Algebra needed for construction of solutions to systems of linear algebraic equations and some applications. We use matrices and vectors as essential elements in obtaining and expressing the solutions.

Cholesky decomposition

Every symmetric, positive definite matrix A can be decomposed into a product of a unique lower triangular matrix L and its transpose:
\[ {\bf A} = {\bf L} \, {\bf L}^{\ast} , \]
where \( {\bf L}^{\ast} = \overline{\bf A}^{\mathrm T} \) denotes the conjugate transpose of L. The matrix L is called the Cholesky factor of A, and can be interpreted as a generalized square root of A, as described in Cholesky decomposition or Cholesky factorization. It was discovered by a French military officer and mathematician André-Louis Cholesky (1875--1918) for real matrices. When it is applicable, the Cholesky decomposition is roughly twice as efficient as the LU decomposition for solving systems of linear equations.

If the matrix A is Hermitian and positive semi-definite, then it still has a decomposition of the form \( {\bf A} = {\bf L}{\bf L}^{\ast} \) if the diagonal entries of L are allowed to be zero.

When A has real entries, L has real entries as well and the factorization may be written \( {\bf A} = {\bf L}{\bf L}^{\mathrm T} . \)

The Cholesky decomposition is unique when A is positive definite; there is only one lower triangular matrix L with strictly positive diagonal entries such that \( {\bf A} = {\bf L}{\bf L}^{\mathrm T} . \) However, the decomposition need not be unique when A is positive semidefinite.

A closely related variant of the classical Cholesky decomposition is the LDL decomposition,

\[ {\bf A} = {\bf L} \, {\bf \Lambda}\,{\bf L}^{\ast} , \]
where L is a lower unit triangular (unitriangular) matrix and \( \Lambda \) is a diagonal matrix.

The LDL variant, if efficiently implemented, requires the same space and computational complexity to construct and use but avoids extracting square roots. Some indefinite matrices for which no Cholesky decomposition exists have an LDL decomposition with negative entries in \( \Lambda . \) For these reasons, the LDL decomposition may be preferred.

Example. Here is the Cholesky decomposition of a symmetric real matrix:

\[ \begin{bmatrix} 4&12&-16 \\ 12&37&-43 \\ -16&-43&93 \end{bmatrix} = \begin{bmatrix} 2&0&0 \\ 6&1&0 \\ -8&5&3 \end{bmatrix} \begin{bmatrix} 2&6&-8 \\ 0&1&5 \\ 0&0&3 \end{bmatrix} . \]
And here is its \( {\bf L}{\bf \Lambda}{\bf L}^{\mathrm T} \) decomposition:
\[ \begin{bmatrix} 4&12&-16 \\ 12&37&-43 \\ -16&-43&93 \end{bmatrix} = \begin{bmatrix} 1&0&0 \\ 3&1&0 \\ -4&5&1 \end{bmatrix} \begin{bmatrix} 4&0&0 \\ 0&1&0 \\ 0&0&9 \end{bmatrix} \begin{bmatrix} 1&3&-4 \\ 0&1&5 \\ 0&0&1 \end{bmatrix} . \]
The Cholesky decomposition is mainly used for the numerical solution of linear equations \( {\bf A}{\bf x} = {\bf b} . \)

The Cholesky algorithm, used to calculate the decomposition matrix L, is a modified version of Gaussian elimination. The Cholesky–Banachiewicz and Cholesky–Crout algorithms are based on rewriting the equation \( {\bf A} = {\bf L}{\bf L}^{\ast} \) as

\[ {\bf A} = {\bf L} {\bf L}^{\mathrm T} = \begin{bmatrix} L_{11} & 0&0 \\ L_{21}&L_{22} &0 \\ L_{31} & L_{32} & L_{33}\end{bmatrix} \begin{bmatrix} L_{11}& L_{21}& L_{31} \\ 0&L_{22} & L_{32} \\ 0&0&L_{33} \end{bmatrix} = \begin{bmatrix} L_{11}^2 &L_{21} L_{11}&L_{31} L_{11} \\ L_{21} L_{11} &L_{21}^2 L_{22}^2 & L_{31} L_{21} +L_{32} L_{22} \\ L_{31} L_{11} & L_{31} L_{21} +L_{32} L_{22} & L_{31}^2 + L_{32}^2 + L_{33}^2 \end{bmatrix} , \]
which leads to
\begin{align*} L_{j,j} &= \sqrt{A_{j,j} - \sum_{k=1}^{j-1} L_{j,k} L_{j,k}^{\ast}} , \\ L_{i,j} &= \frac{1}{L_{j,j}} \left( A_{i,j} - \sum_{k=1}^{j-1} L_{i,k} L_{j,k}^{\ast} \right) , \quad i>j. \end{align*}
So we can compute the (i, j) entry if we know the entries to the left and above. The computation is usually arranged in either of the following orders
  • The Cholesky–Banachiewicz algorithm starts from the upper left corner of the matrix L and proceeds to calculate the matrix row by row.
  • The Cholesky–Crout algorithm starts from the upper left corner of the matrix L and proceeds to calculate the matrix column by column.

Example. Consider the positive definite matrix

\[ {\bf A} = \begin{bmatrix} 2&1&1 \\ 1&2&1 \\ 1&1&2 \end{bmatrix} . \]
Using the Cholesky–Banachiewicz algorithm, we construct the lower triangular matrix L:
\[ L_{1,1} = \sqrt{2}, \quad L_{2,1} = L_{3.1} = \frac{1}{\sqrt{2}} , \quad L_{2,2} = \sqrt{\frac{3}{2}}, \quad L_{3,2} = \frac{1}{\sqrt{6}} , \quad L_{3,3} = \frac{2}{\sqrt{3}} . \]
We verify with Sage that \( {\bf L} \,{\bf L}^{\mathrm T} = {\bf A} : \)

sage: A = 
sage: A.eigenvalues()
sage: A.fcp()
Therefore
\[ {\bf L} = \begin{bmatrix} \sqrt{2} &0&0 \\ 1/\sqrt{2} & \sqrt{3/2} & 0 \\ 1/\sqrt{2} & 1/\sqrt{6} & 2/\sqrt{3} \end{bmatrix} \quad \mbox{and} \quad {\bf L}^{\mathrm T} = \begin{bmatrix} \sqrt{2} & 1/\sqrt{2} & 1/\sqrt{2} \\ 0&\sqrt{3/2} & 1/\sqrt{6} \\ 0&0&2/\sqrt{3} \end{bmatrix} . \]

Table Of Contents

Previous topic

Matrix Algebra

Next topic

Systems of linear ODEs