FilesExpand file tree

 Write-for-current-Python

/

newton_raphson.py

Copy path

Blame

More file actions

Blame

More file actions

Latest commit

 

History

History

History

54 lines (46 loc) · 1.78 KB

 Write-for-current-Python

/

newton_raphson.py

Copy path

Top

File metadata and controls

• Code
• Blame

54 lines (46 loc) · 1.78 KB

Raw

Copy raw file

Download raw file

Open symbols panel

Edit and raw actions

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

"""

Author: P Shreyas Shetty

Implementation of Newton-Raphson method for solving equations of kind

f(x) = 0. It is an iterative method where solution is found by the expression

x[n+1] = x[n] + f(x[n])/f'(x[n])

If no solution exists, then either the solution will not be found when iteration

limit is reached or the gradient f'(x[n]) approaches zero. In both cases, exception

is raised. If iteration limit is reached, try increasing maxiter.

"""

import math as m

def calc_derivative(f, a, h=0.001):

"""

Calculates derivative at point a for function f using finite difference

method

"""

return (f(a + h) - f(a - h)) / (2 * h)

def newton_raphson(f, x0=0, maxiter=100, step=0.0001, maxerror=1e-6, logsteps=False):

a = x0 # set the initial guess

steps = [a]

error = abs(f(a))

f1 = lambda x: calc_derivative(f, x, h=step) # noqa: E731 Derivative of f(x)

for _ in range(maxiter):

if f1(a) == 0:

raise ValueError("No converging solution found")

a = a - f(a) / f1(a) # Calculate the next estimate

if logsteps:

steps.append(a)

if error < maxerror:

break

else:

raise ValueError("Iteration limit reached, no converging solution found")

if logsteps:

# If logstep is true, then log intermediate steps

return a, error, steps

return a, error

if __name__ == "__main__":

from matplotlib import pyplot as plt

f = lambda x: m.tanh(x) ** 2 - m.exp(3 * x) # noqa: E731

solution, error, steps = newton_raphson(

f, x0=10, maxiter=1000, step=1e-6, logsteps=True

)

plt.plot([abs(f(x)) for x in steps])

plt.xlabel("step")

plt.ylabel("error")

plt.show()

print(f"solution = {{{solution:f}}}, error = {{{error:f}}}")