Skip to content
scipython.com

Q8.4.2Failing to find the root with scipy.optimize.newton

Using scipy.optimize.newton to find a root of the following functions (with the given starting point, \(x_0\)) fails. Explain why and find the roots either by modifying the call to newton or by using a different method.

(a)

\[f(x) = x^3 - 5x, \quad x_0 = 1\]

(b)

\[f(x) = x^3 - 3x + 1, \quad x_0 = 1\]

(c)

\[f(x) = 2-x^5, \quad x_0 = 0.01\]

(d)

\[f(x) = x^4 - (4.29)x^2 -5.29, \quad x_0 = 0.8\]
Show solution

Some examples of root-finding for which the Newton-Raphson algorithm fails and how to solve this.

(a)

In [x]: newton(lambda x: x**3 - 5*x, 1, lambda x: 3*x**2 - 5)
...
RuntimeError: Failed to converge after 50 iterations, value is 1.0

The Newton-Raphson algorithm enters an endless cycle of values for \(x\):

\[\begin{align*} x_0 = 1:& \quad x_1 = x_0 - f(x_0)/f'(x_0) = -1\\ x_1 = -1: & \quad x_2 = x_1 - f(x_1)/f'(x_1) = 1\\ x_2 = 1: & \quad x_1 = x_0 - f(x_0)/f'(x_0) = -1\\ \cdots \end{align*}\]

Alternative starting points converge correctly on a root. Even a very small displacement from \(x=0\) ensures convergence:

In [x]: newton(lambda x: x**3 - 5*x, 1.0001, lambda x: 3*x**2 - 5)
Out[x]: 2.23606797749979
In [x]: newton(lambda x: x**3 - 5*x, 1.1, lambda x: 3*x**2 - 5)
Out[x]: -2.23606797749979
In [x]: newton(lambda x: x**3 - 5*x, 0.5, lambda x: 3*x**2 - 5)
Out[x]: 0.0

(b)

In [x]: f, fp = lambda x: x**3 - 3*x+1, lambda x: 3*x**2 - 3
In [x]: newton(f, 1, fp)
Out[x]: 1.0
In [x]: f(1.0)
Out[x]: -1

The algorithm converged, but not on a root! Unfortunately, the gradient of the function is zero at the chosen starting point and because of round-off error this has not led to a ZeroDivisionError. To find the roots, choose different starting points such that \(f'(x_0) \ne 0\), or use a different method after bracketing the roots by inspection of a plot of the function:

In [x]: brentq(f, -0.5, 0.5), brentq(f, -2, -1.5), brentq(f, 1, 2)
Out[x]: (0.34729635533386066, -1.879385241571423, 1.532088886237956)

(c) The function \(f(x) = 2-x^5\) has a flat plateau around \(f(0)=2\) and the small gradient here leads to slow convergence on the root:

In [x]: newton(f, 0.01, fp)
...
RuntimeError: Failed to converge after 50 iterations, value is ...

To find it using newton either move the starting point closer to the root, or increase the maximum number of iterations:

In [x]: newton(f, 0.01, fp, maxiter=100)
Out[x]: 1.148698354997035

(d) This is another example of a function which generates an endless cycle of values from the Newton-Raphson method:

In [x]: f = lambda x: x**4 - 4.29 * x**2 - 5.29
In [x]: fp = lambda x: 4*x**3 - 8.58 * x
In [x]: newton(f, 0.8, fp)
...
RuntimeError: Failed to converge after 50 iterations, value is ...

Unlike the function in (a), the region \(0.6 \le x_0 \le 1.1\) attracts this cyclic behaviour, so one needs to initialize the algorithm outside this range to obtain the roots \(\pm 2.3\). For example,

In [x]: newton(f, 1.2, fp)
Out[x]: -2.3