Programa en python para implementar método de Muller para encontrar raices de una ecuación

A continuación se muestra el código paso a apaso del Programa en python para implementar método de Muller para encontrar raices de una ecuación:



#pmuller.py
#Descripcion del problema: Realizar un programa en python
#Encontrar la raiz de la función: f(x)=x^3-13*x-12, por el método de muller con valores iniciales x0=4.5, x1=5.5, x2=5
#Autor: Alexander Arias
#fecha: 2018-10-06

i=0
x0=4.5
x=x0
f0=x**3-13*x-12
x1=5.5
x=x1
f1=x**3-13*x-12
x2=5
x=x2
f2=x**3-13*x-12
h0=x1-x0
h1=x2-x1
d0=(f1-f0)/h0
d1=(f2-f1)/h1
a=(d1-d0)/(h1+h0)
b=a*h1+d1
c=f2
print("i x0     x1     x2     f0      f1      f2      h0      h1     d0      d1      a       b       c     ")
print(str(i) + " {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f}".format(x0,x1,x2,f0,f1,f2,h0,h1,d0,d1,a,b,c))
if b>0:
    x3=x2-2*c/(b+(b**2-4*c*a)**0.5)
else:
    x3=x2-2*c/(b-(b**2-4*c*a)**0.5)
print("x3 = " + "{:.4f}".format(x3))
print("---------------------------------------")

i=1
x0=x1
x=x0
f0=x**3-13*x-12
x1=x2
x=x1
f1=x**3-13*x-12
x2=x3
x=x2
f2=x**3-13*x-12
h0=x1-x0
h1=x2-x1
d0=(f1-f0)/h0
d1=(f2-f1)/h1
a=(d1-d0)/(h1+h0)
b=a*h1+d1
c=f2
print("i x0     x1     x2     f0      f1      f2      h0      h1     d0      d1      a       b       c     ")
print(str(i) + " {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f}".format(x0,x1,x2,f0,f1,f2,h0,h1,d0,d1,a,b,c))
if b>0:
    x3=x2-2*c/(b+(b**2-4*c*a)**0.5)
else:
    x3=x2-2*c/(b-(b**2-4*c*a)**0.5)
print("x3 = " + "{:.4f}".format(x3))

print("---------------------------------------")

i=2
x0=x1
x=x0
f0=x**3-13*x-12
x1=x2
x=x1
f1=x**3-13*x-12
x2=x3
x=x2
f2=x**3-13*x-12
h0=x1-x0
h1=x2-x1
d0=(f1-f0)/h0
d1=(f2-f1)/h1
a=(d1-d0)/(h1+h0)
b=a*h1+d1
c=f2
print("i x0     x1     x2     f0      f1      f2      h0      h1     d0      d1      a       b       c     ")
print(str(i) + " {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f}".format(x0,x1,x2,f0,f1,f2,h0,h1,d0,d1,a,b,c))
if b>0:
    x3=x2-2*c/(b+(b**2-4*c*a)**0.5)
else:
    x3=x2-2*c/(b-(b**2-4*c*a)**0.5)
print("x3 = " + "{:.4f}".format(x3))
print("---------------------------------------")

i=3
x0=x1
x=x0
f0=x**3-13*x-12
x1=x2
x=x1
f1=x**3-13*x-12
x2=x3
x=x2
f2=x**3-13*x-12
h0=x1-x0
h1=x2-x1
d0=(f1-f0)/h0
d1=(f2-f1)/h1
a=(d1-d0)/(h1+h0)
b=a*h1+d1
c=f2
print("i x0     x1     x2     f0      f1      f2      h0      h1     d0      d1      a       b       c     ")
print(str(i) + " {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f} {:.4f}".format(x0,x1,x2,f0,f1,f2,h0,h1,d0,d1,a,b,c))
if b>0:
    x3=x2-2*c/(b+(b**2-4*c*a)**0.5)
else:
    x3=x2-2*c/(b-(b**2-4*c*a)**0.5)
print("x3 = " + "{:.4f}".format(x3))

El resultado es:




i x0 x1 x2 f0 f1 f2 h0 h1 d0 d1 a b c 0 4.5000 5.5000 5.0000 20.6250 82.8750 48.0000 1.0000 -0.5000 62.2500 69.7500 15.0000 62.2500 48.0000 x3 = 3.9765 --------------------------------------- i x0 x1 x2 f0 f1 f2 h0 h1 d0 d1 a b c 1 5.5000 5.0000 3.9765 82.8750 48.0000 -0.8163 -0.5000 -1.0235 69.7500 47.6949 14.4765 32.8780 -0.8163 x3 = 4.0011 --------------------------------------- i x0 x1 x2 f0 f1 f2 h0 h1 d0 d1 a b c 2 5.0000 3.9765 4.0011 48.0000 -0.8163 0.0368 -1.0235 0.0246 47.6949 34.7310 12.9775 35.0498 0.0368 x3 = 4.0000 --------------------------------------- i x0 x1 x2 f0 f1 f2 h0 h1 d0 d1 a b c 3 3.9765 4.0011 4.0000 -0.8163 0.0368 0.0000 0.0246 -0.0010 34.7310 35.0126 11.9775 35.0000 0.0000 x3 = 4.0000