Método SOR (successive overrelaxation)

Resolver el sistema de ecuaciones $$ \begin{align*} 10 x_{1} - x_{2} + 2 x_{3} &= 6 \\ -x_{1} + 11 x_{2} - x_{3} + 3 x_{4} &= 25 \\ 2 x_{1} - x_{2} + 10 x_{3} - x_{4} &= -11 \\ 3 x_{2} - x_{3} + 8 x_{4} &= 15 \end{align*} $$

Despejando $ x_{i} $ $$ \begin{alignat*}{5} x_{1}^{(k)} &= \cfrac{6}{10} & &+ \cfrac{1}{10} x_{2}^{(k-1)} &- \cfrac{2}{10} x_{3}^{(k-1)} & \\ x_{2}^{(k)} &= \cfrac{25}{11} &+ \cfrac{1}{11} x_{1}^{(k)} & &+ \cfrac{1}{11} x_{3}^{(k-1)} &- \cfrac{3}{11} x_{4}^{(k-1)} \\ x_{3}^{(k)} &= \cfrac{-11}{10} &- \cfrac{2}{10} x_{1}^{(k)} &+ \cfrac{1}{10} x_{2}^{(k)} & &+ \cfrac{1}{10} x_{4}^{(k-1)} \\ x_{4}^{(k)} &= \cfrac{15}{8} & &- \cfrac{3}{8} x_{2}^{(k)} &+ \cfrac{1}{8} x_{3}^{(k)} & \end{alignat*} $$

Usando $ \omega $ para acelerar la convergencia $$ \begin{alignat*}{6} x_{1}^{(k)} &= \omega \biggl[ \cfrac{6}{10} & &+ \cfrac{1}{10} x_{2}^{(k-1)} &- \cfrac{2}{10} x_{3}^{(k-1)} & &\biggr] + (1 - \omega) x_{1}^{(k-1)} \\ x_{2}^{(k)} &= \omega \biggl[ \cfrac{25}{11} &+ \cfrac{1}{11} x_{1}^{(k)} & &+ \cfrac{1}{11} x_{3}^{(k-1)} &- \cfrac{3}{11} x_{4}^{(k-1)} &\biggr] + (1 - \omega) x_{2}^{(k-1)} \\ x_{3}^{(k)} &= \omega \biggl[ \cfrac{-11}{10} &- \cfrac{2}{10} x_{1}^{(k)} &+ \cfrac{1}{10} x_{2}^{(k)} & &+ \cfrac{1}{10} x_{4}^{(k-1)} &\biggr] + (1 - \omega) x_{3}^{(k-1)} \\ x_{4}^{(k)} &= \omega \biggl[ \cfrac{15}{8} & &- \cfrac{3}{8} x_{2}^{(k)} &+ \cfrac{1}{8} x_{3}^{(k)} & &\biggr] + (1 - \omega) x_{4}^{(k-1)} \end{alignat*} $$

Para este ejemplo se usará $ \omega = 1.05 $

Iteración 0

$$ \begin{align*} x_{1}^{(0)} &= \color{blue}{0} \\ x_{2}^{(0)} &= \color{blue}{0} \\ x_{3}^{(0)} &= \color{blue}{0} \\ x_{4}^{(0)} &= \color{blue}{0} \end{align*} $$

Iteración 1

$$ \begin{align*} x_{1}^{(1)} &= 1.05 \biggl[ \frac{6 + x_{2}^{(0)} - 2 x_{3}^{(0)}}{10} \biggr] + (1 - 1.05) x_{1}^{(0)} \\ &= 1.05 \biggl[ \frac{6 + \color{blue}{0} - 2 (\color{blue}{0})}{10} \biggr] + (1 - 1.05) (\color{blue}{0}) = \color{green}{0.63} \\ x_{2}^{(1)} &= 1.05 \biggl[ \frac{25 + x_{1}^{(1)} + x_{3}^{(0)} - 3 x_{4}^{(0)}}{11} \biggr] + (1 - 1.05) x_{2}^{(0)} \\ &= 1.05 \biggl[ \frac{25 + \color{green}{0.63} + \color{blue}{0} - 3 (\color{blue}{0})}{11} \biggr] + (1 - 1.05) (\color{blue}{0}) = \color{green}{2.4465} \\ x_{3}^{(1)} &= 1.05 \biggl[ \frac{-11 - 2 x_{1}^{(1)} + x_{2}^{(1)} + x_{4}^{(0)}}{10} \biggr] + (1 - 1.05) x_{3}^{(0)} \\ &= 1.05 \biggl[ \frac{-11 - 2 (\color{green}{0.63}) + \color{green}{2.4465} + \color{blue}{0}}{10} \biggr] + (1 - 1.05) (\color{blue}{0}) = \color{green}{-1.030418} \\ x_{4}^{(1)} &= 1.05 \biggl[ \frac{15 - 3 x_{2}^{(1)} + x_{3}^{(1)}}{8} \biggr] + (1 - 1.05) x_{3}^{(0)} \\ &= 1.05 \biggl[ \frac{15 - 3 (\color{green}{2.4465}) + (\color{green}{-1.030418})}{8} \biggr] + (1 - 1.05) (\color{blue}{0}) = \color{green}{0.870198} \end{align*} $$

Iteración 2

$$ \begin{align*} x_{1}^{(2)} &= 1.05 \biggl[ \frac{6 + x_{2}^{(1)} - 2 x_{3}^{(1)}}{10} \biggr] + (1 - 1.05) x_{1}^{(1)} \\ &= 1.05 \biggl[ \frac{6 + \color{green}{2.4465} - 2 (\color{green}{-1.030418})}{10} \biggr] + (1 - 1.05) (\color{green}{0.63}) = \color{red}{1.071770} \\ x_{2}^{(2)} &= 1.05 \biggl[ \frac{25 + x_{1}^{(2)} + x_{3}^{(1)} - 3 x_{4}^{(1)}}{11} \biggr] + (1 - 1.05) x_{2}^{(1)} \\ &= 1.05 \biggl[ \frac{25 + \color{red}{1.071770} + (\color{green}{-1.030418}) - 3 (\color{green}{0.870198})}{11} \biggr] + (1 - 1.05) (\color{green}{2.4465}) = \color{red}{2.018793} \\ x_{3}^{(2)} &= 1.05 \biggl[ \frac{-11 - 2 x_{1}^{(2)} + x_{2}^{(2)} + x_{4}^{(1)}}{10} \biggr] + (1 - 1.05) x_{3}^{(1)} \\ &= 1.05 \biggl[ \frac{-11 - 2 (\color{red}{1.071770}) + \color{red}{2.018793} + \color{green}{0.870198}}{10} \biggr] + (1 - 1.05) (\color{green}{-1.030418}) = \color{red}{-1.025207} \\ x_{4}^{(2)} &= 1.05 \biggl[ \frac{15 - 3 x_{2}^{(2)} + x_{3}^{(2)}}{8} \biggr] + (1 - 1.05) x_{4}^{(1)} \\ &= 1.05 \biggl[ \frac{15 - 3 (\color{red}{2.018793}) + (\color{red}{-1.025207})}{8} \biggr] + (1 - 1.05) (\color{green}{0.870198}) = \color{red}{0.995782} \end{align*} $$

Iteración 3

$$ \begin{align*} x_{1}^{(3)} &= 1.05 \biggl[ \frac{6 + x_{2}^{(2)} - 2 x_{3}^{(2)}}{10} \biggr] + (1 - 1.05) x_{1}^{(2)} \\ &= 1.05 \biggl[ \frac{6 + \color{red}{2.018793} - 2 (\color{red}{-1.025207})}{10} \biggr] + (1 - 1.05) (\color{red}{1.071770}) = \color{fuchsia}{1.003678} \\ x_{2}^{(3)} &= 1.05 \biggl[ \frac{25 + x_{1}^{(3)} + x_{3}^{(2)} - 3 x_{4}^{(2)}}{11} \biggr] + (1 - 1.05) x_{2}^{(2)} \\ &= 1.05 \biggl[ \frac{25 + \color{fuchsia}{1.003678} + (\color{red}{-1.025207}) - 3 (\color{red}{0.995782})}{11} \biggr] + (1 - 1.05) (\color{red}{2.018793}) = \color{fuchsia}{1.998213} \\ x_{3}^{(3)} &= 1.05 \biggl[ \frac{-11 - 2 x_{1}^{(3)} + x_{2}^{(3)} + x_{4}^{(2)}}{10} \biggr] + (1 - 1.05) x_{3}^{(2)} \\ &= 1.05 \biggl[ \frac{-11 - 2 (\color{fuchsia}{1.003678}) + \color{fuchsia}{1.998213} + \color{red}{0.995782}}{10} \biggr] + (1 - 1.05) (\color{red}{-1.025207}) = \color{fuchsia}{-1.000142} \\ x_{4}^{(3)} &= 1.05 \biggl[ \frac{15 - 3 x_{2}^{(3)} + x_{3}^{(3)}}{8} \biggr] + (1 - 1.05) x_{4}^{(2)} \\ &= 1.05 \biggl[ \frac{15 - 3 (\color{fuchsia}{1.998213}) + (\color{fuchsia}{-1.000142})}{8} \biggr] + (1 - 1.05) (\color{red}{0.995782}) = \color{fuchsia}{1.000896} \end{align*} $$

Fórmula matemática

$$ \begin{align*} k &= 1, \dots \\ & \quad i = 1, \dots, m \\ & \quad \quad x_{i}^{(k)} = \frac{\omega}{a_{ii}} \biggl( b_{i} - \sum_{j = 1}^{i-1} a_{ij} x_{j}^{(k)} - \sum_{j = i+1}^{n} a_{ij} x_{j}^{(k-1)} \biggr) + (1 - \omega) x_{i}^{(k-1)} \end{align*} $$

Seudocódigo

function metodo_SOR(A, b)
    m, n = tamaño(A)
    x_actual = [0,...,0]
    for k=1 to iteraciones do
        for i=1 to m do
            x_anterior = x_actual
            sumatoria = b(i)
            for j=1 to i-1 do
                sumatoria = sumatoria - A(i,j)*x_actual(j)
            end for
            for j=i+1 to n do
                sumatoria = sumatoria - A(i,j)*x_anterior(j)
            end for
            x_actual(i) = omega * sumatoria/A(i,i) + (1 - omega)*x_anterior(i)
        end for
    end for
    return x
end function

Implementación

import numpy as np

def metodo_SOR(A,b,omega):
    m, n = A.shape
    x_actual = np.zeros(m)
    error_permitido = 0.000001
    print("{0:s} \t {1:9s} \t {2:9s} \t {3:9s} \t {4:9s} \t {5:9s}".format('k', 'x(1)', 'x(2)', 'x(3)', 'x(4)', 'error relativo %'))
    print("{0:d} \t {1:10.9f} \t {2:10.9f} \t {3:10.9f} \t {4:10.9f} \t {5:9s}".format(0, x_actual[0], x_actual[1], x_actual[2], x_actual[3], '??????????????'))
    for k in range(0,30):
        for i in range(0,m):
            x_anterior = np.copy(x_actual)
            sumatoria = b[i]
            for j in range(0,i):
                sumatoria = sumatoria - A[i,j]*x_actual[j]
            for j in range(i+1,n):
                sumatoria = sumatoria - A[i,j]*x_anterior[j]
            x_actual[i] = (omega*sumatoria)/A[i,i] + (1-omega)*x_anterior[i]
        error_relativo = np.linalg.norm(((x_actual - x_anterior)/x_actual)*100)
        print("{0:d} \t {1:10.9f} \t {2:10.9f} \t {3:10.9f} \t {4:10.9f} \t {5:14.10f}".format(k+1, x_actual[0], x_actual[1], x_actual[2], x_actual[3], error_relativo))
        if error_relativo < error_permitido:
            break

A = np.array([[10,-1,2,0],
              [-1,11,-1,3],
              [2,-1,10,-1],
              [0,3,-1,8]])
B = np.array([6,25,-11,15])

metodo_SOR(A,B,1.05)

    k 	 x(1)      	 x(2)      	 x(3)      	 x(4)      	 error relativo %
    0 	 0.000000000 	 0.000000000 	 0.000000000 	 0.000000000 	 ??????????????
    1 	 0.630000000 	 2.446500000 	 -1.030417500 	 0.870198328 	 100.0000000000
    2 	 1.071770175 	 2.018792780 	 -1.025206795 	 0.995782035 	  12.6115658080
    3 	 1.003678160 	 1.998213227 	 -1.000142571 	 1.000895728 	   0.5109116599
    4 	 0.999658421 	 1.999786620 	 -0.999849493 	 1.000058986 	   0.0836692535
    5 	 0.999963068 	 2.000004619 	 -0.999993091 	 0.999996139 	   0.0062847209
    6 	 1.000000881 	 2.000001618 	 -1.000000766 	 0.999999455 	   0.0003316465
    7 	 1.000000287 	 2.000000029 	 -1.000000076 	 1.000000006 	   0.0000550385
    8 	 1.000000005 	 1.999999990 	 -0.999999998 	 1.000000004 	   0.0000001782