Logo Studenta

04_Sistemas_de_Ecuaciones_Algebraicas

¡Este material tiene más páginas!

Vista previa del material en texto

NOTAS DE CLASE
Tema: Sistemas de Ecuaciones Lineales
Métodos Directos
Métodos Iterativos
 Mal condicionamiento
Tema: Sistemas de ecuaciones No Lineales
Métodos de Punto Fijo y Newton-Raphson
14/12/2020
 
Dada una matriz A (m x n) y un vector B (m x 1), 
encontrar el vector X (n x 1) que satisface la ecuación
A X = B (m ecuaciones con n incógnitas)
Condiciones suficientes para la existencia de solución única son 
- m = n 
- determinante de A  0
EXISTENCIA y UNICIDAD DE LA SOLUCIÓN
A es regular <==> det (A)  0 <==> existe A-1 Recordar:
 
TRANSFORMACIONES ELEMENTALES
Cualquiera de las siguientes operaciones aplicadas a un sistema de ecuaciones lineales 
produce un sistema equivalente
- INTERCAMBIO: el orden de las ecuaciones puede alterarse
 - ESCALAMIENTO: multiplicar una ecuación por una constante no nula
 - SUSTITUCIÓN: una ecuación puede ser reemplazada por la suma de ella misma 
más un múltiplo de otra ecuación
A
nxn
 X
nx1
 = B
nx1
 
ELIMINACIÓN GAUSSIANA
sistema de ecuaciones
MATRIZ AMPLIADA
 
ELIMINACIÓN GAUSSIANA
sistema de ecuaciones
MATRIZ AMPLIADA
La primera fila se usa para eliminar los elementos de la primera columna que 
están por debajo de la diagonal principal. El elemento a11 es el pivote. Los valores 
de mk1 = ak1/a11 son los multiplicadores.
pivote
m21= 2
m31= 4
m41= -3
 
Ahora, la segunda fila se usa para eliminar los elementos de la segunda columna 
que están por debajo de la diagonal principal. El elemento a22 es el pivote. Los 
valores de mk2 = ak2/a22 son los multiplicadores.
pivote
m32= 1.5
m42= -1.75
Finalmente, restando la cuarta fila de la tercera multiplicada por m43 = a43/a33 =-1.9 
se obtiene el sistema triangular superior.
pivote
m43= -1.9
 
Utilizando el algoritmo de SUSTITUCIÓN HACIA ATRÁS obtenemos
 
IDEA
ALGORITMO
PASO 1
ELIMINACIÓN GAUSSIANA
Transformar el sistema A X = B en un sistema equivalente 
A* X = B* tal que A* sea una matriz triangular superior.
A matriz invertible de n x n, B vector n x 1.
Almacenar todos los coeficientes en la matriz ampliada (A|B)
PASO 2 Eliminar la incógnita xq en todas las filas desde la q+1 en adelante. 
PASO 3 Sustitución hacia atrás.
 
Código F90 BÁSICO del método de Eliminación de Gauss simple
program sist_ec_lin
use kinds 
implicit none
integer(ikind4) :: n, i
real(rkind8), allocatable :: A(:,:)
real(rkind8), allocatable :: B(:), x(:)
open(unit=10,file="input.in",status="old")
read(10,*) n
allocate(A(n,n),B(n),x(n))
do i=1,n ! Lectura de la matriz de coeficientes y t.i.
 read(10,*) A(i,:), B(i)
end do
close(10)
call elim_gauss_v0(n,A,B,x)
write(*,*)"x=", x
deallocate(A,B,x)
end program sist_ec_lin
Sample input file: input.in
subroutine elim_gauss_v0(n,A,B,x) 
use kinds ! Eliminación de Gauss (versión preliminar)
implicit none
integer(ikind4) :: n, i, j
real(rkind8), dimension(n,n) :: A
real(rkind8), dimension(n) :: B, x
real(rkind8), dimension(n,n+1) :: Amp
real(rkind8) :: mult, sumi
! -----------------------------------------------------------------------------------------------------
do i=1, n ! Armamos la matriz ampliada
 Amp(:,i)=A(:,i)
end do
Amp(:,n+1)=B(:)
do i=1, n-1 ! Reducción de la matriz ampliada
 if ( Amp(i,i) == 0 ) then
 write(*,*)"ERROR: cabecera nula" ; stop
 endif
 do j=i+1, n
 mult=Amp(j,i)/Amp(i,i) ; Amp(j,:)=Amp(j,:)-mult*Amp(i,:)
 end do
end do
x(n)=Amp(n,n+1)/Amp(n,n) ! Sustitución hacia atrás
do i=n-1,1,-1
 sumi=0._rkind8
 do j=i+1, n
 sumi=sumi+Amp(i,j)*x(j)
 end do
 x(i)=(Amp(i,n+1)-sumi)/Amp(i,i)
end do
! -----------------------------------------------------------------------------------------------------
end subroutine elim_gauss_v0
2
1 1 2 
1 -1 0
 
Mejoras simples, necesarias para los códigos anteriores:
- verificar en el programa principal que A.x=B: evaluación de || A.x – B ||
- más mejoras ?
- implementar en la subrutina la permutación de filas cuando se encuentra cabecera nula y
 detectar posible incompatibilidad o indeterminación del sistema de ecuaciones
 
subroutine elim_gauss(n,A,B,x) ! Eliminación de Gauss (SIN pivoteo parcial)
use kinds 
implicit none
integer(ikind4) :: n, i, j
real(rkind8), dimension(n,n) :: A
real(rkind8), dimension(n) :: B, x
real(rkind8), dimension(n,n+1) :: Amp
real(rkind8) :: mult, sumi
real(rkind8), dimension(n+1) :: aux
! ---------------------------------------------------------------------------------------------
do i=1, n !Armarmos la matriz ampliada
 Amp(:,i)=A(:,i)
end do
Amp(:,n+1)=B(:)
do i=1, n-1 ! Reducción de la matriz ampliada
 if ( Amp(i,i) == 0 ) then
 do j=i+1, n
 if ( Amp(j,i) /= 0 ) then
 aux=Amp(i,:) ; Amp(i,:)=Amp(j,:) ; Amp(j,:)=aux
 exit
 elseif ( j == n) then
 write(*,*)"ERROR: no se puede aplicar el método"
 stop
 endif
 end do 
 endif
 do j=i+1, n
 mult=Amp(j,i)/Amp(i,i) ; Amp(j,:)=Amp(j,:)-mult*Amp(i,:)
 end do
end do
x(n)=Amp(n,n+1)/Amp(n,n) ! Sustitución hacia atrás
do i=n-1,1,-1
 sumi=0.
 do j=i+1, n
 sumi=sumi+Amp(i,j)*x(j)
 end do
 x(i)=(Amp(i,n+1)-sumi)/Amp(i,i)
end do
! ---------------------------------------------------------------------------------------------
end subroutine elim_gauss
program sist_ec_lin
use kinds
implicit none
integer(ikind4) :: n, i
real(rkind8), allocatable :: A(:,:)
real(rkind8), allocatable :: B(:), x(:), test(:)
open(unit=10,file="input.in",status="old")
read(10,*) n
allocate(A(n,n),B(n),x(n),test(n))
do i=1,n ! Lectura de la matriz de coeficientes y t.i.
 read(10,*) A(i,:), B(i)
end do
close(10)
call elim_gauss(n,A,B,x)
write(*,*)"x=", x
test=matmul(A,x)-B
write(*,*) "Norma de Ax-B=",norm2(test) 
deallocate(A,B,x,test)
end program sist_ec_lin
Sample input file: input.in
2
0 1 1 
1 1 2
 
ESTRATEGIAS DE PIVOTEO PARA REDUCIR LOS ERRORES DE REDONDEO
1.133 x1+5 .281 x2=6 .414
24 .14 x1−1 .210 x2=22 .03
Los valores
son la solución del sistema
x1=x2=1 .000
Resolvemos el sistema haciendo las operaciones con 4 dígitos significativos, se obtiene:
x1=0 .9956 , x2=1 .001
24 .14 x1−1 .210 x2=22 .03
1.133 x1+5 .281 x2=6 .414
Intercambiamos las ecuaciones del sistema anterior
Resolvemos nuevamente utilizando 4 dígitos significativos, se obtiene: x1=1.000 , x 2=1 .000
Pivoteo trivial y pivoteo parcialPivoteo trivial y pivoteo parcial
 
ELIMINACIÓN GAUSSIANA CON PIVOTEO PARCIAL
IDEA
ALGORITMO
PASO 1
Permanecen sin cambios
PASO 2
 Encontrar la fila k “más grande”, es decir, 
|akq| = max {|aqq|, |aq+1q|, …, |an-1q|, |anq|}. 
Intercambiar la fila q-ésima con la fila k-ésima
Eliminar la incógnita xq en todas las filas desde la q+1 en adelante.
PASO 3 Permanece sin cambios
 
PIVOTEO PARCIAL VS. TOTAL
PIVOTEO TOTAL: En la etapa r-esima del proceso de eliminación, elegimos el elemento 
a
pq
=max{ | a
ij
 | / i,j ≥r } como nuevo pivote. Para ello intercambiamos las filas r y p, y las 
columnas r y q de forma que situamos el elemento a
pq
 en laposicion (r,r). Obviamente 
hemos tomado i,j ≥r para no perturbar los ceros que ya teníamos. Posteriormente 
continuamos la eliminación con el nuevo pivote
PIVOTEO PARCIAL: Se busca solamente en la r-esima columna; es decir, se toma 
a
pr
 = max { | a
ij
 | / i ≥ r } como nuevo pivote. Para ello intercambiamos las filas r y p, 
continuando posteriormente el proceso de eliminación.
En la práctica, el método de Gauss con pivoteo total puede consumir mucho tiempo 
computacional, pues para hallar el máximo en cada paso hay que buscar entre (m - r +1).(n - 
r +1) elementos. En el otro caso, además del ahorro de tiempo, las incógnitas de nuestro sistema 
no cambian de orden en el sistema reducido. Por ello, en general, es suficiente utilizar un pivoteo 
parcial.
 
subroutine elim_gauss_pivoteo_parcial(n,A,B,x)! Eliminación de Gauss CON pivoteo parcial
use kinds
implicit none
integer(ikind4) :: n, i, j, row_max_abs_elem
real(rkind8), dimension(n,n) :: A
real(rkind8), dimension(n) :: B, x
real(rkind8), dimension(n,n+1) :: Amp
real(rkind8) :: mult, sumi, max_abs_elem
real(rkind8), dimension(n+1) :: aux
! ---------------------------------------------------------------------------------------------
do i=1, n ! Armarmos la matriz ampliada
 Amp(:,i)=A(:,i)
end do
Amp(:,n+1)=B(:)
do i=1, n-1 ! Reducción de la matriz ampliada
 max_abs_elem=abs(Amp(i,i)) ; row_max_abs_elem=i
 do j=i+1,n
 if (abs(Amp(j,i)) > max_abs_elem) then
 max_abs_elem=abs(Amp(j,i)) ; row_max_abs_elem=j
 endif
 end do
 aux=Amp(i,:) ; Amp(i,:)=Amp(row_max_abs_elem,:) ; Amp(row_max_abs_elem,:)=aux
 if ( Amp(i,i) == 0 ) then
 write(*,*)"ERROR: no se puede aplicar el método" ; stop
 end if
 do j=i+1, n
 mult=Amp(j,i)/Amp(i,i) ; Amp(j,:)=Amp(j,:)-mult*Amp(i,:)
 end do
end do
x(n)=Amp(n,n+1)/Amp(n,n) ! Sustitución hacia atrás
do i=n-1,1,-1
 sumi=0._rkind8
 do j=i+1, n
 sumi=sumi+Amp(i,j)*x(j)
 end do
 x(i)=(Amp(i,n+1)-sumi)/Amp(i,i)
end do
! ---------------------------------------------------------------------------------------------
end subroutine elim_gauss_pivoteo_parcial
 
 
¿ Qué cambios hay que introducir a los algoritmos 
para implementar eliminación de Gauss-Jordan ?
subroutine elim_gauss_jordan_pivoteo_parcial(n,A,B,x) ! Eliminación de Gauss-Jordan CON pivoteo 
parcial
implicit none
integer(4) :: n, i, j, row_max_abs_elem
real(8), dimension(n,n) :: A
real(8), dimension(n) :: B, x
real(8), dimension(n,n+1) :: Amp
real(8) :: mult, sumi, max_abs_elem
real(8), dimension(n+1) :: aux
do i=1, n !Armarmos la matriz ampliada
 Amp(:,i)=A(:,i)
end do
Amp(:,n+1)=B(:)
do i=1, n ! Reducción de la matriz ampliada
 max_abs_elem=abs(Amp(i,i)) ; row_max_abs_elem=i
 do j=i+1,n
 if (abs(Amp(j,i)) > max_abs_elem) then
 max_abs_elem=abs(Amp(j,i)) ; row_max_abs_elem=j
 end if
 end do
 aux=Amp(i,:) ; Amp(i,:)=Amp(row_max_abs_elem,:)
 Amp(row_max_abs_elem,:)=aux
 if ( Amp(i,i) == 0 ) then
 write(*,*)"ERROR: no se puede aplicar el método" ; stop
 endif
 do j=1, n
 if (j == i) cycle 
 mult=Amp(j,i)/Amp(i,i)
 Amp(j,:)=Amp(j,:)-mult*Amp(i,:)
 end do
 Amp(i,:)=Amp(i,:)/Amp(i,i)
end do
x(:)=Amp(:,n+1)
end subroutine elim_gauss_jordan_pivoteo_parcial
program sist_ec_lin
implicit none
integer(4) :: n, i
real(8), allocatable :: A(:,:)
real(8), allocatable :: B(:), x(:), test(:)
open(unit=10,file="input.in",status="old")
read(10,*) n
allocate(A(n,n),B(n),x(n),test(n))
do i=1,n ! Lectura de la matriz de coeficientes y t.i.
 read(10,*) A(i,:), B(i)
end do
close(10)
call elim_gauss(n,A,B,x) ! Eliminación de Gauss
write(*,*)"Resultado Gauss, x=", x
test=matmul(A,x)-B
write(*,*) "Norma de Ax-B=",norm2(test) 
call elim_gauss_jordan_pivoteo_parcial(n,A,B,x) ! Elim. Gauss-Jordan
write(*,*)"Resultado Gauss-Jordan, x=", x
test=matmul(A,x)-B
write(*,*) "Norma de Ax-B=",norm2(test) 
deallocate(A,B,x,test)
end program sist_ec_lin
 
DEBE SER SIN PIVOTEO PARCIAL !!!!!
 
subroutine elim_gauss_jordan_pivoteo_parcial_msist(n,A,m,B,x)
! Esta rutina resuelve simultaneamente, m sistemas de ecuaciones nxn de la forma: Ax=b_i con i=1,2,...,m; 
! todos con la misma matriz de coeficientes A (nxn) y los m vectores b_i son las columnas de la matriz B (nxm). 
! Las soluciones de los m sistemas de ecuaciones son las m columnas de la matriz x (nxm) que es calculada. 
!
! El método utilizado es Eliminación de Gauss-Jordan CON pivoteo parcial
implicit none
integer(4) :: n, i, j, row_max_abs_elem, m
real(8), dimension(n,n) :: A
real(8), dimension(n,m) :: B, x
real(8), dimension(n,n+m) :: Amp
real(8) :: mult, sumi, max_abs_elem
real(8), dimension(n+m) :: aux
do i=1, n !Armarmos la matriz ampliada
 Amp(:,i)=A(:,i)
end do
Amp(:,n+1:n+m)=B
do i=1, n ! Reducción de la matriz ampliada
 max_abs_elem=abs(Amp(i,i)) ; row_max_abs_elem=i
 do j=i+1,n
 if (abs(Amp(j,i)) > max_abs_elem) then
 max_abs_elem=abs(Amp(j,i)) ; row_max_abs_elem=j
 end if
 end do
 aux=Amp(i,:) ; Amp(i,:)=Amp(row_max_abs_elem,:) ; Amp(row_max_abs_elem,:)=aux
 if ( Amp(i,i) == 0 ) then
 write(*,*)"ERROR: no se puede aplicar el método" ; stop
 end if
 do j=1, n
 if (j == i) cycle 
 mult=Amp(j,i)/Amp(i,i) ; Amp(j,:)=Amp(j,:)-mult*Amp(i,:)
 end do
 Amp(i,:)=Amp(i,:)/Amp(i,i)
end do
x=Amp(:,n+1:n+m)
end subroutine elim_gauss_jordan_pivoteo_parcial_msist
 
subroutine elim_gauss_jordan_msist(n,A,m,B,x)
! Esta rutina resuelve simultaneamente, m sistemas de ecuaciones nxn de la forma: Ax=b_i con i=1,2,...,m; 
! todos con la misma matriz de coeficientes A (nxn) y los m vectores b_i son las columnas de la matriz B (nxm). 
! Las soluciones de los m sistemas de ecuaciones son las m columnas de la matriz x (nxm) que es calculada. 
! Notar que si m=n y B es la matriz identidad nxn, la matriz x calculada es la matriz inversa de A, PERO EN ESTE CASO 
! el método de eliminación de Gauss-Jordan utilizado DEBE SER SIN pivoteo parcial
implicit none
integer(4) :: n, i, j, row_max_abs_elem, m
real(8), dimension(n,n) :: A
real(8), dimension(n,m) :: B, x
real(8), dimension(n,n+m) :: Amp
real(8) :: mult, sumi, max_abs_elem
real(8), dimension(n+m) :: aux
do i=1, n !Armarmos la matriz ampliada
 Amp(:,i)=A(:,i)
end do
Amp(:,n+1:n+m)=B
do i=1, n ! Reducción de la matriz ampliada
 max_abs_elem=abs(Amp(i,i)) ; row_max_abs_elem=i
 do j=i+1,n
 if (abs(Amp(j,i)) > max_abs_elem) then
 max_abs_elem=abs(Amp(j,i)) ; row_max_abs_elem=j
 end if
 end do
 aux=Amp(i,:) ; Amp(i,:)=Amp(row_max_abs_elem,:) ; Amp(row_max_abs_elem,:)=aux
 if ( Amp(i,i) == 0 ) then
 write(*,*)"ERROR: no se puede aplicar el método" ; stop
 end if
 do j=1, n
 if (j == i) cycle 
 mult=Amp(j,i)/Amp(i,i) ; Amp(j,:)=Amp(j,:)-mult*Amp(i,:)
 end do
 Amp(i,:)=Amp(i,:)/Amp(i,i)
end do
x=Amp(:,n+1:n+m)
end subroutine elim_gauss_jordan_msist
SI QUEREMOS USARLO PARA CALCULAR LA MATRIZ INVERSA: SIN PIVOTEO PARCIAL
 
CONSIDERACIONES NUMÉRICAS
 ERROR DE TRUNCAMIENTO: Ninguno. En la ausencia de error de 
redondeo se obtiene la respuesta exacta en un número finito de pasos.
 ERROR DE REDONDEO: Puede ser serio, se recomienda el uso 
Pivoteo Parcial para minimizar el efecto.
 ESFUERZO COMPUTACIONAL: Proceso de eliminación requiere de 
(n3) operaciones. Sustitución hacia atrás (n2) operaciones. Total de 
operaciones (n3). 
 SISTEMAS MAL CONDICIONADOS (pequeños cambios en los 
elementos de A o de B provocan grandes cambios en X): los métodos 
numéricos son proclives a tener más errores de redondeo (veremos más 
adelante).
 
FIN CLASE 1 de SISTEMAS DE ECUACIONES
 
Extienden los métodos iterativos introducidos en el capítulo de “raíces de 
ecuaciones no lineales” a espacios de mayor dimensión. En particular 
consideraremos extensiones del método de punto fijo aplicado a sistemas de 
ecuaciones lineales.
METODOS ITERATIVOS
x = g(x)
x: vector
g: función vectorial de x
 
MÉTODO ITERATIVO DE JACOBI
4 x1− x2+ x3= 7
4 x1− 8 x 2+ x3= −21
−2 x1+ x2+ 5 x3= 15
sistema de ecuaciones
El sistema anterior puede ser escrito en la 
forma siguiente:
x1=
7+x2−x3
4
x2=
21+4 x1+x3
8
x3=
15+2 x1−x2
5
x1
( k )
=
7+x2
(k−1 )
−x3
( k−1)
4
x2
( k )
=
21+4 x1
( k−1 )+x3
(k−1)
8
x3
( k )
=
15+2x1
(k−1 )−x2
( k−1)
5
Lo cual sugiere el siguiente proceso 
iterativo
 
x1
(2 )
=
7+3 .375−3
4
=1 .84375
x2
(2 )
=
21+4∗1.75+3
8
=3 .875
x3
(2 )
=
15+2∗1 .75−3 .375
5
=3 .025
Tomando el vector aproximación inicial x(0) = ( 1 , 2 , 2 ) 
x(1) = ( 1.75 , 3.375 , 3 ) 
x1
(1 )
=
7+2−2
4
=1.75
x2
(1 )
=
21+4∗1+2
8
=3 .375
x3
(1 )
=
15+2∗1−2
5=3
x(2) = ( 1.84375 , 3.875 , 3.025 ) 
MÉTODO ITERATIVO DE JACOBI
 
PASO 1
PASO 2
PASO 3
Sigo? Si ||x(k+1) − x(k)|| < tol, detengo el proceso y retorno el 
valor de x(k+1) como respuesta
Actualizo k ← k+1. Retorno al paso 1
x(k+1 )=g (x(k ) )Fórmula de recurrencia
MÉTODO ITERATIVO DE JACOBI
x j
(k+1 )
=
b j−a j 1 x1
(k )
−.. .−a jj−1 x j−1
(k )
−a jj+1 x j+1
(k )
−.. .−a jN xN
(k )
a jj
Norma
 
PASO 1
PASO 2
PASO 3
MÉTODO ITERATIVO DE GAUSS-SEIDEL
x j
(k+1 )
=
b j−a j 1 x1
(k+1 )
− .. .−a jj−1 x j−1
(k+1)
−a jj+1 x j+ 1
(k )
−. . .−a jN xN
(k )
a jj
Sigo? Si ||x(k+1) − x(k)|| < tol, detengo el proceso y retorno el 
valor de x(k+1) como respuesta
Actualizo k ← k+1. Retorno al paso 1
x(k+1 )=g (x(k ) )Fórmula de recurrencia
Norma
 
CONDICION SUFICIENTE DE CONVERGENCIA
Definición: Matriz de diagonal estrictamente dominante
|aii|>∑
j= 1
n
|aij| ∀ i=1,… ,n
Teorema:
Si la matriz de coeficientes A del sistema de ecuaciones es diagonal 
estrictamente dominante, entonces el proceso iterativo convergerá a la 
solución exacta, independientemente del vector inicial x(0)
 
What does it mean by ill conditioned and well-conditioned system of equations?
- A system of equations is considered to be well conditioned if a small change in the
coefficient matrix or a small change in the right hand side results in a small change in the
solution vector.
- A system of equations is considered to be ill conditioned if a small change in the
coefficient matrix or a small change in the right hand side results in a large change in the
solution vector.
Sistemas de ecuaciones algebraicas 
bien/mal condicionados
 
 
 
 
 
 
 
 
 
 
 
 
 
Some properties of the norm of a matrix
 
 
 
Para poder “confiar” en los 4 primeros dígitos significativos de un número, su error relativo 
debe ser ~ 10-4 (*)
(*) Por ejemplo, el mínimo número que se redondea a 5.043 (manteniendo 4 cifras significativas)
es 5.0425 y el máximo es 5.0434999999999. Es decir, x=5.0430±0.0005.
En una aritmética de punto flotante, se llama épsilon de la máquina (ε
mach
) 
al menor valor de una determinada máquina que cumple lo siguiente: 1.0 + ε
mach
 > 1.0 
ε
mach
 es el número decimal más pequeño que, sumado a 1, la computadora nos arroja 
un valor > 1 es decir, que no es redondeado a 1. Representa la exactitud relativa 
de la aritmética del computador. La existencia del épsilon de la máquina 
es una consecuencia de la precisión finita de la aritmética en coma flotante.
 
 
 
 
 
 
 
Sistemas de ecuaciones 
algebraicas NO lineales
 
MÉTODOS ITERATIVOS PARA SISTEMAS NO LINEALES
f 1 ( x,y )=x
2−2 x− y+0 .5
f 2 ( x,y )=x
2+4 y2−4
Queremos hallar la forma de resolver el sistemas de ecuaciones:
f 1 ( x,y )=0
f 2 ( x,y )=0
Consideremos las funciones
 
MÉTODOS ITERATIVOS PARA SISTEMAS NO LINEALES
f1(x,y) y f2(x,y) definen implícitamente curvas en el plano
f 1 ( x,y )=x
2−2 x− y+0 .5
f 2 ( x,y )=x
2+4 y2−4
Una solución del sistema es 
un punto (p,q) en el que ambas curvas 
se cruzan f1(p,q)=0 y f2(p,q)=0
es una parábola
es una elipse
 
Definición Un Punto Fijo del sistema de 2 ecuaciones
Definición El método de iteración de punto fijo es
es un punto (p,q) tal que p=g1(p,q) y q=g2(p,q).
x=g1 ( x,y )
y=g2 ( x,y )
pk+1=g1( pk ,qk )
qk+1=g2( pk ,qk ) k=0,1,. ..
METODO ITERATIVO DE PUNTO FIJO
 
despejar x de f1 
despejar y de f2 
x=
x2− y+0 .5
2
y=
−x2−4 y2+8 y+4
8
Tomando (p0,q0) generamos una sucesión de valores {(pk+1,qk+1)}
pk+1=g1( pk ,qk )=
pk
2−qk+0 .5
2
qk+1=g2( pk ,qk )=
−pk
2−4 qk
2+8 qk+4
8
(s.m.m. -8y a f2)
f 1 ( x,y )=x
2−2 x− y+0 .5=0
f 2 ( x,y )=x
2+4 y2−4=0
EJEMPLO
Sea el siguiente 
sistema no lineal
 
CASO I Si usamos como punto inicial (p0,q0) = (0,1), entonces
CONVERGE !
CASO II Si usamos como punto inicial (p0,q0) = (2,0), entonces
DIVERGE !
 
Cambiando
las 
generatrices…
(s.m.m. -11y a f2)
(s.m.m. -2x a f1)
CONVERGE !
pk+1=g1( pk ,qk )=
−pk
2+4 pk+qk−0 .5
2
qk+1=g2( pk ,qk )=
−pk
2−4 qk
2+11qk+4
11
 
TEORÍA Cómo determinar las generatrices de forma que la sucesión 
generada sea convergente? 
Sean las funciones g1 y g2 y sus derivadas parciales continuas en una
región del plano llamada D que contiene un punto fijo (p,q)
Si se cumple que
Entonces la sucesión
generada a partir de (p0,q0) ϵ D converge a (p,q)
(|
∂ g i
∂ x
|+|
∂ gi
∂ y
|< 1 ) ∀( x,y )∈D
pk+1=g1( pk ,qk )
qk+1=g2( pk ,qk ) k=0,1,. ..
i=1,2
 
PUNTO FIJO ATRACTIVO
(Caso bidimensional) (Caso tridimensional)
|
∂ g1
∂ x
( p,q ) |+|
∂ g1
∂ y
(p,q ) | < 1
|
∂ g2
∂ x
(p,q ) |+|
∂ g2
∂ y
( p,q )| < 1
|
∂ g1
∂ x
( p,q,r ) |+|
∂ g1
∂ y
(p,q,r ) |+|
∂ g1
∂ z
(p,q,r ) |< 1
|
∂ g2
∂ x
(p,q,r )|+|
∂ g2
∂ y
( p,q,r ) |+|
∂ g2
∂ z
( p,q,r ) | < 1
|
∂ g3
∂ x
(p,q,r )|+|
∂ g3
∂ y
( p,q,r ) |+|
∂ g3
∂ z
( p,q,r ) | < 1
Si alguna de las condiciones anteriores NO SE CUMPLE
 el proceso iterativo puede DIVERGER o CONVERGER.
 
PUNTO FIJO REPULSIVO Si se cumple que:
Entonces la iteración de punto fijo genera una sucesión divergente
(Caso bidimensional) (Caso tridimensional)
|
∂ g1
∂ x
( p,q ) |+|
∂ g1
∂ y
(p,q ) | > 1
|
∂ g2
∂ x
(p,q ) |+|
∂ g2
∂ y
( p,q )| > 1
|
∂ g1
∂ x
( p,q,r ) |+|
∂ g1
∂ y
(p,q,r ) |+|
∂ g1
∂ z
(p,q,r ) |> 1
|
∂ g2
∂ x
(p,q,r )|+|
∂ g2
∂ y
( p,q,r ) |+|
∂ g2
∂ z
( p,q,r ) | > 1
|
∂ g3
∂ x
(p,q,r )|+|
∂ g3
∂ y
( p,q,r ) |+|
∂ g3
∂ z
( p,q,r ) | > 1
 
Método iterativo de Gauss-Seidel
pk+1=g1( pk ,qk )
qk+1=g2( pk+1 ,qk ) k=0,1, . ..
En forma análoga al método iterativo de Gauss-Seidel para sistemas de 
ecuaciones lineales podemos proponer una mejora al método iterativo de Punto 
Fijo para sistemas de ecuaciones no lineales.
 
PASO 1 Evaluar la función Pk+1=g (Pk )=[
g1(Pk )
g2(Pk )
. . .
gn(Pk )
]
PASO 2 Sigo? Si, ‖Pk+1−Pk‖max< tol
‖Pk+1−Pk‖max
‖Pk+1‖max+eps
< tol
detengo el proceso y retorno el valor de Pk+1 como respuesta
ó
PASO 3 k  k + 1; vuelvo a PASO 1
PuntoFijo.mMETODO ITERATIVO DE PUNTO FIJO
IDEA
ALGORITMO
Resolver el sistema no lineal F (x) = 0 mediante una sucesión 
{Pk} que converge a la solución, a partir de un punto inicial P0. 
Evaluar la sucesión Xk+1= G(Xk) 
 
 
 
 
 
 
EJEMPLO
f 1 ( x )=x12−2 x1−x 2+0 .5=0
f 2 ( x )=x12+4 x 22−4=0
Dado el siguiente sistema de 2 ecuaciones no lineales resolver utilizando el método de 
Newton-Raphson
F (X )=[
x
12
−2 x1−x2+0 .5
x
12
+4 x
22
−4 ]; J (X )=[
2 x1−2 −1
2 x1 8 x2 ]
F (2 .00,0 .25 )=[0 .250 .25 ]; J (2.00,0 .25 )=[
2 .0 −1 .0
4 .0 2 .0 ]
Tomando como estimación inicial (2.00, 0.25)
[2 .0 −1 . 04 .0 2 .0 ] [ ΔPP1ΔPP2 ]=−[
0. 25
0 .25 ]
ΔPP=[ ΔPP1ΔPP2 ]=[
−0 .09375
0 .0625 ]
Resolviendo el sistema lineal obtenemos 
donde
Así el punto siguiente es P1=P0+ΔPP=[2 .000 .25 ]+[
−0 .09375
0 .0625 ]=[
1.90625
0 .3125 ]
de manera similar P2=[1 .9006910 .311213 ]
P3=[1.9006770 .311219 ]
 
MÉTODO DE NEWTON-RAPHSON PARA SISTEMAS DE ECUACIONES NO LINEALES
IDEA
ALGORITMO
Resolver el sistema no lineal F (x) = 0 mediante una sucesión 
{x(k)} que esperamos converja a la solución, a partir de un punto 
inicial x(0). La sucesión {x(k)} es obtenida linealizando F(x).
Evaluar la sucesión
PASO 1 Evaluar la función vectorial
PASO 2 Evaluar la matriz Jacobiana
PASO 3
Cálculo el siguiente punto x(k+1) 
resolviendo el SISTEMA LINEAL 
(e.g. para k=0):
PASO 4
detengo el proceso y retorno el valor de x(k+1) como respuesta
De lo contrario, vuelvo al paso 1 hasta que converge o hasta 
alcanzar un cierto número máximo de iteraciones permitidas
Sigo ? Si, || x(k+1)-x(k) ||
max
 < tol
 
FIN NOTAS DE CLASE
Tema: Sistemas de ecuaciones algebraicas
	Diapositiva 1
	Diapositiva 2
	Diapositiva 3
	Diapositiva 4
	Diapositiva 5
	Diapositiva 6
	Diapositiva 7
	Diapositiva 8
	Diapositiva 9
	Diapositiva 10
	Diapositiva 11
	Diapositiva 12
	Diapositiva 13
	Diapositiva 14
	Diapositiva15
	Diapositiva 16
	Diapositiva 17
	Diapositiva 18
	Diapositiva 19
	Diapositiva 20
	Diapositiva 21
	Diapositiva 22
	Diapositiva 45
	Diapositiva 46
	Diapositiva 47
	Diapositiva 48
	Diapositiva 49
	Diapositiva 50
	Diapositiva 51
	Diapositiva 52
	Diapositiva 53
	Diapositiva 54
	Diapositiva 55
	Diapositiva 56
	Diapositiva 57
	Diapositiva 58
	Diapositiva 59
	Diapositiva 60
	Diapositiva 61
	Diapositiva 62
	Diapositiva 63
	Diapositiva 64
	Diapositiva 65
	Diapositiva 66
	Diapositiva 67
	Diapositiva 68
	Diapositiva 69
	Diapositiva 70
	Diapositiva 71
	Diapositiva 72
	Diapositiva 73
	Diapositiva 74
	Diapositiva 75
	Diapositiva 76
	Diapositiva 77
	Diapositiva 78
	Diapositiva 79
	Diapositiva 80
	Diapositiva 81
	Diapositiva 82
	Diapositiva 83
	Diapositiva 87
	Diapositiva 88
	Diapositiva 89
	Diapositiva 90
	Diapositiva 91
	Diapositiva 92
	Diapositiva 93
	Diapositiva 94
	Diapositiva 95
	Diapositiva 96

Continuar navegando

Materiales relacionados