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-1)} & &+ \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-1)} &+ \cfrac{1}{10} x_{2}^{(k-1)} & &+ \cfrac{1}{10} x_{4}^{(k-1)} \\ x_{4}^{(k)} &= \cfrac{15}{8} & &- \cfrac{3}{8} x_{2}^{(k-1)} &+ \cfrac{1}{8} x_{3}^{(k-1)} & \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-1)} & &+ \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-1)} &+ \cfrac{1}{10} x_{2}^{(k-1)} & &+ \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-1)} &+ \cfrac{1}{8} x_{3}^{(k-1)} & &\biggr] + (1 - \omega) x_{4}^{(k-1)} \end{alignat*}Para este ejemplo se usará $\omega = 1.05$
function metodo_gauss_jacobi_ponderado(A,b,omega)
m, n = tamaño(A)
iteraciones = ...
x_actual = [0,...,0]
for k=1 to iteraciones do
x_anterior = x_actual
for i=1 to m do
sumatoria = 0
for j=1 to n do
if (j != i) then
sumatoria = sumatoria + A(i,j)*x_anterior(j)
end if
end for
x_actual(i) = omega/A(i,i) * (b(i) - sumatoria) + (1 - omega)*x_anterior(i)
end for
end for
return x
end function
In [1]:
import numpy as np
def metodo_gauss_jacobi_ponderado(A, b, omega):
m, n = A.shape
x_actual = np.zeros(n)
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):
x_anterior = np.copy(x_actual)
for i in range(0,m):
sumatoria = b[i]
for j in range(0,n):
if j != i:
sumatoria = sumatoria - A[i,j]*x_anterior[j]
x_actual[i] = (omega/A[i,i] * sumatoria) + (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
In [2]:
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])
In [3]:
metodo_gauss_jacobi_ponderado(A,B,1.05)