Método numérico Gauss-Seidel

Las matemáticas son, sin duda alguna, la peor materia para la mayoría de los estudiantes; más sin embargo, y por desgracia para los programadores, los métodos numéricos son algoritmos por excelencia, y por consiguiente son relativamente fáciles de codificar.

El método de Gauss-Seidel no será muy sencillo de resolver debido a la cantidad de operaciones e iteraciones que se realizan, pero para una computadora ésto no es problema alguno.
Trataré de explicarles la lógica que se ha de seguir para hacer un programa que resuelva un sistema de ecuaciones de 3X3 con este método; para ello yo usaré el lenguaje PHP pero el método, como todo buen algoritmo, es mudable a cualquier lenguaje.

Primero que nada es bueno el encasillar todos nuestros programas, o sea dividir lo más posible el problema en funciones pequeñas que se mandan a llamar de manera pertinente. Así que para empezar sería bueno tener funciones que calculen lo que será X1, X2 y X3.

Para ello podemos hacer ésto:

function getX1($b1,$a12,$x2,$a13,$x3,$a11)
{
  return ($b1-($a12*$x2)-($a13*$x3))/$a11;
}


function getX2($b2,$a21,$x1,$a23,$x3,$a22)
{
  return ($b2 - ($a21*$x1) - ($a23 * $x3)) / $a22;
}

function getX3($b3,$a31,$x1,$a32,$x2,$a33)

{
  return ($b3 - ($a31 * $x1) - ($a32 * $x2)) / $a33;
}

Los parámetros que reciben las funciones son los mismo con los cuales se calculan matemáticamente X1, X2 y X3, por ello los nombré de la misma manera que en las funciones tradicionales.

Lo mismo para el error aproximado:

function getEa($x1,$x2)
{
    return abs((($x2-$x1)/$x2)*100);
}

Siendo x1 el valor anterior y x2 el valor actual.

Muy buen, ya tenemos todo lo necesario para que funcione nuestro programa, aun que sera difícil de creer así es, ya que solo es necesario mandar a llamar estos métodos de manera reiterativa hasta que Ea tenga el valor acercado deseado.

¿Cómo hacemos eso? Sencillo, primero necesitamos tener todos los valores necesarios para la función, o sea los que será Es; además de tener los valores iniciales para X1, X2 y X3 y un arreglo de 2 espacios para cada una de éstas variables (o sea 3 arreglos de 2 espacios); ¿para qué? pues por que allí se almacenarán los valores de Ea actual y anterior. Estos arreglos quedarían así (Por motivos propios los nombré como X11, X21 y X31 a cada arreglo):

$x11[0]=$x1;
$x11[1]=0;
$x21[0]=$x2;
$x21[1]=0;
$x31[0]=$x3;
$x31[1]=0;

La posición inicial de los arreglos será la que guarde el valor calculado de las X's anteriores, en ella se guardará en valor inicial de cada X's; en la segunda posición se guardará el nuevo valor de las X's, inicialmente es 0.

Luego, recordando que todos las operaciones se repiten indefinidamente, es obvio que necesitamos un bucle; el "while" es en este caso el más apto. Yo por practico le puse una condición True para que se repita infinitas veces:

while (true)
{
}

Claro, dentro del bucle le pondré un "break" cuando la condición adecuada deba romper el bucle, pero por ahora lo dejamos así.

Inmediatamente, al inicial el bucle, es necesario re-calcular los valores de las X's y depositarlos en la segunda posición de su respectivo arreglo; y bueno, para el cálculo ya tenemos las funciones, nuestro bucle queda, por ahora, así:

while (true)
{
  $x11[1]=getX1($b1,$a12,$x21[0],$a13,$x31[0],$a11);
  $x21[1]=getX2($b2,$a21,$x11[1],$a23,$x31[0],$a22);
  $x31[1]=getX3($b3,$a31,$x11[1],$a32,$x21[1],$a33);
}

Como ven, en algunas ocasiones envío la posición 0 de unos arreglos y luego envío la posición 1; ésto por que recuerden que en el calculo de X1 se usan los antiguos valores de las X's, en X2 se usa el nuevo valor de X1 pero el antiguo valor de X3 (por que aún no se ha calculado el nuevo) y en X3 ya se usan los nuevos valores de X1 y X2.

Después, y solo por que en ocasiones suele ocurrir en la primera interacción, se comprueba si ninguno de los nuevos valores de X's no es igual a 0.

Después se calcula, mediante la función getEa que programamos anteriormente, los Ea para cada uno de los nuevos valores de las X's con los anteriores; y si los 3 son menores o iguales a Es rompemos el ciclo. Entonces nuestro "while" queda así:

while(true)
{
  $x11[1]=getX1($b1,$a12,$x21[0],$a13,$x31[0],$a11);
  $x21[1]=getX2($b2,$a21,$x11[1],$a23,$x31[0],$a22);
  $x31[1]=getX3($b3,$a31,$x11[1],$a32,$x21[1],$a33);
  if($x11[1]!=0 || $x21[1]!=0 || $x31[1]!=0)
  {
    if(getEa($x11[0],$x11[1])<=$es)
    {
      if(getEa($x21[0],$x21[1])<=$es)
      {
        if(getEa($x31[0],$x31[1])<=$es)
        {
          break;
        }
      }
    } 
  }
}

Y para finalizar nuestro bucle, re-posicionamos los nuevos valores de las X's en la posición 0 de sus arreglos y dejamos la posición 1 en 0:

$x11[0]=$x11[1];
$x21[0]=$x21[1];
$x31[0]=$x31[1];

Entonces el programa final queda:

        $x11[0]=$x1;

        $x11[1]=0;

        $x21[0]=$x2;

        $x21[1]=0;

        $x31[0]=$x3;

        $x31[1]=0;

        while(true)

        {

            $x11[1]=getX1($b1,$a12,$x21[0],$a13,$x31[0],$a11);

            $x21[1]=getX2($b2,$a21,$x11[1],$a23,$x31[0],$a22);

            $x31[1]=getX3($b3,$a31,$x11[1],$a32,$x21[1],$a33);

            if($x11[1]!=0 || $x21[1]!=0 || $x31[1]!=0)

            {

                if(getEa($x11[0],$x11[1])<=$es)

                {

                    if(getEa($x21[0],$x21[1])<=$es)

                    {

                        if(getEa($x31[0],$x31[1])<=$es)

                        {

                            break;

                        }

                    }

                }

            }

            $x11[0]=$x11[1];

            $x21[0]=$x21[1];

            $x31[0]=$x31[1];

        }

        $resultados[0]= $x11[1];

        $resultados[1]= $x21[1];

        $resultados[2]= $x31[1];

Las cosas siempre parecen fáciles después de haberlas visto más detenidamente; espero que haya sido un poco claro en éste programa, pero de no haber sido así lo correcto es siempre volver a analizar todo y experimentar con el código.

He aquí como me quedó el programa al subirlo al internet (y ponerle una interfaz con html :P ):
Programa gauss-seidel

Sigan programando.

Edición: Como el código de arriba puede confundir por lo largo y las inexistentes alineaciones en Blogger dejo un enlace de descarga:

Gauss-seidel en PHP