3  Gaussian elimination and reduced echelon form

A general system with m equations and n unknowns has the form \begin{align*} \begin{array}{lllll} a_{11}x_1 &+ a_{12}x_2 &+ \cdots &+ a_{1n}x_n &= b_1 \\ a_{21}x_1 &+ a_{22}x_2 &+ \cdots &+ a_{2n}x_n &= b_2 \\ \vdots & \vdots & \vdots & \vdots & \vdots\\ a_{m1}x_1 &+ a_{m2}x_2 &+ \cdots &+ a_{mn}x_n &= b_m \end{array} \end{align*} In the previous chapter we wrote this more compactly on the form A\mathbf{x}=\mathbf{b}, where the coefficients a_{ij} were assembled in the matrix A, and the right hand sides values b_i in the vector \mathbf{b}. Gaussian elimination solves the system A\mathbf{x}=\mathbf{b} by systematically eliminating variables from left to right. To explain this further, suppose that the variables x_1,...,x_n already have been eliminated. There are now two possibilities:

  1. x_{n+1} does not appear in the remaining equations. Then x_{n+1} has been eliminated as well, and we proceed with the elimination of x_{n+2}.
  2. x_{n+1} appears in at least one of the remaining equations. We do the following:
    1. Take one of the equations where x_{n+1} appears, and place this first among the remaining equations.
    2. Scale this equation so that there is a 1 in front of x_{n+1}. We shall call this a leading one (there should be no contribution left of x_{n+1} in this equation).
    3. Add multiples of this equation to the other remaining equations so that x_{n+1} gets eliminated.
    We have now eliminated x_{n+1}, and proceed with the elimination of x_{n+2}.

Some concrete examples will clarify this

Example 3.1 (A system with a unique solution) Consider the system \begin{array}{rrrr} & -8x_2&+3x_3=& 6 \\ 4x_1& -3x_2&+2x_3=& 8\\ x_1& -2x_2& + x_3=& 3 \end{array} First we eliminate x_1. Since x_1 does not contribute in the first equation we swap the first and third equation and get \begin{array}{rrrr} x_1& -2x_2& + x_3=& 3 \\ \textcolor{red}{4}x_1& -3x_2&+2x_3=& 8\\ & -8x_2&+3x_3=& 6 \end{array} We also need to ensure that there is a leading one in the first equation, but this is already the case. We now eliminate x_1 in the second equation by subtracting 4 times the first equation (we mark in red what we eliminate). Then we get \begin{array}{rrrr} x_1& -2x_2& + x_3=& 3 \\ & 5x_2&-2x_3=& -4\\ & -8x_2&+3x_3=& 6 \end{array} We then make a leading one in the second equation: \begin{array}{rrrr} x_1& -2x_2& + x_3=& 3 \\ & x_2&-\frac{2}{5}x_3=& -\frac{4}{5}\\ & \textcolor{red}{-8}x_2&+3x_3=& 6 \end{array} We eliminate x_2 in the third equation: \begin{array}{rrrr} x_1& -2x_2& + x_3=& 3 \\ & x_2&-\frac{2}{5}x_3=& -\frac{4}{5}\\ & &-\frac{1}{5}x_3=& -\frac{2}{5} \end{array} We make a leading one in the third equation: \begin{array}{rrrr} x_1& -2x_2& + \textcolor{red}{1}x_3=& 3 \\ & x_2&\textcolor{red}{-\frac{2}{5}}x_3=& -\frac{4}{5}\\ & &x_3=& 2 \end{array} We now say that the system is on echelon form, meaning that:

  • each equation is either zero or has a leading one,
  • each equation begins with at least one more 0-coefficient than the equation above.

When we have arrived at an echelon form, we are finished eliminating variables left to right. We proceed by eliminating variables from right to left: x_3 can be eliminated in the first two equations by adding multiples of the third equation. If we perform both these operations we get \begin{array}{rrrr} x_1& \textcolor{red}{-2}x_2& =& 1 \\ & x_2&=& 0\\ & &x_3=& 2 \end{array} Finally we eliminate x_2 in the first equation by adding a multiple of the second equation: \begin{array}{rrrr} x_1& & =& 1 \\ & x_2&=& 0\\ & &x_3=& 2 \end{array} We thus have the unique solution \mathbf{x}=(x_1,x_2,x_3)=(1,0,2). The system is now on what we call reduced echelon form: It is on echelon form, and in addition each variable with a leading one appears in just one equation. These are also called pivot variables. 3 We also call the position of a leading one a pivot element. \clubsuit

We will also assemble both A and \mathbf{b} in the bigger matrix \begin{pmatrix}A & \mathbf{b}\end{pmatrix}, called the augmented matrix of the system. For the concrete example above we have that \begin{aligned} A &= \begin{pmatrix} 0 & -8 & 3 \\ 4 & -3 & 2 \\ 1 & -2 & 1 \end{pmatrix} & \mathbf{b} &= \begin{pmatrix} 6 \\ 8 \\ 3 \end{pmatrix} & \begin{pmatrix}A & \mathbf{b}\end{pmatrix} &= \begin{pmatrix} 0 & -8 & 3 & 6 \\ 4 & -3 & 2 & 8 \\ 1 & -2 & 1 & 3 \end{pmatrix} \end{aligned} An equation can be thought of as a row in the augmented matrix, and a variable as a column in the coefficient matrix A. The columns in A which have a pivot element are called pivot columns. In the same way, the rows in A which have a pivot element are called pivot rows. If we let \sim denote the operations we did above, everything can be written more compactly as:

\begin{aligned} \begin{pmatrix} 0 & -8 & 3 & 6 \\ 4 & -3 & 2 & 8 \\ 1 & -2 & 1 & 3 \end{pmatrix} &\sim \begin{pmatrix} 1 & -2 & 1 & 3 \\ \textcolor{red}{4} & -3 & 2 & 8 \\ 0 & -8 & 3 & 6 \end{pmatrix} \sim \begin{pmatrix} 1 & -2 & 1 & 3 \\ 0 & 5 & -2 & -4 \\ 0 & -8 & 3 & 6 \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & -2 & 1 & 3 \\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} \\ 0 & \textcolor{red}{-8} & 3 & 6 \end{pmatrix} \sim \begin{pmatrix} 1 & -2 & 1 & 3 \\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} \\ 0 & 0 & -\frac{1}{5} & -\frac{2}{5} \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & -2 & \textcolor{red}{1} & 3 \\ 0 & 1 & \textcolor{red}{-\frac{2}{5}} & -\frac{4}{5} \\ 0 & 0 & 1 & 2 \end{pmatrix} \sim \begin{pmatrix} 1 & \textcolor{red}{-2} & 0 & 1 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 2 \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & 0 & 0 & 1 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 2 \end{pmatrix} \end{aligned} \tag{3.1}

After the first five operations the matrix is on echelon form. The last matrix is on reduced echelon form.

We can bring a general matrix to reduced echelon form with the function rref in the class sym.Matrix. If we test this on the matrix above we get

We see that we obtained the same matrix. This is not obvious, since there are many different ways to bring a matrix to reduced echelon form, and each such can in principle lead to different reduced forms. Later we shall see, however, that the reduced echelon form is unique - regardless of the operations, we will always end up with the same reduced echelon form.

In Equation 3.1 we indicated the operations we did. Let us establish some notation for the first operations.

  1. The first operation swapped row 1 and 3. We will write this as I\leftrightarrow III.
  2. The second operation subtracted 4 times row 1 in row 2. We will write this as II-4I.
  3. The third operation multiplied row 2 with 1/5. We will write this as II\cdot(1/5).

With the help of this notation Equation 3.1 can be rewritten as follows: \begin{align*} \begin{pmatrix} 0 & -8 & 3 & 6 \\ 4 & -3 & 2 & 8 \\ 1 & -2 & 1 & 3 \end{pmatrix} &\stackrel{I\leftrightarrow III}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 \\ 4 & -3 & 2 & 8 \\ 0 & -8 & 3 & 6 \end{pmatrix} \stackrel{II-4I}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 \\ 0 & 5 & -2 & -4 \\ 0 & -8 & 3 & 6 \end{pmatrix}\\ &\stackrel{II\cdot(1/5)}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 \\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} \\ 0 & -8 & 3 & 6 \end{pmatrix} \stackrel{III+8II}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 \\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} \\ 0 & 0 & -\frac{1}{5} & -\frac{2}{5} \end{pmatrix}\\ &\stackrel{III\cdot(-5)}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 \\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} \\ 0 & 0 & 1 & 2 \end{pmatrix} \stackrel{II+\frac{2}{5}III,I-III}{\sim} \begin{pmatrix} 1 & -2 & 0 & 1 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 2 \end{pmatrix} \\ &\stackrel{I+2II}{\sim} \begin{pmatrix} 1 & 0 & 0 & 1 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 2 \end{pmatrix} \end{align*} The system was thus brought to reduced echelon form with help of three types of operations:

  1. Add a multiple of one row to another,
  2. multiply a row with a constant \neq 0,
  3. Swap two rows.

These are called elementary row operations, and they are easily implemented as follows:

A[j,:] = A[j,:] + c*A[i,:]
A[i,:] = c*A[i,:]
A[[i,j],:] = A[[j,i], :]

The three first operations from Example 3.1 can in particular be programmed as follows:

A=np.array([[0,-8,3,6],[4,-3,2,8],[1,-2,1,3]])
A[[0,2],:] = A[[2,0],:]
A[1,:] = A[1,:] - 4*A[0,:]
A[1,:] = A[1,:]/5

In practice, the computer finds the steps for us when Gaussian elimination is applied, see Section 3.3.

We say that A and B are row equivalent if we can transform A to B through a series of elementary row operations. From the discussion above we get \begin{pmatrix} 0 & -8 & 3 & 6 \\ 4 & -3 & 2 & 8 \\ 1 & -2 & 1 & 3 \end{pmatrix} \sim \begin{pmatrix} 1 & 0 & 0 & 1 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 2 \end{pmatrix}. Every matrix is row equivalent to a matrix on reduced echelon form.

Elementary row operations can be expressed in terms of matrices. To see this, consider the matrix E=I_m + c\mathbf{e}_i\mathbf{e}_j^T, where i\neq j. E equals the identity matrix, with the exception that E_{ij}=c. Clearly E is obtained by altering the identity matrix by adding row j multiplied with c to row i.

Let \mathbf{x}\in\mathbb{R}^m. It is clear that \mathbf{x} and E\mathbf{x} differ only in component i, where (E\mathbf{x})_i = x_i+cx_j. If A is a matrix, it is thus also clear that A and EA differ only in row i, where (EA)_{i:}=A_{i:}+cA_{j:}.

To summarise, we have shown that the elementary row operation which adds a multiple of one row to another, corresponds to multiplication to the left with the matrix E obtained by applying the same elementary row operation to the identity matrix. It is straightforward to prove that a similar result holds for the other two types of elementary row operations as well. We thus have the following result.

Proposition 3.1 Suppose B is obtained by applying an elementary row operation to the m\times n-matrix A. Then B=EA where E is the matrix obtained by applying the same row operation to the identity matrix I_m.

A matrix is called elementary if it is obtained by applying an elementary row operation to the identity matrix. The following matrices are elementary, and correspond to the row operations where we add 2 times row 2 to row 1, multiply row 3 with 3, and swap rows 2 and 3, respectively: \begin{aligned} E_1 &=\begin{pmatrix} 1 & 2 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix} & E_2 &=\begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 3 \end{pmatrix} & E_3 &=\begin{pmatrix} 1 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 1 & 0 \end{pmatrix} \end{aligned} \tag{3.2} Row operations can be reversed in the sense that one can find a row operation reversing what the first one did. The reverse row operations of Equation 3.2 are adding -2 times row 2 to row 1, multiplying row 3 with 1/3, and swapping rows 2 and 3, respectively. The elementary matrices for these are \begin{align*} F_1 &=\begin{pmatrix} 1 & -2 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix} & F_2 &=\begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1/3 \end{pmatrix} & F_3 &=\begin{pmatrix} 1 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 1 & 0 \end{pmatrix} \end{align*} That F_i reverses E_i means that F_1E_1=I_n, F_2E_2=I_n, F_3E_3=I_n, which also can be verified by multiplying the matrices. From Proposition 3.1 we get the following.

Proposition 3.2 For every matrix A we can write C=E_k\cdots E_2E_1A where C is a matrix on reduced echelon form, and where E_1,E_2,...,E_k are elementary matrices.

If A\mathbf{x}=\mathbf{b} is a system, we can easily obtain other systems by multiplying with the same matrix on each side. We have the following result on solutions of systems obtained in this way.

Proposition 3.3 Let E be an elementary matrix. The system A\mathbf{x}=\mathbf{b} has the same solutions as the system EA\mathbf{x}=E\mathbf{b}.

Proof. Let F be the elementary matrix associated with the reverse row operation.

  • If A\mathbf{x}=\mathbf{b} it is clear that also EA\mathbf{x}=E\mathbf{b}, so that \mathbf{x} also solves the new system.
  • If EA\mathbf{x}=E\mathbf{b} we get by applying F on both sides that A\mathbf{x}=\mathbf{b}, so that \mathbf{x} also solves the original system.

Thus, the original- and the new system have the same solutions.

In Example 3.1 we obtained a unique solution to the system under consideration. Let us also consider systems with no solution, and systems with infinitely many solutions.

Example 3.2 (A system with infinitely many solutions) Consider the system \begin{array}{rrrr} x_1 & + x_2 & + x_3=& 4 \\ 2x_1 & &-x_3=& 2. \end{array} We write down the augmented matrix, and apply elementary row operations: \begin{align*} \begin{pmatrix} 1 & 1 & 1 & 4 \\ 2 & 0 &-1 & 2 \end{pmatrix} &\stackrel{II-2I}{\sim} \begin{pmatrix} 1 & 1 & 1 & 4 \\ 0 & -2 &-3 & -6 \end{pmatrix} \stackrel{II\cdot(-1/2)}{\sim} \begin{pmatrix} 1 & 1 & 1 & 4 \\ 0 & 1 & 3/2 & 3 \end{pmatrix}\\ &\stackrel{I-II}{\sim} \begin{pmatrix} 1 & 0 & -1/2 & 1 \\ 0 & 1 & 3/2 & 3 \end{pmatrix} \end{align*} Again we ended up with a reduced echelon form. In equation form this says \begin{array}{rrrr} x_1 & & -\frac{1}{2}x_3=& 1 \\ & x_2 & +\frac{3}{2}x_3=& 3 \end{array} Which also can be written x_1=1+\frac{1}{2}x_3, x_2 = 3-\frac{3}{2}x_3. We see that x_3 can be chosen freely, and that a general solution can be written as \mathbf{x} =\begin{pmatrix} x_1 \\ x_2 \\ x_3 \end{pmatrix} =\begin{pmatrix} 1+\frac{1}{2}x_3 \\ 3-\frac{3}{2}x_3 \\ x_3 \end{pmatrix} =\begin{pmatrix} 1\\ 3\\ 0 \end{pmatrix} + x_3\begin{pmatrix} \frac{1}{2} \\ -\frac{3}{2} \\ 1 \end{pmatrix} We see that there are infinitely many solutions. x_1 and x_2 are pivot variables, which can be expressed in terms of x_3. x_3 is not a pivot variable since column 3 does not have a leading one. x_3 can be set to anything, and is therefore called a free variable. \clubsuit

Example 3.3 (A system without solutions) Consider the system \begin{array}{rrr} x_1 & +x_2=& 2 \\ 3x_1 & +3x_2=& 5 \end{array} Row reduction of the augmented matrix gives \begin{align*} \begin{pmatrix} 1 & 1 & 2 \\ 3 & 3 & 5 \end{pmatrix} \stackrel{II-3I}{\sim} \begin{pmatrix} 1 & 1 & 2 \\ 0 & 0 & -1 \end{pmatrix} \stackrel{II\cdot(-1)}{\sim} \begin{pmatrix} 1 & 1 & 2 \\ 0 & 0 & 1 \end{pmatrix} \stackrel{I-2II}{\sim} \begin{pmatrix} 1 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix}. \end{align*} Again we arrived at a reduced echelon form. Note here that, when eliminating x_1, x_2 also got eliminated. In equation form this says that \begin{array}{rrr} x_1 & + x_2=& 0 \\ & 0=& 1 \end{array} 0=1 is impossible, so that there are no solutions. \clubsuit

The three examples above illustrate three possibilities: A system can have

  1. a unique solution, as in Example 3.1. Here we had only pivot variables, and no pivot element in the last column.
  2. infinitely many solutions, as in Example 3.2. Here we had a free variable, and no pivot element in the last column.
  3. no solutions, as in Example 3.3. Here we had a pivot element in the last column.

Furthermore, we arrived at a reduced echelon form in all examples. These observations turn out to hold more generally:

Proposition 3.4 Any system can be brought to reduced echelon form, and we have the following possibilities:

  1. The system has no solutions. In this case the last column is a pivot column,
  2. The system has solutions. In this case the last column is not a pivot column and we have that
    1. when there are free variables the system has infinitely many solutions (any value of the free variables give a solution),
    2. when there are no free variables (all columns except the last is a pivot column) the system has a unique solution.

If we restrict to square coefficient matrices, the discussion above gives the following result.

Theorem 3.1 The following are equivalent for an n\times n-matrix A:

  1. A is row equivalent with I_n.
  2. A\mathbf{x}=\mathbf{0} has \mathbf{x}=\mathbf{0} as the only solution.
  3. A\mathbf{x}=\mathbf{b} has a unique solution for all right hand sides \mathbf{b}.
  4. A can be written as the product of elementary matrices.

Later we will extend this result with more equivalences.

Proof. 1.\Rightarrow 3.: Suppose that A can be row reduced to the identity matrix. We can solve A\mathbf{x}=\mathbf{b} by row reducing the augmented matrix \begin{pmatrix} A & \mathbf{b} \end{pmatrix}, using the same row operations. At the end we then arrive at equations of the form x_1=b_1', x_2=b'_2,…,x_n=b'_n, so that the system has a unique solution for all right hand sides \mathbf{b}.

3.\Rightarrow 2. is obvious.

2.\Rightarrow 1.: We solve A\mathbf{x}=\mathbf{0} by row reducing the augmented matrix \begin{pmatrix} A & \mathbf{0} \end{pmatrix}. Since the system only has one solution, there can be no free variables. But then all columns in A are pivot columns. Since A is square, all rows must also be pivot rows. It follows that A is row equivalent with I_n.

1.\Rightarrow 4.: That A\sim I_n means that we can write I_n=E_k\cdots E_1A, where E_i are elementary matrices. Let F_i be the elementary matrix which corresponds to the reverse row operation of E_i. If we multiply with F_k to the left on both sides we get that F_k=E_{k-1}\cdots E_1A. If we then multiply with F_{k-1} to the left on both sides we get that F_{k-1}F_k=E_{k-2}\cdots E_1A. If we continue this way we arrive in the end at A=F_1\cdots F_k, so that A can be written as a product of elementary matrices.

4.\Rightarrow 1.: If A=F_1\cdots F_k we can as above show that I_n=E_k\cdots E_1A (where again F_i represents the reverse row operation of E_i), so that A\sim I_n.

3.1 Gaussian elimination and complexity

Suppose we have a system with n equations and n unknowns. How many operations (additions, multiplications, etc.) are needed to solve the system, i.e., bring it to reduced echelon form? Suppose for simplicity that the coefficient matrix A is row equivalent with the identity matrix. Let us also assume that row interchanges do not cost anything (in terms of operations).

  • In Gaussian elimination we ensure leading ones in the rows by dividing with the first non-zero. This requires (n-1)+(n-2)+\cdots +1=(n-1)n/2 divisions.
  • In the first part of Gaussian elimination we zero out the components below the diagonal, column by column. Under diagonal element k in column k there are n-k rows. For each of these we multiply row k with a number (this corresponds to n-k multiplications, since the row starts with k-1 zeroes), and add this to the row (this gives n-k additions). Together we get (n-k)((n-k)+(n-k))=2(n-k)^2 operations here. Summing over all k we get 2\sum_{k=1}^{n-1} (n-k)^2=2\sum_{k=0}^n k^2=n(n+1)(2n+1)/6\approx 2n^3/3, since n^3/3 is the dominating term when multiplying out n(n+1)(2n+1)/6 (which is the formula for the sum of the n første squares). The interest in the dominating term only comes from that we are usually interested in counting operations when n gets big - we are not interested in the exact number of operations, which typically is an uglier expression.
  • In the second part of Gaussian elimination we zero out the components over the diagonal, column by column. No extra operations are needed for this.

Since (n-1)n/2\approx n^2/2 is dominated by 2n^3/3, we have in totalt about 2n^3/3 operations. We have shown the following.

Proposition 3.5 (Number of operations in Gaussian elimination) Suppose that A\sim I_n. About 2n^3/3 operations are required to bring A to reduced echelon form.

3.2 Simultaneous solution of systems

Suppose we want to solve two systems with the same coefficient matrix, A\mathbf{x}=\mathbf{b}_1 and A\mathbf{x}=\mathbf{b}_2. We bring A to reduced echelon form, C, and apply the same row operations to \mathbf{b}_1 and \mathbf{b}_2, and get \mathbf{b}'_1 and \mathbf{b}'_2. This corresponds to \begin{pmatrix} A & \mathbf{b}_1 & \mathbf{b}_2 \end{pmatrix} \sim \begin{pmatrix} C & \mathbf{b}'_1 & \mathbf{b}'_2 \end{pmatrix}, since row reducing A, \mathbf{b}_1, and \mathbf{b}_2 separately can also be done by row reducing \begin{pmatrix} A & \mathbf{b}_1 & \mathbf{b}_2 \end{pmatrix} as a whole. From \begin{pmatrix} C & \mathbf{b}'_1 & \mathbf{b}'_2 \end{pmatrix} we can read out the solutions to the two systems. Let us consider an example.

Example 3.4 In Example 3.1 we solved the system A\mathbf{x}=\mathbf{b}_1 where \mathbf{b}_1=(6,8,3). By using the row operations from this example we can simultaneously solve the system A\mathbf{x}=\mathbf{b}_2 where \mathbf{b}_2=(1,8,2) as follows:

\begin{align*} \begin{pmatrix}A & \mathbf{b}_1 & \mathbf{b}_2 \end{pmatrix} &= \begin{pmatrix} 0 & -8 & 3 & 6 & 1\\ 4 & -3 & 2 & 8 & 8 \\ 1 & -2 & 1 & 3 & 2\end{pmatrix} \stackrel{I\leftrightarrow III}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 & 2 \\ 4 & -3 & 2 & 8 & 8 \\ 0 & -8 & 3 & 6 & 1\end{pmatrix} \stackrel{II-4I}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 & 2 \\ 0 & 5 & -2 & -4 & 0 \\ 0 & -8 & 3 & 6 & 1\end{pmatrix}\\ &\stackrel{II\cdot(1/5)}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 & 2\\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} & 0 \\ 0 & -8 & 3 & 6 & 1\end{pmatrix} \stackrel{III+8II}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 & 2\\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} & 0 \\ 0 & 0 & -\frac{1}{5} & -\frac{2}{5} & 1 \end{pmatrix}\\ &\stackrel{III\cdot(-5)}{\sim} \begin{pmatrix} 1 & -2 & 1 & 3 & 2 \\ 0 & 1 & -\frac{2}{5} & -\frac{4}{5} & 0 \\ 0 & 0 & 1 & 2 & -5\end{pmatrix} \stackrel{II+\frac{2}{5}III,I-III}{\sim} \begin{pmatrix} 1 & -2 & 0 & 1 & 7 \\ 0 & 1 & 0 & 0 & -2 \\ 0 & 0 & 1 & 2 & -5\end{pmatrix} \\ &\stackrel{I+2II}{\sim} \begin{pmatrix} 1 & 0 & 0 & 1 & 3 \\ 0 & 1 & 0 & 0 & -2 \\ 0 & 0 & 1 & 2 & -5 \end{pmatrix} \end{align*} It follows that A\mathbf{x}_1=\mathbf{b}_1 as before has the solution \mathbf{x}_1=(1,0,2), but also that A\mathbf{x}_2=\mathbf{b}_2 has the solution \mathbf{x}_2=(3,-2,-5). \clubsuit

3.3 Implementation of Gaussian elimination

Let A be an m\times n-matrix. Gaussian elimination is a bit complicated to implement, but consider the following code which is a skeleton for bringing a matrix to echelon form:

piv = np.zeros(m).astype(int) # Stores the positions of the pivots.
i, j = 0, 0
while i < m and j < n: # Start scan for next pivot element
    k = i
    while k < m and abs(A[k,j]) < 10**(-6):
        k += 1
    if k == m: # Need to go to the next column
        j += 1
    else: # Found a nonzero in the current column
        # The next pivot element is in row k. If k>i, 
        # swap row i and k, so that the next pivot element is (i,j).
        # Create a leading one at pos. (i,j)
        # Zero out col. j below row i.
        piv[i] = j # The i'th leading one is in position (i,j)
        i += 1; j += 1

The initial position is (0,0). The outer while-loop scans for a next nonzero entry to be scaled to a leading one/a pivot element. First it looks for this below the current position (i) in the current column (the inner while-loop).

When applying row operations roundoff errors occur. Thus, we may get something nonzero, even if exactly zero should appear. The inner while-loop remedies this by assuming that a_{kj} is zero when |a_{kj}| < 10^{-6}.

If the search in the current column fails (the if-statement succeeds), the search moves to the next column (by continuing in the outer while-loop).

If the search in the current column succeeds, the else-statement is executed. Let us comment on how to implement the next three lines. If the nonzero (which is in row k) was found below the current row (row i), these rows are swapped:

A[[k,i]] = A[[i,k]]

We then create a leading one at position (i,j):

A[i,j:] = A[i,j:]/A[i,j]

(i.e., \mathbf{a}_i = \mathbf{a}_i/a_{ij}, where \mathbf{a}_i is row i in A). Note here that we only modified column j and upwards in row i. We did not need to divide with a_{ij} for columns 1,...,j-1, as the entries here are zero. In column j we can zero out component k by subtracting a multiple of row i:

A[k,j:] = A[k,j:] - A[k,j] * A[i,j:]

(i.e., \mathbf{a}_k = \mathbf{a}_k - a_{kj}\mathbf{a}_i). We need to do this for k=i+1,i+2,...,m, i.e., \begin{aligned} \mathbf{a}_{i+1} &= \mathbf{a}_{i+1} - a_{(i+1)j}\mathbf{a}_i \\ \mathbf{a}_{i+2} &= \mathbf{a}_{i+2} - a_{(i+2)j}\mathbf{a}_i \\ \vdots & \vdots \\ \mathbf{a}_m &= \mathbf{a}_m - a_{mj}\mathbf{a}_i \end{aligned} We here find A_{(i+1):,:} on the left hand side, and on the right hand side the product
\begin{pmatrix} a_{(i+1)j} \\ a_{(i+2)j} \\ \vdots \\ a_{mj}\end{pmatrix} \mathbf{a}_i = A_{(i+1):,j}\mathbf{a}_i This means that we can program this as

A[(i+1):,j:] = A[(i+1):,j:] - A[(i+1):,[j]] @ A[[i],j:]

Once this has been done, the code stores the position of the leading one in row i (piv[i] = j), and moves to the next row/column. Once the outer while-loop finishes, we arrive at echelon form.

After arriving at echelon form, we can similarly zero out all components above the leading ones to arrive at the reduced echelon form. If the pivot in row i is in column j we get \begin{aligned} \mathbf{a}_{1} &= \mathbf{a}_{1} - a_{1j}\mathbf{a}_i \\ \mathbf{a}_{2} &= \mathbf{a}_{2} - a_{2j}\mathbf{a}_i \\ \vdots & \vdots \\ \mathbf{a}_{i-1} &= \mathbf{a}_{i-1} - a_{(i-1)j}\mathbf{a}_i \end{aligned} On the right hand side we recognise the product \begin{pmatrix} a_{1j} \\ a_{2j} \\ \vdots \\ a_{(i-1)j}\end{pmatrix} \mathbf{a}_i = A_{:(i-1),j}\mathbf{a}_i so that this can be programmed as

j = piv[i]
A[:(i-1),j:] = A[:(i-1),j:] - A[:(i-1),[j]] @ A[[i], j:]

We here needed the stored locations of the leading ones.

If we also add some print statements to track the steps in Gaussian elimination (which are included if verbose is set to true), we arrive at the following code.

The following code tests everything on the matrix from Example 3.1.

Note that the matrix A is changed. So, if you also need A for other purposes, you should send a copy of it to gaussian_elimination. Note also that the entries in the matrix have been converted to floating points (float64). Without this the matrix would have integer entries, and all arithmetic would have been integer-based. In particular, division between two matrix entries would give and integer, which clearly is undesirable.

Quiz

--- primary_color: orange secondary_color: lightgray text_color: black shuffle_questions: false shuffle_answers: false --- ## A system with two equations and two unknowns The equations $$\begin{aligned}x-y&=1 \\ x+y&=2\end{aligned}$$ have 1. [x] a unique solution 2. [ ] infinitely many solutions 3. [ ] no solution ## A system with two equations and three unknowns The equations $$\begin{aligned} x+y+z&=2\\x+2y+3z&=3\end{aligned}$$ have 1. [ ] a unique solution 2. [x] infinitely many solutions 3. [ ] no solution ## A system with three equations and three unknowns The equations  $$\begin{aligned} x-y+z&=3\\2x+2y-z&=1\\4y-3z&=0\end{aligned}$$ have 1. [ ] a unique solution 2. [ ] infinitely many solutions 3. [x] no solution ## A system with three equations and three unknowns What can be said about the solutions to the following system? $$ \begin{aligned} 3x_1 &= 7-x_2 \\ 2x_1-3x_2 &= -3-4x_3 \\ 2x_2-4 &= 2x_3 \end{aligned} $$ > $$\begin{aligned} \begin{pmatrix} 3 & 1 & 0 & 7\\ 2 & -3 & 4 & -3 \\ 0 & 2 & -2 & 4 \end{pmatrix} &\sim \begin{pmatrix} 1 & 1/3 & 0 & 7/3\\ 0 & -11/3 & 4 & -23/3 \\ 0 & 1 & -1 & 2 \end{pmatrix} \sim \begin{pmatrix} 1 & 1/3 & 0 & 7/3\\ 0 & 1 & 1 & 2 \\ 0 & 0 & 1/3 & -1/3 \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & 0 & 0 & 2\\ 0 & 1 & 0 & 1 \\ 0 & 0 & 1 & -1 \end{pmatrix} \end{aligned} $$ 1. [ ] The system has no solutions 2. [ ] The system has infinitely many solutions 3. [ ] The system has the unique solution $$\mathbf{x}=\begin{pmatrix} 4\\-5\\-13/2\end{pmatrix}$$ 4. [ ] The system has the unique solution $$\mathbf{x}=\begin{pmatrix} 6\\-11\\-12 \end{pmatrix}$$ 5. [x] The system has the unique solution $$\mathbf{x}=\begin{pmatrix} 2\\1\\-1 \end{pmatrix}$$ ## System with different right hand sides. Let $A$ be a $4\times 2$-matrix. Suppose that the system $$A\mathbf{x} = \begin{pmatrix}1\\2\\4\\-1\end{pmatrix}$$ has the solution $$\mathbf{x}=\begin{pmatrix}1\\0\end{pmatrix}$$, and that the system $$A\mathbf{y} = \begin{pmatrix}2\\0\\-1\\-3\end{pmatrix}$$ has the solution $$\mathbf{y}=\begin{pmatrix} 0\\1\end{pmatrix}.$$ What can you say about the solution to the system $$A\mathbf{z}=\begin{pmatrix} 1\\4\\2\\0\end{pmatrix}$$? 1. [ ] The system has the solution $$\mathbf{z} = \begin{pmatrix} 11 \\ 4\end{pmatrix}.$$ 2. [ ] The system has a unique solution. 3. [x] The systemet has no solutions. 4. [ ] The system has an infinite number of solutions. 5. [ ] The system has the solution $$\mathbf{z} = \begin{pmatrix} 1 \\ 7\end{pmatrix}.$$ ## The reduced echelon form of a $3\times 3$ matrix The reduced echelon form of the matrix $$\begin{pmatrix} 1 & 2 & 3 \\ 1 & 4 & 5 \\ 2 & 6 & 8 \end{pmatrix}$$ is 1. [ ] $$\begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 2 \\ 0 & 0 & 0 \end{pmatrix} $$ 2. [ ] $$\begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix}$$ 3. [x] $$\begin{pmatrix} 1 & 0 & 1 \\ 0 & 1 & 1 \\ 0 & 0 & 0 \end{pmatrix}$$ ## The reduced echelon form of a $3\times 5$ matrix The reduced echelon form of the matrix $$\begin{pmatrix}2 & -1 & 2 & 10 & 11 \\ 5 &1 &1 &23 &12 \\ 4 & -1 & 3 & 20 & 18 \end{pmatrix}$$ is 1. [ ] $$\begin{pmatrix} 1 & 2 & 3 & 4 & 1 \\ 0 & 1 & 3 & -1 & -2 \\ 0 & 0 & 1 & 2 & 2 \end{pmatrix}$$ 2. [ ] $$\begin{pmatrix} 1 & 0 & 0 & 2 & 2 \\ 0 & 1 & 0 & 3 & 1 \\ 0 & 0 & 1 & 3 & 1 \end{pmatrix}$$ 3. [ ] $$\begin{pmatrix} 1 & 0 & 0 & 2 & -1 \\ 0 & 1 & 0 & 3 & 2 \\ 0 & 0 & 1 & 3 & -7 \end{pmatrix}$$ 4. [x] $$\begin{pmatrix} 1 & 0 & 0 & 3 & 2 \\ 0 & 1 & 0 & 4 & -1 \\ 0 & 0 & 1 & 4 & 3 \end{pmatrix}$$ 5. [ ] $$\begin{pmatrix} 1 & 0 & 0 & 0 & 4 \\ 0 & 1 & -3 & 0 & -1 \\ 0 & 0 & 0 & 1 & 2 \end{pmatrix}$$ ## Elementary matrices Which one of the following matrices is not elementary? 1. [ ] $$\begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 3 & 0 & 1 \end{pmatrix}$$ 2. [x] $$\begin{pmatrix} 1 & 0 & 0 \\ 0 & 2 & 0 \\ 3 & 0 & 1 \end{pmatrix}$$ 3. [ ] $$\begin{pmatrix} 1 & 0 & 0 \\ 0 & 2 & 0 \\ 0 & 0 & 1 \end{pmatrix}$$ ## System with parameter $a$ For which values of $a$ does the system $$\begin{array}{lll} x &+y &+z &= 3 \\ x &-y&+z &= 2 \\ x & &+z &= a \end{array}$$ have solution(s)? 1. [ ] 2 2. [x] 5/2 3. [ ] 1 ## System with parameter $a$ Let $$ \begin{aligned} A &= \begin{pmatrix} 3 & -1 & 2 \\ 2 & 4 & 1 \\ -2 & 3 & 2 \\ 1 & 2 & 4 \end{pmatrix} & {\bf b} &= \begin{pmatrix} 1 \\ 4 \\ -1 \\ a \end{pmatrix}. \end{aligned} $$ For which values of $a$ does the system $A{\bf x}={\bf b}$ have a solution? > Row reducing the augmented matrix gives $$\begin{pmatrix} 1 & 2 & 4 & a \\ 0 & -7 & -10 & 1-3a \\ 0 & 0 & -7 & 4-2a \\ 0 & 0 & 0 & -a \end{pmatrix}.$$ It is clear that the last column is a pivot column if and only if $a\neq 0$, so that the system has a solution only when $a=0$. 1. [ ] $a=-1$ 2. [ ] $a=2$ 3. [x] $a=0$ 4. [ ] $a=1$ 5. [ ] No such $a$ exists. ## Pivot columns The pivot columns of the matrix $$\begin{pmatrix} 4 & 2 & 3 & 4 \\ 2 & -2 & 1 & 3 \\ 8 & 10 & 7 & 6 \end{pmatrix}$$ are column 1. [x] 1 and 2 2. [ ] 1 and 3 3. [ ] 1, 2, and 3 4. [ ] 1, 2, 3, and 4 5. [ ] 1, 2, and 4 ## Which claim is true? Suppose $A$ is a matrix. Which claim is true? 1. [ ] $A$ can have different reduced echelon forms - depending on how we row reduce. 2. [x] $A$ can have different echelon forms - depending on how we row reduce. 3. [ ] In the reduced echelon forms of $A$ there can be numbers $\neq 0$ above the leading ones. 4. [ ] In $A$ the number of pivot rows can be different from the number of pivot columns. 5. [ ] $A{\bf x}={\bf b}$ must have a solution if the reduced echelon form of $A$ has a row of zeros. ## Which claim is true? 1. [ ] In a matrix on reduced echelon form there is always only one non-zero in each column. 2. [x] In a matrix on reduced echelon form there is always only one non-zero in each pivot column. 3. [ ] The reduced echelon form is not unique. 4. [ ] An $m\times n$-matrix always has $n$ pivot columns. 5. [ ] An $m\times n$-matrix always has $n$ rows with pivot elements.

Exercises

Exercise 3.1 Find the solution(s) to the following systems.

  1. \begin{array}{rrrr} 2x_1 & -7x_2&+x_3=& -21 \\ x_1& -4x_2&-3x_3=& -23\\ x_1& -x_2& + 2x_3=& 4 \end{array}

We row-reduce the augmented matrix: \begin{aligned} \begin{pmatrix} 2 & -7 &+1 & -21 \\ 1& -4 &-3 & -23 \\ 1& -1& & 2 & 4\end{pmatrix} & \sim \begin{pmatrix} 1& -1 & 2 & 4 \\ 2 & -7 & 1 & -21 \\ 1& -4 &-3 & -23 \end{pmatrix} \sim \begin{pmatrix} 1& -1 & 2 & 4 \\ 0 & -5 & -3 & -29 \\ 0& -3 &-5 & -27 \end{pmatrix} \\ &\sim \begin{pmatrix} 1& -1 & 2 & 4 \\ 0 & 1 & 3/5 & 29/5 \\ 0& -3 &-5 & -27 \end{pmatrix} \sim \begin{pmatrix} 1& -1 & 2 & 4 \\ 0 & 1 & 3/5 & 29/5 \\ 0& 0 & -16/5 & -48/5 \end{pmatrix} \\ &\sim \begin{pmatrix} 1& -1 & 2 & 4 \\ 0 & 1 & 3/5 & 29/5 \\ 0& 0 & 1 & 3 \end{pmatrix} \sim \begin{pmatrix} 1& -1 & 0 & -2 \\ 0 & 1 & 0 & 4 \\ 0& 0 & 1 & 3 \end{pmatrix} \\ &\sim \begin{pmatrix} 1& 0 & 0 & 2 \\ 0 & 1 & 0 & 4 \\ 0& 0 & 1 & 3 \end{pmatrix} \end{aligned} This gives the unique solution \mathbf{x}=(2,4,3).

  1. \begin{array}{rrrr} x_1 & -2x_2 =& 3 \\ x_1& -3x_2=& 4\\ -x_1& +5x_2=& 2 \end{array}

We row-recuce the augmented matrix: \begin{aligned} \begin{pmatrix} 1 & -2 & 3 \\ 1 & -3 & 4 \\ -1 & 5 & 2 \end{pmatrix} \sim \begin{pmatrix} 1 & -2 & 3 \\ 0 & -1 & 1 \\ 0 & 3 & 5 \end{pmatrix} \sim \begin{pmatrix} 1 & -2 & 3 \\ 0 & 1 & -1 \\ 0 & 0 & 8 \end{pmatrix} \sim \begin{pmatrix} 1 & -2 & 3 \\ 0 & 1 & -1 \\ 0 & 0 & 1 \end{pmatrix} \end{aligned} From this echelon form it is clear that the last column is a pivot column, so that the system has no solution.

  1. \begin{array}{rrrr} x_1 & -2x_2&+2x_3=& 2 \\ 2x_1& +3x_2&-4x_3=& 2 \end{array}

We row-reduce the augmented matrix: \begin{aligned} \begin{pmatrix} 1 & -2 & 2 & 2 \\ 2 & 3 & -4 & 2 \end{pmatrix} &\sim \begin{pmatrix} 1 & -2 & 2 & 2 \\ 0 & 7 & -8 & -2 \end{pmatrix} \sim \begin{pmatrix} 1 & -2 & 2 & 2 \\ 0 & 1 & -8/7 & -2/7 \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & 0 & -2/7 & 10/7 \\ 0 & 1 & -8/7 & -2/7 \end{pmatrix} \end{aligned} There are infinitely many solutions since x_3 is a free variable, and we get the general solution where x_1 = \frac{2}{7}x_3 + \frac{10}{7} and x_2 = \frac{8}{7}x_3 - \frac{2}{7}, so that \mathbf{x}=\begin{pmatrix}x_1\\x_2\\x_3\end{pmatrix} = \begin{pmatrix} \frac{2}{7}x_3 + \frac{10}{7} \\ \frac{8}{7}x_3 - \frac{2}{7} \\ x_3 \end{pmatrix} = x_3\begin{pmatrix} \frac{2}{7}\\ \frac{8}{7}\\1\end{pmatrix}+\begin{pmatrix} \frac{10}{7} \\ -\frac{2}{7} \\ 0\end{pmatrix}.

Exercise 3.2 We are given the matrix A=\begin{pmatrix} 1 & -1 & 3 \\ 2 & 1 & 4 \\ -4 & 7 & 2\end{pmatrix}. Find the solution to the systems A\mathbf{x}=\mathbf{b}_1 and A\mathbf{x}=\mathbf{b}_2 where \mathbf{b}_1=(-8,-7,-7) and \mathbf{b}_2=(9,5,-17).

We row-reduce the augmented matrix: \begin{aligned} \begin{pmatrix} 1 & -1 & 3 & -8 & 9 \\ 2 & 1 & 4 & -7 & 5 \\ -4 & 7 & 2 & -7 & -17 \end{pmatrix} &\sim \begin{pmatrix} 1 & -1 & 3 & -8 & 9 \\ 0 & 3 & -2 & 9 & -13 \\ 0 & 3 & 14 & -39 & 19 \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & -1 & 3 & -8 & 9 \\ 0 & 1 & -2/3 & 3 & -13/3 \\ 0 & 0 & 16 & -48 & 32 \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & -1 & 3 & -8 & 9 \\ 0 & 1 & -2/3 & 3 & -13/3 \\ 0 & 0 & 1 & -3 & 2 \end{pmatrix} \sim \begin{pmatrix} 1 & -1 & 0 & 1 & 3 \\ 0 & 1 & 0 & 1 & -3 \\ 0 & 0 & 1 & -3 & 2 \end{pmatrix} \\ &\sim \begin{pmatrix} 1 & 0 & 0 & 2 & 0 \\ 0 & 1 & 0 & 1 & -3 \\ 0 & 0 & 1 & -3 & 2 \end{pmatrix} \end{aligned} The unique solutions for the two systems are thus \mathbf{x}_1=(2,1,-3), and \mathbf{x}_2=(0,-3,2).

Exercise 3.3 Find another set of row operations which brings the matrix from Example 3.1 to reduced echelon form. Verify that you obtain the same reduced echelon form as in the example.

Exercise 3.4 For which value of a does the following system have a solution? \begin{align*} 2x + 3y + 4z &= 5 \\ 2x + 5y + 4z &= 7 \\ 4x + 7y + 8z &= a \end{align*}

We row-reduce the augmented matrix: \begin{pmatrix} 2 & 3 & 4 & 5 \\ 2 & 5 & 4 & 7 \\ 4 & 7 & 8 & a \end{pmatrix} \sim \begin{pmatrix} 2 & 3 & 4 & 5 \\ 0 & 2 & 0 & 2 \\ 0 & 1 & 0 & a-10 \end{pmatrix} \sim \begin{pmatrix} 2 & 3 & 4 & 5 \\ 0 & 1 & 0 & 1 \\ 0 & 0 & 0 & a-11\end{pmatrix} We see that it is only for a=11 that the last column is not a pivot column, so that the system has solution(s) only for a=11.

Exercise 3.5 Find the reduced echelon form of the following matrices

  1. A = \begin{pmatrix} 1 & 2 & -2 & 2 \\ 2 & 4 & -1 & 16 \\ 3 & 6 & -5 & 11 \end{pmatrix}

We get \begin{pmatrix} 1 & 2 & -2 & 2 \\ 2 & 4 & -1 & 16 \\ 3 & 6 & -5 & 11 \end{pmatrix} \sim \begin{pmatrix} 1 & 2 & -2 & 2 \\ 0 & 0 & 3 & 12 \\ 0 & 0 & 1 & 5 \end{pmatrix} \sim \begin{pmatrix} 1 & 2 & -2 & 2 \\ 0 & 0 & 1 & 4 \\ 0 & 0 & 0 & 1 \end{pmatrix} \sim \begin{pmatrix} 1 & 2 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{pmatrix}

  1. \begin{pmatrix} 1 & 2 & 3 & 4 \\ 2 & 5 & 7 & 13 \\ 3 & 7 & 10 & 18\end{pmatrix}

Wew get \begin{pmatrix} 1 & 2 & 3 & 4 \\ 2 & 5 & 7 & 13 \\ 3 & 7 & 10 & 18\end{pmatrix} \sim \begin{pmatrix} 1 & 2 & 3 & 4 \\ 0 & 1 & 1 & 5 \\ 0 & 1 & 1 & 6\end{pmatrix} \sim \begin{pmatrix} 1 & 2 & 3 & 4 \\ 0 & 1 & 1 & 5 \\ 0 & 0 & 0 & 1\end{pmatrix} \sim \begin{pmatrix} 1 & 0 & 1 & 0 \\ 0 & 1 & 1 & 0 \\ 0 & 0 & 0 & 1\end{pmatrix}

Exercise 3.6 What are the pivot columns of the following matrix? \begin{pmatrix}1 & 2 & 4 \\ 2 & 4 & 9 \\ -2 & -4 & 1 \\ -3 & -6 & 0\end{pmatrix}

\begin{pmatrix} 1 & 2 & 4 \\ 2 & 4 & 9 \\ -2 & -4 & 1 \\ -3 & -6 & 0 \end{pmatrix} \sim \begin{pmatrix} 1 & 2 & 4 \\ 0 & 0 & 1 \\ 0 & 0 & 9 \\ 0 & 0 & 12 \end{pmatrix} \sim\begin{pmatrix} 1 & 2 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix} We see that column 1 and 3 er are pivot columns.

Exercise 3.7 (A product of elementary matrices) What is the product \begin{pmatrix} 1 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 1 & 0 \end{pmatrix} \begin{pmatrix} 1 & 0 & 0 \\ 0 & 5 & 0 \\ 0 & 0 & 1 \end{pmatrix} \begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 3 & 1 \end{pmatrix} \begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 4 & 0 & 1 \end{pmatrix} \begin{pmatrix} 1 & 0 & 0 \\ 2 & 1 & 0 \\ 0 & 0 & 1\end{pmatrix}

The product can be written as E_4E_3E_2E_1A, where A is the matrix to the right.

  • E_1 adds 4 times row 1 to row 3. Applying this to A gives E_1A=\begin{pmatrix} 1 & 0 & 0 \\ 2 & 1 & 0 \\ 4 & 0 & 1 \end{pmatrix}.
  • E_2 adds 3 times row 2 to row 3. Applying this to E_1A gives E_2E_1A=\begin{pmatrix} 1 & 0 & 0 \\ 2 & 1 & 0 \\ 10 & 3 & 1\end{pmatrix}.
  • E_3 multiplies row 2 with 5. Applying this to E_2E_1A gives E_3E_2E_1A=\begin{pmatrix} 1 & 0 & 0 \\ 10 & 5 & 0 \\ 10 & 3 & 1\end{pmatrix}.
  • E_4 swaps row 2 and 3. Applying this to E_3E_2E_1A gives E_4E_3E_2E_1A=\begin{pmatrix} 1 & 0 & 0 \\ 10 & 3 & 1 \\ 10 & 5 & 0 \end{pmatrix}.

The product of all five matrices is thus \begin{pmatrix} 1 & 0 & 0 \\ 10 & 3 & 1 \\ 10 & 5 & 0 \end{pmatrix}

Exercise 3.8 Let A\in\mathbb{R}^{3\times 3}. Suppose that the system A\mathbf{x} = \begin{pmatrix}1\\0\\2\end{pmatrix} has the solution \mathbf{x}=\begin{pmatrix} 5\\6\\-3\end{pmatrix}, and that the system A\mathbf{y} = \begin{pmatrix}0\\4\\0\end{pmatrix} has the solution \mathbf{y}=\begin{pmatrix} 3\\-1\\2 \end{pmatrix}. Find a solution to the system A\mathbf{z}=\begin{pmatrix} 2\\8\\4\end{pmatrix}.

We start by expressing A\mathbf{z} in terms of A\mathbf{x} and A\mathbf{y}. This amounts to finding x_1 and x_2 so that x_1A\mathbf{x} + x_2A\mathbf{y}=A\mathbf{z}, i.e., x_1\begin{pmatrix}1\\0\\2\end{pmatrix} + x_2\begin{pmatrix} 0\\4\\0\end{pmatrix} = \begin{pmatrix} 1 & 0 \\ 0 & 4 \\ 2 & 0 \end{pmatrix} \begin{pmatrix} x_1 \\ x_2 \end{pmatrix} = \begin{pmatrix} 2\\8\\4\end{pmatrix} This can be solved by row-reducing the augmented matrix: \begin{pmatrix} 1 & 0 & 2 \\ 0 & 4 & 8 \\ 2 & 0 & 4 \end{pmatrix} \sim \begin{pmatrix} 1 & 0 & 2 \\ 0 & 4 & 8 \\ 0 & 0 & 0 \end{pmatrix} \sim \begin{pmatrix} 1 & 0 & 2 \\ 0 & 1 & 2 \\ 0 & 0 & 0 \end{pmatrix} so that x_1=x_2=2. We thus have that 2A\mathbf{x} + 2A\mathbf{y}=A\mathbf{z}, so that A(2\mathbf{x}+2\mathbf{y})=A\mathbf{z}. A solution of the system is therefore 2\mathbf{x} + 2\mathbf{y} = 2\begin{pmatrix}5\\6\\-3\end{pmatrix} + 2\begin{pmatrix} 3\\-1\\2\end{pmatrix} = \begin{pmatrix} 16 \\10 \\ -2\end{pmatrix}.

Exercise 3.9 Above we programmed the first three row operations in Example 3.1. Write code that applies the last row operations as well.

The last row operations can be programmed as follows

A[2,:] = A[2,:] + 8*A[1,:]
A[2,:]=(-5)*A[2,:]
A[0,:]=A[0,:]-A[2,:]
A[1,:] = A[1,:] +(2/5)*A[2,:]
A[0,:] = A[0,:] + 2*A[1,:]

Exercise 3.10 (Farkas lemma) Show that the system A\mathbf{x}=\mathbf{b} (with m equations and n unknowns) has a solution if and only if there is no \mathbf{c}\in\mathbb{R}^m so that \mathbf{c}^TA=\mathbf{0} and \mathbf{c}^T\mathbf{b}=1. This result is also called Farkas lemma.

We prove 1.\Rightarrow 2. first. Suppose \mathbf{x} is a solution to the system, and assume for contradiction there exists a \mathbf{c} so that \mathbf{c}^TA=\mathbf{0} and \mathbf{c}^T\mathbf{b}=1. Then \mathbf{c}^TA\mathbf{x}=\mathbf{c}^T\mathbf{b} says that 0x_1 + 0x_2 + \cdots + 0x_n = 1, which is impossible. Thus, there can be no such \mathbf{c}, and this proves 1.\Rightarrow 2..

We prove 2.\Rightarrow 1. contrapositively. Assume thus that the system does not have a solution. The augmented matrix must then have the last column as a pivot column, and the last row says that
0x_1 + 0x_2 + \cdots + 0x_n = 1. This equation can be obtained as a linear combination of the rows of \begin{pmatrix} A & \mathbf{b} \end{pmatrix}, i.e., there exists a vector \mathbf{c} so that \mathbf{c}^TA=\mathbf{0} and \mathbf{c}^T\mathbf{b}=1, so that such a \mathbf{c} exists. This is the negation of 2., and this proves 2.\Rightarrow 1..