Metodos Iterativos Jacobi Gauss Seidel 4w1du

5.5 MÉTODOS ITERATIVOS Se sabe que es posible usar un esquema iterativo que mejore la solución hasta que se obtenga convergencia. Esta es la esencia de este tipo de métodos empleados para resolver sistemas de ecuaciones lineales. Los más comunes son el método de jacobi y el de Gauss-Seidel. Ambos parten del mismo principio, pero el esquema difiere en la forma de sustentación. 5.5.1 MÉTODO DE JACOBI El método Jacobi es el método iterativo para resolver sistemas de ecuaciones lineales más simple y se aplica solo a sistemas cuadrados, es decir a sistemas con tantas incógnitas como ecuaciones. 1. Primero se determina la ecuación de recurrencia. Para ello se ordenan las ecuaciones y las incógnitas. De la ecuación i se despeja la incógnita i. En notación matricial se escribirse como: x = c + B x (1) Donde x es el vector de incógnitas. 2. Se toma una aproximación para las soluciones y a ´esta se le designa por xo 3. Se itera en el ciclo que cambia la aproximación Xi +1 = c + Bxi Solución Debemos primeramente despejar de la ecuación la incógnita correspondiente. x = 0.20 + 0.00 x − 0.40 y y = 0.00 + 0.25 x + 0.00 y Escrito en la notación vectorial quedaría: 𝑥 0.20 0.00 − 0.40 𝑥 [𝑦] = [ ]+[ ][ ] 0.00 0.25 0.00 𝑦 Aplicamos la primera iteración partiendo de x0 = 1.00 y y0 = 2.00: x1 = 0.20 + 0.00 (1.00) − 0.40 (2.00) = −0.60 y1 = 0.00 + 0.25 (1.00) + 0.00 (2.00) = 0.25 Aplicamos la segunda iteración partiendo de x1 = −0.60 y y1 = 0.25: x2 = 0.20 + 0.00 (−0.60) − 0.40 (0.25) = 0.10 y2 = 0.00 + 0.25 (−0.60) + 0.00 (0.25) = −0.15 Aplicamos la siguiente iteración partiendo de x2 = 0.10 y y1 = −0.15: x3 = 0.20 + 0.00 (0.10) − 0.40 (−0.15) = 0.26 y3 = 0.00 + 0.25 (0.10) + 0.00 (−0.15) = 0.025 Aplicamos la siguiente iteración partiendo de x3 = 0.26 y y3 = 0.025: x4 = 0.20 + 0.00 (0.26) − 0.40 (0.025) = 0.190 y4 = 0.00 + 0.25 (0.26) + 0.00 (0.025) = 0.065

Aplicamos la siguiente iteración partiendo de x4 = 0.190 y y4 = 0.065: x5 = 0.20 + 0.00 (0.19) − 0.40 (0.065) = 0.174 y5 = 0.00 + 0.25 (0.19) + 0.00 (0.065) = 0.0475 Aplicamos la siguiente iteración partiendo de x5 = 0.174 y y5 = 0.0475: x6 = 0.20 + 0.00 (0.174) − 0.40 (0.0475) = 0.181 y6 = 0.00 + 0.25 (0.174) + 0.00 (0.0475) = 0.0435 Si uno dispone de una hoja de cálculo como EXCEL es fácil realizar los cálculos

anteriores: Donde Di = máx (|xi − xi+1|, |yi − yi+1|) Este Di es utilizado como criterio de paro en las iteraciones: Cuando Di es menos que cierto valor dado (digamos 0.001) uno ya no realiza la siguiente iteración. Si se grafica las aproximaciones obtenidas en el plano x − y se obtendrá algo como:

La forma más conveniente de expresar estos métodos es usar la notación matricial. La matriz A se divide en tres partes, correspondientes a los tres conjuntos de coeficientes. La ecuación AX = B será (L + D + U)X = B Donde L es una matriz triangular inferior con ceros e la diagonal, D es la matriz diagonal que contiene los elementos de la diagonal de A, y U es una matriz triangular superior con ceros en la diagonal. El método de Jacobi para la r-ésima iteración se escribe como sigue 𝑫𝑿(𝒓) = −(𝑳 + 𝑼)𝑿(𝒓−𝟏) + 𝑩, (r = 0, 1, ….) El método de jacobi converge para matrices A que son diagonalmente dominantes, en el sentido que es matemáticamente preciso. Para determinar si el método es convergente, si se tiene que la matriz de iteración − 𝑫−𝟏 (𝑳 + 𝑼) Además de un término aditivo, − 𝑫−𝟏 𝑩, mapea un conjunto de X en el siguiente. Esta matriz de iteración tiene valores propios, cada uno de los cuales refleja el factor por el que se suprime durante una iteración el valor propio particular de algún residual no deseado. Es evidente que estos factores tienen un módulo menor que 1. Con modulo menor que uno. La rapidez de convergencia de este método está dada por la rapidez de disminución del valor propio más lento: es decir, el factor con el modulo más grande. El módulo de este factor esta, por consiguiente, entre 0 y 1, y se le llama radio espectral del operador de relajación (𝒑𝒔 ). El número de iteraciones ® necesarias para reducir el error toral por un factor (𝟏𝟎−𝒑 ) se estima por: 𝑟=

𝑝 ln(10) − ln( 𝑝𝑠 )

En el límite cuando 𝒑𝒔 tiene a 1, el sistema esta en el punto de ruptura. Aquí puede o no converger, dependiendo de las condiciones iniciales. La idea de iterar se puede aplicar fácilmente a ecuaciones simultáneas lineales. Esto se ilustra mediante un ejemplo numérico. EJEMPLO Mediante el método de jacobi resolver el conjunto de ecuaciones lineales dadas por: 10𝑥1 + 𝑥2 + 𝑥3 = 24 −𝑥1 + 20𝑥2 + 𝑥3 = 21 𝑥1 − 2𝑥2 + 100𝑥3 = 300 Despejando de la primera ecuación la primera variable y así sucesivamente, se obtiene

1 ( 24 − 𝑥2 − 𝑥3 ) 10 1 𝑥2 = ( 21 + 𝑥1 − 𝑥3 ) 20 1 𝑥3 = ( 300 − 𝑥1 + 2𝑥2 ) 100 𝑥1 =

Si las primeras aproximaciones son 𝑥1 = 0, 𝑥2 = 0 𝑦 𝑥3 = 0 entonces la siguiente iteración es 𝑥1 = 2.4, 𝑥2 = 1.05 𝑦 𝑥3 = 3 Y las iteraciones sucesivas son, 𝑥1 = 1.995 , 𝑥2 = 1.02 𝑦 𝑥3 = 2.997 Y después 𝑥1 = 1.9983 , 𝑥2 = 0.9999 𝑦 𝑥3 = 2.9995 La sustitución de los valores 𝒙𝟏 = 𝟐 , 𝒙𝟐 = 𝟏 𝒚 𝒙𝟑 = 𝟑 en el sistema inicial demuestra que el proceso está convergiendo rápidamente a la solución verdadera. Si se utiliza las ecuaciones − 𝑫−𝟏 (𝑳 + 𝑼), 𝒓 =

𝒑 𝐥𝐧(𝟏𝟎) − 𝐥𝐧(𝒑𝒔 )

se tiene que el radio espectral y el número de

iteraciones para reducir el error total por un facto (𝟏𝟎−𝟑 ) Son, respectivamente (𝒑𝒔 = 𝟎. 𝟎𝟕𝟓𝟓) y (𝒓 = 𝟑). Es claro que las iteraciones pueden ser un método muy útil de solución. Se demuestra ahora que este también puede ser un método inapropiado en ciertas circunstancias. Las ecuaciones del ejemplo se pueden tomar en orden diferente: 𝒙𝟏 − 𝟐𝒙𝟐 + 𝟏𝟎𝟎𝒙𝟑 = 𝟑𝟎𝟎 𝟏𝟎𝒙𝟏 + 𝒙𝟐 + 𝒙𝟑 = 𝟐𝟒 −𝒙𝟏 + 𝟐𝟎𝒙𝟐 + 𝒙𝟑 = 𝟐𝟏 Con este acomodo, el radio espectral del sistema es (𝒑𝒔 = 𝟐𝟕. 𝟕𝟖𝟎𝟐), lo cual lo hace divergente. Si se calcula el número de iteraciones con la ecuación 𝒓 =

𝒑 𝐥𝐧(𝟏𝟎) , − 𝐥𝐧(𝒑𝒔 )

se

obtiene un numero negativo, lo cual no tiene sentido. A continuación se muestra que esta forma de ordenar las ecuaciones es divergente, resolviendo el sistema numéricamente. Despejando de la primera ecuación la primera variable y así sucesivamente, se obtiene 𝒙𝟏 = 𝟑𝟎𝟎 + 𝒙𝟐 − 𝟏𝟎𝟎𝒙𝟑 𝒙𝟐 = 𝟐𝟒 − 𝟏𝟎𝒙𝟏 − 𝒙𝟑 𝒙𝟑 = 𝟐𝟏 + 𝒙𝟏 − 𝟐𝟎𝒙𝟐 Utilizando los valores 𝒙𝟏 = 𝟎, 𝒙𝟐 = 𝟎 𝒚 𝒙𝟑 = 𝟎 para a primera aproximación, se obtiene las aproximaciones sucesivas 𝒙𝟏 = 𝟑𝟎𝟎, 𝒙𝟐 = 𝟐𝟒 𝒚 𝒙𝟑 = 𝟐𝟏 Y 𝒙𝟏 = −𝟏𝟕𝟓𝟐, 𝒙𝟐 = −𝟐𝟗𝟗𝟕 𝒚 𝒙𝟑 = −𝟑𝟖𝟕𝟗

Y no hay probabilidad de que esta secuencia converja a los valores 𝒙𝟏 = 𝟐, 𝒙𝟐 = 𝟏 𝒚 𝒙𝟑 = 𝟑. Por tanto, es crucial obtener algún conocimiento para que se pueda asegurar que las iteraciones converjan. El método anterior, conocido como el método de jacobi, para vez se usa porque existen varias mejores que se pueden utilizar para elevar la rapidez de convergencia. Sin embargo, si se cuenta con una computadora que efectué los cálculos en paralelo, se recomienda este método por ser altamente paralelizable. 5.5.2 MÉTODO DE GAUSS-SEIDEL El método de Gauss-Seidel es muy semejante al método de Jacobi. Mientras que en el de Jacobi se utiliza el valor de las incógnitas para determinar una nueva aproximación, en el de Gauss-Seidel se va utilizando los valores de las incógnitas recién calculados en la misma iteración, y no en la siguiente. Por ejemplo, en el método de Jacobi se obtiene en el primer cálculo x i+1, pero este valor de x no se utiliza sino hasta la siguiente iteración. Por ejemplo, en el método de Jacobi se obtiene en el primer cálculo x i+1, pero este valor de x no se utiliza sino hasta la siguiente iteración. En el método de Gauss-Seidel en lugar de eso se utiliza de x i+1 en lugar de x i en forma inmediata para calcular el valor de y i+1 de igual manera procede con las siguientes variables; Por ejemplo, en el método de Jacobi se obtiene en el primer cálculo x i+1, pero este valor de x no se utiliza sino hasta la siguiente iteración. En el método de Gauss-Seidel en lugar de eso se utiliza de x i+1 en lugar de x i en forma inmediata para calcular el valor de y i+1 de igual manera procede con las siguientes variables; siempre se utilizan las variables recién calculadas. EJEMPLO 

Ejemplo Partiendo de ( x = 1, y = 2) aplique dos iteraciones del método de GaussSeidel para resolver el sistema: "

5x+2y=1 x−4y=0 Solución Debemos primeramente despejar de la ecuación la incógnita correspondiente. x = 0.20 + 0.00 x − 0.40y y = 0.00 + 0.25 x + 0.00y Aplicamos la primera iteración partiendo de x 0 = 1.00 y y 0 = 2.00: x 1 = 0.20 + 0.00 (+1.000) − 0.40 (2.00) = − 0.600 y 1 = 0.00 + 0.25 ( − 0.600) + 0.00 (2.00) = − 0.15 Aplicamos la segunda iteración partiendo de x 1 = − 0.600 y y 1 = − 0.15:



x 2 = 0.20 + 0.00 ( − 0.600) − 0.40 ( − 0.15) = 0.26 y 2 = 0.00 + 0.25 ( 0.26) + 0.00 ( − 0.15) = 0.065 Ejemplo Partiendo de (x = 1, y = 2, z = 0) aplique dos iteraciones del método de Gauss-Seidel para resolver el sistema: 10 x + 0 y − z = − 1 4 x + 12 y − 4 z = 8

4 x + 4 y + 10 z = 4 Solución Debemos primeramente despejar de la ecuación la incógnita correspondiente. x = − 0.10 + 0.00 x + 0.00 y + 0.10 z y = 0.66 − 0.33 x + 0.00 y + 0.33 z z = 0.40 − 0.40 x − 0.40 y + 0.00 z Aplicamos la primera iteración partiendo de x 0 = 1.00, y 0 = 2.00, y z = 0.00: x 1 = − 0.10 + 0.00(1.00) + 0.00 (2.00) + 0.10 (0.00) = − 0.1 y 1 = 0.66 − 0.33( − 0.10) + 0.00 (2.00) + 0.33 (0.00) = 0.70 z 1 = 0.40 − 0.40( − 0.10) − 0.40 ( 0.70) + 0.00 (0.00) = 0.16 Aplicamos la primera iteración partiendo de x 0 = 1.00, y 0 = 2.00, y z = 0.00: x 1 = − 0.10 + 0.00(1.00) + 0.00 (2.00) + 0.10 (0.00) = − 0.1 y 1 = 0.66 − 0.33( − 0.10) + 0.00 (2.00) + 0.33 (0.00) = 0.70 z 1 = 0.40 − 0.40( − 0.10) − 0.40 ( 0.70) + 0.00 (0.00) = 0.16 Aplicamos la segunda iteración partiendo de x 1 = − 0.10 y y 1 = 0.70 y z 1 = 0.16: x 1 = − 0.10 + 0.00( − 0.10) + 0.00 (0.70) + 0.10 (0.16) = − 0.084 y 1 = 0.66 − 0.33( − 0.084) + 0.00 (0.70) + 0.33 (0.16) = 0.748 z 1 = 0.40 − 0.40( − 0.084) − 0.40 ( 0.748) + 0.00 (0.16) = 0.134 

Utilice el método iterativo de gauss-seidel para encontrar las soluciones aproximadas de 10𝑥1 + 𝑥2 + 2𝑥3 = 6 −𝑥1 + 11𝑥2 − 𝑥3 + 3𝑥4 = 25 2𝑥1 − 𝑥2 + 10𝑥3 − 𝑥4 = −11 3 𝑥2 − 𝑥3 + 8𝑥4 = 15

Comenzando con x=(0,0,0,0)2 e iterando hasta que ||𝑥 (𝑘) − 𝑥 (𝑘−1) ||∞ |||𝑥 (𝑘) ||∞

< 10−3

La solución está por el método de gauss-seidel escribimos el sistema, para cada k = 1, 2,….. Como (𝑘)

𝑥1

(𝑘)

𝑥2

(𝑘)

𝑥3

=

=

1 (𝑘−1) 1 (𝑘−1) 3 𝑥 − 𝑥3 + 10 2 5 5

1 (𝑘) 1 (𝑘−1) 3 (𝑘−1) 25 𝑥1 + 𝑥3 − 𝑥 + 11 11 11 4 11

=−

1 (𝑘) 1 (𝑘) 1 (𝑘−1) 11 𝑥1 + 𝑥2 + 𝑥 − 5 10 10 4 10

(𝑘)

𝑥4

=−

3 (𝑘) 1 (𝑘) 15 𝑥 + 𝑥3 + 8 2 8 8

Cuando 𝑥 (0) = (0,0,0,0)1 , 𝑡𝑒𝑛𝑒𝑚𝑜𝑠 𝑥 (1) =(0.6000, 2.3272, -0.9873, 0.8789)´. Las iteraciones subsecuentes generan los valores de la tabla k (𝑘) 𝑥1 (𝑘) 𝑥2 (𝑘) 𝑥3 (𝑘) 𝑥4

0 0.0000

1 0.6000

2 1.030

3 1.0065

4 1.0009

5 1.0001

0.0000

2.3272

2.037

2.0036

2.0003

2.0000

0.0000

-0.9873

-1.014

-1.0025

-1.0003

-1.0000

0.0000

0.8789

0.9844

0.9983

0.9983

1.0000

Dado que ||𝑥 (5) − 𝑥 (5−1) ||∞ |||𝑥 (5) ||∞

=

0.0008 = 4 ∗ 10−4 2.000

Se acepta 𝑥 (5) como una aproximación razonable a la solución. Note que en el método de Jacobi requirió el doble de iteraciones para alcanzar el mismo grado de exactitud. Si se quiere expresar el método de Gauss-Seidel en forma matricial, multiplicamos ambos lados de la ecuación por 𝑎Ü y reunimos todos los k-ésimos términos de iteración, los que nos da (𝑘)

(𝑘)

(𝑘)

𝑎𝑖1 𝑥3 + 𝑎𝑖2 𝑥2 + ⋯ + 𝑎𝑖𝑖 𝑥𝑖

(𝑘−1)

= −𝑎𝑖,𝑖+1 𝑥𝑖+1

(𝑘−1)

− ⋯ − 𝑎𝑖𝑛 𝑥𝑛

+ 𝑏𝑖

Por cada i= 1, 2, …, n. al escribir todas las n ecuaciones obtenemos (𝑘)

𝑎11 𝑥1

(𝑘−1)

= −𝑎12 𝑥2

(𝑘)

(𝑘)

𝑎21 𝑥1 + 𝑎22 𝑥2

(𝑘−1)

− 𝑎13 𝑥3

(𝑘−1)

= −𝑎23 𝑥3

(𝑘)

(𝑘)

(𝑘−1)

… − 𝑎1𝑛 𝑥𝑛

+ 𝑏1 ,

(𝑘−1)

− ⋯ − 𝑎2𝑛 𝑥𝑛

+ 𝑏2

(𝑘)

𝑎𝑛1 𝑥1 + 𝑎𝑛2 𝑥2 + ⋯ + 𝑎𝑛𝑛 𝑥𝑛 = +𝑏𝑛𝑖

Y con las definiciones de D, L y U, tenemos que el método de Gauss-Seidel esta representado por (𝐷 − 𝐿)𝑥 (𝑘) = 𝑈𝑥 (𝑘−1) + 𝑏 Y 𝑥 (𝑘) = (𝐷 − 𝐿)−1 𝑈𝑥 (𝑘−1) + (𝐷 − 𝐿)−1 𝑏 . 𝑝𝑎𝑟𝑎 𝑐𝑎𝑑𝑎 𝑘 = 1, 2, …. Si usamos 𝑆𝑒𝑖𝑑𝑒𝑙 𝑡𝑖𝑒𝑛𝑒 𝑙𝑎 𝑓𝑜𝑟𝑚𝑎

𝑇𝑔 = (𝐷 − 𝐿)−1 𝑈 𝑦 𝐶𝑔 = (𝐷 − 𝐿)−1 𝑏, 𝑒𝑙 𝑚é𝑡𝑜𝑑𝑜 𝑑𝑒 𝐺𝑎𝑢𝑠𝑠 − 𝑥 (𝑘) = 𝑇𝑔 𝑥 (𝑘−1) + 𝐶𝑔

Para que la matriz triangular inferior D-L sea no singular, es necesario y suficiente que 𝑎𝑖𝑖 ≠ 0 para cada i=1, 2,…..n. Método iterativo de Gauss-Seidel Para resolver Ax = b dada un aproximación inicial 𝑥 (0) :

ENTRADA el número de ecuaciones e incógnitas n; los elementos 𝑎𝑖𝑗 , 1 ≤ 𝑖, 𝑗 ≤ 𝑛 de la matriz A; los elementos 𝑏𝑖 , 1 ≤ 𝑖 ≤ 𝑛 de b; los elementos 𝑋𝑂𝑖 , 1 ≤ 𝑖 ≤ 𝑛 de 𝑋𝑂 = 𝑥 (0) , la tolerancia TOL; el número de iteraciones N SALIDA la solución aproximada 𝑥𝑖 , . . , 𝑥𝑛 o el mensaje de que se rebaso el numero de iteraciones Paso1: haga k=1 Paso2: mientras (k ≤ N) haga los pasos 3-6 Paso3: para i = 1,…, n 1

𝑛 Haga 𝑥1 = 𝑎 [− ∑𝑖−1 𝑗=1 𝑎𝑖𝑗 𝑥𝑗 − ∑𝑗=𝑖+1 𝑎𝑖𝑗 𝑋𝑂𝑗 + 𝑏𝑖 𝑖𝑖

Paso4: Si||x – XO|| < TOL (procedimiento terminado exitosamente)

entonces

PARAR Paso5: Haga k=k+1 Paso6: Para i=1,….., n tome 𝑋𝑂𝑖 = 𝑥𝑖 Paso7: Salida (´Número máximo de iteraciones excedido´); (Procedimiento terminado exitosamente) Parar

SALIDA

(𝑥𝑖 , . . , 𝑥𝑛 ):