GDIDP_AMAT_1314_Practica1

1464 days ago by evalsanjuan

MÉTODOS PARA EL CÁLCULO DE LOS CEROS DE UNA FUNCIÓN

1. Métodos de dos puntos

1.1 Método de la bisección

Comencemos definiendo una función para aplicar el método y un intervalo donde sea continua y tenga un cero. Vamos a aplicar el algoritmo para calcular la raíz de f(x)=-x+cos(x) en el intervalo [0,1].

f(x)=-x+cos(x); a,b=0.0,1.0 
       

Representamos la función para ver que tiene un cero en ese intervalo. Nótese que no haría falta representarla sino comprobar que cambia de signo en los extremos del intervalo.

plot(f(x),(x,a,b),thickness=3) 
       

Cálculo paso a paso

Paso 1

Calculamos el punto medio del intervalo

c=(a+b)/2; c 
       

Vemos cuál es el signo de f en c y comparamos con el signo en los extremos

f(a),f(c),f(b) 
       

Podemos observar que f(c) tiene signo distinto a f(b), así que la raíz estará entre c y b. Tomamos ese intervalo.

a,b=c,b; a,b 
       

Paso 2

Ahora repetimos el proceso.

c=(a+b)/2; c 
       
f(a),f(c),f(b) 
       

Ahora hemos de coger el intervalo [a,c] y volver a repetir el proceso.

a,b=a,c; print[a,b]; 
       

Paso 3

c=(a+b)/2; c f(a), f(c), f(b) 
       
a,b=c,b; print[a,b]; 
       

Automatizando la selección del intervalo

En el proceso anterior, todas las instrucciones se repiten salvo la elección del intervalo siguiente. Vamos a intentar automatizar este paso.

Volvamos a los datos iniciales:

f(x)=-x+cos(x); a,b=0.0,1.0 
       

Ahora hay que calcular $c$ y elegir entre los intervalos $[a,c]$ y $[c,b]$. Podemos automatizar la selección:

c=(a+b)/2; print "c=", c if f(c)*f(a)<0: b=c else: a=c print "[a,b]=",[a,b] 
       

Para hacer un nuevo paso basta volver a ejecutar la celda anterior. Y así sucesivamente hasta realizar tantos pasos como nos indiquen.

Automatizando la repetición

Ejecutando la celda anterior, de manera sucesiva, vamos obteniendo los intervalos del método de bisección. Podemos automatizar este proceso haciendo un bucle.

for i in range(1,6): if f(c)*f(a)<0: b=c else: a=c print "[a,b]=",[a,b] c=(a+b)/2; print "c=", c 
       

 

Análisis del error: podemos estimar el error cometido mediante el cálculo de la longitud de los intervalos.

b-a 
       

Algoritmo completo

Juntándolo todo:

f(x)=-x+cos(x); a,b=0.0,1.0 print "[a,b]=",[a,b], "error=",b-a for i in range(1,6): c=(a+b)/2; print "Iteración",i,": c=", c if f(c)*f(a)<0: b=c else: a=c print "[a,b]=",[a,b],"error=",(b-a) 
       

Respresentación gráfica

Vamos a representar los pasos anteriores:

f(x)=-x+cos(x); a,b=0.0,1.0 P = plot(f, a, b, color='red',thickness=3) print "[a,b]=",[a,b], "error=", b-a n=5 for i in range(1,n+1): c=(a+b)/2; print "Iteración",i,": c=", c if f(c)*f(a)<0: b=c else: a=c print "[a,b]=",[a,b],"error=",b-a h = (P.ymax() - P.ymin())/ (1.5*n) P += line([(a,h*i), (b,h*i)],thickness=3) P += line([(a,h*i-h/4), (a,h*i+h/4)],thickness=3) P += line([(b,h*i-h/4), (b,h*i+h/4)],thickness=3) show(P) 
       

Practica: Aplica 4 pasos del método de la bisección a la función $f(x)=e^{-x}-sen(x)$ para encontrar una raíz entre 0 y 1. Calcula el error cometido.

 
       

1.2 Método de la Regula Falsi

Para aplicar el método de la Regula-Falsi, tengamos en cuenta que el procedimiento es similar al método de la bisección, sólo que el punto c se define de forma diferente. Veamos cómo resuelve el problema anterior.

a=0; b=2.0; print[a,b] 
       

Podemos definir la función que calcula el punto del intervalo [a,b], a partir de la función f.

r(x,y)=(f(y)*x-f(x)*y)/(f(y)-f(x)) 
       

Y aplicarla al intervalo [a,b]

c=r(a,b); c.n() 
       

Podemos ahora utilizar el algoritmo iterativo definido antes para el método de la bisección.

if f(c)*f(a)<0: b=c else: a=c c=r(a,b); print "[a,b]=",[a,b],"c=",c 
       

Practica: Aplica 4 pasos del método de la Regula Falsi a la función $f(x)$ para encontrar una raíz entre 0 y 1. Estima el error cometido.

f(x)= 
       

 

2. Métodos de un punto

2.1 Método de Newton-Raphson

Apliquemos el método para calcular la raíz de la función $f(x)=e^x-4$, partiendo de $x_0=2$

f(x)=e^x-4 plot(f(x),(x,0,4),thickness=3,ymin=-3,ymax=20) 
       
c=2; 
       
c=c-(f(c)/f.diff(x)(c)).n() ; c 
       

Evaluando la celda de manera iterada, vamos obteniendo las sucesivas aproximaciones. 

f(x)=e^x-4 c=2; for i in range(1,5): c=c-(f(c)/f.diff(x)(c)).n() ; c 
       

Como podemos comprobar, llegaremos al valor exacto de la raíz, que no es otra que ln4

ln(4).n() 
       

Podemos representar gráficamente el método:

f(x)=e^x-4 P=plot(f(x),(x,1,2.1),thickness=3) c=2; for i in range(1,5): P+=line([(c,f(c)),(c-(f(c)/f.diff(x)(c)),0)], color='black',thickness=2)+line([(c,0),(c,f(c))], color='red',thickness=2) c=c-(f(c)/f.diff(x)(c)) ; c.n() show(P) 
       
 
       

Análisis del error: para este método, podemos usar la fórmula $E_{n+1}\leq C E_n^2$, donde $$C=\frac{\max |f''(x)|}{2 \min |f'(x)|}$$

Vamos a tomar como intervalo inicial [1.2,1.6]. Entonces, $E_0=1.6-1.2$ pues es la longitud del intervalo inicial donde vamos a buscar la raíz.

e0=1.6-1.2;e0 
       
g(x)=f.diff(x); g 
       

Cálculo exacto del máximo y el mínimo

En este caso, el mínimo de la función $e^x$ en el intervalo $[1.2,1.6]$ es fácil de calcular, pues la función $e^x$ es creciente:

$\min |f'(x)|=e^{1.2}$

 

g(1.2) 
       
h(x)=f.diff(x,2);h 
       

El máximo también es obvio:

$\max |f''(x)| =e^{1.6}$

h(1.6) 
       

Por tanto, $C=\frac{e^{1.6}}{2e^{1.2}}$

C=h(1.6)/(2*g(1.2)); C.n() 
       

Aproximando numéricamente el máximo y el mínimo

También podríamos haber calculado el máximo o el mínimo usando los métodos de Sage para ello. El primer valor de gmin es el mínimo alcanzado por la función y el segundo valor es en qué punto lo alcanza.

gmin=g.find_local_minimum(1.2,1.6);gmin 
       
hmax=h.find_local_maximum(1.2,1.6);hmax 
       
hmax[0] 
       
C=hmax[0]/(2*gmin[0]); C 
       

Estimación del error

Entonces, para calcular estimaciones del error, haremos:

e1=C*e0^2; e1.n() 
       
e2=C*e1^2; e2.n() 
       

Y automatizando el proceso, para 5 iteraciones:

f(x)=e^x-4 x0=2; print "Punto inicial: ",x0 e0=2-1; print "Error inicial:", e0 for i in range(1,6): x0=x0-(f(x0)/f.diff(x)(x0)) ; print "Iteración ",i, ": ",x0.n() e0=C*(e0^2); print "Error:",e0 
       

Practica: Aplica el método de Newton-Raphson a la función $f(x)=e^{-x}-sen(x)$ para encontrar una raíz entre 0 y 1 con 4 cifras decimales de precisión. ¿Cuántos pasos necesitarías para tener 10 cifras decimales?

 
       

2.2 Método de la secante

En este caso, partimos de dos valores iniciales, x0 y x1.

x0=2.3 x1=2.0 
       

De nuevo obtenemos los valores iterando el método:

c=(x0*f(x1)-x1*f(x0))/(f(x1)-f(x0)); x0=x1; x1=c; x1.n() 
       

Y podemos automatizar el proceso usando un bucle:

f(x)=e^x-4 x0=2.3 x1=2.0 for i in range(1,5): c=(x0*f(x1)-x1*f(x0))/(f(x1)-f(x0)); x0=x1; x1=c; x1.n() 
       

Y también representar los valores obtenidos gráficamente:

f(x)=e^x-4 x0=2.3 x1=2.0 P=plot(f(x),(x,1,2.4),thickness=3) for i in range(1,5): c=(x0*f(x1)-x1*f(x0))/(f(x1)-f(x0)); P+=line([(x0,f(x0)),(c,0)], color='red',thickness=2) x0=x1; x1=c; x1.n() show(P) 
       

Practica: Aplica 5 pasos del método de la secante a la función $f(x)=e^{-x}-sen(x)$ para encontrar una raíz entre 0 y 1.

 
       

3. Método del punto fijo

Consideramos el problema f(x)=-3x+e^(-x)=0, partiendo de x0=1. La función para el método de punto fijo será F(x)=e^(-x)/3.

F(x)=e^(-x)/3 
       
x0=1.0; F(x0) 
       

El método consiste en aplicar la función F de forma sucesiva. Por ello, bastará con evaluar sucesivamente la siguiente celda:

x0=F(x0); x0.n() 
       

Lo podemos automatizar mediante un bucle:

F(x)=e^(-x)/3 x0=1.0 for i in range(1,5): x0=F(x0); x0.n() 
       

Y gráficamente:

F(x)=e^(-x)/3 g1=plot(F(x),(x,0,1),thickness=3)+plot(x,(x,0,1),color='red',thickness=2) g1=g1+line([(0,0),(0,1)], color='yellow',linestyle='--',thickness=3)+line([(0,1),(1,1)], color='yellow',linestyle='--',thickness=3)+line([(1,1),(1,0)], color='yellow',linestyle='--',thickness=3)+line([(1,0),(0,0)], color='yellow',linestyle='--',thickness=3) x0=1.0 for i in range(1,5): g1+=line([(x0,x0),(x0,F(x0))], color='green',thickness=2)+line([(x0,F(x0)),(F(x0),F(x0))], color='green',thickness=2) x0=F(x0); x0.n() show(g1) 
       

Podemos observar que $F(x)\in [0,1]$ para todo $x\in[0,1]$ pues la gráfica (azul) está contenida en el cuadrado amarillo. Por otra parte, vamos a ver que la derivada en valor absoluto es menor que $1$ por lo que tenemos asegurada la convergencia para todo punto inicial en $[0,1]$:

plot(F.derivative(x),(x,0,1)) 
       

Practica: Aplica 4 pasos del método del punto fijo a la función $F(x)=e^{-x/4}-sen(x/4)$ para encontrar un punto fijo partiendo de $x_0=1$. 

 
       

Ejercicios

En todos los ejercicios, comenzaremos definiendo una variable llamada dni con parte entera 0 y cifras decimales tu número de dni. Por ejemplo, si mi dni es 123456789A, definimos:

dni=0.123456789 
       

Ejercicio 1: Aplica el método de la bisección a la función $f(x)=e^{-x}-sen(x)$ para encontrar una raíz entre 1 y 2 con 5 cifras decimales. ¿Cuántos pasos son necesarios para obtener la raíz con 10 cifras decimales?

 
       

Ejercicio 2: Aplica el método de Newton-Raphson a la función $f(x)=e^{-x}-sen(x)$ para encontrar una raíz entre 0 y 1 con 4 cifras decimales de precisión. ¿Cuántos pasos necesitarías para tener 10 cifras decimales?

 
       

Ejercicio 3: Aplica el método del punto fijo para encontrar un cero de la función $f(x)=e^{-x/4}-sen(x/4)-(1+dni) x$ partiendo de $x0=0$ con un error menor de $10^-2$

 
       

Problema 1: Un objeto se mueve con ecuación de movimiento

$$(1+(t-1)^2+sin(dni\ 2\ pi  t), 1+(t-1)^3/3+cos(dni\ 2\ pi  t)), \ t\in(0,2)$$

Calcular el instante de tiempo en el que su distancia al origen es mínima (con 10 cifras de precisión). ¿Cuánto vale esa distancia (con 2 cifras de precisión)?

a) Dibujar la curva pedida usando el comando parametric_plot

b) Definir una función que en cada instante de tiempo, proporcione la distancia del objeto al origen (también es válido --incluso mejor-- usar la distancia al cuadrado).

c) Dibujar dicha función y decidir en qué intervalo estará el mínimo ("a ojo").

d) Calcular la derivada de la función anterior (cuando se anule será un mínimo o un máximo).

e) Dibujar la función derivada.

f) Aplicar un método para calcular un cero de la derivada.

g) Calcular el instante en el que la distancia es mínima

h) Calcular la mínima distancia

 
       

Rellena el cuadro de texto, imprime como pdf y envía la tarea al campus virtual.

Alumnos: