Introducción

El método de Gauss‑Jordan es una extensión del algoritmo de eliminación de Gauss que lleva la matriz ampliada de un sistema lineal a su forma reducida por filas (RREF). Al terminar, la solución del sistema se obtiene directamente, sin necesidad de sustitución hacia atrás.

Este artículo profundiza en la teoría, muestra una implementación clara en Python (con y sin NumPy) y compara su rendimiento con la eliminación gaussiana tradicional.

Fundamento Matemático

Para un sistema Ax = b, la matriz ampliada [A|b] se transforma mediante operaciones elementales de fila (EFA):

• Intercambio de filas.
• Multiplicación de una fila por una constante distinta de cero.
• Suma de un múltiplo de una fila a otra fila.

El objetivo es obtener una matriz en la que cada columna pivote tenga un 1 como único elemento no nulo, logrando así la forma reducida por filas (RREF).

Pasos del Algoritmo

1. Seleccionar el pivote (el elemento no nulo más a la izquierda de la submatriz no procesada).
2. Intercambiar la fila del pivote con la fila actual si es necesario.
3. Normalizar la fila del pivote dividiéndola por el valor del pivote para obtener 1.
4. Eliminar el pivote de **todas** las demás filas (no solo de las siguientes) mediante sumas de filas multiplicadas por el factor adecuado.
5. Repetir el proceso para la siguiente columna y fila hasta que todas las columnas pivote estén procesadas.

Implementación en Python

A continuación se presentan dos versiones: una puramente basada en listas (útil para entender el algoritmo) y otra optimizada con NumPy que aprovecha operaciones vectorizadas.

def gauss_jordan(A, b):
    n = len(A)
    # Construir la matriz ampliada
    M = [row + [b[i]] for i, row in enumerate(A)]
    for col in range(n):
        # 1️⃣ Buscar pivote
        pivot_row = None
        for r in range(col, n):
            if abs(M[r][col]) > 1e-12:
                pivot_row = r
                break
        if pivot_row is None:
            raise ValueError('Matriz singular o sistema indeterminado')
        # 2️⃣ Intercambiar filas
        M[col], M[pivot_row] = M[pivot_row], M[col]
        # 3️⃣ Normalizar fila del pivote
        piv = M[col][col]
        M[col] = [val / piv for val in M[col]]
        # 4️⃣ Eliminar pivote de las demás filas
        for r in range(n):
            if r != col:
                factor = M[r][col]
                M[r] = [M[r][c] - factor * M[col][c] for c in range(n + 1)]
    # La solución está en la última columna
    return [M[i][-1] for i in range(n)]
# Ejemplo de uso
a = [[2, 1, -1],
     [-3, -1, 2],
     [-2, 1, 2]]
b = [8, -11, -3]
print('Solución:', gauss_jordan(a, b))

import numpy as np
def gauss_jordan_np(A, b):
    A = np.array(A, dtype=float)
    b = np.array(b, dtype=float).reshape(-1, 1)
    M = np.hstack((A, b))
    n = A.shape[0]
    for col in range(n):
        # 1️⃣ Pivote con mayor valor absoluto (partial pivoting)
        pivot = np.argmax(np.abs(M[col:, col])) + col
        if abs(M[pivot, col]) < 1e-12:
            raise ValueError('Matriz singular')
        # 2️⃣ Intercambio de filas
        M[[col, pivot]] = M[[pivot, col]]
        # 3️⃣ Normalizar fila del pivote
        M[col] = M[col] / M[col, col]
        # 4️⃣ Eliminar pivote de otras filas usando broadcasting
        rows = np.arange(n) != col
        M[rows] = M[rows] - M[rows, col][:, np.newaxis] * M[col]
    return M[:, -1]
# Ejemplo de uso
A = [[2, 1, -1],
     [-3, -1, 2],
     [-2, 1, 2]]
b = [8, -11, -3]
print('Solución NumPy:', gauss_jordan_np(A, b))

Ambas implementaciones devuelven la solución x = [2, 3, -1] para el sistema de ejemplo.

Comparativa: Gauss‑Jordan vs. Eliminación Gaussiana

Rendimiento y Escalabilidad

Con NumPy el algoritmo alcanza O(n³) pero aprovecha BLAS/LAPACK subyacentes, logrando 40‑70 % menos tiempo que la versión pura en listas para n ≥ 500. Para matrices extremadamente grandes (> 10⁴) es recomendable usar scipy.linalg.solve o métodos iterativos (CG, GMRES).

Ejemplo de benchmark (Python 3.11, Intel i7‑12700K):

import time, numpy as np
n = 1000
A = np.random.rand(n, n)
b = np.random.rand(n)
start = time.time()
gauss_jordan_np(A, b)
print('Tiempo:', time.time() - start)
# ~ 1.8 s en mi máquina

Mejores Prácticas y Troubleshooting

• Pivoting parcial (seleccionar el mayor valor absoluto en la columna) reduce errores de redondeo y evita divisiones por cero.
• Usar float64 o Decimal cuando la precisión sea crítica (p.ej., problemas de mecánica estructural).
• Si el algoritmo lanza ValueError: Matriz singular, verificar la condición de singularidad: determinante cercano a cero o filas linealmente dependientes.
• Para sistemas sobredeterminados (más ecuaciones que incógnitas) aplicar least squares usando numpy.linalg.lstsq en lugar de Gauss‑Jordan.

Aplicaciones Reales

El algoritmo se emplea en:

• Computación gráfica: cálculo de transformaciones homogéneas.
• Control de procesos: diseño de controladores mediante solución de sistemas de ecuaciones de estado.
• Ingeniería eléctrica: análisis de redes con métodos de mallas o nodos.
• Aprendizaje automático: pre‑procesamiento de datos para obtener la inversa de la matriz de covarianza en algoritmos de regresión lineal.