Embed Size (px)
Citation preview
![Page 1: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/1.jpg)
Clase No. 9:
Métodos iterativos para resolversistemas de ecuaciones lineales
MAT–251 Dr. Alonso Ramírez ManzanaresCIMAT A.C.e-mail: [email protected]: http://www.cimat.mx/salram/met_num/
Dr. Joaquín Peña AcevedoCIMAT A.C.e-mail: [email protected]
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 1 / 21
![Page 2: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/2.jpg)
Introducción a los métodos iterativos
Tenemos que la solución x� de Ax = b satisface
x� =Mx� + c (2)
Restando (2) de (1) se obtiene
xt � x� =M(xt�1 � x�)
=) kxt � x�k kMk kxt�1 � x�k
Proposición
Si kMk < 1, entonces la sucesión {xt}, con xt = Mxt�1 + c, converge paracualquier x0.
kxt � x�k kMk kxt�1 � x�k < kxt�1 � x�k
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 5 / 21
![Page 3: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/3.jpg)
Un comentario sobre las normas en Rn
En general, no conocemos x�. Por ello, un criterio para terminar de iterar elalgoritmo es cuando se cumpla
kxt � xt�1kp�m + kxtk < tol o
kAxt � bkp�m + kbk < tol
Proposición
Si k · ka y k · kb son dos normas en Rn, entonces existen constantes �,� 2 Rtales que para todo x 2 Rn
�kxka kxkb �kxka
En general, cuando dos normas cumplen las desigualdades anteriores sedice que son equivalentes.
Desde el punto de vista computacional, esto nos da la libertad de usar enlos algoritmos una norma que no sea costosa de calcular y que nointroduzca demasiados errores al calcularla.
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 6 / 21
![Page 4: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/4.jpg)
Un comentario sobre las normas en Rn
En general, no conocemos x�. Por ello, un criterio para terminar de iterar elalgoritmo es cuando se cumpla
kxt � xt�1kp�m + kxtk < tol o
kAxt � bkp�m + kbk < tol
Proposición
Si k · ka y k · kb son dos normas en Rn, entonces existen constantes �,� 2 Rtales que para todo x 2 Rn
�kxka kxkb �kxka
En general, cuando dos normas cumplen las desigualdades anteriores sedice que son equivalentes.
Desde el punto de vista computacional, esto nos da la libertad de usar enlos algoritmos una norma que no sea costosa de calcular y que nointroduzca demasiados errores al calcularla.
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 6 / 21
![Page 5: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/5.jpg)
Un comentario sobre las normas en Rn
En general, no conocemos x�. Por ello, un criterio para terminar de iterar elalgoritmo es cuando se cumpla
kxt � xt�1kp�m + kxtk < tol o
kAxt � bkp�m + kbk < tol
Proposición
Si k · ka y k · kb son dos normas en Rn, entonces existen constantes �,� 2 Rtales que para todo x 2 Rn
�kxka kxkb �kxka
En general, cuando dos normas cumplen las desigualdades anteriores sedice que son equivalentes.
Desde el punto de vista computacional, esto nos da la libertad de usar enlos algoritmos una norma que no sea costosa de calcular y que nointroduzca demasiados errores al calcularla.
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 6 / 21
![Page 6: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/6.jpg)
Método de Jacobi
Queremos resolver Ax = b haciendo una separación de la matrizA = [aij] =Q� P. Dado que x =Q�1(Px+b), conviene elegir Q de modo quesu inversa sea fácil de calcular.
En el método de Jacobi se elige
Q =
2664
a11 0 · · · 00 a22 · · · 0...
.... . .
...0 0 · · · ann
3775 , P =
2664
0 �a22 · · · �a1n�a21 0 · · · �a2n...
.... . .
...�an1 �an2 · · · 0
3775
Si b = (b1, ...,bn)>, xt�1 =Äxt�11 ,xt�12 , ...,xt�1n
ä>, entonces para
i = 1,2, ...,n, la componente i-ésima de xt está dada por
xti =1
aii
0B@bi �
nX
j=1j 6=i
aijxt�1j
1CA
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 7 / 21
![Page 7: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/7.jpg)
Método de Jacobi
Queremos resolver Ax = b haciendo una separación de la matrizA = [aij] =Q� P. Dado que x =Q�1(Px+b), conviene elegir Q de modo quesu inversa sea fácil de calcular.
En el método de Jacobi se elige
Q =
2664
a11 0 · · · 00 a22 · · · 0...
.... . .
...0 0 · · · ann
3775 , P =
2664
0 �a22 · · · �a1n�a21 0 · · · �a2n...
.... . .
...�an1 �an2 · · · 0
3775
Si b = (b1, ...,bn)>, xt�1 =Äxt�11 ,xt�12 , ...,xt�1n
ä>, entonces para
i = 1,2, ...,n, la componente i-ésima de xt está dada por
xti =1
aii
0B@bi �
nX
j=1j 6=i
aijxt�1j
1CA
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 7 / 21
![Page 8: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/8.jpg)
Método de Jacobi
Queremos resolver Ax = b haciendo una separación de la matrizA = [aij] =Q� P. Dado que x =Q�1(Px+b), conviene elegir Q de modo quesu inversa sea fácil de calcular.
En el método de Jacobi se elige
Q =
2664
a11 0 · · · 00 a22 · · · 0...
.... . .
...0 0 · · · ann
3775 , P =
2664
0 �a22 · · · �a1n�a21 0 · · · �a2n...
.... . .
...�an1 �an2 · · · 0
3775
Si b = (b1, ...,bn)>, xt�1 =Äxt�11 ,xt�12 , ...,xt�1n
ä>, entonces para
i = 1,2, ...,n, la componente i-ésima de xt está dada por
xti =1
aii
0B@bi �
nX
j=1j 6=i
aijxt�1j
1CA
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 7 / 21
![Page 9: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/9.jpg)
Alonso Ramírez Manzanares Métodos Numéricos 09.09.2015
Método Iterativo de Jacobi
• Ejemplo para el SEL con solución exacta
• Para tener un sistema
![Page 10: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/10.jpg)
Convergencia del método de Jacobi
Note que
M =Q�1P = I� Q�1Adonde I es la matriz identidad de tamaño n. Entonces
kMk� = max1in
nX
j=1j 6=i
|aij||aii|
Si para todo i tenemos que
|aii| >nX
j=1j 6=i
|aij| =) kMk� < 1
y por tanto la sucesión es convergente.
Proposición
Si la matriz A es estrictamente diagonal dominante, entonces el método deJacobi converge para cualquier vector inicial x0.
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 8 / 21
![Page 11: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/11.jpg)
Ejemplo 1.
A =
26664
15 5 0 2 65 14 0 1 �14 2 �18 1 63 6 2 12 05 �4 4 1 15
37775
Elegimos
x� =
0BBB@
100020003000�2000�1000
1CCCA , =) b = Ax� =
0BBB@
1500032000�54000�3000�8000
1CCCA
Iniciamos con x0 =
0BBB@
00000
1CCCA.
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 9 / 21
![Page 12: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/12.jpg)
Ejemplo 1.
t=1 t=3 t=5 t=7 t=9 t=11x1 1000.000 1065.794 996.807 996.451 998.268 999.303x2 2285.714 2191.327 2050.923 2012.788 2003.548 2001.011x3 3000.000 2853.889 2986.866 2997.922 2999.722 3000.026x4 -250.000 -1872.778 -1969.974 -1995.201 -1999.320 -1999.953x5 -533.333 -919.048 -982.239 -992.414 -997.533 -999.193
t=13 t=15 t=17 t=19 t=21 t=23x1 999.747 999.912 999.970 999.990 999.997 999.999x2 2000.297 2000.090 2000.028 2000.009 2000.003 2000.001x3 3000.031 3000.015 3000.006 3000.002 3000.001 3000.000x4 -2000.023 -2000.016 -2000.007 -2000.003 -2000.001 -2000.000x5 -999.736 -999.914 -999.972 -999.991 -999.997 -999.999
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 10 / 21
![Page 13: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/13.jpg)
Otra condición suficiente.
Decimos que v es un eigenvector de la matriz A si existe � 2 R tal que
Av = �v,
y llamamos a � es un eigenvalor asociado a v.
El radio espectral �(A) de una matriz A se define como
�(A) =max{ |�| : � es eigenvalor de A}
Proposición
Sea A una matriz n⇥ n. Entonces �(A) kAk para cualquier norma natural.Además, Ak ! 0 si y sólo si �(A) < 1.
Proposición
La sucesión xt =Mxt�1 + c converge para cualquier x0 si �(M) < 1.
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 11 / 21
![Page 14: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/14.jpg)
Ejemplo 2.Consideremos las matrices de la forma
A=
2666666664
1 0 0 · · · 0 0 0�1 2 �1 · · · 0 0 00 �1 2 · · · 0 0 0...
......
. . ....
......
0 0 0 · · · 2 �1 00 0 0 · · · �1 2 �10 0 0 · · · 0 0 1
3777777775
50 100 150 200
3.90
3.92
3.94
3.96
3.98
4.00
Dimension
Rad
io e
spec
tral
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 12 / 21
![Page 15: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/15.jpg)
Ejemplo 2.
Fijamos n = 100. Sea b = (0.5,0,0, ...,0,0,0.5)>. Entonces, inicializamos elmétodo de Jacobi con el siguiente vector x0:
0 20 40 60 80 100
0.0
0.5
1.0
1.5
3
Componente
Res
ulta
do p
arci
al
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 13 / 21
![Page 16: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/16.jpg)
Ejemplo 2.
0 20 40 60 80 100
0.0
0.5
1.0
1.5
4
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.5
1.0
1.5
13
ComponenteR
esul
tado
par
cial
Iteración 1 Iteración 10
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 14 / 21
![Page 17: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/17.jpg)
Ejemplo 2.
0 20 40 60 80 100
0.0
0.5
1.0
1.5
23
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.5
1.0
1.5
53
ComponenteR
esul
tado
par
cial
Iteración 20 Iteración 50
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 15 / 21
![Page 18: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/18.jpg)
Ejemplo 2: Iteración 350
0 20 40 60 80 100
0.0
0.5
1.0
1.5
350
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
A*x^t
0 50 100 150 200 250 300 350
0.1
0.2
0.3
0.4
0.5
Iteracion
Error
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 16 / 21
![Page 19: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/19.jpg)
Ejemplo 2: Iteración 950
0 20 40 60 80 100
0.0
0.5
1.0
1.5
950
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
A*x^t
0 200 400 600 800
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
Error
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 17 / 21
![Page 20: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/20.jpg)
Ejemplo 2: Iteración 1250
0 20 40 60 80 100
0.0
0.5
1.0
1.5
1250
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
A*x^t
0 200 400 600 800 1000 1200
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
Error
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 18 / 21
![Page 21: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/21.jpg)
Ejemplo 2: Iteración 1550
0 20 40 60 80 100
0.0
0.5
1.0
1.5
1550
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
A*x^t
0 500 1000 1500
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
Error
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 19 / 21
![Page 22: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/22.jpg)
Ejemplo 2: Iteración 1850
0 20 40 60 80 100
0.0
0.5
1.0
1.5
1850
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
A*x^t
0 500 1000 1500
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
Error
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 20 / 21
![Page 23: Métodos iterativos para resolver SELalram/met_num/clases/clase08a.pdfMétodo de Jacobi Queremos resolver Ax = b haciendo una separación de la matriz A =[aij]=Q P. Dado que x = Q1(Px+b),](https://reader033.document.onl/reader033/viewer/2022041807/5e5505095d53750726522189/html5/thumbnails/23.jpg)
Ejemplo 2: Iteración 2150
0 20 40 60 80 100
0.0
0.5
1.0
1.5
2150
Componente
Res
ulta
do p
arci
al
0 20 40 60 80 100
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
A*x^t
1000 1200 1400 1600 1800 2000
0.005
0.006
0.007
0.008
Iteracion
Error
0 500 1000 1500 2000
0.0
0.1
0.2
0.3
0.4
0.5
Iteracion
Error
Joaquín Peña (CIMAT) Métodos Numéricos (MAT–251) 07.09.2015 21 / 21
