MAT 281 - Real Newton Basins

3060 days ago by kcrisman

So you want to know what happens with Newton's method? 

In particular, what guesses will lead to which roots?

Here are a number of useful things to engage with that over the real numbers.


First, evaluate the following cell of preliminaries.  Then scroll ahead to some interactive stuff.

def newt(func): vari = func.args()[0] fprime = diff(func,vari) return vari-func/fprime class RootFinder: def __init__(self, roots, newtf): self.newtf = newtf self.roots = roots self.cutoff = .1 def which_root(self,Varia): varia = Varia for counter in [0..19]: varia = self.newtf(varia) for root in self.roots: if (varia-root).abs()<self.cutoff: return root ls = [] for root in self.roots: ls.append((varia-root).abs()) return self.roots[ls.index(min(ls))] 

This first cell allows you to experiment with different functions and what zeros they converge to.  

  • $f$ is the function; you might want to have a good idea of what its zeros look like ahead of time.
  • "iterations" is the number of times you want Newton's method to go. 
  • "guess" is your initial guess.
@interact def _(f=x^3-x,iterations=(3,[1..20]),guess=1.5): f(x) = f L = [(guess,0)] for i in range(iterations): L.append((guess,f(guess))) guess = newt(f)(guess) L.append((guess,0)) myxmax = max([l[0] for l in L])+.5 myxmin = min([l[0] for l in L])-.5 myymax = max([l[1] for l in L])+.5 myymin = min([l[1] for l in L])-.5 P = plot(f,(x,myxmin,myxmax),ymax=myymax,ymin=myymin) Q = line(L,color='red') html("Root approximation is $%s$"%L[-1][0]) show(P+Q) 

Click to the left again to hide and once more to show the dynamic interactive window

Try to figure out what starting values lead to which zeros.  I start with $x^3-x$ because it's a fun one to explore, but you should try for others as well.

In particular, try to see what intervals of starting guesses will all lead to the same zeros.  Write down some of your explorations!


Our next interact is different in style.  

@interact def _(roots = '-1,0,1',detail=.01,xmin=-2,xmax=2,auto_update=False): roots = eval('['+roots+']') colors = rainbow(len(roots)) f(x) = prod([x-r for r in roots]) html("The function is $f(x) = %s$"%latex(f(x).expand())) html("And the roots are ") for i in range(len(roots)): html("$%s$ in <font color=%s>this color</font>"%(roots[i],colors[i])) P = plot(f,(x,xmin,xmax),color='black') Finder = RootFinder(roots,newt(f)) f1 = Finder.which_root for t in srange(xmin,xmax,detail): P += point((t,0),color=colors[roots.index(f1(t))]) show(P) 

Click to the left again to hide and once more to show the dynamic interactive window

Here, we are using color to identify which initial guesses will end up near which roots.  The parameters are

  • The roots of the polynomial, as a list
  • The detail - how far apart we should look at guesses
  • The window we want to look in on on the $x$-axis

Notice that we don't want to ask Sage to guess what the roots are; instead, we create the polynomial by telling Sage what its roots are ahead of time.


Questions you may enjoy asking:

  • Can we guarantee anything about what root is found?
  • Are there any polynomials which have strange behavior when you zoom in?
  • Are there any polynomials which do not have strange behavior?
  • Does the behavior depend on the number of roots?

Finally, I'd ask what happens if you add in a couple of complex roots.  Try the list "-1,i,-i" for kicks.  What function does that give?  Is the result surprising?