D.

C D
D. R P
Análisis numérico

Ricardo Prato
Universidad de Antioquia

Semestre I de 2024 Semana 01 May 24, 2025

2 / 31
Métodos iterativos

D. C D

D. R P
2 / 31
 3 / 31

D. C D

D. R P
Teorema
Sea 𝐴 ∈ ℝ𝑛×𝑛 y ‖ ⋅ ‖ una norma sobre ℝ𝑛×𝑛 entonces

𝜌(𝐴) ≤ ‖𝐴‖.

Para cada matriz 𝐴 y para cada 𝜀 > 0 existe una norma sobre ℂ𝑛 tal
que
‖𝐴‖ ≤ 𝜌(𝐴) + 𝜀

3 / 31
 4 / 31
Lema
Sea 𝐵 ∶ 𝕂 → 𝕂 una transformación linear sobre 𝕂 (𝕂 = ℝ ó 𝕂 = ℂ) con

D. C D

‖𝐵‖ < 1 , sea 𝐼 la trasformación identidad . Entonces

D. R P

𝐼 −𝐵

es invertible, y para cada z ∈ 𝕂 la ecuación

(𝐼 − 𝐵)x = x − 𝐵x = z

tiene solución única.

Para un z fijo pero arbitrario, se define 𝐴x ∶= 𝐵x + z, x ∈ 𝕂 entonces

‖𝐴x − 𝐴y‖ = ‖𝐵x − 𝐵y‖ ≤ ‖𝐵‖‖x − y‖,

por tanto, 𝐴 es una contracción y tiene un único punto fijo x̃

𝐴x̃ = x̃ = 𝐵x̃ + z = x̃ ⇒ x − 𝐵x = (𝐼 − 𝐵)x = z

4 / 31
 5 / 31
Métodos iterativos

D. C D

D. R P
Recordar:
Sea 𝐴 ∈ ℝ𝑛×𝑛 y ‖ ⋅ ‖ una norma sobre ℝ𝑛×𝑛 entonces

radio espectral de 𝐴 = 𝜌(𝐴) ≤ ‖𝐴‖

Teorema
Sea 𝐵 ∈ ℂ𝑛×𝑛 . Las sucesivas aproximaciones

x𝜈+1 ∶= 𝐵 x𝜈 + z, 𝜈 ∈ ℕ0

converge para cada z ∈ ℂ𝑛 y cada x0 ∈ ℂ𝑛 , si y solo si,

𝜌(𝐵) < 1

5 / 31
 6 / 31
Métodos iterativos

D. C D

Teorema

D. R P

Sea 𝐵 ∈ ℂ𝑛×𝑛 . Las sucesivas aproximaciones

x𝜈+1 ∶= 𝐵 x𝜈 + z, 𝜈 ∈ ℕ0

converge para cada z ∈ ℂ𝑛 y cada x0 ∈ ℂ𝑛 , si y solo si,

𝜌(𝐵) < 1

Sketch de la prueba
1. Si 𝜌(𝐵) < 1 Existe ‖ ⋅ ‖ en ℂ𝑛 tal que ‖𝐵‖ < 1 , por tanto 𝐴𝑥 ∶= 𝐵𝑥 + 𝑧 define
una contracción. Es decir, Las sucesivas aproximaciones x𝜈+1 ∶= 𝐵 x𝜈 + z, 𝜈 ∈ ℕ0
converge para cada z ∈ ℂ𝑛 y cada x0 ∈ ℂ𝑛
2. Por RAA. Si 𝜌(𝐵) ≥ 1 ∃𝜆 ∈ 𝜎(𝐵) tal que |𝜆| > 1 Con (𝜆, 𝑤) el par
autovalor-autovector respectivo. Tome 𝑧 = 𝑥0 = 𝑤 ¿Qué concluye?

6 / 31
 3 / 54
Métodos iterativos

D. C D

Sea 𝐵 ∈ ℂ𝑛×𝑛 . Las sucesivas aproximaciones

D. R P

x𝜈+1 ∶= 𝐵 x𝜈 + z, 𝜈 ∈ ℕ0

converge para cada z ∈ ℂ𝑛 y cada x0 ∈ ℂ𝑛 , si y solo si, 𝜌(𝐵) < 1

Resolución sistemas lineales 𝐴x = 𝑏

• 𝐴 = 𝑃 − 𝑁 donde 𝑃 es no singular
• Resolver 𝐴x = 𝑏 es equivalente a resolver el problema de punto
fijo
x = 𝑃 −1 𝑁 x + 𝑃 −1 𝑏

x(𝑘+1) = 𝑃 −1 𝑁 x(𝑘) + 𝑃 −1 𝑏
x(𝑘+1) = x(𝑘) + 𝑃 −1 𝑟(𝑘) , 𝑟(𝑘) = 𝑏 − 𝐴x(𝑘)

3 / 54
 4 / 54
Métodos iterativos estacionarios

D. C D

0

D. R P

⎛
⎜ ⎞
⎟
⎜ 𝑎21 0 ⎟
𝐴𝐿 ⎜
∶= ⎜ ⎟
⎟
⎜
⎜ ⋮ ⋱ ⋱ ⎟
⎟ 𝐴 = 𝐷 + 𝐴𝑅 + 𝐴 𝐿
⎜ ⎟
⎝
𝑎𝑛1 ⋯ 𝑎𝑛,𝑛−1 0
⎠ 𝐴x = y
(𝐷 + 𝐴𝑅 + 𝐴𝐿 )x = y
𝑎11 𝐷x = −(𝐴𝑅 + 𝐴𝐿 )x + y
⎛ ⎞
𝐷 ∶= ⎜
⎜ ⋱
⎟
⎟
⎜ ⎟ x = −𝐷−1 (𝐴𝑅 + 𝐴𝐿 )x + 𝐷−1 y
⎝ 𝑎𝑛𝑛 ⎠

0 𝑎12 ⋯ 𝑎1𝑛 Método de Jacobi

⎛
⎜ ⎞
⎟
⎜ 0 ⋱ ⋮ ⎟ x𝑛+1 ∶= −𝐷−1 (𝐴𝑅 + 𝐴𝐿 )x𝑛 + 𝐷−1 y
𝐴𝑅 ∶= ⎜
⎜
⎜
⎟
⎟
⎜ ⋱ 𝑎𝑛−1,𝑛 ⎟
⎟
⎜ ⎟
⎝ 0 ⎠

4 / 54
 5 / 54
Ejemplo del Sistema 4x4

D. C D

Consideremos el siguiente sistema de ecuaciones:

D. R P

4𝑥1 + 𝑥2 − 𝑥3 + 0𝑥4 = 7
𝑥1 − 3𝑥2 + 𝑥3 + 𝑥4 = −6
−𝑥1 + 𝑥2 + 5𝑥3 − 2𝑥4 = 8
0𝑥1 + 𝑥2 − 2𝑥3 + 4𝑥4 = 9

En forma matricial:

4 1 −1 0 ⎛𝑥1 ⎞ 7
⎛
⎜ ⎞
⎟ ⎜ ⎟ ⎛⎜ ⎞
⎜
⎜ 1 −3 1 1 ⎟⎟
⎜
⎜𝑥2⎟⎟ ⎜ −6⎟
⎟
⎜ ⎟ ⎜
⎜ ⎟
⎟ =⎜ ⎟
⎜ ⎟
⎜
⎜−1 1 5 −2⎟ ⎟⎜ 𝑥3 ⎟ ⎜⎜ 8⎟⎟
⎜ ⎟
⎝ 0 1 −2 4 ⎠ ⎝𝑥4 ⎠ ⎝ ⎠9

5 / 54
 6 / 54

D. C D

Despejando cada variable, el método de Jacobi equivale a un problema de

D. R P

punto fijo

(𝑘+1) 1 (𝑘) (𝑘)

𝑥1 = (7 − 𝑥2 + 𝑥3 )
4
(𝑘+1) 1 (𝑘) (𝑘) (𝑘)
𝑥2 = (−6 − 𝑥1 − 𝑥3 − 𝑥4 )
−3
(𝑘+1) 1 (𝑘) (𝑘) (𝑘)
𝑥3 = (8 + 𝑥1 − 𝑥2 + 2𝑥4 )
5
(𝑘+1) 1 (𝑘) (𝑘)
𝑥4 = (9 − 𝑥2 + 2𝑥3 )
4
(𝑘)
Donde 𝑥𝑖 nos indica la variable 𝑥𝑖 en la 𝑘-ésima iteración de Jacobi.

6 / 54
 7 / 54
La Matriz de Iteración de Jacobi

D. C D

D. R P
Podemos expresar el método de Jacobi en forma matricial como:

𝑥(𝑘+1) = 𝐵𝐽 𝑥(𝑘) + 𝑐𝐽

donde
• 𝐵𝐽 es la matriz de iteración de Jacobi
• 𝑐𝐽 es el vector constante
Descomponemos la matriz 𝐴 en 𝐴 = 𝐷 + 𝐴𝐿 + 𝐴𝑅
4 0 00 0 0 00 0 −1 1 0
⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎜0 −3 0 0⎟ ⎜−1 0 0 0 ⎟ ⎜0 0 −1 −1⎟⎟
𝐷=⎜
⎜ ⎟
⎟ , −𝐴 𝐿 = ⎜
⎜ ⎟
⎟ , −𝐴 𝑅 = ⎜
⎜ ⎟
⎟
⎜
⎜0 0 5 0⎟
⎟ ⎜
⎜ 1 −1 0 0 ⎟
⎟ ⎜
⎜0 0 0 2 ⎟
⎟
⎝0 0 0 4⎠ ⎝ 0 −1 2 0⎠ ⎝0 0 0 0 ⎠

7 / 54
 8 / 54
Calculando
La matriz delaiteración
Matriz de
de Jacobi
Iteración
está𝐵dada
𝐽 por:

D. C D

𝐵𝐽 = −𝐷−1 (𝐴𝐿 + 𝐴𝑅 )

D. R P

Calculando 𝐷−1 :
1/4 0 0 0
⎛
⎜ ⎞
⎟
⎜
⎜ 0 −1/3 0 0 ⎟
⎟
𝐷−1 =⎜
⎜ ⎟
⎟
⎜
⎜ 0 0 1/5 0 ⎟ ⎟
⎜ ⎟
⎝ 0 0 0 1/4 ⎠
Calculando 𝐴𝐿 + 𝐴𝑅 :

0 −1 1 0
⎛
⎜ ⎞
⎟
⎜−1 0 −1 −1 ⎟
−(𝐴𝐿 + 𝐴𝑅 ) = ⎜
⎜ ⎟
⎟
⎜
⎜ 1 −1 0 2 ⎟⎟
⎝ 0 −1 2 0 ⎠

8 / 54
 9 / 54
La Matriz de Iteración 𝐵𝐽

D. C D

Finalmente, la matriz de iteración de Jacobi para este sistema es:

D. R P

1/4 0 0 0 0 −1 1 0
⎛
⎜ ⎞
⎟ ⎛ ⎞
⎜
⎜ 0 −1/3 0 0 ⎟
⎟ ⎜
⎜−1 0 −1 −1 ⎟
⎟
𝐵𝐽 = −𝐷−1 (𝐴𝐿 + 𝐴𝑅 ) = ⎜
⎜ ⎟
⎟ ⎜
⎜ ⎟
⎟
⎜
⎜ 0 0 1/5 0 ⎟
⎟ ⎜
⎜ 1 −1 0 2 ⎟
⎟
⎜ ⎟
⎝ 0 0 0 1/4⎠ ⎝ 0 −1 2 0 ⎠

0 −1/4 1/4 0
⎛
⎜ ⎞
⎟
⎜
⎜1/3 0 1/3 1/3 ⎟
⎟
=⎜
⎜ ⎟
⎟
⎜
⎜1/5 −1/5 0 2/5⎟ ⎟
⎜ ⎟
⎝ 0 −1/4 1/2 0 ⎠

9 / 54
 10 / 54

D. C D

D. R P
El vector constante 𝑐𝐽 está dado por:

1/4 0 0 0 7 7/4
⎛
⎜ ⎞
⎟ ⎛ ⎞ ⎛ ⎞
⎜ ⎜ ⎟
⎜ 0 −1/3 0 0 ⎟ ⎟ ⎜
⎜−6 ⎟
⎟ ⎜
⎜ 2 ⎟
⎟
−1 ⎜
𝑐𝐽 = 𝐷 𝑏 = ⎜ ⎟
⎟ ⎜
⎜ ⎟
⎟ =⎜ ⎟
⎜
⎜ 0 0 1/5 0 ⎟ ⎜ 8 ⎟ ⎜8/5⎟
⎟ ⎜ ⎟ ⎜ ⎟
⎜ ⎟ ⎜ ⎟
0 0 0 1/4 ⎝ 9⎠ 9/4
⎝ ⎠ ⎝ ⎠

10
/ 54
 11 / 54
Paso 1 de la Iteración 0
⎛
⎜ ⎞
0⎟

D. C D

⎜ ⎟
Sea el vector inicial 𝑥(0) =⎜ ⎟
⎜ ⎟. Entonces, la primera iteración es:

D. R P

⎜
⎜0⎟⎟
⎝0⎠
𝑥(1) = 𝐵𝐽 𝑥(0) + 𝑐𝐽
0 −1/4 1/4 0 0 7/4
⎛
⎜ ⎞
⎟ ⎛ ⎞ ⎛ ⎞
⎜ ⎜ ⎟
⎜ 1/3 0 1/3 1/3⎟ ⎟ ⎜
⎜ 0⎟
⎟ ⎜
⎜ 2 ⎟
⎟
=⎜
⎜ ⎟
⎟ ⎜
⎜ ⎟
⎟ +⎜
⎜ ⎟
⎟
⎜
⎜ ⎟ ⎜ 0⎟
1/5 −1/5 0 2/5⎟ ⎜ ⎟ ⎜ 8/5 ⎟
⎜ ⎟ ⎜ ⎟
0 −1/4 1/2 0 ⎝ 0⎠ 9/4
⎝ ⎠ ⎝ ⎠
0 1.75 1.75
⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎜
⎜ 0⎟
⎟ ⎜
⎜ 2.00⎟⎟ ⎜
⎜ 2.00⎟⎟
⎟
=⎜ ⎟+⎜ ⎟ =⎜ ⎟
⎜ 0⎟
⎜ ⎟ ⎜⎜ 1.60 ⎟
⎟ ⎜⎜ 1.60 ⎟
⎟
0
⎝ ⎠ ⎝ 2.25 ⎠ ⎝ 2.25 ⎠
11
/ 54
 12 / 54
Paso 2 de la Iteración
1.75

D. C D

⎛
⎜ ⎞
⎜ 2.00⎟⎟

D. R P

Usando 𝑥(1) ⎜
=⎜ ⎟
⎟, la segunda iteración es:
⎜1.60⎟
⎜ ⎟
⎝2.25⎠

𝑥(2) = 𝐵𝐽 𝑥(1) + 𝑐𝐽
0 −1/4 1/4 0 1.75 1.75
⎛
⎜ ⎞
⎟ ⎛ ⎞ ⎛ ⎞
⎜
⎜ 1/3 0 1/3 1/3⎟
⎟ ⎜
⎜ 2.00⎟⎟ ⎜
⎜ 2.00⎟⎟
=⎜
⎜ ⎟
⎟ ⎜
⎜ ⎟
⎟ + ⎜
⎜ ⎟
⎟
⎜
⎜ 1/5 −1/5 0 2/5⎟
⎟ ⎜
⎜ 1.60 ⎟
⎟ ⎜
⎜ 1.60 ⎟
⎟
⎜ ⎟
0 −1/4 1/2 0 ⎝ 2.25 ⎠ ⎝ 2.25 ⎠
⎝ ⎠
0(1.75) − 0.25(2.00) + 0.25(1.60) + 0(2.25) 1.75
⎛
⎜ ⎞
⎟ ⎛ ⎞
⎜
⎜ 1/3(1.75) + 0(2.00) + 1/3(1.60) + 1/3(2.25)⎟ ⎟ ⎜
⎜ 2.00⎟⎟
=⎜
⎜ ⎟
⎟ + ⎜
⎜ ⎟
⎟
⎜
⎜ 1/5(1.75) − 1/5(2.00) + 0(1.60) + 2/5(2.25) ⎟
⎟ ⎜
⎜ 1.60 ⎟
⎟
⎜ ⎟
0(1.75) − 0.25(2.00) + 0.5(1.60) + 0(2.25) ⎝ 2.25⎠
⎝ ⎠

12
/ 54
 13 / 54
Paso 2 de la Iteración
1.75

D. C D

⎛
⎜ ⎞
⎟

D. R P

⎜ 2.00 ⎟
Usando 𝑥(1) =⎜
⎜ ⎟
⎟, la segunda iteración es:
⎜
⎜ 1.60 ⎟
⎟
⎝2.25⎠

𝑥(2) = 𝐵𝐽 𝑥(1) + 𝑐𝐽
0 − 0.5 + 0.4 + 0 1.75
⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎜
⎜ 0.5833 + 0 + 0.5333 + 0.75⎟
⎟ ⎜
⎜ 2.00⎟⎟
⎟
=⎜ ⎟ +⎜ ⎟
⎜
⎜ 0.35 − 0.4 + 0 + 0.9 ⎟ ⎜1.60⎟
⎟ ⎜ ⎟
⎝ 0 − 0.5 + 0.8 + 0 ⎠ ⎝ 2.25 ⎠
−0.1 1.75 1.65
⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎜
⎜ 1.8666⎟
⎟ ⎜
⎜ 2.00⎟
⎟ ⎜
⎜ 3.8666⎟
⎟
⎟
=⎜ ⎟ +⎜ ⎟ =⎜ ⎟
⎜
⎜ 0.85 ⎟ ⎜1.60⎟ ⎜ 2.45 ⎟
⎟ ⎜ ⎟ ⎜ ⎟
⎝ 0.3 ⎠ ⎝2.25⎠ ⎝ 2.55 ⎠
13
/ 54
 14 / 54
Paso 3 de la Iteración
1.65

D. C D

⎛
⎜ ⎞
⎟
⎜ 3.8666⎟

D. R P

Usando 𝑥(2) ⎜
=⎜ ⎟
⎟, la tercera iteración es:
⎜
⎜ 2.45 ⎟⎟
⎝ 2.55 ⎠

𝑥(3) = 𝐵𝐽 𝑥(2) + 𝑐𝐽
0 −1/4 1/4 0 1.65 1.75
⎛
⎜ ⎞
⎟ ⎛ ⎞ ⎛ ⎞
⎜
⎜1/3 0 1/3 1/3 ⎟
⎟ ⎜
⎜ ⎟
3.8666⎟ ⎜
⎜ 2.00⎟⎟
=⎜
⎜ ⎟
⎟ ⎜
⎜ ⎟
⎟ + ⎜
⎜ ⎟
⎟
⎜
⎜1/5 −1/5 0 2/5 ⎟
⎟ ⎜
⎜ 2.45 ⎟
⎟ ⎜
⎜ 1.60 ⎟
⎟
⎜ ⎟
0 −1/4 1/2 0 ⎝ 2.55 ⎠ ⎝ 2.25 ⎠
⎝ ⎠
0(1.65) − 0.25(3.8666) + 0.25(2.45) + 0(2.55) 1.75
⎛
⎜ ⎞
⎟ ⎛ ⎞
⎜
⎜1/3(1.65) + 0(3.8666) + 1/3(2.45) + 1/3(2.55)⎟ ⎟ ⎜
⎜ 2.00⎟⎟
=⎜
⎜ ⎟
⎟ + ⎜
⎜ ⎟
⎟
⎜
⎜1/5(1.65) − 1/5(3.8666) + 0(2.45) + 2/5(2.55) ⎟
⎟ ⎜
⎜ 1.60 ⎟
⎟
⎜ ⎟
0(1.65) − 0.25(3.8666) + 0.5(2.45) + 0(2.55) ⎝ 2.25⎠
⎝ ⎠

14
/ 54
 15 / 54
Paso 3 de la Iteración
1.65

D. C D

⎛
⎜ ⎞
⎟

D. R P

⎜ 3.8666⎟
Usando 𝑥(2) =⎜
⎜ ⎟
⎟ , la tercera iteración es:
⎜
⎜ 2.45 ⎟
⎟
⎝ 2.55 ⎠

𝑥(3) = 𝐵𝐽 𝑥(2) + 𝑐𝐽
0 − 0.96665 + 0.6125 + 0 1.75
⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎜
⎜ 0.55 + 0 + 0.81665 + 0.85⎟
⎟ ⎜
⎜ 2.00⎟⎟
⎟
=⎜ ⎟ +⎜ ⎟
⎜0.33 − 0.77332 + 0 + 1.02⎟ ⎜1.60⎟
⎜ ⎟ ⎜ ⎟
⎝ 0 − 0.96665 + 1.225 + 0 ⎠ ⎝ 2.25 ⎠
−0.35415 1.75 1.39585
⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎜
⎜ 2.21665 ⎟
⎟ ⎜
⎜ 2.00⎟⎟ ⎜
⎜ 4.21665⎟⎟
⎟
=⎜ ⎟ +⎜ ⎟ =⎜ ⎟
⎜
⎜ 0.57668 ⎟ ⎜1.60⎟ ⎜2.17668⎟
⎟ ⎜ ⎟ ⎜ ⎟
⎝ 0.25835 ⎠ ⎝2.25⎠ ⎝2.50835⎠
15
/ 54
 16 / 54

1
0 𝑎12 ⋯ 𝑎1𝑛
⎛ ⎞

D. C D

⎛ 𝑎11 ⎞ ⎜
⎜ ⎟
⎟
⎜ ⎟ ⎜ 𝑎 21 0 … ⋮ ⎟
𝐷−1 (𝐴𝑅 + 𝐴𝐿 ) = ⎜ ⎟

D. R P

⎜ ⋱ ⎟ ⎜
⎜ ⎟
⎟
⎜
⎜ ⎟
⎟ ⎜ ⋯ ⋱ 𝑎 𝑛−1,𝑛 ⎟
1 ⎜
⎜ ⎟
⎟
⎝ 𝑎𝑛𝑛 ⎠
𝑎𝑛1 𝑎𝑛,𝑛−1 0
⎝ ⎠

𝑎1𝑛
0 𝑎𝑎12 ⋯
⎛
⎜
11 𝑎11 ⎞
⎟
⎜
⎜ 𝑎21 ⎟
⎟
⎜ 𝑎22 0 ⋱ ⋮ ⎟
⎜
= ⎜ ⎟
⎜ 𝑎𝑛−1,𝑛 ⎟⎟
⎜ ⋯ ⋱ 𝑎𝑛−1,𝑛−1 ⎟
⎜
⎜ ⎟
⎟
⎜ 𝑎𝑛1 𝑎𝑛,𝑛−1 ⎟
0
⎝ 𝑎𝑛𝑛 𝑎𝑛𝑛 ⎠

x𝑛+1 ∶= −𝐷−1 (𝐴𝑅 + 𝐴𝐿 )x𝑛 + 𝐷−1 y

converge si y solo si
16
/ 54 𝜌(−𝐷−1 (𝐴𝑅 + 𝐴𝐿 )) < 1
 17 / 54
Teorema: Condiciones suficientes método de Jacobi
Sea 𝐴 = (𝑎𝑗𝑘 ) que satisface

D. C D

𝑛 𝑛

D. R P

𝑎𝑗𝑘 𝑎𝑗𝑘
𝑞∞ ∶= max ∑ ∣ ∣ < 1, ó, 𝑞1 ∶= max ∑ ∣ ∣<1 ó,
𝑗=1,⋯,𝑛
𝑘=1,𝑘≠𝑗
𝑎𝑗𝑗 𝑘=1,⋯,𝑛
𝑗=1,𝑗≠𝑘
𝑎𝑗𝑗
2 1/2
𝑛
𝑎𝑗𝑘
𝑞𝐹 ∶= ( ∑ ∣ ∣ ) <1
𝑗,𝑘=1,𝑗≠𝑘
𝑎𝑗𝑗

entonces el método iterativo de Jacobi,

𝑛
𝑎𝑗𝑘 𝑦𝑗
𝑥𝑛+1,𝑗 ∶= − ∑ 𝑥𝑛,𝑘 +
𝑎
𝑘=1,𝑘≠𝑗 𝑗𝑗
𝑎𝑗𝑗

converge para cada y ∈ ℂ𝑛 y para cada x0 ∈ ℂ𝑛 a la solución única de 𝐴x = y (en

cualquier norma de ℂ𝑛 ). Para 𝜇 = 1, 2 ó ∞ si 𝑞𝜇 < 1 entonces

𝑞𝜇𝑛 𝑞𝜇
‖x𝑛 − x‖𝜇 ≤ ‖x − x0 ‖ , ‖x𝑛 − x‖𝜇 ≤ ‖x − x𝑛−1 ‖
1 − 𝑞𝜇 1 1 − 𝑞𝜇 𝑛

17
/ 54
 18 / 54

Demostración:

D. C D

En el método de Jacobi, la matriz de iteración está dada por

D. R P

𝐵 = −𝐷−1 (𝐴𝐿 + 𝐴𝑅 ),

así,

𝑞1 = ‖−𝐷(𝐴𝐿 + 𝐴𝑅 )‖1 = ‖𝐵‖1

𝑞∞ = ‖−𝐷(𝐴𝐿 + 𝐴𝑅 )‖∞ = ‖𝐵‖2
𝑞2 = ‖−𝐷(𝐴𝐿 + 𝐴𝑅 )‖𝐹 = ‖𝐵‖𝐹

por tanto la sucesión

x𝑛+1 = 𝐵x𝑛 + z ∶= −𝐷−1 (𝐴𝐿 + 𝐴𝑅 )x𝑛 + 𝐷−1 y

converge al punto fijo x, el cual es solución del sistema 𝐴x = y.

18
/ 54
 20 / 54
Ejemplo del Sistema 4x4 (Recordatorio)
Consideramos el sistema con matriz de iteración:

D. C D

D. R P
0 −1/4 1/4 0 7/4
⎛
⎜ ⎞
⎟ ⎛ ⎞
⎜ ⎟ ⎜
⎜ ⎟
⎟
⎜ 1/3 0 1/3 1/3⎟ 2
𝐵𝐽 = ⎜
⎜ ⎟
⎟ , 𝑐𝐽 = ⎜
⎜
⎜
⎟
⎟
⎟
⎜
⎜ 1/5 −1/5 0 2/5⎟ ⎟ ⎜ 8/5 ⎟
⎜ ⎟ ⎜ ⎟
⎝ 0 −1/4 1/2 0 ⎠ ⎝9/4⎠
8 7 13 11 13
‖𝐵𝐽 ‖1 = max { , , , }= ≈ 1.08333
15 10 12 15 12
1 4 3
‖𝐵𝐽 ‖∞ = max { , 1, , } = 1
2 5 4
1 1 1 1 1 1 1 4 1 1
‖𝐵𝐽 ‖∞ = √0 + + +0+ +0+ + + + +0+ +0+ + +0
16 16 9 9 9 25 25 25 16 4
≈ 1.0054
‖𝐵𝐽 ‖2 = 0.653837
𝜌(𝐵𝐽 ) = 0.488639
20
/ 54
 4𝑥 − 𝑦 + 𝑧 = 7
Resolver 4𝑥 − 8𝑦 + 𝑧 = −21

D. C D

D. R P
−2𝑥 + 𝑦 + 5𝑧 = 15 1
0 −1 1 ⎛
⎜
4 ⎞
⎟
⎛
⎜ ⎞
⎟, ⎜ ⎟
𝐴𝐿 + 𝐴 𝑅 = ⎜
⎜ 4 0 1⎟ ⎟ 𝐷 −1
=⎜
⎜ − 18 ⎟
⎟ ⇒
⎜
⎜ ⎟
⎟
⎝−2 1 0 ⎠ 1
⎝ 5⎠

1 1
⎛ 0 4 −4⎞
⎜
⎜ ⎟
⎟
⎜
⎜ ⎟
⎟
⎜ 1 ⎟
𝐵 = −𝐷(𝐴𝐿 + 𝐴𝑅 ) = ⎜
⎜
1
2 0 ⎟,
8 ⎟
⎜
⎜ ⎟
⎟
⎜
⎜ ⎟
⎟
⎜ ⎟
2 1
⎝5 5− 0 ⎠
1 5 3 5
𝑞∞ = ‖𝐵‖∞ = max { , , } =
2 8 5 8
9 9 3 9 ⇒ el método iterativo de Jacobi
𝑞1 = ‖𝐵‖1 = max { , , } =
10 20 8 10 converge.
𝑞 = ‖𝐵‖ ≈ 0.7685
 22 / 54
Método iterativo de Jacobi

D. C D

7+𝑦−𝑧

D. R P

𝑥=
4𝑥 − 𝑦 + 𝑧 = 7 4
21 + 4𝑥 + 𝑧
4𝑥 − 8𝑦 + 𝑧 = −21 ⇒ 𝑦=
8
−2𝑥 + 𝑦 + 5𝑧 = 15 15 + 2𝑥 − 𝑦
𝑧=
5
lo que sugiere el siguiente proceso iterativo

7 + 𝑦 𝑛 − 𝑧𝑛
𝑥𝑛+1 =
4
21 + 4𝑥𝑛 + 𝑧𝑛
𝑦𝑛+1 =
8
15 + 2𝑥𝑛 − 𝑦𝑛
𝑧𝑛+1 =
5

22
/ 54
 23 / 54
Condiciones suficientes - Método de Jacobi

D. C D

𝑛
𝑎𝑗𝑘
𝑞∞ ∶= max ∑ ∣ ∣ < 1, si y solo si

D. R P

𝑗=1,⋯,𝑛
𝑘=1,𝑘≠𝑗
𝑎𝑗𝑗
𝑛
∑ |𝑎𝑗𝑘 | < |𝑎𝑗𝑗 | 𝑗 = 1, … , 𝑛, es decir, si 𝐴 es estrictamente diagonal
𝑘=1, 𝑘≠𝑗
dominante por filas.

Ejemplo
4 −1 1
⎛
⎜ ⎞
⎟
𝐴=⎜
⎜ 4 −8 1⎟
⎟
⎝−2 1 5⎠
es estrictamente diagonalmente dominante por filas, pues

|1| + |1| < 4 |4| + |1| < | − 8|, | − 2| + |1| < 5

23
/ 54
 24 / 54

import numpy as np

D. C D

D. R P
def jacobi(A,b,x,tol,imax):
n = len(A[0])
D = np.diag(A) # vector diagonal
R = A - np.diag(D) # matriz R = A_R + A_L
r = tol + 1
k = 0
while (r > tol) and (k <=imax):
x =(b-np.matmul(R,x))/D.reshape((n,1))
r = np.linalg.norm(b-np.matmul(A,x),'fro')
k = k +1

if k == imax+1:
print('numero maximo de iteraciones alcanzado, residuo=',r)
return [x, k, r]
24
/ 54
import numpy as np 25 / 54
def jacobiBasico(A,b,x,tol,imax):
n = len(A[0])

D. C D

r = tol + 1

D. R P

k = 0
while (r > tol) and (k <=imax):
for j in range(0,n):
xsum = 0.0
for k in range(0,n):
if k!=j:
xsum = xsum + A[j][k]*x[k]/A[j][j]

x[j] = -xsum + b[j]/A[j][j]

r = np.linalg.norm(b-np.matmul(A,x),'fro')
k = k +1

if k == imax+1:
print('numero maximo de iteraciones alcanzado, residuo=',
return
25 [x, k, r]
/ 54
 26 / 54
Métodos iterativos

D. C D

𝐴 = 𝐷 + 𝐴𝑅 + 𝐴 𝐿
0

D. R P

⎛
⎜ ⎞
⎟ (𝐷 + 𝐴𝑅 + 𝐴𝐿 )x = y
⎜ 𝑎21 0 ⎟
𝐴𝐿 ⎜
∶= ⎜ ⎟
⎟
⎜ ⎟ (𝐷 + 𝐴𝐿 )x = −𝐴𝑅 x + y
⎜ ⋮ ⋱ ⋱ ⎟
⎜ ⎟
𝑎𝑛1 ⋯ 𝑎𝑛,𝑛−1 0 x = −(𝐷 + 𝐴𝐿 )−1 𝐴𝑅 x + (𝐷 + 𝐴𝐿 )−1 y
⎝ ⎠

𝑎11
⎛ ⎞
𝐷 ∶= ⎜
⎜
⎜ ⋱
⎟
⎟
⎟ x𝑛+1 = − (𝐷 + 𝐴𝐿 )−1 𝐴𝑅 x𝑛
⎝ 𝑎𝑛𝑛 ⎠ + (𝐷 + 𝐴𝐿 )−1 y

0 𝑎12 ⋯ 𝑎1𝑛 Gauss-Seidel

⎛
⎜ ⎞
⎟
⎜ 0 ⋱ ⋮ ⎟
𝐴𝑅 ∶= ⎜
⎜
⎜
⎟
⎟
⎜ ⋱ 𝑎𝑛−1,𝑛 ⎟
⎟
⎜ ⎟ (𝐷 + 𝐴𝐿 )x𝑛+1 = −𝐴𝑅 x𝑛 + y
⎝ 0 ⎠

26
/ 54
 27 / 54
Forma matricial del método

D. C D

D. R P
El método de Gauss-Seidel se escribe como:

(𝐷 + 𝐴𝐿 )x(𝑘+1) = −𝐴𝑅 x(𝑘) + b

x(𝑘+1) = −(𝐷 + 𝐴𝐿 )−1 𝐴𝑅 x(𝑘) + (𝐷 + 𝐴𝐿 )−1 b

Sin embargo, no calculamos la inversa de 𝐷 + 𝐴𝐿

27
/ 54
 28 / 54
Ecuación por componentes

D. C D

(𝐷 + 𝐴𝐿 )𝑖 x(𝑘+1) = −(𝐴𝑅 )𝑖 x(𝑘) + b𝑖

D. R P

Aplicado en iteraciones, tenemos:
𝑛
(𝑘+1) (𝑘+1)
• 𝐷𝑖 x(𝑘+1) = ∑ 𝑑𝑖𝑗 𝑥𝑗 = 𝑎𝑖𝑖 𝑥𝑖
𝑗=1
𝑖−1
(𝑘+1)
• (𝐴𝐿 )𝑖 x(𝑘+1) = ∑ 𝑎𝑖𝑗 𝑥𝑗
𝑗=1
𝑛
(𝑘) (𝑘)
• (𝐴𝑅 )𝑖 x = ∑ 𝑎𝑖𝑗 𝑥𝑗
𝑗=𝑖+1
Entonces:
𝑖−1 𝑛
(𝑘+1) (𝑘+1) (𝑘)
𝑎𝑖𝑖 𝑥𝑖 + ∑ 𝑎𝑖𝑗 𝑥𝑗 + ∑ 𝑎𝑖𝑗 𝑥𝑗 = 𝑏𝑖
𝑗=1 𝑗=𝑖+1

𝑖−1 𝑛
(𝑘+1) 1 (𝑘+1) (𝑘)
𝑥𝑖 = (𝑏𝑖 − ∑ 𝑎𝑖𝑗 𝑥𝑗 − ∑ 𝑎𝑖𝑗 𝑥𝑗 )
𝑎𝑖𝑖 𝑗=1 𝑗=𝑖+1
28
/ 54
 29 / 54
Método Iterativo de Gauss-Seidel

D. C D

D. R P
La matriz de iteración de Gauss-Seidel está dada por:

𝐵 ∶= −(𝐷 + 𝐴𝐿 )−1 𝐴𝑅

por lo que el método toma la forma:

𝑗−1 𝑛
(𝜈+1) 𝑎𝑗𝑘 (𝜈+1) 𝑎𝑗𝑘 (𝜈) 𝑦𝑗
𝑥𝑗 = −∑ 𝑥𝑘 − ∑ 𝑥𝑘 + , 𝑗 = 1, … , 𝑛, 𝜈 ∈ ℕ0
𝑘=1
𝑎𝑗𝑗 𝑎
𝑘=𝑗+1 𝑗𝑗
𝑎𝑗𝑗

𝑗−1 𝑛
1 (𝜈+1) (𝜈)
= ( − ∑ 𝑎𝑗𝑘 𝑥𝑘 − ∑ 𝑎𝑗𝑘 𝑥𝑘 + 𝑦𝑗 )
𝑎𝑗𝑗 𝑘=1 𝑘=𝑗+1

29
/ 54
 30 / 54
Ejemplo 4x4

D. C D

D. R P
Consideremos nuevamente el sistema de ecuaciones lineales 𝐴𝑥 = 𝑏:

4 1 −1 0 ⎛𝑥1 ⎞ 7
⎛
⎜ ⎞
⎟ ⎜ ⎟ ⎛
⎜ ⎞
⎜
⎜ 1 −3 1 1 ⎟
⎟
⎜
⎜𝑥2 ⎟
⎟ ⎜ −6⎟
⎟
⎜ ⎟ ⎜ ⎟ =⎜ ⎟
⎜
⎜
⎜−1 1 5 −2⎟⎟⎜
⎜𝑥 ⎟ ⎜8⎟ ⎟
⎜ 3⎟⎟ ⎜ ⎟
⎝ 0 1 −2 4 ⎠ ⎝𝑥4 ⎠ ⎝ 9 ⎠

30
/ 54
 31 / 54
El Método Iterativo de Gauss-Seidel

D. C D

El método de Gauss-Seidel es una mejora del método de Jacobi. En cada

D. R P

iteración, utiliza los valores ya calculados en la misma iteración para acelerar
la convergencia:

(𝑘+1) 1 (𝑘) (𝑘) (𝑘)

𝑥1 = (𝑏1 − 𝑎12 𝑥2 − 𝑎13 𝑥3 − 𝑎14 𝑥4 )
𝑎11
(𝑘+1) 1 (𝑘+1) (𝑘) (𝑘)
𝑥2 = (𝑏2 − 𝑎21 𝑥1 − 𝑎23 𝑥3 − 𝑎24 𝑥4 )
𝑎22
(𝑘+1) 1 (𝑘+1) (𝑘+1) (𝑘)
𝑥3 = (𝑏3 − 𝑎31 𝑥1 − 𝑎32 𝑥2 − 𝑎34 𝑥4 )
𝑎33
(𝑘+1) 1 (𝑘+1) (𝑘+1) (𝑘+1)
𝑥4 = (𝑏4 − 𝑎41 𝑥1 − 𝑎42 𝑥2 − 𝑎43 𝑥3 )
𝑎44

31
/ 54
 32 / 54
Reescribiendo las Ecuaciones para Gauss-Seidel

D. C D

D. R P
Para nuestro sistema:
(𝑘+1) 1 (𝑘) (𝑘)
𝑥1 = (7 − 𝑥2 + 𝑥3 )
4
(𝑘+1) 1 (𝑘+1) (𝑘) (𝑘)
𝑥2 = (−6 − 𝑥1 − 𝑥3 − 𝑥4 )
−3
(𝑘+1) 1 (𝑘+1) (𝑘+1) (𝑘)
𝑥3 = (8 + 𝑥1 − 𝑥2 + 2𝑥4 )
5
(𝑘+1) 1 (𝑘+1) (𝑘+1)
𝑥4 = (9 − 𝑥2 + 2𝑥3 )
4

32
/ 54
La Matriz de Iteración de Gauss-Seidel
Recordemos la descomposición 𝐴 = 𝐷 + 𝐴𝑅 + 𝐴𝐿 . La matriz de iteración de

D. C D

Gauss-Seidel, 𝐵𝐺𝑆 , y el vector constante 𝑐𝐺𝑆 están dados por:

D. R P

𝐵𝐺𝑆 = −(𝐷 + 𝐴𝐿 )−1 𝐴𝑅 , 𝑐𝐺𝑆 = (𝐷 + 𝐴𝐿 )−1 𝑏

4 0 00 0 −1 1 0 7
⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟
⎜ −1 −3 0 0⎟ ⎜0 0 −1 −1 ⎟ ⎜−6 ⎟
𝐷 + 𝐴𝐿 = ⎜
⎜ ⎟
⎟, −𝐴𝑅 = ⎜
⎜ ⎟
⎟, 𝑏=⎜
⎜ ⎟
⎟
⎜
⎜ 1 −1 5 0⎟⎟ ⎜
⎜0 0 0 2 ⎟ ⎟ ⎜
⎜8⎟ ⎟
⎝ 0 −1 2 4⎠ ⎝0 0 0 0 ⎠ ⎝9⎠

0 − 14 14 0 ⎞ ‖𝐵𝐺𝑆 ‖1 = 0.820833
⎛
⎜ ⎟
⎜ 1 5 1 ⎟ ‖𝐵𝐺𝑆 ‖∞ = 0.833333
⎜
⎜ 0 − 12 12
⎟
3 ⎟
𝐵𝐺𝑆 =⎜
⎜ ⎟
⎟ ‖𝐵𝐺𝑆 ‖𝐹 = 0.742673
⎜
⎜
1
0 − 30 − 301 1 ⎟
⎜ 3 ⎟
⎟ ‖𝐵𝐺𝑆 ‖2 = 0.623721
⎜ ⎟
1 29 1
⎝ 0 240 − 240 12 ⎠ 𝜌(𝐵𝐺𝑆 ) = 0.214087
 34 / 54
Primera Iteración (k=0) ( Sol. Exacta ( 79 , 27
7
𝑇
, 2, 16
7
) ≈ (1.28571 3.85714 2. 2.28571) )
𝑇

D. C D

𝑇
• Asumimos un vector inicial: x(0) = (0 0 0 0)

D. R P

▶ Calculamos 𝐵𝐺𝑆 = −(𝐷 + 𝐴𝐿 )−1 𝐴𝑅 y
▶ Calculamos (𝐷 + 𝐴𝐿 )−1 :
𝑐𝐺𝑆 = (𝐷 + 𝐴𝐿 )−1 𝑏:
1
0 0 0⎞
⎛
⎜
4
⎟ ⎛ 0 − 14 14 0 ⎞ ⎛
7
4 ⎞
⎜
⎜ 1 1 ⎟
⎟ ⎜ ⎟ ⎜ ⎟
⎜ 12 − 3 0 0 ⎟ ⎜
⎜ 0 − 1 5 1 ⎟
⎟ ⎜
⎜ 31 ⎟
⎟
(𝐷+𝐴𝐿 )−1 ⎜
=⎜ ⎟
⎟ ⎜ 12 12 3 ⎟ ⎜ 12 ⎟
⎜ 𝐵𝐺𝑆 ⎜
=⎜ ⎟ ⎜
⎟ , 𝑐𝐺𝑆 = ⎜ ⎟
⎜
1
30
1 1
15 5 0⎟⎟ ⎜ 1 1 1 ⎟ ⎜ ⎟
43 ⎟
⎜
⎜ ⎟
⎟ ⎜ 0 − 30 − 30 3 ⎟ ⎜ 30 ⎟
1 7 1 1
⎜
⎜ ⎟
⎟ ⎜
⎜ ⎟ ⎟
− 240 1 29 1 557
⎝ 60 10 4 ⎠ 0 −
⎝ 240 240 12 ⎠ ⎝ 240 ⎠

• Calculamos
7 7
0 ⎛ 4 ⎞
⎜ ⎟ ⎛⎜
4 ⎞
⎟ ⎛ 1.75 ⎞
⎛
⎜ ⎞
⎟ ⎜
⎜ 31 ⎟
⎟ ⎜
⎜ 31 ⎟
⎟ ⎜
⎜ 0 ⎟ ⎜ 12 ⎟ ⎜ 12 ⎟ ⎜2.58333⎟⎟
x(1) = 𝐵𝐺𝑆 𝑥(0) + 𝑐𝐺𝑆 = 𝐵𝐺𝑆 ⋅ ⎜ ⎜ ⎟
⎟ +⎜⎜ ⎟
⎟ =⎜
⎜ ⎟
⎟ ≈⎜
⎜ ⎟
⎟
⎜ 0 ⎟
⎜ ⎟ ⎜ ⎜ 43 ⎟
⎟ ⎜
⎜
43 ⎟
⎟ ⎜
⎜ 1.43333⎟
⎟
⎜ 30 ⎟
⎜ ⎟ ⎜ ⎟ ⎜ 30 ⎟
0
⎝ ⎠ 557 557 ⎝ 2.32083⎠
⎝ 240 ⎠ ⎝ 240 ⎠
𝑇
• Calcular el residuo 𝑟 = 𝐴x(1) − 𝑏 = (1.15 3.75417 −4.64167 0)
• Calcular la norma del residuo ‖𝑟‖2 ≈ 6.07958360461928
34
/ 54
 35 / 54
Segunda Iteración (k=1) ( Sol. Exacta ( 9 27
, , 2, 16
7 7 7
)
𝑇
≈ (1.28571 3.85714 2. 2.28571) )
𝑇

D. C D

𝑇

D. R P

• Usamos x(1) = ( 74 31
12
43
30
557
240 )
(2) (1)
• Calculamos x = 𝐵𝐺𝑆 𝑥 + 𝑐𝐺𝑆 :

⎛ 0 − 14 14 0 ⎞ ⎛ 74 ⎞ ⎛ 74 ⎞
⎜ ⎟⎜ ⎟ ⎜ ⎟
⎜
⎜
⎜
1
0 − 12 5
12
⎟⎜
1 ⎟ ⎜ 31 ⎟
3 ⎟ ⎜ 12 ⎟
⎟ ⎜⎜
⎜
31 ⎟
⎟
12 ⎟
x(2) =⎜
⎜ ⎟
⎟ ⎜
⎜ ⎟
⎟ +⎜
⎜ ⎟
⎟
⎜
⎜ 0 − 30 − 30 3 ⎟
1 1 1
⎟ ⎜
⎜
43 ⎟
⎟ ⎜
⎜
43 ⎟
⎟
⎜
⎜ ⎟
⎟⎜⎜ ⎟
30
⎟ ⎜⎜ ⎟
30
⎟
1 29 1 557 557
0 −
⎝ 240 240 12 ⎠ ⎝ 240 ⎠ ⎝ 240 ⎠
117
⎛
⎜
80 ⎞
⎟ 1.4625
⎜
⎜ 673 ⎟
⎟ ⎛
⎜ ⎞
⎟
⎜ 180 ⎟ ⎜3.7389 ⎟
=⎜
⎜ ⎟
⎟ ⎜
≈⎜ ⎟
⎟
⎜
⎜
7463 ⎟ ⎜
⎜2.0731 ⎟
⎟
⎜ 3600 ⎟
⎟
⎜ ⎟ 2.3518
16933 ⎝ ⎠
⎝ 7200 ⎠

• Calcular el residuo
𝑇
𝑟 = 𝐴x(2) − 𝑏 = (0.515833 0.670694 −0.0619444 0)
• Calcular la norma del residuo ‖𝑟‖2 ≈ 0.848382095393314
35
/ 54
 36 / 54
Tercera Iteración (k=2) ( Sol. Exacta ( 79 , 27
7
𝑇
, 2, 16
7
) ≈ (1.28571 3.85714 2. 2.28571) )
𝑇

D. C D

D. R P
𝑇
• Usamos x(2) = (1.4625 3.7389 2.0731 2.3518)
• Calculamos x(3) = 𝐵𝐺𝑆 𝑥(2) + 𝑐𝐺𝑆 :

⎛ 0 − 41 14 0 ⎞ ⎛
7
4 ⎞
⎜ ⎟ 1.4625 ⎜ ⎟
⎜
⎜ 1 5 1 ⎟ ⎛ ⎞ ⎜
⎜ 31 ⎟
⎟
⎜0 − 12 12 3 ⎟
⎜ ⎟
⎟
⎜
⎜3.7389⎟ ⎟ ⎜ 12 ⎟
x(3) =⎜ ⎟ ⎜
⎜ ⎟
⎟ ⎜
+⎜ ⎟ ⎟
⎜ 1 1 1 ⎟ ⎜2.0731⎟ ⎜ 43 ⎟
⎜ 0 − 30 − 30 3 ⎟ ⎜ ⎟ ⎜ 30 ⎟
⎜
⎜ ⎟
⎟ 2.3518 ⎜
⎜ ⎟ ⎟
0 1
− 29 1 ⎝ ⎠ 557
⎝ 240 240 12 ⎠ ⎝ 240 ⎠
117
⎛
⎜
80 ⎞
⎟ 1.3335
⎜
⎜ 673 ⎟
⎟ ⎛
⎜ ⎞
⎜ 180 ⎟ ⎜3.9195⎟⎟
=⎜
⎜ ⎟ ⎜
⎟ ≈⎜ ⎟
⎟
⎜
⎜
7463 ⎟ ⎜
⎜2.0235 ⎟
⎟
⎜ 3600 ⎟
⎟
⎜ ⎟
16933 ⎝2.2819⎠
⎝ 7200 ⎠
𝑇
• Calcular el residuo 𝑟 = 𝐴x(3) − 𝑏 = (0.230095 −0.1194 0.139801 0)
• Calcular la norma del residuo ‖𝑟‖2 ≈ 0.294524192810059
36
/ 54
 (𝑘) (𝑘) (𝑘) (𝑘)
𝑘 𝑥1 𝑥2 𝑥3 𝑥4 𝑟 = 𝐴x(𝑘) 37
−/𝑏54
0 0.00000000 0.00000000 0.00000000 0.00000000 15.16575089

D. C D

1 1.75000000 2.58333333 1.43333333 2.32083333 6.07958360

D. R P

2 1.46250000 3.73888889 2.07305556 2.35180556 0.84838210
3 1.33354167 3.91946759 2.02353704 2.28190162 0.29453810
4 1.27601736 3.86048534 1.99586705 2.28281219 0.04122913
5 1.28384543 3.85417489 1.99905898 2.28598577 0.01308057
6 1.28622102 3.85708859 2.00022079 2.28583825 0.00204572
7 1.28578305 3.85728070 2.00003577 2.28569771 0.00057203
8 1.28568877 3.85714075 1.99998869 2.28570916 0.00010207
9 1.28571198 3.85713661 1.99999874 2.28571522 0.00002465
10 1.28571553 3.85714316 2.00000056 2.28571449 0.00000507
11 1.28571435 3.85714313 2.00000004 2.28571424 0.00000105
12 1.28571423 3.85714283 1.99999997 2.28571428 0.00000025
13 1.28571428 3.85714285 2.00000000 2.28571429 0.00000004
14 1.28571429 3.85714286 2.00000000 2.28571429 0.00000001
1.28571429 3.85714286 2.00000000 2.28571429
37
/ 54
function X = GaussSeidel(A,y,P,delta, max1) 38 / 54

N = length(y);

D. C D

D. R P
for k=1:max1
for j=1:N
if j==1
X(1) = (y(1)-A(1,2:N)*P(2:N))/A(1,1);
elseif j == N
X(N)=(B(N) - A(N,1:N-1)*(X(1:N-1))')/A(N,N);
else
% X: contiene las aproximaciones k-esimas y P las (k-1)-esimas
X(j)=(y(j)-A(j,1:j-1)*X(1:j-1)'-A(j,j+1:N)*P(j+1:N))/A(j,j);
end
end
err=abs(norm(X'-P));
relerr=err/(norm(X)+eps);
P=X';
if (err<delta)|(relerr<delta)
break
end
end

X=X';
38
/ 54
En python ... (𝐷 + 𝐴𝐿 )x𝑛+1 = −𝐴𝑅 x𝑛 + 𝑏

D. C D

D. R P
import numpy as np

def GaussSeidel(A, b, x0, tol, nmax):

n = len(A[0])
x = np.zeros((n, 1), 'float')

R = np.tril(A)
U = np.triu(A, 1)
r = tol + 1
i = 0
while (r > tol) and (i < nmax):
y = -np.matmul(U, x0) + b
x = np.linalg.solve(R, y)
r = np.linalg.norm(b - np.matmul(A, x))
x0=x
i = i + 1
if i == nmax:
print('El esquema alcanzo el numero max. de iteraciones, r=', r)
break
return [x, r, i]
 3 / 17
Resolución sistemas lineales 𝐴x = y

D. C D

Si 𝐴 = 𝑃 − 𝑁 , x = 𝑃 −1 𝑁 x + 𝑃 −1 y

D. R P

Esquema iterativo 1ra forma: x𝑛+1 = 𝑃 −1 𝑁 x𝑛 + 𝑃 −1 y,
Esquema iterativo 2da forma: x𝑛+1 = x𝑛 + 𝑃 −1 𝑟𝑛 , 𝑟𝑛 = y − 𝐴x𝑛

Modificación método de Jacobi

La nueva aproximación x𝑛+1 es obtenida al corregir la anterior:

x𝑛+1 = x𝑛 + 𝐷−1 𝑟𝑛 , 𝑟𝑛 = y − 𝐴x𝑛

Si se multiplica la corrección por un factor 𝜔, la velocidad de conver-

gencia podría mejorar:

x𝑛+1 = x𝑛 + 𝜔𝐷−1 𝑟𝑛 , 𝑟𝑛 = y − 𝐴x𝑛

3 / 17
 4 / 17

JOR (Jacobi over-relaxation)

D. C D

D. R P
x𝑛+1 = x𝑛 + 𝜔𝐷−1 𝑟𝑛 , 𝑟𝑛 = y − 𝐴x𝑛

JOR: Matriz de iteración

𝐵𝜔 ∶ = 𝐼 + 𝜔𝐷−1 (−𝐴) = (1 − 𝜔)𝐼 − 𝜔𝐷−1 (𝐴𝑅 + 𝐴𝐿 )

= (1 − 𝜔)𝐼 + 𝜔𝐵
con 𝐵 ∶= −𝐷−1 (𝐴𝑅 + 𝐴𝐿 )

Si 𝐵 tiene autovalores reales, JOR es minimal para el parámetro de

relajación
2
𝜔opt =
2 − 𝜆max − 𝜆min

4 / 17
Modificación Gauss-Seidel

D. C D

(𝐷 + 𝐴𝐿 )x𝑛+1 = y − 𝐴𝑅 x𝑛 ,

D. R P

𝐷x𝑛+1 = y − 𝐴𝐿 x𝑛+1 − 𝐴𝑅 x𝑛 ,
x𝑛+1 = 𝐷−1 y − 𝐷−1 𝐴𝐿 x𝑛+1 − 𝐷−1 𝐴𝑅 x𝑛 + x𝑛 − x𝑛
= x𝑛 + 𝐷−1 (y − 𝐴𝐿 x𝑛+1 − (𝐷 + 𝐴𝑅 )x𝑛 )

x𝑛+1 = x𝑛 + 𝜔𝐷−1 (y − 𝐴𝐿 x𝑛+1 − (𝐷 + 𝐴𝑅 )x𝑛 )

SOR

x𝑛+1 = 𝜔𝐷−1 (y − 𝐴𝐿 x𝑛+1 − 𝐴𝑅 x𝑛 ) + (1 − 𝜔)x𝑛

𝑗−1 𝑛
(𝑛+1) 𝜔 (𝑛+1) (𝑛) (𝑛)
x𝑗 = (y𝑗 − ∑ 𝑎𝑗,𝑘 x𝑘 − ∑ 𝑎𝑗,𝑘 x𝑘 ) + (1 − 𝜔)x𝑗
𝑎𝑗,𝑗 𝑘=1 𝑘=𝑗+1
 6 / 17
x𝑛+1 = x𝑛 + 𝜔𝐷−1 (y − 𝐴𝐿 x𝑛+1 − (𝐷 + 𝐴𝑅 )x𝑛 )
x𝑛+1 + 𝜔𝐷−1 𝐴𝐿 x𝑛+1 = x𝑛 + 𝜔𝐷−1 y − 𝜔𝐷−1 (𝐷 + 𝐴𝑅 )x𝑛

D. C D

(𝐼 + 𝜔𝐷−1 𝐴𝐿 ) x𝑛+1 = x𝑛 + 𝜔𝐷−1 y − 𝜔𝐷−1 𝐷x𝑛 − 𝜔𝐷−1 𝐴𝑅 x𝑛

D. R P

(𝐼 + 𝜔𝐷−1 𝐴𝐿 ) x𝑛+1 = x𝑛 − 𝜔x𝑛 − 𝜔𝐷−1 𝐴𝑅 x𝑛 + 𝜔𝐷−1 y
(𝐼 + 𝜔𝐷−1 𝐴𝐿 ) x𝑛+1 = ((1 − 𝜔)𝐼 − 𝜔𝐷−1 𝐴𝑅 ) x𝑛 + 𝜔𝐷−1 y
−1
x𝑛+1 = 𝐵SOR x𝑛 + 𝜔 (𝐼 + 𝜔𝐷−1 𝐴𝐿 ) 𝐷−1 y
−1
x𝑛+1 = 𝐵SOR x𝑛 + 𝜔 (𝐷 + 𝜔𝐴𝐿 ) y

donde
−1
𝐵SOR = (𝐼 + 𝜔𝐷−1 𝐴𝐿 ) ((1 − 𝜔)𝐼 − 𝜔𝐷−1 𝐴𝑅 )
−1
= (𝐼 + 𝜔𝐷−1 𝐴𝐿 ) ((1 − 𝜔)𝐷−1 𝐷 − 𝜔𝐷−1 𝐴𝑅 )
−1
= (𝐼 + 𝜔𝐷−1 𝐴𝐿 ) 𝐷−1 ((1 − 𝜔)𝐷 − 𝜔𝐴𝑅 )
−1
= [𝐷 (𝐼 + 𝜔𝐷−1 𝐴𝐿 )] ((1 − 𝜔)𝐷 − 𝜔𝐴𝑅 )
−1
= (𝐷 + 𝜔𝐴𝐿 ) ((1 − 𝜔)𝐷 − 𝜔𝐴𝑅 )

−1 −1
x𝑛+1 = (𝐷 + 𝜔𝐴𝐿 ) ((1 − 𝜔)𝐷 − 𝜔𝐴𝑅 ) x𝑛 + 𝜔 (𝐷 + 𝜔𝐴𝐿 ) y
6 / 17
 4𝑥 − 𝑦 + 𝑧 = 7 1
⎛ ⎞
Resolver 4𝑥 − 8𝑦 + 𝑧 = −21 tomando 𝑥0 = ⎜
⎜−1⎟

D. C D

⎟

D. R P

−2𝑥 + 𝑦 + 5𝑧 = 15 ⎝0⎠
−1 −1
x𝑛+1 = (𝐷 + 𝜔𝐴𝐿 ) ((1 − 𝜔)𝐷 − 𝜔𝐴𝑅 ) x𝑛 + 𝜔 (𝐷 + 𝜔𝐴𝐿 ) y

1
Tomando 𝜔 = 2

4 4
⎛ ⎞ ⎛ ⎞
(𝐷 + 𝜔𝐴𝐿 ) = ⎜
⎜ ⎟ +1⎛
⎜4 ⎞ ⎟ ⎜
⎜ 2 −8 ⎟⎟
⎜ −8 ⎟
⎟ 2⎜⎜ ⎟
⎟= ⎜
⎜ ⎟
⎟
1
⎝ 5⎠ ⎝−2 1 ⎠ ⎝−1 2 5⎠
−40
−1 1 ⎛⎜ ⎞
⎟
(𝐷 + 𝜔𝐴𝐿 ) = − ⎜ −10 20 ⎟
160 ⎜ ⎟
⎝ −7 −2 −32 ⎠
4 1 −1
1⎛⎜ ⎞
((1 − 𝜔)𝐷 − 𝜔𝐴𝑅 ) = ⎜
⎜ 0 −8 −1⎟⎟
⎟
2
⎝ 0 0 5 ⎠
 8 / 17
−40 4 1 −1
1 ⎛⎜ ⎞
⎟ 1⎛⎜ ⎞
𝐵SOR =− ⎜
⎜ −10 20 ⎟
⎟ ⋅ ⎜
⎜ −8 −1⎟⎟
⎟
160 2

D. C D

⎝ −7 −2 −32⎠ ⎝ 5⎠

D. R P

−160 −40 40
1 ⎛⎜ ⎞
=− ⎜ −40 −170 −10 ⎟
⎟
320 ⎜ ⎟
⎝ −28 9 −151 ⎠
−160 −40 40 1 −40 7
1 ⎛⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟ 1 ⎛⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟
x1 = − ⎜
⎜ −40 −170 −10 ⎟
⎟ ⎜
⎜ −1 ⎟
⎟ − ⎜
⎜−10 20 ⎟
⎟ ⎜
⎜ −21⎟
⎟
320 320
⎝ −28 9 −151⎠ ⎝ 0 ⎠ ⎝ −7 −2 −32⎠ ⎝ 15 ⎠
1.2500
⎛ ⎞
=⎜
⎜
⎜1.1250⎟
⎟
⎟
⎝1.6375⎠
−160 −40 40 1.2500 −40 7
1 ⎛⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟ 1 ⎛⎜ ⎞
⎟ ⎛
⎜ ⎞
⎟
x2 = − ⎜ −40 −170 −10 ⎟
⎟⎜⎜1.1250⎟ ⎟ − 320 ⎜
⎜−10 20 ⎟
⎟⎜⎜−21⎟
320 ⎜ ⎟
⎝ −28 9 −151⎠ ⎝1.6375⎠ ⎝ −7 −2 −32⎠ ⎝ 15 ⎠
𝑇
= (1.43594 2.33633 2.37230)
8 / 17
 9 / 17
Ejemplo 4x4

D. C D

D. R P
Consideremos nuevamente el sistema de ecuaciones lineales 𝐴𝑥 = 𝑏:

4 1 −1 0 ⎛𝑥1 ⎞ 7
⎛
⎜ ⎞
⎟ ⎜ ⎟ ⎛
⎜ ⎞
⎜
⎜ 1 −3 1 1 ⎟
⎟
⎜
⎜𝑥2 ⎟
⎟ ⎜ −6⎟
⎟
⎜ ⎟ ⎜ ⎟ =⎜ ⎟
⎜
⎜
⎜−1 1 5 −2⎟⎟⎜
⎜𝑥 ⎟ ⎜8⎟ ⎟
⎜ 3⎟⎟ ⎜ ⎟
⎝ 0 1 −2 4 ⎠ ⎝𝑥4 ⎠ ⎝ 9 ⎠

9 / 17
 10 / 17
𝑇 𝑇
Sol. Exacta ( 79 , 27
7
, 2, 16
7
) ≈ (1.28571 3.85714 2. 2.28571)

D. C D

D. R P
(10) (10) (10) (10)
𝜔 ‖𝐴x(10) − 𝑏‖2 𝑥1 𝑥2 𝑥3 𝑥4
0.27 1.31921967 1.40251619 3.50913059 1.96895893 2.27450286
0.34 0.53507918 1.36965991 3.74318256 2.01887573 2.30913989
0.41 0.20398501 1.33173615 3.83876438 2.02269130 2.30601781
0.48 0.08539991 1.30511488 3.86506127 2.01351385 2.29575259
0.55 0.03723608 1.29138438 3.86533853 2.00531246 2.28883237
0.62 0.01292444 1.28637080 3.86066772 2.00121273 2.28610488
0.69 0.00288954 1.28541934 3.85785074 2.00000641 2.28558993
0.76 0.00057408 1.28557714 3.85709937 1.99992128 2.28566268
0.83 0.00012591 1.28571270 3.85710852 1.99999190 2.28571537
0.90 0.00001991 1.28571788 3.85714681 2.00000212 2.28571479
0.97 0.00000532 1.28571363 3.85714157 1.99999962 2.28571438
1.04 0.00000727 1.28571343 3.85714451 1.99999953 2.28571359
1.11 0.00008102 1.28570946 3.85712266 1.99999903 2.28572127
1.18 0.00076212 1.28586651 3.85725062 2.00009980 2.28572210
1.25 0.00556657 1.28636781 3.85830121 2.00040481 2.28534738
1.32 0.01638024 1.28545241 3.86052473 1.99980985 2.28328740
1.39 0.06935489 1.26921144 3.85241619 1.99014811 2.27845105
1.46 0.60864958 1.17059579 3.75269211 1.93416081 2.28663656
1.5310 3.51675540 0.73336014 3.18376911 1.69765198 2.43301399
/ 17
 11 / 17

Gráfico de Error vs ω

D. C D

D. R P
100

10−1
Error=‖x10 − x * ‖

10−2

10−3

10−4

10−5

11 0.4 0.6 0.8 1.0 1.2 1.4

/ 17
ω
 12 / 17

Gráfico de Número de iteraciones vs ω

45

D. C D

D. R P
40

35
Número de iteraciones

30

25

20

15

10

12 0.4 0.6 0.8 1.0 1.2 1.4

/ 17
ω
 13 / 17
Condiciones necesarias: SOR
Kahan [Kress, Th. 4.11]

D. C D

D. R P
Si el esquema iterativo

x(𝜈+1) ∶= x(𝜈) + 𝜔𝐷−1 (y − 𝐴𝐿 x(𝜈+1) − (𝐷 + 𝐴𝑅 )x(𝜈) )

es convergente, entonces 0 < 𝜔 < 2

Sketch de la prueba
1.
𝑛
∏ 𝜆𝑖 = det 𝐵(𝜔) = det(𝐷 + 𝜔𝐴𝐿 )−1 [(1 − 𝜔)𝐷 − 𝜔𝐴𝑅 ]
𝑗=1

= det(𝐼 + 𝜔𝐷−1 𝐴𝐿 )−1 [(1 − 𝜔)𝐼 − 𝜔𝐷−1 𝐴𝑅 ]

2. (𝐼 + 𝜔𝐷−1 𝐴𝐿 )−1 triangular inferior

3. (1 − 𝜔)𝐼 − 𝜔𝐷−1 𝐴𝑅 triangular superior
4. 𝜌(𝐵(𝜔)) ≥ | 1 − 𝜔 |

13
/ 17
 14 / 17
function [x,iter]= JOR(A,b,x0,nmax,tol,omega)
% Resuelve un sistema Ax=b usando el metodo JOR

D. C D

D. R P
iter = 0;
r = b - A * x0;
r0 = norm(r);
err = norm (r);
x = x0;
P = diag(diag(A).^(-1));

while err > tol && iter < nmax

iter = iter +1;
xnew = x + omega*P*r;
x = xnew;
r = b-A*x;
err = norm(r)/r0;
end
end
14
/ 17
function [x,iter] = SOR(A,b,x0,nmax,tol,omega) 15 / 17
% Resuelve un sistema Ax=b usando el metodo SOR

[n, m]=size(A);

D. C D

if n ~=m, error('A debe ser una matriz cuadrada'); end

D. R P

iter = 0;
r = b-A*x0;
r0 = norm(r);
err = norm(r);
xold = x0;

while err>tol && iter < nmax

iter = iter +1;
for i = 1:n
s = 0;
for j = 1:i-1
s = s+A(i,j)*x(j);
end
for j=i+1:n
s = s +A(i,j)*xold(j);
end
x(i,1) = omega*(b(i)-s)/A(i,i) + (1-omega)*xold(i);
end
xold = x;
r = b-A*x;
err = norm(r)/r0;
end
return
15
/ 17
 16 / 17
Ejercicio

D. C D

2 −1
⎛ ⎞

D. R P

⎜
⎜ −1 2 −1 ⎟
⎟
⎜
⎜ ⎟
⎟
⎜ −1 2 −1 ⎟
𝐴=⎜
⎜ ⎟
⎟
⎜
⎜ −1 2 1 ⎟
⎟
⎜
⎜ ⎟
⎟
⎜ −1 2 −1 ⎟
⎝ −1 2 ⎠
1. Converge el método de Jacobi? En caso afirmativo, determine el numero
de iteraciones.
2. Converge el método de Gauss? En caso afirmativo, determine el numero
de iteraciones.
3. Converge el método de JOR? En caso afirmativo, determine el numero
de iteraciones.
4. Realice una gráfica de 𝜔 vs. iteraciones.
16
/ 17
Resultados de convergencia

1. 𝜌(𝐷−1 (𝐴𝐿 + 𝐴𝑅 )) < 1 ⇔ Jacobi converge.

D. C D

2. 𝜌((𝐷 + 𝐴𝐿 )−1 𝐴𝑅 ) < 1 ⇔ Gauss-Seidel converge.

D. R P

3. 𝜌(𝐼 − 𝜔𝐷−1 𝐴) < 1 ⇔ JOR converge.
4. Si 𝐴 es diagonalmente dominante por filas entonces los métodos
de Jacobi y Gauss Seidel convergen.
5. Si 𝐴 es diagonalmente dominante por columnas entonces el
métodos de Jacobi converge.
6. Si el método de Jacobi converge entonces JOR converge para
0 < 𝜔 ≤ 1.
Kress, Rainer 17 / 17
Numerical analysis,
Graduate Texts in Mathematics

D. C D

181,

D. R P

Springer-Verlag,
1998,
Quarteroni, Alfio M. and Saleri, Fausto and Sacco, Riccardo
Numerical Mathematics,
Text in applied mathematics 37,
Springer Verlag, New York,
2000,
Trefethen L; Bau David
Numerical Linear Algebra,
Society for Industrial and Applied Mathematics (SIAM),
1997,
Mathews, John H. Fink, Kurtis D.
Métodos numéricos con MATLAB 3ED
Prentice Hall
2000
17
/ 17