Elusive Zeros under Newton’s Method

Though  well-known for its simplicity and efficiency, Newton’s method applied to a  complex polynomial can fail quite miserably, even on a relatively large open  set of initial guesses. In this work, we present some analytic and numerical  results for Newton’s method applied to the complex quartic family  where  is a parameter. The  symmetric location of the roots of allows for some easy  reductions. In particular, when λ is either real or  purely imaginary, standard techniques from real dynamical systems theory can be  employed for rigorous analysis. Classifying those λ-values where Newton’s method fails on an open set leads to  complex and aesthetically intriguing geometry in the λ-parameter plane, complete with  fractal-like figures such as Mandelbrot-like sets, tricorns and swallows.

section_content: One of the most common iterative algorithms for finding solutions to an equation is Newton's method. Given an equation ( ) 0 f x = and an initial guess 0 x, Newton's method attempts to locate a better approximation, 1 x, given by ( ) ( ) ( )

Here, the numerical technique uses information about the first derivative of f at 0 x to obtain an improved approximation 1

x. The process then begins anew from 1 x, generating a sequence of numbers 0 1 2,,, x x x  intended to converge to a solution of the equation. 

Interestingly, though perhaps of little surprise to those familiar with iterative algorithms, numerical methods such as this have a tendency for failing in unpredictable manners. For example, applying Newton's method to find the roots of ( ) 4 2 6 11 p z z z = -leads to an entire region of the complex plane for which initial seeds eventually bounce back and forth between 1 and 1 -, neither of which is a solution to ( ) 0 p z =. In this case, 1 z = and 1 z =lie on a super-attracting cycle of period two for the map p N. It is this failure of Newton's method to converge on an open set of initial guesses that is investigated in this work. 

The notion of studying the failure of Newton's method applied to a complex polynomial dates back to the pioneering work of Curry, Garnett and Sullivan in 1983 [1]. Focusing on Newton's method applied to a particular family of cubic polynomials, the authors amazingly discover one of the most famous fractals of all, the Mandelbrot set, lurking throughout the parameter plane. Each parameter value in these special sets corresponds to a "bad" polynomial where Newton's method fails on an open region in the complex plane due to the existence of an extraneous attracting cycle distinct from the roots. Similar results in the case of a complex cubic were later obtained by Blanchard [2], Head [3], Lei [4], Roberts and Horgan-Kobelski [5] (also verifying the phenomenon for Halley's method), and Haeseler and Kriete [6] who applied quasiconformal surgery to prove the existence of rogue attractors for relaxed Newton's method. The existence of Mandelbrot-like sets in the parameter plane and the connection with quadratic-like dynamics was thoroughly explained by Douady and Hubbard using their theory of polynomial-like mappings [7]. 

In this work, we study Newton's method applied to the complex quartic family ( ) ( ) ( ) ( ) (

)

where λ ∈  is a complex parameter. The symmetric placement of the four roots is motivated by the nearest-root principle, the notion that initial seeds "typically" converge to the closest root. This is precisely the case for Newton's method applied to a complex quadratic map with two distinct roots 1 r and 2 r. The only points that fail to converge to either root under iteration lie on the perpendicular bisector l of the segment joining 1 r and 2 r. These initial seeds cannot decide which root to converge toward and consequently, remain on l for all time. For any other initial guess 0 z, Newton's method converges to the root located on the same side of l as 0 z. In general, the invariance of l occurs whenever there is a line of symmetry amid the configuration of roots. For the family p λ, the real axis is a line of symmetry and is therefore invariant under Newton's method. This will play a key role in studying the dynamics as it allows for a reduction to a map of one real variable. The primary difference between the cubic and quartic cases is the number of "free" critical points. In addition to the roots of a polynomial, its inflection points also turn out to be critical points of the Newton map. Unlike the fixed roots, these points can iterate freely around the complex plane. There is only one such point for a cubic map, but two free critical points for the quartic case since the second derivative is quadratic. Not surprisingly, due to the identical numbers of critical points, our problem turns out to have many similarities with the general complex cubic. In [8], Milnor classifies the types of hyperbolic components possible in the parameter plane for the general cubic, obtaining fractal-like figures such as Mandelbrot-like sets and tricorns as well as swallow and product configurations. The distinguishing element between these fractal-like sets is the type and behavior of the critical points. Following Milnor's work, we are able to locate and explain the existence of these same types of fractals in the parameter plane for Newton's method applied to p λ. 

The main tool in our analysis is the reduction to a map in one real variable. This occurs in the cases where λ is real or purely imaginary due to the invariance and symmetry described above. We apply standard arguments from real dynamical systems theory to prove that there are no extraneous attracting cycles in the case λ ∈ . Using a bisection algorithm, we numerically locate an abundance of values β such that Newton's method applied to i p β contains a super-attracting n -cycle. Bifurcations are explored as β varies. Whether the free critical points lie on the same or distinct periodic cycles for these special parameter values has important consequences for the resulting figures in the λ -parameter plane. f z z =. We refer to 0 z as the initial seed of the orbit.
section_name: Introduction

section_content: A point 0 z is said to be a periodic point if ( ) 0 0 n f z z = for some n ∈ , and the smallest such n is known as the period of 0 z. In this case, we say that 0 z lies on an n -cycle or periodic orbit of period n. A periodic point of Period 1 is known as a fixed point. 

A periodic point 0 z with period n is said to be attracting if ( ) ( )

It is straight-forward to show that seeds in a sufficiently small neighborhood of an attracting periodic orbit are attracted to that orbit under iteration. A periodic point 0 z with period n is said to be

title that corresponds to the rate at which nearby points converge. 

Definition 1 Suppose that O is an attracting periodic orbit of period n We now restrict our attention to studying Newton's method applied to a complex polynomial p. For a fixed p, this produces a rational map of the extended complex plane  denoted p N and given by ( ) ( ) ( )

Studying the convergence of Newton's method is equivalent to investigating the orbits of initial seeds under iteration of the map p N, placing our study squarely in the field of complex dynamical systems. There are two complimentary sets used to describe the dynamics of a map in complex dynamics, the Julia set, where the interesting and chaotic behavior occurs, and its tame cousin, the Fatou set, where attracting periodic cycles and their basins of attraction lie. The Julia set for a rational map is the closure of the set of repelling periodic points [9]. This is an invariant, perfect and fractal-like set displaying sensitive dependence on initial conditions. Any neighborhood of a point in the Julia set is mapped under iteration to cover all of the extended complex plane except at most two points. For Newton's method, this implies that arbitrarily close to any point in the Julia set are seeds that will iterate to each root of p. Thus, choosing an initial seed in the Julia set of p N is not a major difficulty since a small perturbation will ensure convergence toward a root. 

The roots of ( ) p z and their basins of attraction are in the Fatou set. It is easy to see that if ( )

so that the root 0 z is a fixed point of the dynamical system p N. Moreover, a short calculation gives

If 0 z is a simple root of p, then it follows that ( )

and thus 0 z is a super-attracting fixed point. This is certainly a good property for a numerical root-finding algorithm to possess as it implies that nearby points are strongly attracted towards a root of the polynomial. However, equation (1. 1) is significant for another reason, for it indicates that inflection points of p are also critical points of p N. Since these inflection points are typically not roots of p (and thus not fixed), we will refer to them as the free critical points of p N. In complex dynamics, it is the orbit of the critical point that governs the behavior of the underlying dynamics. In particular, we have the following important theorem of Fatou and Julia (see [9] [10]) for rational maps on  : Theorem 1 Every attracting cycle of a rational map attracts at least one critical point. In the case of Newton's method, the critical points corresponding to roots are themselves attractors. Generically, we also expect the free critical points to be attracted towards the roots. However, an intriguing situation occurs when one or more of the free critical points converges to an extraneous attracting cycle, that is, a cycle distinct from one of the roots of the polynomial. The basin of attraction for such a cycle is an entire region in C for which initial seeds never converge to a root. In this case, the roots elude detection, for a small perturbation of a failing initial seed may not lead to convergence to one of the roots. The dynamical plane for such an example is displayed in Figure 1. These "bad'' polynomials, with elusive zeros and extraneous attracting cycles, are the main focus of this paper.
section_name: Newton's Method and Complex Dynamics

section_content: We now restrict our investigation to applying Newton's method to the family of fourth degree polynomials defined by

)

where λ ∈  is a complex parameter. This family was briefly considered by Sutherland at the end of his doctoral thesis [11]. By expanding p λ, we find that

and thus p λ always has real coefficients. For the remainder of this work, we will denote N λ as the complex rational map obtained by applying New- ton's method to p λ. The free critical points for N λ, arising as the inflection points of p λ, are given by

Note that by Theorem 1, there can be at most two distinct extraneous attracting cycles. If we let

The equation 0 δ = defines a hyperbola in the λ -parameter plane and serves as an important boundary to the possible types of free critical points c ± described below and in Within the curves, the two free critical points are real, while outside, they form a complex conjugate pair. 

• 0

• 0 δ < ⇒ two free critical points as a complex conjugate pair with real part ( )

While this rational map may appear imposing, it possesses certain symmetries under conjugacy that significantly reduce the size of the λ -parameter plane. 

, and iii) N λ is topologically conjugate to N λ via the conjugacy ( )

Proof Since p λ is a polynomial with real coefficients, so is its derivative. This leads to all coefficients in N λ being real and the invariance of the real axis. To prove item (ii), note that p

Since attracting cycles are preserved under conjugacy, the symmetries described in items (ii) and (iii) above enable us to restrict the parameter plane to the region

Being only concerned with those λ -values for which λ N possesses an extraneous attracting orbit, it suffices to follow the orbit of the free critical points for S λ ∈. There is an important additional symmetry when λ is purely ima- ginary, as indicated by the following lemma.
section_name: A Symmetric Fourth-Degree Family of Polynomials

section_content: i N β is conjugate to

Proof It is straight-forward to show that two systems arising from Newton's method applied to polynomials whose roots differ by an affine map, are topologically conjugate (see [5] for example). Since h is an affine map and a bijection between the roots of i p β and the roots of

, the result follows. Alternatively, direct substitution into Equation (1. 4) shows that
section_name: Lemma 2

section_content: Restricting to the case λ ∈  leads to some nice reductions that allow for a complete classification of the orbits of the free critical points. We will show that the orbits of c ± converge to one of the roots of p λ and con- sequently, there are no extraneous attracting cycles for N λ when λ is real. For λ ∈ , N λ now contains the repeated root λ in addition to the roots 1 and 1 -. Because λ is a repeated root, it is not super-attracting in this case. The Newton map simplifies to ( ) ( )

Since the free critical points are real for λ ∈ , we restrict the domain of the Newton map to the real axis, hence the use of the variable x as opposed to z. The free critical points are

Theorem 2 For all λ ∈ , the orbits of the two free critical points under N λ always converge to roots of p λ. Specifically, we distinguish the following two cases:

1) For ( )

2) For ( )

( )

Proof Due to item (iii) of Lemma 1, it suffices to consider the case 0 λ ≥. We begin with the case 0 1 λ ≤ <, where a series of straight-forward estimates show that the critical points, poles and roots of p λ are arranged in the following manner:

and the ordering described in inequality (1. 6), we see that

<. Moreover, it follows from equation (1. 5) that N λ is strictly increasing on the open interval ( )

, we have ( )

It follows that for

, the orbit of x under N λ is a strictly increasing sequence bounded above by the fixed point λ and for ( ]

, the orbit of x is a strictly decreasing sequence bounded below by the fixed point λ. Since the only other fixed points for N λ are 1

and 1, standard arguments using the con-

show that all points in [ ]

converge to the attracting fixed point at λ. For the case 1 λ =, the Newton map reduces even further since p λ has the root 1 with multiplicity three in addition to the simple root 1 -. In this case, N λ becomes ( ) ( )

+ + + and 1 c + = while 0 c -=. The pole p + vanishes and 1 2 p -= -. Here, the free critical point + c has been transformed into an attracting fixed point with derivative 2 3. A simple web-diagram of 1 N shows that the orbit of 0 c -= converges to the repeated root at 1 λ =. To prove the rest of item (ii), we compute that for 1 5

λ =, we have 1 c -= (super-attracting fixed point). The major change between these two inequalities and inequality (1. 6 ) is that the two roots λ and 1 have interchanged positions. For 1 λ >, the super-attract- ing fixed point 1 lies between the two poles, while for 0 1 λ ≤ <, it is the super-attracting fixed point at λ that is located between the two poles. In addition, the critical point c + jumps to the other side of the pole + p, be- coming a local min instead of a local max (see Figure 4 ). However, the dynamical behavior of the critical points is unchanged. Similar arguments as with case (i) show that the orbit of + c monotonically increases as it converges to λ for 1 λ >. The orbit of c -monotonically increases as it converges to 1 for 1 5 λ < < and monotonically decreases to 1 for 5 λ >. As before, the continuity of N λ on the appropriate intervals in addition to the number and location of the fixed points allows us to prove these convergence results. □
section_name: The Case λ Real

section_content: Due to Theorem 2, the interesting Newton maps occur when λ becomes complex. However, the reduction to a real map in the previous subsection suggests a practical approach to investigating the case where λ is purely imaginary. 

Consider the case i λ β = with β ∈ . By Lemma 2, we can restrict to the regime 0 1 β < ≤. For ease of notation, we let N β correspond to Newton's method applied to i p β. In this case, the roots lie on the vertices of a rhombus in the complex plane. Based on previous work applying Newton's method to cubics [1] [2] [5], the symmetry of this configuration suggests that Newton's method may fail badly for certain β -values. We shall see that this is indeed the case. The driving force behind the interesting dynamics is that initial seeds on the real- axis are equidistant from the two imaginary roots. If these seeds are far enough away from the roots at 1 and 1 -, their orbits bounce around the real axis, unable to converge to a root. For the special case 1 β =, the roots of i p are the four roots of unity and as a result of this symmetry, the only "free" critical point is 0 of multiplicity two. However, 0 is also a pole for the Newton map and therefore iterates to ∞ which is a repelling fixed point. Thus, there can be no extraneous attracting cycles for the case 1 β =. For 0 1 β < <, the critical points of N β are real and given by Due to the invariance of the real axis under our Newton map, this implies that any extraneous attracting cycle must lie on the real axis. As with the case λ ∈ , this allows us to restrict our study to a rational map N β of one real variable,

Of particular importance is the fact that N β is an odd function whose critical points are symmetric about the origin. If the orbit of c + lies on or converges to an attracting n -cycle, then so will its partner c -. However, it may be the case that + c and - c lie on or converge to the same attracting orbit, a phenomenon that will be critical to understanding the larger parameter plane picture. 

The poles of N β are 0 and p ±, where 3 p c

Using some basic estimates for the regime 0 1

Moreover, for 0 x p + < <, N β is concave down with a local maximum at c + and ( )
section_name: The Case λ Pure Imaginary

section_content: β <. For x p + >, N β is concave up with a local minimum at the super-attracting fixed point 1. It follows that each ( )

p + ∈ ∞ converges to the root 1 under iteration of N β. Since N β is odd, we also have that each ( )

converges to the root 1 -. Thus, the interesting dynamics of N β occurs between the poles p -and p +. 

The following lemma captures some key information about the behavior of the critical points as β varies. 

Lemma 3 For 1 3 1 β ≤ <, the orbits of the free critical points ± c converge to the roots 1 or 1 -. Spe- cifically, we have the following: ) ( )

which is clearly negative, proving that ( )

is decreasing with respect to β. 

To prove item (ii), letting 1 3 β =, we see that ( )

Item (iii) follows easily from (ii) and the fact that all points in ( ), p + ∞ converge to 1 under iteration while all points in ( )

Our goal is to search for attracting cycles of N β distinct from roots when ( )

. While it is difficult to make precise calculations locating specific examples, there are three values worthy of mention. 

• For ( ) A β = agrees with the value given in Exercise 10, Chapter 13 of Devaney's undergraduate textbook [12]. The second value, 1 β, is found by solving ( )
section_name: N x x

section_name: 2

section_content: β + + = for β. This leads to an even eighth-degree polynomial that can be solved using Ferrari's formula for a quartic equation (calculated using Maple [13] ). It is the only new solution (other than its reciprocal 1 1 β. ) Thus, 1 β and 2 β are the only parameter values in ( ) 0,1 3 for which the critical points are pe- riodic of period two. The value 0 β is significant in that it defines a regime 0 0 β β < < with more complicated dynamics for the free critical points. This will be apparent in some of our numerical investigations. 

Due to the challenges of studying higher iterations of N β analytically, a computer program was written in C++ to locate super-attracting periodic points numerically. For a fixed period where the free critical points lie on super-attracting periodic cycles. Thus, Newton's method for these cases will fail on an open set of initial conditions in the complex plane. As the period n increases, the number of periodic orbits increases exponentially. Letting the program run exhaustively on period 16 yields 2525 periodic points. It should also be noted that error checking devices were used to ensure that poles were not mistaken for periodic points. These values were also confirmed using Maple. 

From this numerical evidence, we observe some noteworthy behavior. For instance, a period-doubling cascade to chaos is readily apparent from our data, as shown in Table 2. These values correspond well to the universal rate of doubling given by Feigenbaum's constant. The period doubling can also be seen in the bifurcation diagram obtained by following the orbits of both free critical points (see Figure 6 ). >, the signs of each element in the orbit of either critical point must alternate, thus making only even periodic cycles possible. This is indicated clearly in Table 1 for 3 n = and 5 n =. Moreover, for β -values close to 2 3 -, the image of both free critical points is close to the pole at zero and consequently, further iteration leads to convergence to either 1 or 1 -. As long as the image of the critical points is trapped close to the pole at zero, the long-term behavior will be convergence to roots. This explains the gap in the center of Figure 6. 

It is important to note that in most cases, the free critical points lie on disjoint periodic orbits of the same period. Except for ( ), all the parameter values in Table 1 are for two disjoint periodic cycles. This is always the case when the period is odd.
section_name: N c c

section_content: For β -values corresponding to super-attracting periodic cycles of odd period, c + and c -do not lie on the same orbit. 

Proof Without loss of generality, suppose

for some primitive odd period n. By contradiction, suppose that c -and c + lie on the same periodic orbit. This implies that ( )

for some k n <. However, since N β is an odd function and c c

Then, by a simple substitution we have ( )

. 

However, we also have that

Therefore, we see that c -is periodic with even period 2k. Again, given the symmetry of N β, this also implies that c + is periodic with even period 2k. Since k n < and n is the least period of the orbit, it follows that 2k n =. But this contradicts the fact that n is odd. □

The birth of two disjoint super-attracting period n -cycles usually arises from a pitch-fork bifurcation. In the rare case when there is a single cycle containing both free critical points on the orbit, a saddle-node bifurcation occurs (see Figure 7 ). These two types of bifurcations transpire remarkably close together in parameter space, with the saddle-node bifurcation taking place first. For example, when 2 n =, the successive values are 2 0. 4439 β ≈ (same cycle) followed by 1 0. 3836 β ≈ (different cycles). For period 4, the parameter values are only 0. 0045 apart.
section_name: Lemma 4

section_content: Unfortunately, once we leave the real and imaginary axes in the λ -parameter plane, we are no longer able to work with the simplified versions of our Newton map. Using a computer, we follow the orbits of both free critical points under N λ as λ varies, producing an intriguing picture of the parameter plane (indicated on the left in Figure 8 ). The deeper red colors indicate faster convergence to a root, while the lighter colors indicate a slower convergence. Black represents λ -values where one or both critical points fails to converge to within 6 10 -of a root after 100 iterations. The same color scheme is used for all figures shown in the parameter plane. 

Much of the analytic work we have done up to this point will be useful in explaining the general structure and interesting figures in the parameter plane for N λ. Although quite intricate, there are a few phenomena that can be explained via connections with the of Milnor on the general case of a complex cubic polynomial [8]. We find a striking similarity between the dynamical behavior of N λ and that of a general complex cubic- polynomial. This is not entirely surprising, since a generic cubic map has two free critical points, exactly the case for N λ. In [8], Milnor classifies four types of hyperbolic components in the parameter plane for the generic cubic based on orbits of the critical points. These four cases are referred to as adjacent critical points, bitransitive, capture and disjoint periodic sinks. Only the second and fourth cases are relevant to our study. We adapt Milnor's definition to our problem. Definition 2 (Milnor [8] ) Suppose that both free critical points c ± converge toward attracting periodic orbits O ± distinct from the roots of p λ. Let U ⊂  be the open set of all points in the basin of attraction of O ±. 

Bitransitive: The two free critical points belong to different components 0 U and 1 U of U, but there exist natural numbers p and q such that ( )

( )

We assume that p and q are primitive, so that both 0 U and 1 U have period p q +.
section_name: The General Case for N λ

section_content: The two free critical points belong to different components 0 U and 1 U, where no forward image of 1 U is equal to 0 U and vice versa. In this case, there exist natural numbers p and q with ( )

( )

Milnor characterizes the types of fractal-like figures we should expect to see in the parameter plane for each case, along with a prototype map. Not surprisingly, there are different figures depending on the type and behavior of the critical points. In the bitransitive case, one finds either a swallow configuration (indicated on the right in Figure 8 ) or a three-pointed configuration Milnor refers to as a tricorn (indicated on the left in Figure 9 ). For the swallow configuration, the prototypical dynamical system is the family of real maps ( )

with parameters 1 2, c c ∈  (so the model swallow lives in 2  ). For our Newton map N λ, the center point of each swallow configuration appears to correspond to the situation where the two real free critical points lie on the same periodic cycle. As discussed in Section 3. 2, this occurs for the parameter values For the tricorn, the model is the complex map ( )

. The tricorn actually contains three embedded copies of the Mandelbrot set, where the cusp of each has been stretched out over a triangular region, joining them in a peculiar fashion. The difference between the tricorn and swallow configuration is the type of critical points. For the swallow, we have real distinct critical points, while for the tricorn, we have a complex conjugate pair. Due to Lemma 2, inverting the key parameter values 2 β and 4 β on the imaginary axis gives a conjugate dynamical system. However, the real, symmetric free critical points are mapped to a pure-imaginary conjugate pair under the conjugacy. This explains why we see tricorns in the parameter plane centered at the λ - values ( ) In the case of disjoint periodic sinks, there is either a product configuration (two long, thin, black strips stretched across each other) or an actual copy of the Mandelbrot set itself (indicated on the right in Figure 9 ). The defining map for the product configuration of the two disjoint real functions with c ∈ , the usual defining map for the Mandelbrot set. The disjoint periodic sink cases yield ``simpler'' dynamical phenomenon because the orbits of the critical points stay away from each other, leading to the decoupled prototype maps. 

Throughout our discussion, we see that the different possible orbits of the free critical points determine the types of fractal-like figures found in the parameter plane. Our detailed analysis of the special case restricting λ to the imaginary axis has provided a useful guide explaining the locations of these figures. However, there are some exceptions worth mentioning. For example, there are disjoint periodic sink values very close to bitransitive values that are contained inside a swallow configuration. Given our limited resolution and computational resources, we are unable to numerically substantiate the claim that every parameter value corresponding to a disjoint periodic sink lies at the center of a product configuration or a Mandelbrot-like set. This discrepancy is most likely attributable to the fact that N λ is a cubic-like map, as introduced by Douady and Hubbard [7], and not an actual cubic polynomial. Nevertheless, after plotting the long-term behavior of each critical point separately and investigating the location of several swallow configurations, we have obtained strong numerical evidence in support of the following conjecture:

Conjecture 1 Each bitransitive λ -value corresponding to the two free critical points sharing the same super-attracting n -cycle lies at the center of a swallow configuration in the parameter plane.
section_name: Disjoint Periodic Sinks:

section_content: In conclusion, we point out one more feature of the parameter plane that can be analytically derived. The yellow diamond shaped boundary that encompasses the interesting dynamical behavior in the parameter plane (indicated on the left in Figure 8 ) is defined by those λ -values where both p λ ′ and p λ ′′ simultaneously vanish. Letting a bi λ = +, this occurs on the algebraic curve Γ defined by the equation For any parameter value on Γ, one or both of the free critical points coincide with poles of N λ and thus map to the fixed point at ∞. Taking successive pre-images of these curves appears to define the sequence of intertwining yellow "leaves" that approach the real axis. 

We have analytically and numerically investigated Newton's method N λ applied to a highly symmetric family of fourth degree complex polynomials. Using techniques from both real and complex dynamical systems theory, we were able to study some reductions of our Newton map that shed significant light on the varied dynamical behavior of this system. Remarkably, we find quite intricate and complicated fractal-like figures throughout the parameter plane for this simple system. Milnor's work on the complex cubic polynomial combined with our understanding of the dynamics for N λ with λ purely imaginary provided some rationale for the location and types of fractals encountered.
section_name: Conclusions

section_name: Acknowledgements

section_content: TOB would like to thank the Department of Mathematics and Computer Science at the College of the Holy Cross for providing an enjoyable and supportive atmosphere throughout his undergraduate career. Funding for TOB was provided by a Council on Undergraduate Research Student Summer Research Fellowship in Mathematics and Science. GR would like to thank Paul Blanchard for pointing him in the direction of this problem and for several interesting discussions concerning this work. He would also like to thank Bruce Peckham for ongoing and insightful discussions on Newton's method. 

T. M. O'Brien, G. E. Roberts
section_name: Acknowledgements