stringtranslate.com

Método de Gauss-Seidel

En álgebra lineal numérica , el método de Gauss-Seidel , también conocido como método de Liebmann o método de desplazamiento sucesivo , es un método iterativo utilizado para resolver un sistema de ecuaciones lineales . Lleva el nombre de los matemáticos alemanes Carl Friedrich Gauss y Philipp Ludwig von Seidel , y es similar al método de Jacobi . Aunque se puede aplicar a cualquier matriz con elementos distintos de cero en las diagonales, la convergencia solo está garantizada si la matriz es estrictamente diagonalmente dominante , [1] o simétrica y definida positiva . Sólo se mencionó en una carta privada de Gauss a su alumno Gerling en 1823. [2] Seidel no entregó una publicación antes de 1874. [3]

Descripción

Sea un sistema cuadrado de n ecuaciones lineales, donde:

Cuando y se conocen y se desconocen, podemos utilizar el método de Gauss-Seidel para aproximar . El vector denota nuestra estimación inicial de (a menudo, de ). Lo denotamos como la k -ésima aproximación o iteración de , y es la siguiente (o k +1) iteración de .

Fórmula basada en matrices

La solución se obtiene iterativamente mediante

triangular inferiortriangular estrictamente superior[4]

Por qué funciona la fórmula matricial

El sistema de ecuaciones lineales se puede reescribir como:

El método Gauss-Seidel ahora resuelve el lado izquierdo de esta expresión para , utilizando el valor anterior para en el lado derecho. Analíticamente, esto puede escribirse como:

Fórmula basada en elementos

Sin embargo, aprovechando la forma triangular de , los elementos de se pueden calcular secuencialmente para cada fila mediante sustitución directa : [5]

Observe que la fórmula utiliza dos sumas por iteración que se pueden expresar como una suma que utiliza la iteración calculada más recientemente de . El procedimiento generalmente continúa hasta que los cambios realizados por una iteración estén por debajo de cierta tolerancia, como un residuo suficientemente pequeño .

Discusión

La fórmula de elementos del método de Gauss-Seidel es similar a la del método de Jacobi .

El cálculo de utiliza los elementos que ya se han calculado, y solo los elementos que no se han calculado en la ( k +1) -ésima iteración. Esto significa que, a diferencia del método de Jacobi, sólo se requiere un vector de almacenamiento ya que los elementos se pueden sobrescribir a medida que se calculan, lo que puede resultar ventajoso para problemas muy grandes.

Sin embargo, a diferencia del método de Jacobi, los cálculos para cada elemento son generalmente mucho más difíciles de implementar en paralelo , ya que pueden tener una ruta crítica muy larga y, por lo tanto, son más factibles para matrices dispersas . Además, los valores en cada iteración dependen del orden de las ecuaciones originales.

Gauss-Seidel es lo mismo que una relajación excesiva sucesiva con .

Convergencia

Las propiedades de convergencia del método de Gauss-Seidel dependen de la matriz A. Es decir, se sabe que el procedimiento converge si:

El método de Gauss-Seidel a veces converge incluso si no se cumplen estas condiciones.

Golub y Van Loan dan un teorema para un algoritmo que se divide en dos partes. Supongamos que no es singular. Sea el radio espectral de . Entonces las iteraciones definidas por convergen para cualquier vector inicial si es no singular y . [8]

Algoritmo

Dado que los elementos se pueden sobrescribir a medida que se calculan en este algoritmo, solo se necesita un vector de almacenamiento y se omite la indexación de vectores. El algoritmo es el siguiente:

algoritmo El método de Gauss-Seidel es  entradas:  A , b  salida:  φ Elija una suposición inicial φ para la solución  repita hasta la convergencia para  i  desde 1 hasta  n  haga  σ ← 0  para  j  desde 1 hasta  n  haga  si  ji  entonces  σσ + a ij φ j  final si  final ( j -bucle) φ i ← ( b iσ ) / a ii  final ( i -bucle) comprobar si se alcanza la convergencia terminar (repetir)

Ejemplos

Un ejemplo para la versión matricial.

Un sistema lineal mostrado como está dado por:

queremos usar la ecuacion

Debemos descomponer en la suma de una componente triangular inferior y una componente triangular superior estricta :

La inversa de es:

Ahora podemos encontrar:

Ahora tenemos y y podemos usarlos para obtener los vectores de forma iterativa.

En primer lugar, tenemos que elegir : sólo podemos adivinar. Cuanto mejor sea la suposición, más rápido funcionará el algoritmo.

Elegimos un punto de partida:

Entonces podemos calcular:

Como era de esperar, el algoritmo converge a la solución exacta:

De hecho, la matriz A es estrictamente diagonalmente dominante (pero no definida positiva).

Otro ejemplo para la versión matricial.

Otro sistema lineal mostrado está dado por:

queremos usar la ecuacion

Debemos descomponer en la suma de una componente triangular inferior y una componente triangular superior estricta :

La inversa de es:

Ahora podemos encontrar:

Ahora tenemos y y podemos usarlos para obtener los vectores de forma iterativa.

En primer lugar, tenemos que elegir : sólo podemos adivinar. Cuanto mejor sea la suposición, más rápido ejecutará el algoritmo.

Suponemos:

Entonces podemos calcular:

Si probamos la convergencia encontraremos que el algoritmo diverge. De hecho, la matriz A no es diagonalmente dominante ni definida positiva. Entonces, la convergencia a la solución exacta.

Un ejemplo para la versión de la ecuación.

Supongamos que se dan k ecuaciones donde x n son vectores de estas ecuaciones y el punto de partida x 0 . De la primera ecuación resuelve x 1 en términos de Para las siguientes ecuaciones sustituye los valores anteriores de  x s.

Para que quede claro, consideremos un ejemplo.

Resolviendo para y da:

Supongamos que elegimos (0, 0, 0, 0) como aproximación inicial, entonces la primera solución aproximada viene dada por

Utilizando las aproximaciones obtenidas, se repite el procedimiento iterativo hasta alcanzar la precisión deseada. Las siguientes son las soluciones aproximadas después de cuatro iteraciones.

La solución exacta del sistema es (1, 2, −1, 1) .

Un ejemplo usando Python y NumPy

El siguiente procedimiento numérico simplemente se itera para producir el vector solución.

importar  numpy  como  npITERACIÓN_LIMIT  =  1000# inicializamos la matriz A  =  np . matriz (  [  [ 10.0 ,  - 1.0 ,  2.0 ,  0.0 ],  [ -1.0 , 11.0 , - 1.0 , 3.0 ], [2.0 , - 1.0 , 10.0 , - 1.0 ] , [ 0.0 , 3.0 , - 1.0 , 8.0 ] , ] ) # inicializar el vector RHS b = np . matriz ([ 6.0 , 25.0 , - 11.0 , 15.0 ])                 print ( "Sistema de ecuaciones:" ) para  i  en  el rango ( A . forma [ 0 ]):  fila  =  [ f " { A [ i , j ] : 3g } *x { j + 1 } "  para  j  en  el rango ( A. forma [ 1 ])] print ( "[ {0} ] = [ {1:3g} ] " . formato ( " + " . join ( fila ), b [ i ] ))  x  =  np . zeros_like ( b ,  np . float_ ) para  it_count  en  el rango ( 1 ,  ITERATION_LIMIT ):  x_new  =  np . zeros_like ( x ,  dtype = np . float_ )  print ( f "Iteración { it_count } : { x } " )  para  i  en  el rango ( A . forma [ 0 ]):  s1  =  np . punto ( A [ i ,:  i ], x_new [: i ] ) s2 = np . punto ( A [ i , i + 1 :], x [ i + 1 :]) x_new [ i ] = ( b [ i ] - s1 - s2 ) / A [ i , i ] si np . allclose ( x , x_new , rtol = 1e-8 ): romper x = x_new                              print ( f "Solución: { x } " ) error  =  np . punto ( A ,  x )  -  b imprimir ( f "Error: { error } " )

Produce la salida:

Sistema de ecuaciones: [ 10*x1 + -1*x2 + 2*x3 + 0*x4] = [ 6] [ -1*x1 + 11*x2 + -1*x3 + 3*x4] = [ 25] [ 2*x1 + -1*x2 + 10*x3 + -1*x4] = [-11] [ 0*x1 + 3*x2 + -1*x3 + 8*x4] = [ 15] Iteración 1: [ 0 . 0. 0. 0.] Iteración 2: [ 0.6 2.32727273 -0.98727273 0.87886364] Iteración 3: [ 1.03018182 2.03693802 -1.0144562 0.98434122] Iteración 4: [ 1.006585 04 2.00355502 -1.00252738 0.99835095] Iteración 5 : [ 1.00086098 2.00029825 -1.00030728 0.99984975] Iteración 6 : [ 1.00009128 2.00002134 -1.00003115 0.9999881 ] Iteración 7: [ 1.00000836 2.00000117 -1.00000275 0.99999922] Iteración 8: [ 1.00000067 2.0000 0002 -1.00000021 0.99999996] Iteración 9: [ 1.00000004 1.99999999 -1.00000001 1. ] Iteración 10: [ 1. 2. -1. 1.] Solución: [ 1. 2. -1. 1.] Error: [2.06480930e-08 -1.25551054e-08 3.61417563e-11 0.00000000e+00]

Programa para resolver arbitrario no. de ecuaciones usando Matlab

El siguiente código utiliza la fórmula

función  x = gauss_seidel ( A, b, x, iteros ) para i = 1 : iteros para j = 1 : tamaño ( A , 1 ) x ( j ) = ( b ( j ) - suma ( A ( j ,:) ' .* x ) + A ( j , j ) * x ( j )) / A ( j , j ); fin fin fin                     

Ver también

Notas

  1. ^ Sauer, Timoteo (2006). Análisis numérico (2ª ed.). Pearson Education, Inc. pág. 109.ISBN​ 978-0-321-78367-7.
  2. ^ Gauss 1903, pag. 279; enlace directo.
  3. ^ Seidel, Ludwig (1874). "Über ein Verfahren, die Gleichungen, auf welche die Methode der kleinsten Quadrate führt, sowie lineäre Gleichungen überhaupt, durch sucesive Annäherung aufzulösen" [Sobre un proceso para resolver por aproximaciones sucesivas las ecuaciones a las que conduce el método de mínimos cuadrados y el método lineal ecuaciones en general]. Abhandlungen der Mathematisch-Physikalischen Klasse der Königlich Bayerischen Akademie der Wissenschaften (en alemán). 11 (3): 81–108.
  4. ^ Préstamo Golub y Van 1996, pág. 511.
  5. ^ Golub & Van Loan 1996, ecuación (10.1.3)
  6. ^ Golub & Van Loan 1996, Thm 10.1.2.
  7. ^ Bagnara, Roberto (marzo de 1995). "Una prueba unificada de la convergencia de los métodos de Jacobi y Gauss-Seidel". Revisión SIAM . 37 (1): 93–97. CiteSeerX 10.1.1.26.5207 . doi :10.1137/1037008. JSTOR  2132758. 
  8. ^ Golub & Van Loan 1996, Thm 10.1.2

Referencias

Este artículo incorpora texto del artículo Gauss-Seidel_method en CFD-Wiki que se encuentra bajo la licencia GFDL .


enlaces externos