Mostrando entradas con la etiqueta método de la secante. Mostrar todas las entradas
Mostrando entradas con la etiqueta método de la secante. Mostrar todas las entradas

lunes, 28 de octubre de 2013

Métodos iterativos para resolver ecuaciones (V): el método de "regula falsi"

Esta entrada pertenece a la serie "Métodos iterativos para resolver ecuaciones" (parte I: introducción, parte II: método de Newton, parte III: método de bisección, parte IV: método de la secante).

Existe otro método iterativo para resolver numéricamente ecuaciones denominado "regula falsi". Una vez entendidos el método de bisección y el método de la secante, el método de "regula falsi" es sencillo de comprender porque se basa en ambos. Por una parte, realiza el mismo proceso de calcular el punto de corte de la secante a la función en dos puntos con el eje X. Una vez hallado el punto de corte, utiliza el criterio de los signos opuestos de las imágenes de las abscisas (del método de bisección) para decidir cuál de las dos abscisas anteriores se mantiene en la nueva iteración junto con la nueva abscisa proveniente del punto de corte calculado.


Conceptualmente no tengo mucho más que añadir, sin embargo sí que quiero compartir algunos detalles de la algoritmia del método. En un principio pensé que escribir el bucle del procedimiento sería poco más que un "copiar y pegar" de los dos procedimientos que ya están implementados en entradas anteriores (bisección y secante). Sin embargo solamente añadir el criterio de los signos opuestos de las imágenes de las abscisas al final del bucle del método de la secante no funciona y a continuación explico la razón.

La razón principal es que con este método, a diferencia del método de bisección, el intervalo en el que trabajamos en las sucesivas iteraciones no tiene por qué estrecharse tanto como se quiera para alcanzar una precisión predefinida. En la figura anterior puede observarse un caso en el que el intervalo se estrecha sólo por uno de los dos extremos. ¿Qué consecuencia tiene esto? Pues que el criterio que estábamos utilizando en métodos anteriores sobre la cercanía de los extremos de los intervalos para decidir cuándo se había alcanzado la precisión deseada no sirve.

Así que introduciendo la comprobación inicial de que las abscisas con las que empezamos tengan imágenes de signo opuesto, entramos en un bucle que se repite mientras los extremos de los intervalos disten más que la precisión definida. Dentro del bucle se calcula el punto de corte de la secante con el eje X y se determina el nuevo intervalo en función de dónde quede la alternancia de signo de las imágenes de las abscisas. Aquí se añade una comprobación más que en caso de cumplirse se sale del bucle: si dos iteraciones consecutivas producen abscisas tan cercanas como se haya definido en la precisión (¿podría esto parar el bucle en algunos casos antes de estar cerca de la solución?).

Las instrucciones en wxMaxima quedan de la siguiente manera:
x[0]:0; x[1]:1;
x[3]:x[0];
prec:10^(-4);
if p(x[0])*p(x[1]) < 0 then while abs(x[1]-x[0])>prec do (x[2]: float(x[0] - p(x[0])*(x[1]-x[0])/(p(x[1])-p(x[0]))), if p(x[2])=0 then return(x[2]), if p(x[0])*p(x[2]) < 0 then x[1]:x[2] else x[0]:x[2], if abs(x[3]-x[2])
float(x[2]);

Que en ejemplo que estamos resolviendo en esta serie

devuelve como solución .6629131443657456

jueves, 24 de octubre de 2013

Métodos iterativos para resolver ecuaciones (IV): el método de la secante

Esta entrada pertenece a la serie "Métodos iterativos para resolver ecuaciones" (parte I: introducción, parte II: método de Newton, parte III: método de bisección).

En el método de Newton que vimos en una entrada anterior calculábamos la recta tangente a la función en los sucesivos puntos de las iteraciones. Para ello necesitábamos que la función fuera derivable.

Existe un método parecido que en lugar de utilizar la recta tangente a la función en un punto calcula la recta secante a la función en dos puntos. Una ventaja es que no necesitamos que la función sea derivable.
http://www.webtag.com.ar/matcpp/images/derivada/derivada1.gif

Dicho método se conoce con el nombre de método de la secante. Entonces partiendo de dos abscisas, a y b, calculamos la recta secante a la función en dichos puntos. Es decir, buscamos la recta que pasa por los puntos (a,f(a)) y (b,f(b)).

Hay varias formas de calcular dicha recta secante. Por ejemplo, podemos hallar la pendiente de la recta a partir de los puntos por los que pasa: (f(b)-f(a))/(b-a) y utilizar la ecuación punto-pendiente:
y-f(a) = [(f(b)-f(a))/(b-a)]·(x-a)
[Nota: también podíamos haber escogido el punto (b,f(b)) y la ecuación hubiera quedado y-f(b) = [(f(b)-f(a))/(b-a)]·(x-b), pero la recta sigue siendo la misma]

Si ahora buscamos el punto de corte de la recta secante con el eje X (y=0) y despejamos la x tenemos:
-f(a) = [(f(b)-f(a))/(b-a)]·(x-a)
-f(a)/[(f(b)-f(a))/(b-a)] = (x-a)
a - f(a)/[(f(b)-f(a))/(b-a)] = x
a - f(a)·(b-a)/(f(b)-f(a)) = x

Esta x es la nueva abscisa, la llamamos c, que substituye en la próxima iteración a la abscisa a. Es decir, ahora se repetiría el proceso con los puntos (b,f(b)) y (c,f(c)).

Escrito en forma de sucesión, la iteración general queda definida, empezando con x[0]=a y x[1]=b, de la siguiente manera:
x[n]= x[n-2] - f(x[n-2])·(x[n-1]-x[n-2])/(f(x[n-1])-f(x[n-2]))
[Nota: si cogemos la otra ecuación de la recta secante la expresión general de la iteración queda x[n]= x[n-1] - f(x[n-1])·(x[n-1]-x[n-2])/(f(x[n-1])-f(x[n-2])), ambas expresiones son igual de válidas]

Vamos a utilizar wxMaxima para aplicar el método con la ecuación de la primera entrada de la serie:

Definimos el polinomio y los puntos iniciales (por ejemplo 1 y 2):
p(x):=x^5+3*x^4−2*x^3−x^2+5*x−3;
x[0]:1; x[1]:2;

Realizamos la primera iteración y mostramos la aproximación decimal del resultado:
x[2]: x[0] - p(x[0])*(x[1]-x[0])/(p(x[1])-p(x[0]));
float(%);
0.953125

Gráficamente lo que hemos hecho en la iteración anterior es:

Realizamos unas cuantas iteraciones más:
x[3]: x[1] - p(x[1])*(x[2]-x[1])/(p(x[2])-p(x[1]));
float(%);
.9144356310802643




x[4]: x[2] - p(x[2])*(x[3]-x[2])/(p(x[3])-p(x[2]));
float(%);
.7451112080395356

x[5]: x[3] - p(x[3])*(x[4]-x[3])/(p(x[4])-p(x[3]));
float(%);
.6868819480267324

x[6]: x[4] - p(x[4])*(x[5]-x[4])/(p(x[5])-p(x[4]));
float(%);
.6651550727881196

x[7]: x[5] - p(x[5])*(x[6]-x[5])/(p(x[6])-p(x[5]));
float(%);
.6629972438289449

x[8]: x[6] - p(x[6])*(x[7]-x[6])/(p(x[7])-p(x[6]));
float(%);
.6629400060916464

Podemos observar que las 3 últimas iteraciones mantienen dos cifras decimales (0.66) por lo que damos por válida dicha aproximación (por truncamiento) de la solución de la ecuación.

El proceso anterior puede automatizarse en wxMaxima de la siguiente manera:

- Definimos las abscisas iniciales (en nuestro ejemplo a=1 y b=2)
x[0]:1; x[1]:2;

- Definimos la precisión deseada en la aproximación decimal de la solución:
prec:10^(-4);

- Mientras las abscisas disten más que la precisión definida se van haciendo iteraciones, salvo que en alguna iteración se encuentre la solución exacta:
while abs(x[1]-x[0])>prec do (x[2]: float(x[0] - p(x[0])*(x[1]-x[0])/(p(x[1])-p(x[0]))), if p(x[2])=0 then return(x[2]), x[0]:x[1], x[1]:x[2]);

- Si el bucle no para con la solución de la ecuación, pedimos a wxMaxima que nos muestre la aproximación decimal de la última iteración:
float(x[2]);
.6629400060916464

Al igual que en los otros métodos, te propongo que pruebes con otros valores iniciales para encontrar las otras soluciones de la ecuación.