Gauss-Seidel Method for Solving Linear Systems

Problem

Implement the Gauss-Seidel iterative method for solving a system of linear equations Ax = b. The method updates each variable sequentially using the latest available values.

Solution

import numpy as np

def gauss_seidel(A, b, x0=None, n_iterations=100, tolerance=1e-10):
    A = np.array(A, dtype=np.float64)
    b = np.array(b, dtype=np.float64).flatten()
    n = len(b)

    if x0 is None:
        x = np.zeros(n)
    else:
        x = np.array(x0, dtype=np.float64).flatten()

    for _ in range(n_iterations):
        x_old = x.copy()
        for i in range(n):
            sum_val = b[i]
            for j in range(n):
                if j != i:
                    sum_val -= A[i][j] * x[j]
            x[i] = sum_val / A[i][i]

        if np.linalg.norm(x - x_old) < tolerance:
            break

    return x.tolist()

Explanation

1. Start with an initial guess (default is zeros).
2. For each variable x_i, compute its new value using the current values of all other variables (which may include already-updated values from this iteration).
3. The update formula is: x_i = (b_i - sum_{j != i} A_{ij} * x_j) / A_{ii}.
4. The key difference from Jacobi iteration is that Gauss-Seidel uses the latest values immediately, often converging faster.
5. Stop when the change between iterations falls below the tolerance.

Complexity

• Time: O(n_iterations * n^2) where n is the number of variables
• Space: O(n) for the solution vector