1+ '''
2+ Author: P Shreyas Shetty
3+ Implementation of Newton-Raphson method for solving equations of kind
4+ f(x) = 0. It is an iterative method where solution is found by the expression
5+ x[n+1] = x[n] + f(x[n])/f'(x[n])
6+ If no solution exists, then either the solution will not be found when iteration
7+ limit is reached or the gradient f'(x[n]) approaches zero. In both cases, exception
8+ is raised. If iteration limit is reached, try increasing maxiter.
9+ '''
10+
11+ import math as m
12+
13+ def calc_derivative (f , a , h = 0.001 ):
14+ '''
15+ Calculates derivative at point a for function f using finite difference
16+ method
17+ '''
18+ return (f (a + h )- f (a - h ))/ (2 * h )
19+
20+ def newton_raphson (f , x0 = 0 , maxiter = 100 , step = 0.0001 , maxerror = 1e-6 ,logsteps = False ):
21+
22+ a = x0 #set the initial guess
23+ steps = [a ]
24+ error = abs (f (a ))
25+ f1 = lambda x :calc_derivative (f , x , h = step ) #Derivative of f(x)
26+ for _ in range (maxiter ):
27+ if f1 (a ) == 0 :
28+ raise ValueError ("No converging solution found" )
29+ a = a - f (a )/ f1 (a ) #Calculate the next estimate
30+ if logsteps :
31+ steps .append (a )
32+ error = abs (f (a ))
33+ if error < maxerror :
34+ break
35+ else :
36+ raise ValueError ("Itheration limit reached, no converging solution found" )
37+ if logsteps :
38+ #If logstep is true, then log intermediate steps
39+ return a , error , steps
40+ return a , error
41+
42+ if __name__ == '__main__' :
43+ import matplotlib .pyplot as plt
44+ f = lambda x :m .tanh (x )** 2 - m .exp (3 * x )
45+ solution , error , steps = newton_raphson (f , x0 = 10 , maxiter = 1000 , step = 1e-6 , logsteps = True )
46+ plt .plot ([abs (f (x )) for x in steps ])
47+ plt .xlabel ("step" )
48+ plt .ylabel ("error" )
49+ plt .show ()
50+ print ("solution = {%f}, error = {%f}" % (solution , error ))
0 commit comments