Ejercicio: 1Eva_IIT2007_T2 Aplicar Gauss-Seidel
Probando solución con Jacobi, enviado como tarea
Desarrrollo Analítico
– Verificar que la matriz es diagonal dominante
No es necesario realizar el pivoteo por filas, ya la matriz tiene la diagonal dominante.
A = [[7.63, 0.30, 0.15, 0.50, 0.34, 0.84],
[0.38, 6.40, 0.70, 0.90, 0.29, 0.57],
[0.83, 0.19, 8.33, 0.82, 0.34, 0.37],
[0.50, 0.68, 0.86, 10.21, 0.53, 0.70],
[0.71, 0.30, 0.85, 0.82, 5.95, 0.55],
[0.43, 0.54, 0.59, 0.66, 0.31, 9.25]]
B = [ -9.44, 25.27, -48.01, 19.76, -23.63, 62.59]
– revisar el número de condición
cond(A) = || A||.||A-1||
El número de condición no es «muy alto», los valores de la diagonal son los mayores en toda la fila, por lo que el sistema converge.
>>> np.linalg.cond(A) 2.0451853966291011
– realizar iteraciones
Dado que no se establece en el enunciado el vector inicial, se usará el vector cero. La tolerancia requerida es 10-5
X0 = [ 0. 0. 0. 0. 0. 0.]
iteración 1
x0 = (-9.44 -0.30(0) -0.15(0)
-0.50(0) -0.34(0) -0.84(0))/7.63
= -9.44/7.63 = -1.23722149
x1 = (25.27 -0.38(0) -0.70(0)
-0.90(0) -0.29(0) -0.57(0))/6.40
= 25.27/6.40 = 3.9484375
x2 = (-48.01 -0.83(0) -0.19(0) -0.82(0)
-0.34(0) -0.37(0))/8.33
= -48.01/8.33 = -5.7635054
x3 = (19.76 -0.50(0) -0.68(0) -0.86(0)
-0.53(0) -0.70(0))/10.21
= 19.76/10.21 = 1.93535749
x4 = (-23.63 - 0.71(0) -0.30(0) -0.85(0)
-0.82(0) -0.55(0))/5.95
= -23.63/5.95 = -3.97142857
x5 = (62.59 - 0.43(0) -0.54(0) -0.59(0)
-0.66(0) -0.31(0))/9.25
= 62.59/9.25 = 6.76648649
X1 = [-1.23722149 3.9484375 -5.7635054 1.93535749 -3.97142857 6.76648649]
diferencia = X1 - X0 = X1
errado = max(|diferencia|) = 6.76648649
iteración 2
x0 = (-9.44 -0.30(3.9484375) -0.15(-5.7635054)
-0.50(1.93535749) -0.34(-3.97142857) -0.84(6.76648649))/7.63
= -9.44/7.63 = -1.23722149
x1 = (25.27 -0.38(-1.23722149) -0.70(-5.7635054)
-0.90(1.93535749) -0.29(-3.97142857) -0.57(6.76648649))/6.40
= 25.27/6.40 = 3.9484375
x2 = (-48.01 -0.83(-1.23722149) -0.19(3.9484375)
-0.82(1.93535749) -0.34(-3.97142857) -0.37(6.76648649))/8.33
= -48.01/8.33 = -5.7635054
x3 = (19.76 -0.50(-1.23722149) -0.68(3.9484375)
-0.86(-5.7635054) -0.53(-3.97142857) -0.70(6.76648649))/10.21
= 19.76/10.21 = 1.93535749
x4 = (-23.63 - 0.71(-1.23722149) -0.30(3.9484375)
-0.85(-5.7635054) -0.82(1.93535749) -0.55(6.76648649))/5.95
= -23.63/5.95 = -3.97142857
x5 = (62.59 - 0.43(-1.23722149) -0.54(3.9484375)
-0.59(-5.7635054) -0.66(1.93535749) -0.31(-3.97142857))/9.25
= 62.59/9.25 = 6.76648649
X1 = [-1.97395113 3.95743644 -6.05925771 1.96068604 -4.09171178 6.95612152]
diferencia = X1 - X0 = [-0.73672964, 0.00899894, -0.29575231, 0.02532855, -0.12028321, 0.18963504]
errado = max(|diferencia|) = 0.736729635697
iteración 3
x0 = (-9.44 -0.30(3.95743644) -0.15(-6.05925771)
-0.50(1.96068604) -0.34(-4.09171178) -0.84(6.95612152))/7.63
= -9.44/7.63 = -1.23722149
x1 = (25.27 -0.38(-1.97395113) -0.70(-6.05925771)
-0.90(1.96068604) -0.29(-4.09171178) -0.57(6.95612152))/6.40
= 25.27/6.40 = 3.9484375
x2 = (-48.01 -0.83(-1.97395113) -0.19(3.95743644)
-0.82(1.96068604) -0.34(-4.09171178) -0.37(6.95612152))/8.33
= -48.01/8.33 = -5.7635054
x3 = (19.76 -0.50(-1.97395113) -0.68(3.95743644)
-0.86(-6.05925771) -0.53(-4.09171178) -0.70(6.95612152))/10.21
= 19.76/10.21 = 1.93535749
x4 = (-23.63 - 0.71(-1.97395113) -0.30(3.95743644)
-0.85(-6.05925771) -0.82(1.96068604) -0.55(6.95612152))/5.95
= -23.63/5.95 = -3.97142857
x5 = (62.59 - 0.43(-1.97395113) -0.54(3.95743644)
-0.59(-6.059257710) -0.66(1.96068604) -0.31(-4.09171178))/9.25
= 62.59/9.25 = 6.76648649
X1 = [-1.98566781 4.0185268 -5.9920623 2.01431955 -3.98302285 7.01093224]
diferencia = X1 - X0 = [-0.01171668, 0.06109037, 0.0671954 , 0.05363351, 0.10868893, 0.05481072]
errado = max(|diferencia|) = 0.108688931048
Desarrollo numérico con Python
Se verifica el resultado obtenido realizando A.Xi y comparando con el vecto B
en la tabla se usa el signo de errado para la gráfica.
X0: [ 0. 0. 0. 0. 0. 0.] Xi: [-1.23722149 3.9484375 -5.7635054 1.93535749 -3.97142857 6.76648649] errado: 6.76648648649 Xi: [-1.97395113 3.95743644 -6.05925771 1.96068604 -4.09171178 6.95612152] errado: -0.736729635697 Xi: [-1.98566781 4.0185268 -5.9920623 2.01431955 -3.98302285 7.01093224] errado: 0.108688931048 Xi: [-2.00378293 3.99452422 -6.00443878 1.99576482 -4.0067623 6.9961552 ] errado: -0.0240025781157 Xi: [-1.99869529 4.00195452 -5.99863447 2.00153846 -3.99769932 7.00130746] errado: 0.00906298532211 Xi: [-2.00045097 3.99933614 -6.00047801 1.99948184 -4.00078219 6.99955127] errado: -0.00308287453133 Xi: [-1.99984629 4.00022733 -5.99983706 2.00017793 -3.99973154 7.00015339] errado: 0.0010506528253 Xi: [-2.00005265 3.9999222 -6.00005579 1.99993915 -4.00009178 6.9999475 ] errado: -0.0003602430897 Xi: [-1.99998199 4.00002662 -5.99998091 2.00002082 -3.99996859 7.00001796] errado: 0.000123195668206 Xi: [-2.00000616 3.99999089 -6.00000653 1.99999287 -4.00001075 6.99999385] errado: -4.21623029112e-05 respuesta de A.X=B : [-2.00000616 3.99999089 -6.00000653 1.99999287 -4.00001075 6.99999385] iteraciones: 10 A.Xi: [ -9.44006312 25.26992175 -48.01007303 19.75990236 -23.63008584 62.58992368] B: [ -9.44 25.27 -48.01 19.76 -23.63 62.59]
el gráfico de los errores vs iteraciones es:
Aplicando Gauss-Seidel
la tabla de aproximaciones sucesivas para el vector X es:
Tabla de iteraciones con AX=B: [[ 0. 0. 0. 0. 0. 0. ] [-1.23722149 4.02189753 -5.73196479 2.21089228 -3.51242077 6.91477852] [-2.0322951 3.92844105 -6.03202587 1.98957766 -3.97864919 7.00774962] [-1.9976784 4.00317297 -6.00049341 1.99807691 -4.00081785 6.99990294] [-1.99994191 4.00036665 -5.99978715 2.00000392 -4.00004739 6.99996363] [-2.00001274 3.99998231 -5.99999516 2.00000635 -3.99999579 7.00000072] [-2.00000008 3.99999833 -6.00000078 1.99999991 -3.99999985 7.00000015] [-1.99999994 4.00000007 -6.00000001 1.99999997 -4.00000002 7. ]] >>>
que se obtiene aplicando la función de Gauss-Seidel, tomando como vector inicial el vector de ceros.
Tarea: X=TX+C
# Algoritmo Gauss-Seidel # solución de matrices # métodos iterativos # Referencia: Chapra 11.2, p.310, pdf.334 # Rodriguez 5.2 p.162 import numpy as np def gauss_seidel_tabla(A,B,tolera,X = [0], itermax=100): tamano = np.shape(A) n = tamano[0] m = tamano[1] if (len(X)==n): X = np.zeros(n, dtype=float) diferencia = np.ones(n, dtype=float) errado = np.max(diferencia) tabla = [np.copy(X)] itera = 0 while (errado>tolera or itera>itermax): for i in range(0,n,1): xi = B[i] for j in range(0,m,1): if (i!=j): xi = xi-A[i,j]*X[j] xi = xi/A[i,i] diferencia[i] = np.abs(xi-X[i]) X[i] = xi errado = np.max(diferencia) itera = itera + 1 tabla.append(np.copy(X)) tabla = np.array(tabla) # No converge if (itera>itermax): X=0 return(tabla) # Programa de prueba ####### # INGRESO A = np.array([[7.63, 0.30, 0.15, 0.50, 0.34, 0.84], [0.38, 6.40, 0.70, 0.90, 0.29, 0.57], [0.83, 0.19, 8.33, 0.82, 0.34, 0.37], [0.50, 0.68, 0.86, 10.21, 0.53, 0.70], [0.71, 0.30, 0.85, 0.82, 5.95, 0.55], [0.43, 0.54, 0.59, 0.66, 0.31, 9.25]]) B = np.array([[ -9.44], [ 25.27], [-48.01], [ 19.76], [-23.63], [ 62.59]]) tolera = 0.00001 # PROCEDIMIENTO n = len(A) X = np.zeros(n, dtype=float) respuesta = gauss_seidel_tabla(A,B, tolera, X) # SALIDA print('Tabla de iteraciones con AX=B: ') print(respuesta)
En el caso de la norma infinito, para la matriz A, se puede usar el algoritmo desarrollado en clase.
Como valor para verificar su algoritmo, se obtuvo:
>>> np.linalg.norm(A, np.inf) 13.479999999999999
Tarea: incluir la norma infinito para T