Tomemos o exemplo de encontrar a solução de $x^2-5=0$. Existem muitas formas de transformar este problema num problema de ponto fixo:
Vamos retomar o algoritmo do ponto fixo da aula passada para ver a convergência destes sistemas
def PontoFixo(g,x0,epsilon,N):
xold,xnew = x0,g(x0)
iteracao=1
while abs(xnew-xold)>epsilon and iteracao <= N:
print("{0:<4} ** {1:6.7f} ** {2} ".format(iteracao,xnew,abs(xnew-xold)))
xold,xnew = xnew,g(xnew)
iteracao+=1
return xnew
Acima temos quatro funções para achar o ponto fixo:
g1 = lambda x : 5+x-x**2 # função 1.
g2 = lambda x : 5/x # função 2
g3 = lambda x : 1+x-(x**2)/5 # função 3
g4 = lambda x : 0.5*(x+5/x) # função 4
#1
PontoFixo(g1,2.5,0.0001,10)
1 ** 1.2500000 ** 1.25 2 ** 4.6875000 ** 3.4375 3 ** -12.2851562 ** 16.97265625 4 ** -158.2102203 ** 145.92506408691406 5 ** -25183.6840394 ** 25025.473819054896 6 ** -634243120.4799571 ** 634217936.7959177 7 ** -402264336510396544.0000000 ** 4.022643358761534e+17 8 ** -161816596428149550818404783450226688.0000000 ** 1.6181659642814955e+35 9 ** -26184610879590621838544711141932033442920917644155845961599972955652096.0000000 ** 2.6184610879590622e+70 10 ** -685633846915575549629699855410106399266881391154855323356420383329712846863964030212182519807694982522014869944060618863533103659069871226880.0000000 ** 6.8563384691557555e+140
-4.700937720362509e+281
PontoFixo(g2,2.5,0.0001,10)
1 ** 2.0000000 ** 0.5 2 ** 2.5000000 ** 0.5 3 ** 2.0000000 ** 0.5 4 ** 2.5000000 ** 0.5 5 ** 2.0000000 ** 0.5 6 ** 2.5000000 ** 0.5 7 ** 2.0000000 ** 0.5 8 ** 2.5000000 ** 0.5 9 ** 2.0000000 ** 0.5 10 ** 2.5000000 ** 0.5
2.0
PontoFixo(g3,2.5,0.0001,10)
1 ** 2.2500000 ** 0.25 2 ** 2.2375000 ** 0.012500000000000178 3 ** 2.2362187 ** 0.0012812499999999005 4 ** 2.2360839 ** 0.00013485957031234008
2.2360696574218544
PontoFixo(g4,2.5,0.0001,10)
1 ** 2.2500000 ** 0.25 2 ** 2.2361111 ** 0.01388888888888884
2.2360679779158037
Lembremos que: $$|x_{n+1} - x_f| = \Phi^\prime(\xi_n)|x_n-x_f)|$$ daí $$\lim_{n\to\infty} \frac{|x_{n+1}-x_f|}{|x_n-x_f|} = \Phi^\prime(x_f)$$ Então se $0<|\Phi^\prime(x_f)|<1 $, a ordem da convergência é linear.
Para uma primeira estimativa do erro, lembramos que, por hipótese, $\max|\Phi^\prime(x)|=\lambda \lt 1 $, e $$|x_0 - x_f| = |x0 -x_1 +x_1 -x_f|\leq |x_1-x_0| + |x_1-x_f|\leq |x_1-x_0| + \lambda|x0-x_f| \implies |x_0-x_f|\leq \frac{1}{1-\lambda}|x_1-x_0|$$ Da última aula temos que $$ |x_n-x_f|\leq \frac{\lambda^n}{1-\lambda}|x_1-x_0|$$ Na estimativa de Aitken temos o erro estimado como: $$ |x_n-x_f| \cong \frac{\lambda_n}{1-\lambda_n}|x_n-x_{n-1}|$$ com $$ \lambda_n = \frac{|x_n-x_{n-1}|}{|x_{n-1}-x_{n-2}|}$$
%matplotlib inline
from matplotlib.pyplot import plot
import numpy as np
def Coweb(g,x0,N,a=0,b=1):
""" Traça o diagrama de Verhulst ou Coweb de g no intervalo [a,b] """
x = arange(a,b,(b-a)/100)
plot(x,x,x,g(x))
grid()
# construção dos vetores de web
y0=[]
y1=[]
xnew = x0
for i in range(N):
y0.append(xnew)
y1.append(g(xnew))
y0.append(g(xnew))
y1.append(g(xnew))
xnew = g(xnew)
if xnew > b or xnew < a:
raise ValueError(" Valor de g(x) fora do intervalo original ")
plot(y0,y1)
return
Coweb(g3,2.1,5,2,2.4)
Coweb(g2,2.1,5,1,5)
Coweb(g1,2.1,5,2,2.4)
--------------------------------------------------------------------------- ValueError Traceback (most recent call last) <ipython-input-26-f8ff9de74235> in <module>() ----> 1 Coweb(g1,2.1,5,2,2.4) <ipython-input-22-fa557f4d58ca> in Coweb(g, x0, N, a, b) 19 xnew = g(xnew) 20 if xnew > b or xnew < a: ---> 21 raise ValueError(" Valor de g(x) fora do intervalo original ") 22 plot(y0,y1) 23 return ValueError: Valor de g(x) fora do intervalo original
Coweb(g4,2.1,5,2,2.4)
%pastbin '02