Gaussian Elimination Algorithm

linear-algebra
Published

April 2, 2026

Introduction to Gaussian Elimination

Now you’re ready to learn a very famous and classic algorithm called Gaussian elimination. You will see that this is just the elimination method you studied before, but restructured to solve a system of equations and formalized so that it can be followed by hand or implemented in code.

Recall that when you learned about the singularity of a matrix, you ignored the constant values on the right-hand side of your equations. To actually solve a system of equations, you will now need to pay attention to them.

The Augmented Matrix

Let’s start with a system of equations:

\[\begin{cases} 2a - b + c = 1 \\ 2a + 2b + 4c = -2 \\ 4a + b = -1 \end{cases}\]

To begin, make a matrix from the coefficients in your system of equations, just as you always have:

\[\begin{bmatrix} 2 & -1 & 1 \\ 2 & 2 & 4 \\ 4 & 1 & 0 \end{bmatrix}\]

Now add another column to the right side of your matrix which holds the constant values \(1\), \(-2\), and \(-1\). This is called the augmented matrix:

\[\left[\begin{array}{ccc|c} 2 & -1 & 1 & 1 \\ 2 & 2 & 4 & -2 \\ 4 & 1 & 0 & -1 \end{array}\right]\]

The vertical line is used to separate the constants, so you remember they’re not part of the variables.

The Gaussian Elimination Process

Now if you proceed with the elimination method as normal, you can use the augmented matrix to solve your system of equations. The process involves:

  1. Find a pivot on the diagonal of the matrix
  2. Set the pivot to 1 using row operations
  3. Set all values below the pivot to 0 using row operations
  4. Repeat for each row
  5. Back-substitute to find the solution

Whatever row operations you perform on the matrix will also be applied to the column of constants.

Step-by-Step Example

Step 1: Make the First Pivot Equal to 1

To start, select the top-left cell as your pivot. Your first task is to use row operations to set your pivot to 1.

Since the pivot is currently 2, multiply the first row by \(\frac{1}{2}\):

\[R_1 \leftarrow \frac{1}{2}R_1\]

\[\left[\begin{array}{ccc|c} 1 & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \\ 2 & 2 & 4 & -2 \\ 4 & 1 & 0 & -1 \end{array}\right]\]

Step 2: Eliminate Below the First Pivot

Next, set all the values below the pivot to 0 using row operations.

For the second row, subtract 2 times the first row: \[R_2 \leftarrow R_2 - 2R_1\]

\[\left[\begin{array}{ccc|c} 1 & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \\ 0 & 3 & 3 & -3 \\ 4 & 1 & 0 & -1 \end{array}\right]\]

For the third row, subtract 4 times the first row: \[R_3 \leftarrow R_3 - 4R_1\]

\[\left[\begin{array}{ccc|c} 1 & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \\ 0 & 3 & 3 & -3 \\ 0 & 3 & -2 & -3 \end{array}\right]\]

Step 3: Make the Second Pivot Equal to 1

Moving along the diagonal, your new pivot is 3. Multiply the second row by \(\frac{1}{3}\):

\[R_2 \leftarrow \frac{1}{3}R_2\]

\[\left[\begin{array}{ccc|c} 1 & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \\ 0 & 1 & 1 & -1 \\ 0 & 3 & -2 & -3 \end{array}\right]\]

Step 4: Eliminate Below the Second Pivot

Subtract 3 times the second row from the third row: \[R_3 \leftarrow R_3 - 3R_2\]

\[\left[\begin{array}{ccc|c} 1 & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \\ 0 & 1 & 1 & -1 \\ 0 & 0 & -5 & 0 \end{array}\right]\]

Step 5: Make the Third Pivot Equal to 1

The third pivot is \(-5\). Multiply the third row by \(-\frac{1}{5}\):

\[R_3 \leftarrow -\frac{1}{5}R_3\]

\[\left[\begin{array}{ccc|c} 1 & -\frac{1}{2} & \frac{1}{2} & \frac{1}{2} \\ 0 & 1 & 1 & -1 \\ 0 & 0 & 1 & 0 \end{array}\right]\]

The matrix is now in row echelon form: the diagonal contains all ones and below the diagonal there are only zeros.

Back-Substitution

Now we perform back-substitution, starting from the bottom row and working upward. We use the pivot from each row to cancel the values in the cells above it.

Eliminate Above the Third Pivot

Cancel the \(\frac{1}{2}\) in the first row, third column: \[R_1 \leftarrow R_1 - \frac{1}{2}R_3\]

Cancel the \(1\) in the second row, third column: \[R_2 \leftarrow R_2 - R_3\]

\[\left[\begin{array}{ccc|c} 1 & -\frac{1}{2} & 0 & \frac{1}{2} \\ 0 & 1 & 0 & -1 \\ 0 & 0 & 1 & 0 \end{array}\right]\]

Eliminate Above the Second Pivot

Cancel the \(-\frac{1}{2}\) in the first row, second column: \[R_1 \leftarrow R_1 + \frac{1}{2}R_2\]

\[\left[\begin{array}{ccc|c} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & -1 \\ 0 & 0 & 1 & 0 \end{array}\right]\]

The Solution

The final matrix has the identity matrix on the left side, giving us the solution:

\[\begin{cases} a = 0 \\ b = -1 \\ c = 0 \end{cases}\]

Singular Cases

What happens if the matrix is singular? During Gaussian elimination, you will encounter a row that becomes all zeros. When this happens:

Case 1: Infinitely Many Solutions

If the constant value in the row of zeros is also zero:

\[\left[\begin{array}{ccc|c} 1 & 2 & 3 & 4 \\ 0 & 1 & 2 & 1 \\ 0 & 0 & 0 & 0 \end{array}\right]\]

The last row says \(0a + 0b + 0c = 0\), which is always true. The system has infinitely many solutions.

Case 2: No Solutions

If the constant value in the row of zeros is non-zero:

\[\left[\begin{array}{ccc|c} 1 & 2 & 3 & 4 \\ 0 & 1 & 2 & 1 \\ 0 & 0 & 0 & 4 \end{array}\right]\]

The last row says \(0a + 0b + 0c = 4\), which is impossible. The system has no solutions.

Summary: The Gaussian Elimination Algorithm

Here’s the complete process:

  1. Create the augmented matrix by adding the constants to a new column on the right
  2. Get the matrix in reduced row echelon form using row operations:
    • Make each pivot equal to 1
    • Make all entries below each pivot equal to 0
    • Make all entries above each pivot equal to 0
  3. Read the solution from the final column

If you encounter a row of zeros:

  • If the constant is also zero: infinitely many solutions
  • If the constant is non-zero: no solutions

Quiz: Practice the Algorithm

Solve the following system using Gaussian elimination:

\[\begin{cases} x + 2y + z = 6 \\ 2x + y - z = 1 \\ x - y + 2z = 5 \end{cases}\]

Step 1: Set up the augmented matrix: \[\left[\begin{array}{ccc|c} 1 & 2 & 1 & 6 \\ 2 & 1 & -1 & 1 \\ 1 & -1 & 2 & 5 \end{array}\right]\]

Step 2: Eliminate below the first pivot: \[R_2 \leftarrow R_2 - 2R_1, \quad R_3 \leftarrow R_3 - R_1\] \[\left[\begin{array}{ccc|c} 1 & 2 & 1 & 6 \\ 0 & -3 & -3 & -11 \\ 0 & -3 & 1 & -1 \end{array}\right]\]

Step 3: Make the second pivot equal to 1: \[R_2 \leftarrow -\frac{1}{3}R_2\] \[\left[\begin{array}{ccc|c} 1 & 2 & 1 & 6 \\ 0 & 1 & 1 & \frac{11}{3} \\ 0 & -3 & 1 & -1 \end{array}\right]\]

Step 4: Eliminate below the second pivot: \[R_3 \leftarrow R_3 + 3R_2\] \[\left[\begin{array}{ccc|c} 1 & 2 & 1 & 6 \\ 0 & 1 & 1 & \frac{11}{3} \\ 0 & 0 & 4 & 10 \end{array}\right]\]

Step 5: Make the third pivot equal to 1: \[R_3 \leftarrow \frac{1}{4}R_3\] \[\left[\begin{array}{ccc|c} 1 & 2 & 1 & 6 \\ 0 & 1 & 1 & \frac{11}{3} \\ 0 & 0 & 1 & \frac{5}{2} \end{array}\right]\]

Step 6: Back-substitution: \[R_1 \leftarrow R_1 - R_3, \quad R_2 \leftarrow R_2 - R_3\] \[\left[\begin{array}{ccc|c} 1 & 2 & 0 & \frac{7}{2} \\ 0 & 1 & 0 & \frac{7}{6} \\ 0 & 0 & 1 & \frac{5}{2} \end{array}\right]\]

\[R_1 \leftarrow R_1 - 2R_2\] \[\left[\begin{array}{ccc|c} 1 & 0 & 0 & \frac{5}{3} \\ 0 & 1 & 0 & \frac{7}{6} \\ 0 & 0 & 1 & \frac{5}{2} \end{array}\right]\]

Answer: \(x = \frac{5}{3}\), \(y = \frac{7}{6}\), \(z = \frac{5}{2}\)

Lab: Solving Systems with NumPy

You’ve been doing Gaussian elimination by hand. NumPy does it for you with np.linalg.solve. Let’s verify our examples.

import numpy as np

Solve the system from the example above:

\[\begin{cases} 2a - b + c = 1 \\ 2a + 2b + 4c = -2 \\ 4a + b + 0c = -1 \end{cases}\]

A = np.array([[2, -1, 1],
              [2,  2, 4],
              [4,  1, 0]])
b = np.array([1, -2, -1])

solution = np.linalg.solve(A, b)
print("Solution:", solution)
print("Verify A @ x =", A @ solution)
Solution: [ 0. -1.  0.]
Verify A @ x = [ 1. -2. -1.]

The quiz system:

\[\begin{cases} x + 2y + z = 6 \\ 2x + y - z = 1 \\ x - y + 2z = 5 \end{cases}\]

A2 = np.array([[1,  2,  1],
               [2,  1, -1],
               [1, -1,  2]])
b2 = np.array([6, 1, 5])

print("Solution:", np.linalg.solve(A2, b2))
Solution: [1.16666667 1.16666667 2.5       ]

What about a singular system?

A_singular = np.array([[1, 1, 1],
                       [2, 2, 2],
                       [3, 3, 3]])
b_singular = np.array([1, 2, 3])

try:
    np.linalg.solve(A_singular, b_singular)
except np.linalg.LinAlgError as e:
    print(f"Cannot solve: {e}")

You can check singularity beforehand with the determinant:

print(f"det(A)  = {np.linalg.det(A):.2f} -> solvable")
print(f"det(A2) = {np.linalg.det(A2):.2f} -> solvable")
print(f"det(A_singular) = {np.linalg.det(A_singular):.2f} -> singular!")
det(A)  = -30.00 -> solvable
det(A2) = -12.00 -> solvable
det(A_singular) = -0.00 -> singular!

np.linalg.solve does Gaussian elimination under the hood. Now you know what it’s actually doing.

Next: Vector Algebra