Universidad Tecnológica de Panamá

Facultad de Ingeniería Mecánica

Métodos Numéricos

Raíces de una función real

Ing. Fernando Castillo Balboa

Raíces de una función real
Métodos que utilizan intervalo
Bisección.
Regula Falsi (Regla Falsa).
Muller (para raíces complejas)
Métodos de intervalo abierto.
Punto fijo.
Newton-Raphson.
Sistema de ecuaciones. Jacobiano
Secante.
Raíces de una función real
El objeto del cálculo de las raíces de una ecuación es determinar los valores de x
para los que se cumple:
f(x) = 0 (1)
La determinación de las raíces de una ecuación es uno de los problemas más
antiguos en matemáticas y se han realizado un gran número de esfuerzos en este
sentido. Su importancia radica en que si podemos determinar las raíces de una
ecuación también podemos determinar máximos y mínimos, valores propios de
matrices, resolver sistemas de ecuaciones lineales y diferenciales, etc..

La determinación de las soluciones de la ecuación (1) puede llegar a ser un

problema muy difícil. Si f(x) es una función polinómica de grado 1 ó 2, se
conocen expresiones simples que permitirán determinar sus raíces. Para
polinomios de grado 3 ó 4 es necesario emplear métodos complejos y
laboriosos. Sin embargo, si f(x) es de grado mayor de cuatro o bien no es
polinómica, no hay ninguna fórmula conocida que permita determinar los ceros
de la ecuación (excepto en casos muy particulares).
Raíces de una función real
Existen una serie de reglas que pueden ayudar a determinar las raíces de una
ecuación:
– El teorema de Bolzano, que establece que si una función continua, f(x), toma
en los extremos del intervalo [a,b] valores de signo opuesto, entonces la
función admite, al menos, una raíz en dicho intervalo.
– En el caso en que f(x) sea una función algebraica (polinómica) de grado n y
coeficientes reales, podemos afirmar que tendrá n raíces reales o complejas.
– La propiedad más importante que verifican las raíces racionales de una
ecuación algebraica establece que si p/q es una raíz racional de la ecuación de
coeficientes enteros:

a0 + a1x + a2x2 + … + anxn = 0 ( ai ɛ £ )

entonces el denominador q divide al coeficientes an y el numerador p divide

al término independiente a0.
Raíces de una función real
Ejemplo: Pretendemos calcular las raíces racionales de la ecuación:
3x3 + 3x2 - x - 1 = 0
Primero es necesario efectuar un cambio de variable x = y/3:

y después multiplicamos por 32 :

y3 + 3y2 -3y -9 = 0
con lo que los candidatos a raíz del polinomio son:

Sustituyendo en la ecuación, obtenemos que la única raíz real es y = -3, es

decir, x = -3 / 3 = -1 (que es además, la única raíz racional de la ecuación).
Lógicamente, este método es muy poco potente, por lo que sólo nos puede
servir a modo de orientación.
Raíces de una función real
Método de la bisección
Es el método más elemental y antiguo para determinar las raíces de una
ecuación. Está basado directamente en el teorema de Bolzano explicado
con anterioridad. Consiste en partir de un intervalo [x0,x1] tal
que f(x0).f(x1) < 0, por lo que sabemos que existe, al menos, una raíz real.
A partir de este punto se va reduciendo el intervalo sucesivamente hasta
hacerlo tan pequeño como exija la precisión que hayamos decidido
emplear.

En la figura siguiente, se presenta el algoritmo y el diagrama de flujo del

método de bisección.
Dada la ecuación f(x) = 0, el indicador de precisión
ɛ y dos puntos a y b en los que f(a) . f(b) < 0,
 Programa del método de la Bisección: Matlab y Scilab
disp(' METODO DE LA BISECCION '); function y=f(x)
disp(' ---------------------- '); y=exp(-x^2)-x;
f=input('INGRESE FUNCION: ','s'); endfunction
xai=input('INGRESE LIMITE INFERIOR DEL INTERVALO: ');
xbi=input('INGRESE LIMITE SUPERIOR DEL INTERVALO: function xr=biseccion(xai,xbi,tol)
'); i=1;
tol=input('INGRESE PORCENTAJE DE ERROR: '); ea(1)=100;
f=inline(f); if f(xai)*f(xbi) < 0 xa(1)=xai;
i=1; xb(1)=xbi;
ea(1)=100; xr(1)=(xa(1)+xb(1))/2; p
if f(xai)*f(xbi) < 0 rintf('It.\t\t Xa\t\t Xr\t\t Xb\t Error \n');
xa(1)=xai; printf('%2d \t %11.7f \t %11.7f \t %11.7f
xb(1)=xbi; \n',i,xa(i),xr(i),xb(i));
xr(1)=(xa(1)+xb(1))/2; while abs(ea(i)) >= tol
fprintf('It. Xa Xr Xb Error aprox \n'); if f(xa(i))*f(xr(i))< 0 xa(i+1)=xa(i);
fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \n',i,xa(i),xr(i),xb(i)); xb(i+1)=xr(i);
while abs(ea(i)) >= tol, end
if f(xa(i))*f(xr(i))< 0 if f(xa(i))*f(xr(i))> 0 xa(i+1)=xr(i);
xa(i+1)=xa(i); xb(i+1)=xb(i);
xb(i+1)=xr(i); end
end xr(i+1)=(xa(i+1)+xb(i+1))/2;
if f(xa(i))*f(xr(i))> 0 ea(i+1)=abs((xr(i+1)-xr(i))/(xr(i+1))*100);
xa(i+1)=xr(i); printf('%2d \t %11.7f \t %11.7f \t %11.7f \t %7.3f
xb(i+1)=xb(i); \n',i+1,xa(i+1),xr(i+1),xb(i+1),ea(i+1)); i
end =i+1;
xr(i+1)=(xa(i+1)+xb(i+1))/2; end
ea(i+1)=abs((xr(i+1)-xr(i))/(xr(i+1))*100); else
fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \t %7.3f \n',... printf('No existe una raíz en ese intervalo');
i+1,xa(i+1),xr(i+1),xb(i+1),ea(i+1)); end
i=i+1; endfunction
end
else
fprintf('No existe una raíz en ese intervalo');
end
Raíces de una función real
Método de la bisección

El algoritmo empleado se esquematiza en la figura anterior. Inicialmente, es necesario

suministrar al programa el número máximo de iteraciones MaxIter, la tolerancia d,
que representa las cifras significativas con las que queremos obtener la solución y dos
valores de la variable independiente, x0 y x1, tales que cumplan la relación f(x0) f(x1) <
0. Una vez que se comprueba que el intervalo de partida es adecuado, lo dividimos en
dos subintervalos tales que y , y determinamos en qué subintervalo
se encuentra la raíz (comprobando de nuevo el producto de las funciones). Repetimos
el proceso hasta alcanzar la convergencia (hasta que D  d) o bien hasta que
se excede el número de iteraciones permitidas (Iter > MaxIter), en cuyo caso es
necesario imprimir un mensaje de error indicando que el método no converge.
Raíces de una función real
Método de las aproximaciones sucesivas
Dada la ecuación f(x) = 0, el método de las aproximaciones sucesivas reemplaza esta
ecuación por una equivalente, x=g(x), definida en la forma g(x)=f(x)+x. Para
encontrar la solución, partimos de un valor inicial x0 y calculamos una nueva
aproximación x1=g(x0). Reemplazamos el nuevo valor obtenido y repetimos el
proceso. Esto da lugar a una sucesión de valores , que si converge,
tendrá como límite la solución del problema.

Figura: Interpretación
geométrica del método
de las aproximaciones
sucesivas.
Raíces de una función real
Método de las aproximaciones sucesivas
En la figura anterior se representa la interpretación geométrica del método.
Partimos de un punto inicial x0 y calculamos y = g(x0). La intersección de
esta solución con la recta y=x nos dará un nuevo valor x1 más próximo a la
solución final.
Sin embargo, el método puede divergir fácilmente. Es fácil comprobar que
el método sólo podrá converger si la derivada g'(x) es menor en valor
absoluto que la unidad (que es la pendiente de la recta definida por y=x).
Un ejemplo de este caso se muestra en la figura. Esta condición, que a priori
puede considerarse una severa restricción del método, puede obviarse
fácilmente. Para ello basta elegir la función g(x) del siguiente modo:
g(x) = x + α f(x)
de forma que tomando un valor de α adecuado, siempre podemos hacer
que g(x) cumpla la condición de la derivada.
Programa de aproximaciones sucesivas o punto fijo: (Scilab)
function y=g(x)
y=(x^2-exp(x))/5;
endfunction

function x=puntofijo(x0,tol)
i=1;
ea(1)=100; x(1)=x0;
while abs(ea(i))>=tol,
x(i+1) = g(x(i));
ea(i+1) = abs((x(i+1)-x(i))/x(i+1))*100;
i=i+1;
end
printf(' i \t X(i) Error aprox (i) \n');
for j=1:i;
printf('%2d \t %11.7f \t %7.3f \n',j-1,x(j),ea(j));
end
endfunction
Ejemplo:
Aproxima la solución de cos(x) − x = 0 mediante la iteración de punto fijo para la
forma;
a) x = (x + cos(x))/2
b) x = cos(x) , a partir del valor inicial x0 = 1.
¿Cuántas iteraciones son necesarias para obtener 5 decimales exactos?
Solución: En a), la fórmula de recurrencia es:

Resulta en:
x(1)=1;
for i=1:7
x(i+1)=(x(i)+cos(x(i)))/2;
end
x'
Solución: En b), la fórmula de recurrencia es:
i xi
1
for i=1:17 2
x0 = 1, x(i+1)=cos(x(i)); 3
end 4
Xi = cos(Xi) 5
x' 6
7
1.2 8
9
1.1 10
11
1
Como se puede observar, 12
13
0.9
la solución converge de 14
0.8 forma basculante 15
16
0.7 17
18
0.6
19
20
0.5
0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 21
22
23
24
Raíces de una función real
Método de Newton-Raphson
Este método parte de una aproximación inicial xi, y obtiene una
aproximación mejor, xi+1, dada por la fórmula:

La expresión anterior puede derivarse a partir de un desarrollo en serie de

Taylor. Efectivamente, sea r un cero de f y sea x una aproximación a r tal que
r = x + h. Si f '' existe y es continua, por el teorema de Taylor tenemos:
0 = f(r) = f(x+h) = f(x) + hf'(x) + O(h2)
en donde h = r - x. Si x está próximo a r (es decir h es pequeña), es razonable
ignorar el término O(h2):
0 = f(x) + hf'(x)
Raíces de una función real
Método de Newton-Raphson (cont.)

por lo que obtenemos la siguiente expresión para h:

A partir de la ecuación anterior y teniendo en cuenta que a = x + h es

fácil derivar la ecuación

Figura: Interpretación geométrica

del método de Newton.
Raíces de una función real
Método de Newton-Raphson (cont.)
El método de Newton tiene una interpretación geométrica sencilla, como se
puede apreciar del análisis de la figura anterior. De hecho, el método de
Newton consiste en una linealización de la función, es decir, f se reemplaza
por una recta tal que contiene al punto (xi ,f(xi)) y cuya pendiente coincide
con la derivada de la función en el punto, f'(xi). La nueva aproximación a la
raíz, xi+1, se obtiene de la intersección de la función linear con el eje X de
ordenadas.
Veamos como podemos obtener la ecuación a partir de lo dicho en el párrafo
anterior. La ecuación de la recta que pasa por el punto (xi,f(xi)) y de
pendiente f'(xi) es:
y - f(xi) = f'(xi)(x-xi)
haciendo y = 0 y despejando x obtenemos la ecuación de Newton-Raphson
Raíces de una función real
Método de Newton-Raphson (cont.)
El método de Newton es muy rápido y eficiente ya que la convergencia es
de tipo cuadrático (el número de cifras significativas se duplica en cada
iteración). Sin embargo, la convergencia depende en gran medida de la
forma que adopta la función en las proximidades del punto de iteración.
En la figura se muestran dos situaciones en las que este método no es
capaz de alcanzar la convergencia (figura (a)) o bien converge hacia un
punto que no es un cero de la ecuación (figura (b)).
Para ilustrar, el método de Newton, se
resolverá el problema de mecanismos de
cuatro eslabones. La ecuación es:
R1 cos(a) – R2 cos(f) + R3 – cos(a – f)
=0
donde

haciendo

su derivada

entonces

haciendo
Sustituyendo se tiene:

Para la primera iteración se hace

Sustituyendo

Y así sucesivamente
disp(' METODO DE NEWTON-RAPHSON ');
disp(' ---------------------- ');
f=input('INGRESE FUNCION: ','s');
fd=input('INGRESE LA DERIVADA: ','s');
xai=input('INGRESE VALOR INICIAL: ');
tol=input('INGRESE EL ERROR: ');
f=inline(f);
fd=inline(fd);
i=1;
ea(1)=100;
fprintf('It. Xa f(xa) fd(xa) \n');
fprintf('%2d \t %11.7f \t %11.7f \t %11.7f
\n',i,xai,f(xai),fd(xai));
while abs(f(xai))> tol
xai=xai-f(xai)/fd(xai);
i=i+1;
fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \n'
,i,xai,f(xai),fd(xai));
if f(xai)==0
fprintf('No hay raiz')
end
end
Sistemas de Ecuaciones No Lineales:
Al igual que para una sola ecuación no lineal, el método de Newton para
sistemas de ecuaciones no lineales se basa en tomar los dos primeros
términos en la serie de Taylor de las funciones y resolver. Sea el sistema de
dos ecuaciones no lineales,
f1(x1; x2) = 0
f2(x1; x2) = 0
Sea (x1 (0) ; x2 (0)) una solución inicial aproximada. Tomamos los dos primeros
términos de la serie de Taylor de las funciones en el punto (x1 (1) ; x2 (1)) e
igualando a cero se tiene,

La ecuación anterior puede ponerse bajo forma matricial,

donde J es la matriz Jacobiana de las funciones f1 y f2,

Evidentemente es necesario que la matriz J sea no-singular.

Entonces, el método de Newton para dos ecuaciones no lineales en dos

variables consiste en comenzar a partir de una solución inicial aproximada
(x1 (0) ; x2 (0)) y aplicar la iteración siguiente:
Ejemplo de Sistemas de Ecuaciones No Lineales: Sea el sistema de
ecuaciones no lineales,

Calcular una solución del sistema tomando como valores iniciales x1 = 6 y

x2 = 5.
Solución: La matriz Jacobiana del sistema es,

Aplicando las iteraciones de Newton,

con las condiciones iniciales x1(0) = 6 y x2(0)= 5, se obtiene que
converge hacia la solución x1 = 5 y x2 = 4.

// Programa que realiza los pasos para encontrar

// las raíces de un sistema de ecuaciones
// en Scilab
x=[6; 5];
deff('[y]=f(x)','y=[ x(1)-x(2)^3+5*x(2)^2-2*x(2)-13,
x(1)+x(2)^3+x(2)^2-14*x(2)-29]');
for i=1:4,
J(1,1)=1; J(1,2)=-3*x(2)^2+10*x(2)-2;
J(2,1)=1; J(2,2)=3*x(2)^2+2*x(2)-14;
x = x - inv(J)*f(x) *
end; x
x1 x2
x = 2.6122449 4.244898
x = 4.8262832 4.0202435
Solución: x = 4.9987174 4.0001558
x = 4.9999999 4.
x = 4.9999999 4.
Resolver el problema siguiente mediante el método de
Newton para sistemas de ecuaciones

Solución: La matriz Jacobiana es:

Aplicando el método de iteraciones de Newton:

Programa en Scilab para resolver el problema raíces de
sistemas de ecuaciones no lineales
x=[1;1];
deff('[y]=f(x)','y=[0.5*sin(x(1)*x(2))-x(2)/(4*3.1415926)-0.5*x(1);
(1-(1/(4*3.1415926)))*(exp(2*x(1))-exp(1))+
exp(1)*x(2)/3.1415926 - 2*exp(1)*x(1)]');
for i=1:4
j(1,1)=x(2)*cos(x(1)*x(2))/2 - 1/2;
j(1,2)=x(1)*cos(x(1)*x(2))/2 - 1/(4*3.1415926);
j(2,1)=(1 - 2/(4*3.1415916))*2*exp(2*x(1)) - 2*exp(1);
j(2,2)=exp(1)/3.1415926;
x = x - inv(j)*f(x);
end
x

La ejecución para x=[ 0.4 ; 3.0] es:

Método de la Secante:
Un problema del método de Newton-Raphson es el de la evaluación de la
derivada. Al aproximar dicha derivada por medio de una línea secante
(derivada), en donde:

Resulta en:
 Método de la Secante: ejemplo

Usar el método de la Secante para calcular la raíz aproximada de la función

f (x) = x2 − 4 . Comenzando con x0 = 4 y x1 = 3 y hasta que |∈r |≤1% .
Aplicando para la primera iteración con la fórmula

se tendría un valor para x2 = 2.2857 . Si se calcula el error relativo con x2

como valor real y x1 como valor aproximado, se tendrá:

x=-4.5:0.1:4.5 x=-4.5:0.1:4.5
f=inline('x.^2-4','x') f = @(x) x.^2-4;
plot(x,f(x),) plot(x,f(x))
grid on grid on
 Gráfica del problema
20

f(x) = x 2- 4
15

10
ordenada

5
Raiz

-5
-5 -4 -3 -2 -1 0 1 2 3 4 5
abcisa
Raíces complejas de polinomios con constantes complejas
El Método de Muller
Suponga que se tiene un polinomio de grado n

Donde a los coeficientes ai se les conoce como los coeficientes del

polinomio f, siendo an ≠ 0. La función f(x) = 0, se asume que es un
polinomio sin grado.
• Para la ecuación de orden n, hay n raíces reales o complejas. Se
debe notar que esas raíces no son necesariamente distintas.
• Si n es impar, hay al menos una raíz real.
• Si las raíces complejas existen, existe un par conjugado.
Un predecesor del método de Muller, es el método de la secante, el cual
obtiene raíces, estimando una proyección de una línea recta en el eje x, a
través de dos valores de la función (Figura 1). El método de Muller toma un
punto de vista similar, pero proyecta una parábola a través de tres puntos
(Figura 2).
El método consiste en obtener los coeficientes de los tres puntos, sustituirlos
en la fórmula cuadrática y obtener el punto donde la parábola intercepta el eje
x. La aproximación es fácil de escribir, en forma conveniente esta sería:

f(x) Línea recta

Raíz estimada

X1 X0 X

Raíz
Figura 1
Así, se busca esta parábola para intersectar los tres puntos [x0, f(x0)], [x1, f(x1)] y [x2,
f(x2)]. Los coeficientes de la ecuación anterior se evalúan al sustituir uno de esos tres
puntos para dar:

f(x) Parábola

Raíz 0

x x

X2 X1 X0 X

Raíz estimada

Figura 2
La última ecuación genera que, f(x) = C, de esta forma, se puede tener un sistema de
dos ecuaciones con dos incógnitas:

Definiendo de esta forma:

Sustituyendo en el sistema:
Teniendo como resultado los coeficientes:

Hallando la raíz, se implementa la solución convencional, pero debido al error de

redondeo potencial, se usará una formulación alternativa:

despejando

La gran ventaja de este método es que se pueden localizar tanto las raíces reales
como las imaginarias. Hallando el error este será:

Al ser un método de aproximación, este se realiza de forma secuencial e

iterativamente, donde x1, x2, x3 reemplazan los puntos x0, x1, x2 llevando el error a
un valor cercano a cero.
Ejemplo:
 Programa Muller para obtener raíces complejas
function [res,fval,it] = muller (f,Z0,itmax,ztol,ftol,option)
z0 = Z0(1); z1 = Z0(2); z2 = Z0(3);
y0 = feval ( f, z0); y1 = feval ( f, z1);
y2 = feval ( f, z2);
for it = 1:itmax Comandos en Matlab
q = (z2 - z1)/(z1 - z0);
A = q*y2 - q*(1+q)*y1 + q^2*y0;
x=[0 0.1 1]; % Intervalo donde hay una raíz
B = (2*q + 1)*y2 - (1 + q)^2*y1 + q^2*y0; f=@(x) x.^2-x+12; % Definición de la función
C = (1 + q)*y2;
if ( A ~= 0 ) [resultado,fval,iteraciones] = muller(f,x)
disc = B^2 - 4*A*C;
den1 = ( B + sqrt ( disc ) );
den2 = ( B - sqrt ( disc ) ); Resultado
if ( abs ( den1 ) < abs ( den2 ) ) resultado = 0.5000 + 3.4278i
z3 = z2 - (z2 - z1)*(2*C/den2);
else fval = -1.2967e-013 -3.8636e-014i
z3 = z2 - (z2 - z1)*(2*C/den1); iteraciones =1
end
elseif ( B ~= 0 )
z3 = z2 - (z2 - z1)*(2*C/B);
else
warning('Muller Method failed to find a root. Last iteration result used as an output. Result may not be accurate')
res = z2;
fval = y2;
return
end
y3 = feval ( f, z3);
if opt == 1
if ( abs (z3 - z2) < ztol || abs ( y3 ) < ftol )
res = z3;
fval = y3;
return
end
else
if ( abs (z3 - z2) < ztol && abs ( y3 ) < ftol )
res = z3;
fval = y3;
return
end
end
z0 = z1; z1 = z2; z2 = z3; y0 = y1; y1 = y2; y2 = y3;
end
res = z2; fval = y2;
warning('Maximum number of iterations reached. Result may not be accurate')
Método de Bairstow