Solving Linear Equations with Gauss Elimination: A Complete Guide with C Code

Have you ever wondered how engineering software solves systems of 50 equations with 50 unknowns in milliseconds? The answer at the heart of most linear solvers is the Gauss Elimination Method — one of the most important algorithms in numerical computing.

Today I'll break it down from scratch: the concept, the algorithm, a working C implementation, and a real example with output.

The Problem We're Solving

Suppose you have a system of linear equations like this:

2x₁ +  x₂ − x₃ =  8
−3x₁ − x₂ + 2x₃ = −11
−2x₁ + x₂ + 2x₃ = −3

You need to find the values of x₁, x₂, and x₃ that satisfy all three equations simultaneously.

For 2 or 3 variables, you could solve by hand using substitution. But what about 10 variables? Or 100? That's where Gauss Elimination comes in.

The Core Idea

Gauss Elimination works in two phases:

Phase 1: Forward Elimination

Transform the system so that each equation has one fewer variable than the one above it — an upper triangular form.

Original:               After elimination:
[ 2   1  -1 |  8 ]     [ 2   1  -1 |  8  ]
[-3  -1   2 | -11]  →  [ 0  0.5 0.5|  1  ]
[-2   1   2 | -3 ]     [ 0   0   1 | -1  ]

Phase 2: Back Substitution

Once the matrix is triangular, solve from the bottom up:

• Last row gives x₃ directly
• Second-to-last gives x₂ using known x₃
• And so on upward

The Algorithm (Step by Step)

INPUT:  n equations, augmented matrix [A|b]
OUTPUT: solution x[0], x[1], ..., x[n-1]

FORWARD ELIMINATION:
For k = 0 to n-2:
    If A[k][k] == 0 → zero pivot, STOP
    For i = k+1 to n-1:
        factor = A[i][k] / A[k][k]
        For j = k to n:
            A[i][j] = A[i][j] - factor × A[k][j]

BACK SUBSTITUTION:
x[n-1] = A[n-1][n] / A[n-1][n-1]
For i = n-2 down to 0:
    sum = A[i][n] - Σ (A[i][j] × x[j])
    x[i] = sum / A[i][i]

The key operation is computing a factor (the multiplier) for each row below the pivot, then using it to zero out the entry in the pivot column.

C Implementation

#include <stdio.h>
#include <math.h>

#define MAX 10

int main() {
    float a[MAX][MAX+1], x[MAX];
    int n, i, j, k;
    float factor, sum;

    printf("Gauss Elimination Method\n");
    printf("========================\n");
    printf("Enter number of equations: ");
    scanf("%d", &n);

    /* Input the augmented matrix */
    printf("Enter augmented matrix [A|b] row by row:\n");
    for (i = 0; i < n; i++) {
        printf("Row %d: ", i+1);
        for (j = 0; j <= n; j++)
            scanf("%f", &a[i][j]);
    }

    /* ---- Forward Elimination ---- */
    for (k = 0; k < n-1; k++) {
        if (fabs(a[k][k]) < 1e-9) {
            printf("Zero pivot encountered. Method fails!\n");
            return 1;
        }
        for (i = k+1; i < n; i++) {
            factor = a[i][k] / a[k][k];
            for (j = k; j <= n; j++)
                a[i][j] -= factor * a[k][j];
        }
    }

    /* ---- Back Substitution ---- */
    x[n-1] = a[n-1][n] / a[n-1][n-1];
    for (i = n-2; i >= 0; i--) {
        sum = a[i][n];
        for (j = i+1; j < n; j++)
            sum -= a[i][j] * x[j];
        x[i] = sum / a[i][i];
    }

    /* ---- Display solution ---- */
    printf("\nSolution:\n");
    for (i = 0; i < n; i++)
        printf("x[%d] = %.4f\n", i+1, x[i]);

    return 0;
}

Running the Example

Compile and run:

gcc gauss_elimination.c -o gauss_elimination -lm
./gauss_elimination

Enter the augmented matrix for our system:

Enter number of equations: 3
Row 1: 2  1  -1  8
Row 2: -3  -1  2  -11
Row 3: -2  1  2  -3

Output:

Solution:
x[1] =  2.0000
x[2] =  3.0000
x[3] = -1.0000

Let's verify: 2(2) + (3) − (−1) = 4 + 3 + 1 = 8 ✓

What is the Augmented Matrix?

Instead of writing the equations separately, we combine the coefficients and constants into one matrix:

2x₁ + x₂ - x₃ = 8   becomes row → [ 2   1  -1 |  8 ]

The vertical bar separates the left-hand side coefficients from the right-hand side constants. This compact form makes it easy to apply row operations systematically.

Why Check for Zero Pivot?

The factor is computed as A[i][k] / A[k][k]. If A[k][k] is zero, we get a division-by-zero error. This is called a zero pivot and means the method fails in its basic form.

The fix is partial pivoting — before each elimination step, swap the current row with the row below that has the largest absolute value in the pivot column. This both avoids division by zero and reduces numerical error.

Complexity

The algorithm runs in O(n³) time — three nested loops. This is fine for small-to-medium systems (up to a few hundred equations). For very large systems (thousands of variables), iterative methods or sparse solvers are preferred.

Advantages

• Produces an exact solution (not iterative approximation)
• Works for any n×n non-singular system
• Foundation for LU decomposition, which is used in almost all professional linear algebra libraries

Limitations

• Fails with zero pivot without partial pivoting
• Numerical instability can accumulate for large n
• O(n³) makes it slow for very large systems

Key Takeaway

Gauss Elimination is the foundational algorithm behind solving linear systems. Every time you use NumPy's linalg.solve(), MATLAB's backslash operator, or any finite element solver — there's a descendant of Gauss Elimination working under the hood.

Understanding it at this level gives you a solid foundation for numerical analysis, scientific computing, and engineering simulations.