Reputation: 57
I was trying to use the newton raphson method to compute the derivative of a function and I got the following error:
import numpy as np
import matplotlib.pyplot as plt
import sympy as sym
acc = 10**-4
x = sym.Symbol('x')
def p(x): #define the function
return 924*x**6 - 2772*x**5 + 3150*x**4 - 1680*x**3 + 420*x**2 - 42*x + 1
p_prime = sym.diff(p(x))
def newton(p,p_prime, acc, start_val):
x= start_val
delta = p(x)/p_prime(x)
while abs(delta) > acc:
delta = p(x)/p_prime(x)
print(abs(delta))
x =x- delta;
return round(x, acc);
a = newton(p, p_prime,-4, 0)
The error was:
delta = p(x)/p_prime(x)
TypeError: 'Add' object is not callable
Upvotes: 1
Views: 649
Reputation: 23139
There are a few mistakes in your code, correcting them as pointed out in the following modified code snippet, it works fine:
def p_prime(xval): # define as a function and substitute the value of x
return sym.diff(p(x)).subs(x, xval)
def newton(p, p_prime, prec, start_val): # rename the output precision variable
x = start_val # otherwise it's shadowing acc tolerance
delta = p(x) / p_prime(x) # call the function p_prime()
while abs(delta) > acc:
delta = p(x) / p_prime(x)
x = x - delta
return round(x, prec) # return when the while loop terminates
a = newton(p, p_prime,4, 0.1)
a
# 0.0338
The following animation shows how Newton-Raphson converges to a root of the polynomial p(x)
.
Upvotes: 1
Reputation: 80574
In sympy functions usually are represented as expressions. Therefore, sym.diff()
returns an expression and not a callable function. It makes sense to also represent p
as an expression involving x
.
To "call" a function p
on x
, p.subs(x, value)
is used. To get a numeric answer, use p.subs(x, value).evalf()
.
import sympy as sym
acc = 10 ** -4
x = sym.Symbol('x')
p = 924 * x ** 6 - 2772 * x ** 5 + 3150 * x ** 4 - 1680 * x ** 3 + 420 * x ** 2 - 42 * x + 1
p_prime = sym.diff(p)
def newton(p, p_prime, acc, start_val):
xi = start_val
delta = sym.oo
while abs(delta) > acc:
delta = (p / p_prime).subs(x, xi).evalf()
print(xi, delta)
xi -= delta
return xi
a = newton(p, p_prime, 10**-4, 0)
Intermediate values:
0 -0.0238095238095238
0.0238095238095238 -0.00876459050881560
0.0325741143183394 -0.00117130331903862
0.0337454176373780 -1.98196463560775e-5
0.0337652372837341
PS: To plot a function represented as an expression:
sym.plot(p, (x, -0.03, 1.03), ylim=(-0.5, 1.5))
Upvotes: 0
Reputation: 472
Your main problem is to call
something which is not callable
.
import sympy as sym
def f(x):
return x**2
x = sym.Symbol("x")
g = 2*x
print(callable(f)) # True
print(callable(g)) # False
print(f(0)) # print(0)
print(g(0)) # error
So, in your case
def p(x):
return 924*x**6 - 2772*x**5 + 3150*x**4 - 1680*x**3 + 420*x**2 - 42*x + 1
print(p(0)) # p is callable, gives the result 1
print(p(1)) # p is callable, gives the result 1
print(p(2)) # p is callable, gives the result 8989
print(callable(p)) # True
And now, if you use a symbolic variable
from sympy
you get:
x = sym.Symbol("x")
myp = p(x)
print(callable(myp)) # False
print(myp(0)) # gives an error, because myp is an expression, which is not callable
So, if you use diff
to get the derivative of the function, you will get a sympy expression
. You must transform this expression to a callable
function. To do it, you can use the lambda
function:
def p(x):
return 924*x**6 - 2772*x**5 + 3150*x**4 - 1680*x**3 + 420*x**2 - 42*x + 1
x = sym.Symbol("x")
myp = p(x)
mydpdx = sym.diff(myp, x) # get the derivative, once myp is an expression
dpdx = lambda x0: mydpdx.subs(x, x0) # dpdx is callable
print(dpdx(0)) # gives '-42'
Another way to make a sympy expression callable is to use lambdify from sympy:
dpdx = sym.lambdify(x, mydpdx)
Upvotes: 0