Wednesday, 28 January 2015 23:25

Gauss elimination

Written by
Rate this item
(0 votes)

The purpose of this note is to clarify questions and statements in class, and a step-by-step way to find an LU decomposition by hand.

 

We start with notations. Let \(A\) be an \(n\)-by-\(n\) non-singular matrix, \(b\) is a vector in \(\mathbb{R}^n\). The goal is to solve the linear system \(Ax=b\), where \(x\) is the unique solution of the system.

Define \(L_{ij}:=\mathbb{I}+l_{ij}E_{ij}\) be the matrix associate to the row operation between \(i\)-th and \(j\)-th row when performing Gauss elimination. Here, \(E_{ij}\) is an \(n\)-by-\(n\) matrix whose entries are all zeros except the \((i,j)\)-entry is 1. The Gauss elimination procedure (without pivoting) can be described as

\[U=L_{n,(n-1)}\cdots L_{32}L_{31}L_{21}A,\]

where \(U\) is an upper triangular matrix. The system is equivalent to

\[Ux=L_{n,(n-1)}\cdots L_{32}L_{31}L_{21}b,\]

and can be solved easily.

The corresponding LU decomposition of \(A\) is given by

\[A=LU,\quad\text{where}~~L=L_{21}^{-1}L_{31}^{-1}L_{32}^{-1}\cdots L_{n,(n-1)}^{-1}.\]

We get \(U\) and \(L_{ij}\) directly from the Gauss elimination procedure. Here, we derive a way to calculate \(L\) by hand.

First, we compute the inverse of \(L_{ij}\).

Proposition 1. For \(i>j\), \(L_{ij}^{-1}=\mathbb{I}-l_{ij}E_{ij}\).

Proof. We check the following identity

\[L_{ij}L_{ij}^{-1}=(\mathbb{I}-l_{ij}E_{ij})(\mathbb{I}-l_{ij}E_{ij})=\mathbb{I}-l_{ij}^2E_{ij}^2=\mathbb{I}.\]

Here, we use the fact that \(E_{ij}^2=0\) if \(i\neq j\), which is a direct consequence of proposition 2.

Example 1. For  \(L_{32}=\begin{pmatrix}1&0&0\\0&1&0\\0&5&1\end{pmatrix}\), its inverse \(L_{32}^{-1}=\begin{pmatrix}1&0&0\\0&1&0\\0&-5&1\end{pmatrix}\).

Proposition 2. The matrix \(E_{ij}\) satisfies the following identity.

\[E_{i_1j_1}E_{i_2j_2}=\begin{cases}0&\text{if~} j_1\neq i_2\\E_{i_1j_2}&\text{if~} j_1=i_2 \end{cases}\]

This in particular implies \(L_{i_1j_1}\) and \(L_{i_2j_2}\) are commutable if \(i_2\neq j_1\) and \(i_1\neq j_2\). {\it Note: they do not always commute.}

Next, we compute the product of \(L_{ij}^{-1}\).

Proposition 3. If \(j_1\neq i_2\), then \(L_{i_1j_1}^{-1}L_{i_2j_2}^{-1}=I-l_{i_1j_1}E_{i_1j_1}-l_{i_2j_2}E_{i_2j_2}\).

Example 2. If \(L_{21}^{-1}=\begin{pmatrix}1&0&0\\-1&1&0\\0&0&1\end{pmatrix}\) and \(L_{31}^{-1}=\begin{pmatrix}1&0&0\\0&1&0\\2&0&1\end{pmatrix}\), then \(L_{21}^{-1}L_{31}^{-1}=\begin{pmatrix}1&0&0\\-1&1&0\\2&0&1\end{pmatrix}\).

As \(L=L_{21}^{-1}L_{31}^{-1}L_{32}^{-1}\cdots L_{n,(n-1)}^{-1}\), we are always in the case that \(j_1<i_2\). Hence, we can simply assemble each term together to obtain \(L\).

Example 3 (Complete procedure of LU decomposition (without pivoting)) Find the LU decomposition of the following matrix  \(A=\begin{pmatrix}20&30&-10\\6&49&-13\\4&26&3\end{pmatrix}\).

Solution. We proceed with Gauss elimination.

\begin{align*}\begin{matrix}\textbf{1}\\ \textbf{2}\\ \textbf{3}\end{matrix}\begin{pmatrix}20&30&-10\\6&49&-13\\4&26&3\end{pmatrix}\longrightarrow\begin{matrix}\textbf{1}\\ \textbf{2}+(-\frac{6}{20})\times\textbf{1}\\ \textbf{3}\end{matrix}\begin{pmatrix}20&30&-10\\0&40&-10\\4&26&3\end{pmatrix}&\qquad L_{21}=\begin{pmatrix}1&0&0\\-.3&1&0\\0&0&1\end{pmatrix}\\ \longrightarrow\begin{matrix}\textbf{1}\\ \textbf{2}+(-\frac{6}{20})\times\textbf{1}\\ \textbf{3}+(-\frac{4}{20})\times\textbf{1}\end{matrix}\begin{pmatrix}20&30&-10\\0&40&-10\\0&20&5\end{pmatrix}&\qquad L_{31}=\begin{pmatrix}1&0&0\\0&1&0\\-.2&0&1\end{pmatrix}\\ \begin{matrix}\textbf{1}\\ \textbf{2}\\ \textbf{3}\end{matrix}\begin{pmatrix}20&30&-10\\0&40&-10\\0&20&5\end{pmatrix}\longrightarrow\begin{matrix}\textbf{1}\\ \textbf{2}\\ \textbf{3}+(-\frac{20}{40})\times\textbf{2}\end{matrix}\begin{pmatrix}20&30&-10\\0&40&-10\\0&0&10\end{pmatrix}&\qquad L_{32}=\begin{pmatrix}1&0&0\\0&1&0\\0&-.5&1\end{pmatrix}\end{align*}

Therefore, we conclude \(U=\begin{pmatrix}20&30&-10\\0&40&-10\\0&0&10\end{pmatrix}\), and

\[L=L_{21}^{-1}L_{31}^{-1}L_{32}^{-1}=\begin{pmatrix}1&0&0\\.3&1&0\\0&0&1\end{pmatrix} \begin{pmatrix}1&0&0\\0&1&0\\.2&0&1\end{pmatrix}\begin{pmatrix}1&0&0\\0&1&0\\0&.5&1\end{pmatrix} =\begin{pmatrix}1&0&0\\.3&1&0\\.2&.5&1\end{pmatrix}.\]

Read 1910 times Last modified on Saturday, 28 March 2015 16:46

Related items