Matrices dispersas

Sistemas de Ecuaciones Lineales V
´ ´ • Metodos Iterativos: Matrices dispersas. Esquema general. Metodos de Jacobi y de Gauss-Seidel.

521230

-1-

´ DIM – Universidad de Concepcion

Matrices dispersas
´ • Cuando la matriz A del sistema a resolver es dispersa, pero no banda, los metodos ´ (directos) estudiados hasta ahora (eliminacion de Gauss o Cholesky) presentan el defectodenominado llenado (fill-in). ´ • El llenado consiste en que, a medida que el proceso de eliminacion avanza, se van creando elementos no nulos en posiciones de L y U en donde la matriz A tiene ceros.

• Como consecuencia del llenado se tiene, por una parte, el aumento del numero de flop y ´
con ello el aumento del error de redondeo. Por otra parte se tiene el aumento en las necesidades de memoria paraalmacenar las matrices L y U , lo que puede llegar a hacer ´ ˜ imposible aplicar estos metodos cuando A es de gran tamano. ´ • Los metodos que estudieremos en seguida, llamados iterativos, evitan el llenado y sus consecuencias, al trabajar resolviendo reiteradamente sistemas con matriz diagonal o triangular-dispersa.

521230

-2-

´ DIM – Universidad de Concepcion

Llenado de matricesdispersas

521230

-3-

´ DIM – Universidad de Concepcion

Llenado de matrices dispersas (cont.)
>> load data.0125.mat >> [L,U]=lu(A); >> whos Name Size A L U b 10821x10821 10821x10821 10821x10821 10821x1

Bytes 951580 151325524 151325524 129860

Class sparse sparse sparse sparse array array array array

Grand total is 25300219 elements using 303732496 bytes

521230

-4-

´ DIM –Universidad de Concepcion

Esquema general
• Considere el sistema de ecuaciones Ax = b,
con A

∈ Rn×n no singular y b ∈ Rn .

´ Un metodo iterativo para resolver el sistema construye, a partir de un vector inicial x(0) , ´ una sucesion de vectores x(1) , x(2) , . . . , x(k) , . . . la que, bajo condiciones apropiadas, ´ resultara convergente a x.

• Si suponemos A = N − P , donde N debeser invertible, entonces Ax = b ⇐⇒ Nx = P x + b ⇐⇒ x = N −1 P x + N −1 b.

521230

-5-

´ DIM – Universidad de Concepcion

Esquema general (cont.)
• Se usa la igualdad N x = P x + b para definir un esquema general para construir la ´ sucesion {x(k) }. • Algoritmo del esquema general:
Dado el vector inicial x(0) , para k

= 1, 2, . . . resolver:

N x(k) = P x(k−1) + b,
´ hasta que sesatisfaga un criterio de detencion. ´ • Definiendo M := N −1 P (matriz de iteracion) y e(k) := x − x(k) (error de x(k) ), para cada k = 1, 2, . . . se tiene

e(k) = M k e(0) ,

k = 1, 2, . . .

521230

-6-

´ DIM – Universidad de Concepcion

´ Convergencia de metodos iterativos.
´ ´ • Teorema. (Convergencia) La sucesion x(k) converge a la solucion x de Ax = b, si y ´ solo si, ρ(M ) <1. ´ ´ ´ • Observacion. Si la sucesion x(k) converge, necesariamente lo hace a la solucion x de Ax = b.

• Lema. (Cota para el radio espectral) Sea A una matriz cuadrada. Para cualquier norma
matricial se tiene que

ρ(A) ≤ A .

´ ´ • Corolario. (Condicion suficiente de convergencia) Una condicion suficiente para que la ´ ´ sucesion x(k) sea convergente a la solucion x de Ax = b es que

M <1,
´ ´ donde M es la matriz de iteracion del metodo que genera a
521230 -7-

x(k) .
´ DIM – Universidad de Concepcion

´ Criterio de detencion
´ ´ • Detencion del proceso. Cuando el proceso iterativo es convergente, este se debe detener para un x(k+1) tal que e(k+1) = x − x(k+1) ≤ tol, donde tol indica un nivel de tolerancia prefijado para el error.

• Lema. Para M < 1, se tiene que x−x(k+1)

M ≤ 1− M

x(k+1) − x(k) .

´ • Un criterio de detencion usual consiste en detener el proceso cuando

x(k+1) − x(k) ≤ tol.
Sin embargo, este criterio resulta muchas veces inadecuado!

M En efecto, si 1− M

≫ 1, usualmente, M 1− M x(k+1) − x(k) > tol.
´ DIM – Universidad de Concepcion

521230

-8-

´ Criterio de detencion (cont.)
´ • Observacion. Al graficar

M , F...