Newton’s Method
Contents
Newton’s Method#
The goal of this notebook is to show some examples of using Python to implement Newton’s Method, which attempts to approximate a root of a given function.
Mathematical Introduction to Newton’s Method#
Newton’s Method is a beautiful application of calculus, but for purposes of time, it is a topic which is often skipped over in calculus courses. Here is a video introducing Newton’s Method, by Yasmeen Baki.
Applying at a good point#
Newton’s Method is what’s called an iterative method, in that we start with a value, and then iteratively generate new values from the previous value. How well Newton’s Method works depends greatly on the value we start with. In this section, we will consider an example where Newton’s Method works well.
Let \(f(x) = x^3 - 2x + 2\). Starting at \(x_0 = -1.5\), apply Newton’s method three times to estimate a root of \(f\).
import matplotlib.pyplot as plt
import numpy as np
The goal of Newton’s Method is to approximate a root of a given function \(f\). To apply Newton’s Method, we need both \(f\) and the derivative of \(f\). In this case, we are just going to explicitly write out the derivative of \(f\).
f = lambda x: x**3 - 2*x + 2
df = lambda x: 3*x**2 - 2
The iterative formula for Newton’s Method is \(z \leadsto z - f(z)/f'(z)\). At each stage of Newton’s Method, we have some approximate root \(z\), and for the next stage we will use \(z - f(z)/f'(z)\) as our approximate root. We define this function in the next cell.
newt = lambda z: z - f(z)/df(z)
We want to apply Newton’s method three times, starting at \(x_0 = -1.5\). We will hold the successive values in a length-4 NumPy array arr. The zeroth element in arr will be -1.5.
z = -1.5
reps = 3
arr = np.zeros(reps+1)
arr[0] = z
The zeros in arr are just place-holders.
arr
array([-1.5, 0. , 0. , 0. ])
Because we have already defined newt above, it is very easy to fill in the subsequent values produced by Newton’s Method. (I don’t see a way to find these numbers using, for example, list comprehension in place of the for loop.)
for i in range(reps):
arr[i+1] = newt(arr[i])
In good situations, the values produced by Newton’s Method will be getting gradually closer to an actual root of \(f\).
arr
array([-1.5 , -1.84210526, -1.77282692, -1.76930129])
As a reality check, let’s see if newt(arr[2]) is actually arr[3].
newt(-1.77282692)
-1.7693012925504912
To get a sense for what Newton’s Method is doing, let’s plot \(y = f(x)\). We will use one of the pre-defined styles mentioned in a previous section.
plt.style.use('seaborn-darkgrid')
Looking at the graph, it is clear that this cubic polynomial has only one root, and it looks to be at approximately \(x \approx -1.8\). You should imagine starting at \((-1.5, f(-1.5))\), drawing the tangent line to \(y = f(x)\) at that point, and finding the root of the tangent line (i.e., finding where the tangent line crosses the \(x\)-axis). According to our implementation above, the tangent line should cross the \(x\)-axis at approximately \(x = -1.842\).
fig, ax = plt.subplots()
x = np.linspace(-2,2,1000)
ax.plot(x, f(x))
[<matplotlib.lines.Line2D at 0x1118ca2c0>]
Let’s actually apply f to the elements in arr and see if the outputs are getting closer and closer to 0. (One of the reasons we made arr as a NumPy array, as opposed to a list, is so that we could make this evaluation using f(arr). That syntax would not work if arr were a list.) In this case, f(arr[3]) is \(-0.00006605...\), so by the third iteration of Newton’s Method, we have indeed gotten very close to an actual root of \(f\).
f(arr)
array([ 1.62500000e+00, -5.66700685e-01, -2.61909903e-02, -6.60651520e-05])
Applying at a worse point#
As mentioned above, the success of Newton’s Method depends heavily on the starting point. In this section, we will see some of what can go wrong if we do not begin at a good point.
Let \(f(x) = x^3 - 2x + 2\). Starting at \(x_0 = -0.5\), apply Newton’s method three times to estimate a root of \(f\).
import numpy as np
import matplotlib.pyplot as plt
plt.style.use('seaborn-darkgrid')
We will use the same setup as in the previous ssection.
f = lambda x: x**3 - 2*x + 2
df = lambda x: 3*x**2 - 2
newt = lambda z: z - f(z)/df(z)
x = np.linspace(-2,2,1000)
fig, ax = plt.subplots()
ax.plot(x, f(x))
[<matplotlib.lines.Line2D at 0x1114a84c0>]
Again, this code is all the same as in the previous section, with the exception of the fact that we start here at -0.5 instead of at -1.5.
z = -0.5
reps = 3
arr = np.zeros(reps+1)
arr[0] = z
for i in range(reps):
arr[i+1] = newt(arr[i])
arr
array([-0.5 , 1.8 , 1.25181347, 0.71203147])
As a reality check, we make sure that newt(arr[1]) really is arr[2].
newt(1.8)
1.2518134715025906
We will get a much better sense of Newton’s Method if we can draw tangent lines. The following function takes as input an \(x\)-value x0 and as output returns the function corresponding to the tangent line of f at x0.
def tang(x0):
return lambda x: df(x0)*(x-x0)+f(x0)
For each of our three iterations of Newton’s method, we draw the corresponding tangent line using ax.plot(x, tang(arr[i])(x)) and we draw a bullet point at the corresponding root estimate (i.e., where the tangent line crosses the \(x\)-axis) using ax.plot(arr[i+1],0,'ko'). The syntax 'ko', standing for “black circle” marker, might look familiar from Matlab. Instead of letting Matplotlib automatically choose the ranges for x and y, we set them explicitly using ax.set(xlim=(-1,2), ylim=(-2,3.5)).
It’s worth staring at the resulting plot and understanding how what it shows relates to the values in arr.
x = np.linspace(-2,2,1000)
fig, ax = plt.subplots()
ax.plot(x, f(x))
z = -0.5
reps = 3
arr = np.zeros(reps+1)
arr[0] = z
for i in range(reps):
ax.plot(x, tang(arr[i])(x))
arr[i+1] = newt(arr[i])
ax.plot(arr[i+1],0,'ko')
ax.set(xlim=(-1,2), ylim=(-2,3.5))
[(-1.0, 2.0), (-2.0, 3.5)]
Let’s use the exact same code but with 10 iterations of Newton’s Method. Notice how we don’t seem to be converging on the an actual root, even after 10 iterations.
x = np.linspace(-2,2,1000)
fig, ax = plt.subplots()
ax.plot(x, f(x))
z = -0.5
reps = 10
arr = np.zeros(reps+1)
arr[0] = z
for i in range(reps):
ax.plot(x, tang(arr[i])(x))
arr[i+1] = newt(arr[i])
ax.plot(arr[i+1],0,'ko')
ax.set(xlim=(-1,2), ylim=(-2,3.5))
[(-1.0, 2.0), (-2.0, 3.5)]
On the other hand, if we choose a better starting point like -1.5, then we very quickly converge on the actual root. The following chart includes 10 iterations of Newton’s Method, but most of these iterations are so close to each other that we can’t distinguish between them.
x = np.linspace(-2,2,1000)
fig, ax = plt.subplots()
ax.plot(x, f(x))
z = -1.5
reps = 10
arr = np.zeros(reps+1)
arr[0] = z
for i in range(reps):
ax.plot(x, tang(arr[i])(x))
arr[i+1] = newt(arr[i])
ax.plot(arr[i+1],0,'ko')
ax.set(xlim=(-2,-1), ylim=(-2,3.5))
[(-2.0, -1.0), (-2.0, 3.5)]
If we instead start at 0 for this function \(f\), we find another example of what can go wrong when applying Newton’s Method. In this case, Newton’s Method oscillates forever between 0 and 1, neither of which is anywhere near the actual root.
x = np.linspace(-2,2,1000)
fig, ax = plt.subplots()
ax.plot(x, f(x))
z = 0
reps = 10
arr = np.zeros(reps+1)
arr[0] = z
for i in range(reps):
ax.plot(x, tang(arr[i])(x))
arr[i+1] = newt(arr[i])
ax.plot(arr[i+1],0,'ko')
ax.set(xlim=(-2,2), ylim=(-2,3.5))
[(-2.0, 2.0), (-2.0, 3.5)]
arr
array([0., 1., 0., 1., 0., 1., 0., 1., 0., 1., 0.])
f(arr)
array([2., 1., 2., 1., 2., 1., 2., 1., 2., 1., 2.])
Analyzing Newton’s method#
We’ve seen one starting value where Newton’s Method works well, and we’ve seen two starting values where Newton’s Method does not work well. Here we are going to investigate Newton’s Method more systematically.
Let \(f(x) = x^3 - 2x + 2\). For each value of \(z\) in
np.linspace(-2,2,100), apply Newton’s method 7 times, and call the result \(z_{7}\). For how many of these 100 values do we have \(-0.001 \leq f(z_{7}) \leq 0.001\)?
The initial setup is the same as in the previous sections.
import numpy as np
import matplotlib.pyplot as plt
plt.style.use('seaborn-darkgrid')
f = lambda x: x**3 - 2*x + 2
df = lambda x: 3*x**2 - 2
newt = lambda z: z - f(z)/df(z)
Next we make the array that will hold the 100 starting values we will consider.
z = np.linspace(-2,2,100)
We will make a \(100 \times 8\) NumPy array to hold the results. We think of each row as holding a starting value as well as the first 7 applications of Newton’s Method to that starting value. We initialize all the values to 0.
reps = 7
arr = np.zeros((len(z), reps+1))
We fill in the zeroth column with the starting values, and then we one-by-one fill in the next columns using a for loop. This code is very similar to what we did above, except instead of applying newt to a single entry, here we are applying newt to an entire column.
arr[:, 0] = z
for i in range(reps):
arr[:, i+1] = newt(arr[:, i])
Here is what arr looks like. For example, the top row of arr is saying that, if we start at \(z_0 = -2\) and apply Newton’s Method 7 times, we reach approximately \(z_7 = -1.76929\).
arr
array([[-2.00000000e+00, -1.80000000e+00, -1.76994819e+00,
-1.76929266e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.95959596e+00, -1.79093206e+00, -1.76962110e+00,
-1.76929243e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.91919192e+00, -1.78321559e+00, -1.76942955e+00,
-1.76929237e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.87878788e+00, -1.77700676e+00, -1.76933475e+00,
-1.76929236e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.83838384e+00, -1.77248656e+00, -1.76929966e+00,
-1.76929235e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.79797980e+00, -1.76986592e+00, -1.76929259e+00,
-1.76929235e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.75757576e+00, -1.76939218e+00, -1.76929236e+00,
-1.76929235e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.71717172e+00, -1.77135720e+00, -1.76929541e+00,
-1.76929235e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.67676768e+00, -1.77610789e+00, -1.76932547e+00,
-1.76929236e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.63636364e+00, -1.78405978e+00, -1.76944655e+00,
-1.76929237e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.59595960e+00, -1.79571479e+00, -1.76978006e+00,
-1.76929252e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.55555556e+00, -1.81168492e+00, -1.77052745e+00,
-1.76929345e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.51515152e+00, -1.83272408e+00, -1.77199981e+00,
-1.76929760e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.47474747e+00, -1.85977167e+00, -1.77465685e+00,
-1.76931290e+00, -1.76929235e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.43434343e+00, -1.89401326e+00, -1.77915850e+00,
-1.76936154e+00, -1.76929236e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.39393939e+00, -1.93696679e+00, -1.78643431e+00,
-1.76949961e+00, -1.76929239e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.35353535e+00, -1.99060765e+00, -1.79777892e+00,
-1.76985803e+00, -1.76929258e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.31313131e+00, -2.05755477e+00, -1.81498754e+00,
-1.77072260e+00, -1.76929382e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.27272727e+00, -2.14135575e+00, -1.84055561e+00,
-1.77268324e+00, -1.76930058e+00, -1.76929235e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.23232323e+00, -2.24693804e+00, -1.87798643e+00,
-1.77690004e+00, -1.76933359e+00, -1.76929236e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.19191919e+00, -2.38135163e+00, -1.93228936e+00,
-1.78555984e+00, -1.76947918e+00, -1.76929238e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.15151515e+00, -2.55504910e+00, -2.01083181e+00,
-1.80264288e+00, -1.77006384e+00, -1.76929278e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.11111111e+00, -2.78421900e+00, -2.12488859e+00,
-1.83521470e+00, -1.77220938e+00, -1.76929845e+00,
-1.76929235e+00, -1.76929235e+00],
[-1.07070707e+00, -3.09534448e+00, -2.29267033e+00,
-1.89571317e+00, -1.77941316e+00, -1.76936513e+00,
-1.76929236e+00, -1.76929235e+00],
[-1.03030303e+00, -3.53493070e+00, -2.54579403e+00,
-2.00645071e+00, -1.80156357e+00, -1.77001551e+00,
-1.76929273e+00, -1.76929235e+00],
[-9.89898990e-01, -4.19283168e+00, -2.94481740e+00,
-2.20998294e+00, -1.86429691e+00, -1.77518108e+00,
-1.76931710e+00, -1.76929235e+00],
[-9.49494949e-01, -5.26809655e+00, -3.62311891e+00,
-2.59814796e+00, -2.03148836e+00, -1.80791637e+00,
-1.77032155e+00, -1.76929311e+00],
[-9.09090909e-01, -7.30721003e+00, -4.94570843e+00,
-3.41754056e+00, -2.47681584e+00, -1.97445304e+00,
-1.79411694e+00, -1.76972357e+00],
[-8.68686869e-01, -1.25489622e+01, -8.40579361e+00,
-5.66676469e+00, -3.87913658e+00, -2.75233274e+00,
-2.10844708e+00, -1.83003047e+00],
[-8.28282828e-01, -5.39311359e+01, -3.59625627e+01,
-2.39879225e+01, -1.60116576e+01, -1.07048754e+01,
-7.18419619e+00, -4.86522360e+00],
[-7.87878788e-01, 2.16214141e+01, 1.44334331e+01,
9.64996988e+00, 6.47249122e+00, 4.36860038e+00,
2.98162239e+00, 2.06782443e+00],
[-7.47474747e-01, 8.75499163e+00, 5.87909731e+00,
3.97681505e+00, 2.72387826e+00, 1.89646945e+00,
1.32445388e+00, 8.11227019e-01],
[-7.07070707e-01, 5.41234046e+00, 3.66896797e+00,
2.52132147e+00, 1.76065060e+00, 1.22137662e+00,
6.64168077e-01, 2.08979749e+00],
[-6.66666667e-01, 3.88888889e+00, 2.66603410e+00,
1.85781432e+00, 1.29564866e+00, 7.74023344e-01,
5.29225028e+00, 3.58981161e+00],
[-6.26262626e-01, 3.02561426e+00, 2.09696269e+00,
1.46909445e+00, 9.70187358e-01, -2.10728426e-01,
1.08138869e+00, 3.50851023e-01],
[-5.85858586e-01, 2.47567367e+00, 1.72983568e+00,
1.19714636e+00, 6.22490415e-01, 1.81199535e+00,
1.26098987e+00, 7.25624810e-01],
[-5.45454545e-01, 2.09905020e+00, 1.47056765e+00,
9.71633849e-01, -1.98763933e-01, 1.07134095e+00,
3.18232191e-01, 1.14111631e+00],
[-5.05050505e-01, 1.82839634e+00, 1.27346525e+00,
7.43555739e-01, 3.45020712e+00, 2.37727063e+00,
1.66306528e+00, 1.14323583e+00],
[-4.64646465e-01, 1.62731129e+00, 1.11343008e+00,
4.42477131e-01, 1.29313611e+00, 7.70657255e-01,
4.96922487e+00, 3.37699055e+00],
[-4.24242424e-01, 1.47440442e+00, 9.75388643e-01,
-1.68663359e-01, 1.04958485e+00, 2.39488421e-01,
1.07910154e+00, 3.43610662e-01],
[-3.83838384e-01, 1.35628849e+00, 8.49733800e-01,
-4.65204983e+00, -3.23172447e+00, -2.36957020e+00,
-1.92727627e+00, -1.78464424e+00],
[-3.43434343e-01, 1.26416389e+00, 7.30245184e-01,
3.05123116e+00, 2.11391907e+00, 1.48104525e+00,
9.81845963e-01, -1.19900269e-01],
[-3.03030303e-01, 1.19201601e+00, 6.13193205e-01,
1.76479672e+00, 1.22460920e+00, 6.69472450e-01,
2.13587312e+00, 1.49646726e+00],
[-2.62626263e-01, 1.13560207e+00, 4.97077616e-01,
1.39373967e+00, 8.92143634e-01, -1.49537969e+00,
-1.84514226e+00, -1.77311649e+00],
[-2.22222222e-01, 1.09185185e+00, 3.82690537e-01,
1.20969840e+00, 6.44519621e-01, 1.94290031e+00,
1.35859771e+00, 8.52433243e-01],
[-1.81818182e-01, 1.05849802e+00, 2.73219078e-01,
1.10312474e+00, 4.14836015e-01, 1.25172274e+00,
7.11895268e-01, 2.66553223e+00],
[-1.41414141e-01, 1.03384003e+00, 1.74051236e-01,
1.04208023e+00, 2.09298936e-01, 1.06051703e+00,
2.80563528e-01, 1.10884021e+00],
[-1.01010101e-01, 1.01658906e+00, 9.19650586e-02,
1.01206157e+00, 6.82750213e-02, 1.00672095e+00,
3.90185922e-02, 1.00222936e+00],
[-6.06060606e-02, 1.00576401e+00, 3.36178027e-02,
1.00166006e+00, 9.87840470e-03, 1.00014543e+00,
8.71955853e-04, 1.00000114e+00],
[-2.02020202e-02, 1.00062081e+00, 3.71332121e-03,
1.00002063e+00, 1.23781368e-04, 1.00000002e+00,
1.37885051e-07, 1.00000000e+00],
[ 2.02020202e-02, 1.00060431e+00, 3.61492542e-03,
1.00001955e+00, 1.17316569e-04, 1.00000002e+00,
1.23858900e-07, 1.00000000e+00],
[ 6.06060606e-02, 1.00531632e+00, 3.10739766e-02,
1.00142044e+00, 8.46258044e-03, 1.00010683e+00,
6.40627796e-04, 1.00000062e+00],
[ 1.01010101e-01, 1.01449580e+00, 8.11339231e-02,
1.00943313e+00, 5.40602026e-02, 1.00424437e+00,
2.49380548e-02, 1.00091821e+00],
[ 1.41414141e-01, 1.02800913e+00, 1.47645749e-01,
1.03047689e+00, 1.58977231e-01, 1.03522821e+00,
1.80153278e-01, 1.04502798e+00],
[ 1.81818182e-01, 1.04584980e+00, 2.24678939e-01,
1.06965321e+00, 3.12539390e-01, 1.13590519e+00,
4.97782800e-01, 1.39524085e+00],
[ 2.22222222e-01, 1.06814815e+00, 3.07408149e-01,
1.13131327e+00, 4.86987776e-01, 1.37289512e+00,
8.68889872e-01, -2.59722977e+00],
[ 2.62626263e-01, 1.09519343e+00, 3.92441041e-01,
1.22181852e+00, 6.64896525e-01, 2.09594177e+00,
1.46837377e+00, 9.69478766e-01],
[ 3.03030303e-01, 1.12747281e+00, 4.77769474e-01,
1.35483044e+00, 8.48023211e-01, -4.95648203e+00,
-3.42438586e+00, -2.48081406e+00],
[ 3.43434343e-01, 1.16573570e+00, 5.62556971e-01,
1.56477448e+00, 1.05933958e+00, 2.76293132e-01,
1.10549512e+00, 4.21334452e-01],
[ 3.83838384e-01, 1.21109854e+00, 6.46916878e-01,
1.95908288e+00, 1.37039301e+00, 8.66040962e-01,
-2.80266122e+00, -2.13447158e+00],
[ 4.24242424e-01, 1.26521822e+00, 7.31770980e-01,
3.09068759e+00, 2.14002146e+00, 1.49937511e+00,
9.99407793e-01, -3.56379859e-03],
[ 4.64646465e-01, 1.33058821e+00, 8.18845144e-01,
-7.82770391e+01, -5.21904800e+01, -3.48024160e+01,
-2.32149389e+01, -1.54970330e+01],
[ 5.05050505e-01, 1.41106821e+00, 9.10869311e-01,
-9.98948796e-01, -4.01903504e+00, -2.83775191e+00,
-2.15284808e+00, -1.84436813e+00],
[ 5.45454545e-01, 1.51289009e+00, 1.01212504e+00,
6.86139803e-02, 1.00678672e+00, 3.93879565e-02,
1.00227130e+00, 1.34749066e-02],
[ 5.85858586e-01, 1.64672381e+00, 1.12970318e+00,
4.83144541e-01, 1.36525448e+00, 8.60149412e-01,
-3.31202551e+00, -2.41559829e+00],
[ 6.26262626e-01, 1.83237831e+00, 1.27648215e+00,
7.47802127e-01, 3.60959356e+00, 2.48223775e+00,
1.73427248e+00, 1.20065819e+00],
[ 6.66666667e-01, 2.11111111e+00, 1.47906865e+00,
9.79929417e-01, -1.33997171e-01, 1.03015086e+00,
1.57493532e-01, 1.03458666e+00],
[ 7.07070707e-01, 2.58521156e+00, 1.80363793e+00,
1.25460126e+00, 7.16198954e-01, 2.74355396e+00,
1.90959817e+00, 1.33415600e+00],
[ 7.47474747e-01, 3.59661518e+00, 2.47369348e+00,
1.72849684e+00, 1.19608506e+00, 6.20581591e-01,
1.80196241e+00, 1.25331781e+00],
[ 7.87878788e-01, 7.41858586e+00, 4.99410616e+00,
3.39337826e+00, 2.33982171e+00, 1.63751143e+00,
1.12200826e+00, 4.64340933e-01],
[ 8.28282828e-01, -1.48478114e+01, -9.93159828e+00,
-6.67292539e+00, -4.53143285e+00, -3.15588288e+00,
-2.32659456e+00, -1.90938038e+00],
[ 8.68686869e-01, -2.61113064e+00, -2.03778989e+00,
-1.80958417e+00, -1.77041046e+00, -1.76929325e+00,
-1.76929235e+00, -1.76929235e+00],
[ 9.09090909e-01, -1.03761755e+00, -3.44266052e+00,
-2.49150308e+00, -1.98116568e+00, -1.79561389e+00,
-1.76977639e+00, -1.76929252e+00],
[ 9.49494949e-01, -4.08706234e-01, 1.42542737e+00,
9.26007184e-01, -7.19546860e-01, 6.14447620e+00,
4.15197471e+00, 2.83910547e+00],
[ 9.89898990e-01, -6.38458417e-02, 1.00641391e+00,
3.72911176e-02, 1.00203833e+00, 1.21067392e-02,
1.00021813e+00, 1.30737305e-03],
[ 1.03030303e+00, 1.58186516e-01, 1.03488558e+00,
1.78653530e-01, 1.04429414e+00, 2.18385786e-01,
1.06583278e+00, 2.99411544e-01],
[ 1.07070707e+00, 3.16101609e-01, 1.13915147e+00,
5.05270632e-01, 1.41155797e+00, 9.11391114e-01,
-9.87872052e-01, -4.23436798e+00],
[ 1.11111111e+00, 4.36392915e-01, 1.28355057e+00,
7.57620974e-01, 4.06524968e+00, 2.78205426e+00,
1.93526117e+00, 1.35301155e+00],
[ 1.15151515e+00, 5.32764976e-01, 1.47808842e+00,
9.78977278e-01, -1.41115928e-01, 1.03368695e+00,
1.73374096e-01, 1.04175933e+00],
[ 1.19191919e+00, 6.13016042e-01, 1.76393659e+00,
1.22393911e+00, 6.68377506e-01, 2.12610270e+00,
1.48961078e+00, 9.90099232e-01],
[ 1.23232323e+00, 6.81908103e-01, 2.25755078e+00,
1.58103819e+00, 1.07367449e+00, 3.25998976e-01,
1.14842883e+00, 5.26047693e-01],
[ 1.27272727e+00, 7.42511823e-01, 3.41379640e+00,
2.35327848e+00, 1.64670393e+00, 1.12968662e+00,
4.83104863e-01, 1.36517643e+00],
[ 1.31313131e+00, 7.96894925e-01, 1.04123609e+01,
6.97833523e+00, 4.70291660e+00, 3.20163989e+00,
2.21333913e+00, 1.55047066e+00],
[ 1.35353535e+00, 8.46499785e-01, -5.25675646e+00,
-3.61586348e+00, -2.59382463e+00, -2.02939786e+00,
-1.80736909e+00, -1.77029315e+00],
[ 1.39393939e+00, 8.92362474e-01, -1.48818753e+00,
-1.85004035e+00, -1.77360559e+00, -1.76930565e+00,
-1.76929235e+00, -1.76929235e+00],
[ 1.43434343e+00, 9.35245831e-01, -5.83137473e-01,
2.44587023e+00, 1.70966541e+00, 1.18107689e+00,
5.92755621e-01, 1.67398490e+00],
[ 1.47474747e+00, 9.75723500e-01, -1.66043710e-01,
1.04791525e+00, 2.32919980e-01, 1.07483082e+00,
3.29803648e-01, 1.15209851e+00],
[ 1.51515152e+00, 1.01423479e+00, 7.97688608e-02,
1.00912412e+00, 5.23658902e-02, 1.00398608e+00,
2.34499722e-02, 1.00081263e+00],
[ 1.55555556e+00, 1.05112154e+00, 2.45462063e-01,
1.08309815e+00, 3.56194707e-01, 1.17922943e+00,
5.89217371e-01, 1.65980849e+00],
[ 1.59595960e+00, 1.08665348e+00, 3.67127699e-01,
1.19138467e+00, 6.12036749e-01, 1.75920727e+00,
1.22024983e+00, 6.62305782e-01],
[ 1.63636364e+00, 1.12104608e+00, 4.61937263e-01,
1.32578507e+00, 8.12888734e-01, 5.24903769e+01,
3.50018118e+01, 2.33467014e+01],
[ 1.67676768e+00, 1.15447341e+00, 5.39116520e-01,
1.49514654e+00, 9.95389365e-01, -2.83180399e-02,
1.00122705e+00, 7.31744098e-03],
[ 1.71717172e+00, 1.18707749e+00, 6.04072885e-01,
1.72226184e+00, 1.19113269e+00, 6.11574414e-01,
1.75698919e+00, 1.21851669e+00],
[ 1.75757576e+00, 1.21897520e+00, 6.60190718e-01,
2.05721771e+00, 1.44093685e+00, 9.42002486e-01,
-4.95680390e-01, 1.77652335e+00],
[ 1.79797980e+00, 1.25026354e+00, 7.09699829e-01,
2.62810251e+00, 1.83241422e+00, 1.27650934e+00,
7.47840244e-01, 3.61111261e+00],
[ 1.83838384e+00, 1.28102360e+00, 7.54131332e-01,
3.88701543e+00, 2.66480155e+00, 1.85698965e+00,
1.29503001e+00, 7.73196448e-01],
[ 1.87878788e+00, 1.31132363e+00, 7.94572930e-01,
9.40622918e+00, 6.31083609e+00, 4.26182442e+00,
2.91137196e+00, 2.02123769e+00],
[ 1.91919192e+00, 1.34122141e+00, 8.31819615e-01,
-1.12032491e+01, -7.51405549e+00, -5.08117412e+00,
-3.50374256e+00, -2.52738502e+00],
[ 1.95959596e+00, 1.37076616e+00, 8.66466648e-01,
-2.77048628e+00, -2.11778748e+00, -1.83295710e+00,
-1.77201911e+00, -1.76929768e+00],
[ 2.00000000e+00, 1.40000000e+00, 8.98969072e-01,
-1.28877933e+00, -2.10576730e+00, -1.82919995e+00,
-1.77171581e+00, -1.76929656e+00]])
We are interested in the results after the 7 applications of Newton’s Method. This corresponds to the last column of arr. We store this last column with the variable name results.
results = arr[:, -1]
Let’s get ready to plot the results. First we plot \(f\) over a slightly larger domain.
fig, ax = plt.subplots()
x = np.linspace(-2.5,2.5,1000)
ax.plot(x, f(x))
[<matplotlib.lines.Line2D at 0x1133eb8b0>]
Here is the Boolean array corresponding to our condition \(-0.001 \leq f(z_{7}) \leq 0.001\). For example, the zeroth entry being True indicates that \(z_0 = -2\) is a starting value for which this condition is satisfied.
good_bool = np.abs(f(results)) <= 0.001
good_bool
array([ True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, True,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, True,
True, False, False, False, False, False, False, False, False,
False, False, False, True, False, False, False, False, False,
False, False, False, False, False, False, False, False, True,
True])
We are getting ready to plot the good starting values. As a first step, we use Boolean indexing to select the starting values where the condition is satisfied: pts = z[good_bool].
fig, ax = plt.subplots()
x = np.linspace(-2.5,2.5,1000)
ax.plot(x, f(x))
pts = z[good_bool]
This variable pts represents a one-dimensional NumPy array of length 32.
pts.shape
(32,)
Here we try to add in black bullet points along the \(x\)-axis at all of the good points. The error message tells us that the error is in the last line ax.plot(pts, 0, 'ko');, and that the specific error is “ValueError: x and y must have same first dimension, but have shapes (32,) and (1,)”.
fig, ax = plt.subplots()
x = np.linspace(-2.5,2.5,1000)
ax.plot(x, f(x))
pts = z[good_bool]
ax.plot(pts, 0, 'ko');
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
Input In [13], in <cell line: 7>()
3 ax.plot(x, f(x))
5 pts = z[good_bool]
----> 7 ax.plot(pts, 0, 'ko')
File ~/miniconda3/envs/math9f22/lib/python3.10/site-packages/matplotlib/axes/_axes.py:1632, in Axes.plot(self, scalex, scaley, data, *args, **kwargs)
1390 """
1391 Plot y versus x as lines and/or markers.
1392
(...)
1629 (``'green'``) or hex strings (``'#008000'``).
1630 """
1631 kwargs = cbook.normalize_kwargs(kwargs, mlines.Line2D)
-> 1632 lines = [*self._get_lines(*args, data=data, **kwargs)]
1633 for line in lines:
1634 self.add_line(line)
File ~/miniconda3/envs/math9f22/lib/python3.10/site-packages/matplotlib/axes/_base.py:312, in _process_plot_var_args.__call__(self, data, *args, **kwargs)
310 this += args[0],
311 args = args[1:]
--> 312 yield from self._plot_args(this, kwargs)
File ~/miniconda3/envs/math9f22/lib/python3.10/site-packages/matplotlib/axes/_base.py:498, in _process_plot_var_args._plot_args(self, tup, kwargs, return_kwargs)
495 self.axes.yaxis.update_units(y)
497 if x.shape[0] != y.shape[0]:
--> 498 raise ValueError(f"x and y must have same first dimension, but "
499 f"have shapes {x.shape} and {y.shape}")
500 if x.ndim > 2 or y.ndim > 2:
501 raise ValueError(f"x and y can be no greater than 2D, but have "
502 f"shapes {x.shape} and {y.shape}")
ValueError: x and y must have same first dimension, but have shapes (32,) and (1,)
We need to replace the single number 0 with a NumPy array of 32 zeros. (A length 32 list would also work.)
fig, ax = plt.subplots()
x = np.linspace(-2.5,2.5,1000)
ax.plot(x, f(x))
pts = z[good_bool]
ax.plot(pts, np.zeros(len(pts)), 'ko');
We make the bullet points a little smaller by setting a markersize. (I just guessed-and-checked values to see what size looked good.)
fig, ax = plt.subplots()
x = np.linspace(-2.5,2.5,1000)
ax.plot(x, f(x))
pts = z[good_bool]
ax.plot(pts, np.zeros(len(pts)), 'ko', markersize=2);
Let’s make a similar plot for the bad points. Here is a reminder of how good_bool looked.
good_bool
array([ True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, True,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, True,
True, False, False, False, False, False, False, False, False,
False, False, False, True, False, False, False, False, False,
False, False, False, False, False, False, False, False, True,
True])
To take the elementwise negation of this Boolean array in Numpy, we use the tilde symbol ~.
~good_bool
array([False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, True,
True, True, True, True, True, True, True, True, False,
False, True, True, True, True, True, True, True, True,
True, True, True, False, True, True, True, True, True,
True, True, True, True, True, True, True, True, False,
False])
Here we use Boolean indexing with ~good_bool to keep only the bad values from z. We also change the color from black to red, by changing k to r.
fig, ax = plt.subplots()
x = np.linspace(-2.5,2.5,1000)
ax.plot(x, f(x))
pts = z[~good_bool]
ax.plot(pts, np.zeros(len(pts)), 'ro', markersize=2);
Here is the same plot from above with the good values.
fig, ax = plt.subplots()
x = np.linspace(-2.5,2.5,1000)
ax.plot(x, f(x))
pts = z[good_bool]
ax.plot(pts, np.zeros(len(pts)), 'ko', markersize=2);
The plots say that, initially, all the starting values near -2 work well. Once we pass the local maximum of \(f\), after which Newton’s Method starts moving us away from the root, we have mostly bad values. There are a few outliers (good points surrounded by bad points). It’s worth looking at these outliers and thinking about the corresponding tangent lines. Can you see why Newton’s Method could work well for these particular values? If we used more than 7 repetitions, we would see many more “good points” scattered throughout the bad points.