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.

QR-decomposition

The QR decomposition expresses a matrix as the product of an orthogonal matrix and an upper triangular matrix. The letter Q is a substitute for the letter O from "orthogonal" and the letter R is from "right", an alternative for "upper".

Suppose that A is an \( m \times n \) matrix with linearly independent columns (full column rank), and Q is the matrix that results by applying the Gram--Schmidt process to the column vectors of A. Then A can be factored as

\[ {\bf A} = {\bf Q} \, {\bf R} = \left[ {\bf q}_1 \ {\bf q}_2 \ \cdots \ {\bf q}_n \right] \begin{bmatrix} R_{1,1} & R_{1,2} & \cdots & R_{1,n} \\ 0& R_{2,2} & \cdots & R_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ 0&0& \cdots & R_{n,n} \end{bmatrix} , \]

where Q is an \( m \times n \) matrix with orthonormal column vectors (\( \| {\bf q}_i \| =1, \quad \mbox{and} \quad {\bf q}^{\mathrm T}_i {\bf q}_j =0 \) if \( i \ne j \) ), and R is an \( n \times n \) invertible upper triangular matrix. This factorization is unique if we require that the diagonal elements of R be positive. If \( R_{i,i} <0 , \) we can switch the signs of \( R_{i,i} , \ldots , R_{i,n} , \) and the vector \( {\bf q}_{i} \) to make the diagonal term \( R_{i,i} \) positive.

In the above definition, \( {\bf Q}^{\mathrm T} {\bf Q} = {\bf I} ; \) If A is a square matrix, then Q is orthogonal: \( {\bf Q}^{\mathrm T} {\bf Q} = {\bf Q} {\bf Q}^{\mathrm T} = {\bf I} . \) If instead A is a complex square matrix, then there is a decomposition \( {\bf Q}\, {\bf R} = {\bf A} , \) where Q is a unitary matrix (so \( {\bf Q}^{\ast} {\bf Q} = {\bf I} \) ).

There are several methods for actually computing the QR decomposition, such as by means of the Gram--Schmidt process ( \( 2mn^2 \) flops, sensitive to rounding errors), Householder transformations (\( 2mn^2 - (2/3) n^3 \) flops), or Givens rotations. Each has a number of advantages and disadvantages.

Example. There is an example of QR decomposition:

\begin{align*} {\bf A} &= \begin{bmatrix} -1&-1&1 \\ 1&3&3 \\ -1&-1&5 \\ 1&3&7 \end{bmatrix} = \frac{1}{2} \begin{bmatrix} -1&1&-1 \\ 1&1&-1 \\ -1&1&1 \\ 1&1&1 \end{bmatrix} \begin{bmatrix} 2&-4&2 \\ 0&2&8 \\ 0&0&4 \end{bmatrix} \\ &= \left[ {\bf q}_1 \ {\bf q}_2 \ {\bf q}_3 \right] \begin{bmatrix} R_{1,1}&R_{1,2}&R_{1,3} \\ 0&R_{2,2}&R_{2,3} \\ 0&0&R_{3,3} \end{bmatrix} = {\bf Q}\, {\bf R} , \end{align*}

where column vectors \( {\bf q}_1 , {\bf q}_2 , {\bf q}_3 \) are identified from the above equation. ■

Recall that a matrix is full row rank when each of the rows of the matrix are linearly independent and full column rank when each of the columns of the matrix are linearly independent. Using QR decomposition for a full column rank matrix A, we can express its pseudo-inverse through matrices Q and R:

\begin{align*} {\bf A}^{\dagger} &= \left( {\bf A}^{\mathrm T} {\bf A} \right)^{-1} {\bf A}^{\mathrm T} = \left( \left( {\bf Q}\,{\bf R} \right)^{\mathrm T} {\bf Q}\,{\bf R} \right)^{-1} \left( {\bf Q}\,{\bf R} \right)^{\mathrm T} \\ &= \left( {\bf R}^{\mathrm T} {\bf Q}^{\mathrm T} {\bf Q} {\bf R} \right)^{-1} {\bf R}^{\mathrm T} {\bf Q}^{\mathrm T} \\ &= \left( {\bf R}^{\mathrm T} {\bf R} \right)^{-1} {\bf R}^{\mathrm T} {\bf Q}^{\mathrm T} = {\bf R}^{-1} \left( {\bf R}^{\mathrm T} \right)^{-1} {\bf R}^{\mathrm T} {\bf Q}^{\mathrm T} \\ &= {\bf R}^{-1} {\bf Q}^{\mathrm T} . \end{align*}

For a square nonsingular A this is the inverse: \( {\bf A}^{-1} = \left( {\bf Q}\,{\bf R} \right)^{-1} = {\bf R}^{-1} {\bf Q}^{\mathrm T} . \)

Let us emphasize the main property of QR decomposition: Q has the same range as A.

Spectral decomposition

One feature of Sage is that we can easily extend its capabilities by defining new commands. Here will create a function that checks if an augmented matrix represents a consistent system of equations. The syntax is just a bit complicated. lambda is the word that indicates we are making a new function, the input is temporarily named A (think A ugmented), and the name of the function is consistent. Everything following the colon will be evaluated and reported back as the output.

sage: consistent = lambda A: not(A.ncols()-1 in A.pivots())

Execute this block above. There will not be any output, but now the consistent function will be defined and available. Now give it a try (after making sure to have run the code above that defines aug). Note that the output of consistent() will be either True or False.

sage: consistent(aug)
False

The consistent() command works by simply checking to see if the last column of A is not in the list of pivots. We can now test many different augmented matrices, such as perhaps changing the vector of constants while keeping the coefficient matrix fixed. Again, make sure you execute the code above that defines coeff and const.

sage: consistent(coeff.augment(const))
False
sage: w = vector(QQ, [3,1,2])
sage: consistent(coeff.augment(w))
True
sage: u = vector(QQ, [1,3,1])
sage: consistent(coeff.augment(u))
False

Why do some vectors of constants lead to a consistent system with this coefficient matrix, while others do not ? This is a fundamental question, which we will come to understand in several different ways.