Newton method for complex roots

Consider the function $$ f(x) = x^5 + 1 $$ The roots lie on the unit circle in the complex plane.

In [2]:
from numpy import abs
In [3]:
def f(x):
    return x**5 + 1.0

def df(x):
    return 5.0*x**4
In [4]:
def newton(f,df,x0,M=100,eps=1.0e-15):
    x = x0
    for i in range(M):
        dx = - f(x)/df(x)
        x  = x + dx
        print i,x,abs(f(x))
        if abs(dx) < eps * abs(x):
            return x
In [5]:
x0 = 1.0 + 1.0j
x = newton(f,df,x0)
0 (0.85+0.8j) 1.48449260687
1 (0.786945048682+0.653022872974j) 0.358786849698
2 (0.800009825693+0.58871686636j) 0.044669919092
3 (0.8091271136+0.587660739882j) 0.000831133887574
4 (0.809016956688+0.587785211881j) 2.76287012166e-07
5 (0.809016994375+0.587785252292j) 3.01410727401e-14
6 (0.809016994375+0.587785252292j) 2.22044604925e-16
7 (0.809016994375+0.587785252292j) 2.22044604925e-16