Métodos numéricos para
la resolución de Sistemas
de Ecuaciones no
Lineales
Contenido
Planteamiento del problema
Método de Punto Fijo
Método de Newton
Variantes del método de Newton
• Evaluación diferida del jacobiano
• Aproximación por diferencias finitas
• Newton unidimensional
Métodos cuasi-Newton (Broyden)
Notación
Escalar
f1 ( x1 , x 2 ,..., x n ) 0
f2 ( x1 , x 2 ,..., x n ) 0 fi : IR n IR
( x1 ,..., x n ) fi ( x1 ,..., x n )
f n ( x1 , x 2 ,..., x n ) 0
Vectorial
F : IR n IR n
F( x ) 0
x ( x 1 ,..., x n ) ( f 1 ( x ),... f n ( x ))
Resolución iterativa
x(0) estimación inicial de la solución
Iteraciones: x(1), x(2), …, x(k)
Criterio de convergencia
• | x(k+1) - x(k) | < tol
Criterio de parada
• k > maxiter
Esquema del algoritmo
Entrada: f, x0, tol, maxiter
Proceso
• Inicializar incr, iter
• Mientras incr > tol & iter < maxiter
– Obtener x
– incr = norm(x - x0)
– Actualizar x0, iter
Salida: x, iter, incr
• Si incr > tol no converge
Método de Punto Fijo
Punto fijo
F( x) 0 x G( x)
Estimación inicial
x( 0) ( x1( 0) ,..., x(n0) )
Iteraciones
x ( k 1) G( x ( k ) )
Criterio de parada
( k 1)
x -x (k)
tol
Algoritmo de Punto Fijo
function [x,iter,incr] = pfijo(g,x0,tol,
maxiter)
iter = 0;
incr = tol + 1;
while incr > tol & iter < maxiter
x = feval(g,x0);
incr = norm(x - x0);
iter = iter + 1;
x0 = x;
end
if incr > tol, disp(‘No converge’), end
Ejemplo
Sistema no lineal
3x 1 - cos( x 2 x 3 ) - 1 2 0
x 1 - 81( x 2 01
2
. ) sen( x 3 106
2
. 0
- x1 x 2
e 20x 3 10 / 3 - 1 0
Problema de Punto Fijo
x 1 cos( x 2 x 3 ) / 3 1 6
x 2 9 x 1 sen x 3 1.06 - 0.1
1 2
- x1 x 2
x 3 20 (1 - e
1
)-/6
Punto Fijo con desplazamientos simultáneos
x 1( k 1) cos( x (2k ) x (3k ) ) / 3 1 6
x (2k 1) 1
9 x1 sen x 3 1.06 - 0.1
(k) 2 (k)
x (3k 1) 1
20 (k) (k)
1 - exp - x 1 x 2 - / 6
Punto Fijo con desplazamientos sucesivos
x 1( k 1) cos( x (2k ) x (3k ) ) / 3 1
6
x (2k 1) 1
9 x 1 sen x 3 1.06 - 0.1
( k 1) 2 (k)
x (3k 1) 1
20
1 - exp - x 1 x 2 - / 6
( k 1) ( k 1)
Código de la función
function y=f(x)
% Función para el método de punto
% fijo con desplazamientos simultáneos
y(1) = cos(x(2)*x(3))/3 + 1/6;
y(2) = sqrt(x(1)^2+sin(x(3))+1.06)/9-0.1;
y(3) = (1-exp(-x(1)*x(2)))/20 - pi/6;
Ejemplo 1: Desp. simultáneos
Iter x1(k) x2(k) x3(k)
0 0.10000000 0.10000000 -0.10000000
1 0.49998333 0.00944115 -0.52310127
2 0.49999593 0.00002557 -0.52336331
3 0.5 0.00001234 -0.52359814
4 0.5 3.41679E-8 -0.52359847
5 0.5 1.64870 E-8 -0.52359877
Código de la función
function y=f(x)
% Función para el método de punto
% fijo con desplazamientos sucesivos
y(1) = cos(x(2)*x(3))/3 + 1/6;
y(2) = sqrt(y(1)^2+sin(x(3))+1.06)/9-0.1;
y(3) = (1-exp(-y(1)*y(2)))/20 - pi/6;
Ejemplo 1: Desp. sucesivos
Iter x1(k) x2(k) x3(k)
0 0.10000000 0.10000000 -0.10000000
1 0.49998333 0.02222979 -0.52304613
2 0.49997747 0.00002815 -0.52359807
3 0.5 3.762202E-8 -0.52359877
4 0.5 5.028E-11 -0.5235987756
Método de Newton
Sistema de ecuaciones
F : IR n IR n
F( x ) 0
x ( x 1 ,..., x n ) ( f 1 ( x ),... f n ( x ))
Aproximación por el plano tangente
F( x) F( x ( 0) ) DF( x ( 0) ) ( x - x ( 0) )
Paso de Newton
x (1) x ( 0) - DF( x ( 0) ) -1 F( x ( 0) )
Algoritmo de Newton
function [x,iter,incr] = newton(f,x,tol,
maxiter)
iter = 0; incr = tol+1;
while incr > tol & iter < maxiter
[fx,dfx] = feval(f,x);
delta = - dfx \ fx;
El archivo f.m
incr = norm(delta);
evalúa la función
iter = iter+1; y el jacobiano
x = x + delta;
end
if incr>tol, disp(‘No converge’), end
Método de Newton. Ejemplo 2
Sistema
x2 y2 - 1 0
x - y 2 0
2 2 1
Sol: x 1
2 , y 3
4
Estimación inicial x 0 1, y0 3
Primera x x
0 -1
2 2
-1 x y
iteración - 2x0 2y0
1 0 0
2x0 - 2y0
y1 y0 1
x 0 - y 0 2
2 2
Resultados Newton Ejemplo 2
k x y
0 1 3
1 0.62500000000000 1.62500000000000
2 0.51250000000000 1.04326923076923
3 0.50015243902439 0.88108161999291
4 0.50000002323057 0.86615404660332
5 0.50000000000000 0.86602541333757
6 0.50000000000000 0.86602540378444
Método de Newton. Ejemplo 3
Sistema no lineal
3x 1 - cos( x 2 x 3 ) - 1 2 0
x 1 - 81( x 2 01
2
. ) sen( x 3 106
2
. 0
- x1 x 2
e 20x 3 10 / 3 - 1 0
Jacobiana
3 x 3 sen( x 2 x 3 ) x 2 sen( x 2 x 3 )
DF( x) 2 x 1 - 162( x 2 0.1) cos( x 3 )
- x1 x 2
- x 2 e - x 1 e - x1 x 2 20
Resultados Newton. Ejemplo 3
k x1 x2 x3
0 0.10000000 0.10000000 -0.10000000
1 0.49986967 0.01946685 -0.52152047
2 0.50001423 0.00160764 -0.52313166
3 0.50000012 1.48294E-5 -0.52355872
4 0.50000000 2.08910E-8 -0.52359840
5 0.50000000 2.792E-11 -0.52359878
6 0.50000000 4.E-14 -0.52359878
Variantes de Newton (Ejercicio...)
Actualización periódica del
Jacobiano
Aproximación del Jacobiano por
diferencias divididas
Newton con desplazamiento
unidimensional
Métodos casi-Newton
Idea de la secante f ( x1 ) - f ( x 0 )
f ' ( x1 ) a 1
• No usa las x1 - x 0
derivadas (1)
f ( x )
parciales x ( 2) x -
(1)
• Convergencia a1
superlineal
Formulación DF( x (1) ) A1
matricial -1
x ( 2 ) x (1) A1 F( x (1) )
Método de Broyden
(k 1) -1
Iterar x x - Ak F(x )
(k) (k)
(yk - Ak -1sk )
Ak Ak -1
T
sk
siendo sk
2
(k -1)
yk F(x ) - F(x (k)
)
(k -1)
sk x - x (k)
Actualización de la inversa
-1
( y k - A k -1sk ) T
Ak -1
A k -1 sk
2
s k
-1 -1
-1
(s k - A y )s A
k -1 k k
T
k -1
A k -1 T -1 k 1,2,...
ks A y k -1 k
Algoritmo de Broyden
Entrada
• x0 ,tol, maxiter
Inicio
• M: Inversa del Jacobiano en x0
• x1 = x0 - M*F(x0)
• incr, iter
Iteraciones: k = 1, 2, ...
• Actualizar M % Ak-1-1 Ak-1
• xk+1 = xk - M*F(xk)
Actualización de M
w = v; % F(xk-1)
v = F(x); % F del iterado actual
y = v - w; % F(xk) - F(xk-1)
z = - M*y; % - Ak-1-1 * yk
p = - s' *z; % (sk - xk-1)T * Ak-1-1 * yk
q = s' *M; % sk T * Ak-1-1
R = (s+z)*q/p; % Transformación rango 1
M = M+R; % Inversa nueva: Ak-1
s = - M*v; % Paso de Broyden: sk+1
Algoritmo while incr > tol
w = v; % F(x(k-1))
de Broyden v = F(x);
y = v-w; % F(x(k)) - F(x(k-1))
% Inicio z = - M*y; % -inv(A(k-1))*y(k)
v = F(x0) p = - s' *z;
M = inv(DF(x0)) q = s' *M; % s(k)'*inv(A(k-1)
% Inversa Jacobiano R = (s+z)*q/p;
s = - M*v; M = M+R; % inversa de A(k)
x = x0+s; s = - M*v;
% Paso de Newton x = x+s; % Paso de Broyden
incr = norm(s); incr = norm(s);
end
Resultados de Broyden. Ejemplo 3
k x1 x2 x3
0 0.10000000 0.10000000 -0.10000000
1 0.49986967 0.01946684 -0.52152047
2 0.49998637 0.00873783 -0.52317457
3 0.50000660 0.00086727 -0.52357234
4 0.50000032 0.00003953 -0.52359768
5 0.50000000 0.00000019 -0.52359877
Alternativas al primer paso
Estimar el Jacobiano por diferencias
divididas
Estimación unidimensional del
Jacobiano
A0 diag ((F(x1 ) - F(x 0 )). /( x1 - x 0 ))
Fin
