"""
Ubiquitous Newton-Raphson algorithm for solving
f(x) = 0
where a root is repeatedly estimated by
x = x - f(x)/f'(x)
until |dx|/(1+|x|) < TOL is achieved. This termination condition is a
compromise between
|dx| < TOL, if x is small
|dx|/|x| < TOL, if x is large
"""
from scipy import exp
def f(x):
return (3.0-x)*exp(x) - 3.0
def fprime(x):
return (2.0-x)*exp(x)
def newton(func, funcd, x, TOL=1e-6):
f, fd = func(x), funcd(x)
count = 0
while 1:
dx = f / float(fd)
if abs(dx) < TOL * (1 + abs(x)): return x - dx
x = x - dx
f, fd = func(x), funcd(x)
count = count + 1
print "newton(%d): x=%s, f(x)=%s" % (count, x, f)
print newton(f,fprime,5,TOL=1e-3)
newton(1): x=4.32659538633, f(x)=-103.404917649
newton(2): x=3.73937134331, f(x)=-34.1064860902
newton(3): x=3.27329578708, f(x)=-10.2145165658
newton(4): x=2.969407321, f(x)=-2.40404327646
newton(5): x=2.84210428177, f(x)=-0.291801085581
newton(6): x=2.82190151635, f(x)=-0.00638132870494
2.82143960878