Newton Raphson method equation solver algorithm

Question

In the code below, when I choose for example "max_n_iterations" to be equal to 1, the list "approximations", when printed, displays two elements where it should only display one (the initial x).

What is the reason for this?

#This exercise shows an immediate way to find the root of a real valued funciton, using   successive better approximations
#This method is known as Newton Raphson method

print 'Find the root of a given function - NEWTON RAPHSONS METHOD'
print 'The function is the following: ...'
x=input('Choose an initial estimate:') #An initial estimate is necessary to be choosen   to start the iteration process
x=float(x) #The number inserted is transformed into a floating point number
max_n_iterations=input('Choose max n. iterations:') #The user decides the maximum number   of iterations to run
approximations = [] #Vector collecting all the intermediate solutions; i.e. the  approximated roots evaluated before reaching the final solution 
iterations= []

def f(x):  #The given function has to be inserted manually in the code
return 2*x**3+45*x+1/x**2+16

def D(f):  #Evaluates the first derivative of the given function, using the definition of derivative 
def df(x, h): #x is the initial estimate, h is the increment
    return (f(x+h) - f(x))/h #Difference quotient
return df #First derivative

def newtons_method(f, x, h=0.000001, epsilon=0.000001,): #This is the main process: f is  the given function, x and h are the same as above, epsilon is the tolerance level. Epsilon  and h have to be choosen sufficiently small
df = D(f) #df is just the first derivative, as before
for i in range(max_n_iterations): #The code runs until the maximum number of iterations  is reached
    x1 = x - f(x)/df(x, h) #Essence of Newton Raphson method: the iteration process
    approximations.append(x) #Every intermediate solution is collected into a vector, as  defined above
    iterations.append(1)
    if abs(x1 - x) < epsilon: #When the absolute difference between two successive roots  is less than the tolerance level (i.e. the sequence of solutions strongly converges), the  program exists the cycle
        break
    x = x1 #The next solution becomes the starting solution in the new cycle
return x #Final solution

def final_solution(): #The final solution from the Newton Raphson method
return newtons_method(f,x) 

df=D(f) #These values have to be inserted again to allow the execution of the final step
h=0.000001
epsilon=0.000001
x=newtons_method(f,x)

if abs((x-f(x)/df(x,h))-x) < epsilon: #If (strong) convergence has been reached
   print 'Solution is:', final_solution() #Prints the final solution
   print 'Approximations before reaching convergence:', approximations #Prints the vector of intermediate solutions
   print 'Convergence has been reached after', len(iterations), 'iterations'
   print 'Newton Raphson method was successful!'
elif abs((x-f(x)/df(x,h))-x) >= epsilon: #If (strong) convergence has not been reached
   print 'Approximated solution is:', final_solution()
   print 'Approximations evaluated:', approximations
   print 'Convergence has not been reached after', max_n_iterations, 'iterations'
   print 'Newton Raphson method was not successful'

A good place to start would be Carmack's fast inverse square root, which does 2 (hard wired) iterations over an initial estimation for a function. Its as short a NR implementation can get, and should give you the general idea to expand into a more precise solution. Link: en.wikipedia.org/wiki/Fast_inverse_square_root — bconstanzo
Thanks but my intention was to use the Newton Raphson method. I know it's possible to use this method but I'm just stuck on that loop — user3957049
Ah you posted some code, let me read it and then I can give a more proper response! — bconstanzo
I'm just leaving my workplace and couldn't get a clear read of this, I'm getting an itchy feeling from the rounding (why dot it?!?!). I'll check this later at home, but I think someone else might help you in between. — bconstanzo

Alex Riley Alex Riley · Accepted Answer · 2014-08-19T16:42:27

It looks like before newtons_method() can return a value, it has to call itself again. For example:

def newtons_method(f, x, h=0.000001, epsilon=0.000001,):
    ...
    if abs((x - f(x)/df(x, h))-x)< epsilon:
       print 'Solution is:', round(newtons_method(f,x),6) # function called again here
       ...
       return x

    else:
       print 'Approximated solution is:', round(newtons_method(f,x),4) # and again here
       ...
       return x

So the first call to newtons_method() never returns because it has to call itself before return, and then that function call has to call itself before return, and then...

Can you modify your code so that newtons_method() isn't called recursively in this way?

Newton Raphson method equation solver algorithm

2 Answers