BIB-VERSION:: CS-TR-v2.0
ID:: CORNELLCS//TR93-1392
ENTRY:: 1994-03-17
ORGANIZATION:: Cornell University, Computer Science Department
LANGUAGE:: English
TITLE:: Computing the Newtonian Graph (Extended abstract)
AUTHOR:: Kozen, Dexter
AUTHOR:: Stefansson, Kjartan
DATE:: October 1993
PAGES:: 17
ABSTRACT::
 A polynomial $f\in \complex[z]$ defines a vector field $N_f(z) =
  -f(z)/f'(z)$ on $\complex$.  Certain degenerate curves of flow in
  $N_f$ give the edges of the Newtonian graph, as defined by
  \cite{Sma85}.  These give a relation between the roots of $f$ and
  $f'$, much similar to the linear order, when $f$ has real roots
  only.

  We give a purely algebraic algorithm to compute the Newtonian
  graph and the basins of attraction in the Newtonian field.  The
  resulting structure can be used to query whether two points in
  $\complex$ are within the same basin of attraction.  This gives
  us an algebraic approach to root-finding using Newton's method.
  This method extends to rational functions and more generally to
  any functions on $\complex$ whose flow is algebraic over
  $\complex(e)$.

END:: CORNELLCS//TR93-1392
BODY::
Computing the Newtonian Graph
(Extended Abstract)
Dexter Kozen
Kjartan Stefansson
TR 93-1392
October 1993
Department of Computer Science
Cornell University
Ithaca, NY 14853-7501
Computing the Newtonian Graph
(Extended abstract)
Dexter Kozen*			Kiartan Stefansson*
kozen?cs cornell. edu			stefan?cs . cornell, edu
October 13,1993
Abstract
A polynomial f ? C[z] defines a vector field Nf(z) =
on c. Certain degenerate curves of flow in NJ give the edges of the
Newtonian graph, as defined by [6]. These give a relation between the
roots of f and fi, much similar to the linear order, when f has real
roots only.
We give a purely algebraic algoritlim to compute the Newtonian
graph and the basins of attraction in the Newtonian field. The result-
ing structure can be used to query whether two points in C are within
the same basin of attraction. This gives us an algebraic approach to
root-finding using Newton's method. This method extends to ratio-
nal functions and more generally to any functions on C whose flow is
algebraic over
1 Introduction
We follow the definitions of Smale [6] and define the Newtonian vector field
of a polynomial f ? C[z] by Nf(z) = The name is derived from the
fact that Xk+1 & Xk + NJ(xk) is Newton's method for root approximation.
The vector field NJ defines a flow on ? where the flow comes (almost
everywhere) from a pole (in the case of f a polynomial, from infinity) and
*Computer Science Department, Cornell University, Ithaca, NY 14853
converges (almost everywhere) to a root of f. There are degenerate curves of
flow connecting roots of f and f', and the basins of flow split ? into finitely
many regions. These connecting curves will be the edges of our Newtonian
graph (to be defined more formally later). This graph has been studied and
the possible graphs that can arise have been classified for polynomials [7].
We will give a symbolic algorithm to compute the graph from the poly-
nomial. The output of the algorithm is a structure that can act as an oracle
to answer simple questions such as
1. Given a, b ? C, are a and b in the same basin?
2. Given a, b ? C, are a and b on the same curve of flow?
3. Given a ? ? is a on a basin boundary?
4. Given a E ? is a on a graph edge?
So not only do we get the topology of the graph, we get a method for mem-
bership testing for the interesting regions in the field. Such a structure can
for instance be used in a "guaranteed" Newton's method, modifying the step
size at every point to ensure that we stay within a basin.
We then sketch how to extend the definition of a Newtonian graph for
rational functions. We also observe that the resulting fields on ? satisfy
certain algebraic conditions. Given such conditions we can define the graph
and compute it.
2 The Newtonian Graph
We have defined the Newtonian field of a polynomial. A vector field such as
NJ on C defines a flow on ? Given z ? ? the flow through z is a function
: I H ?, where I C ? is an interval containing zero, ? differentiable with
dt
= Nj(?z(t))
That is, ? parameterizes the flow starting at z and at every point the speed
and direction agrees with the field. An example of flow of f of degree 4 is
given in Figure 1
2
mm
Figure 1: Flow in the Newtonian field of a degree 4 polynomial. Every curve
of flow is directed to a root, except the basin boundaries (dotted lines). There
is a root of f' on every basin boundary, and a curve of flow from there to
"adjacent" roots (also dotted lines).
The flow exists on all of C --H Vft (where Vft = fz ? : f'(z) = 01). The
existence and uniqueness follows from the theory of differential equations and
the fact that N1 is a C1 function on --H Vji (see e.g. [3], 8.2 and 8.5).
The following lemma [7] gives us a very important property of the flow:
Lemma 1 Let f ? ?[z], and ?? be flow thwugh z in the Newtonian field
N1. Then f maps the curve t?z(t) : t e 11 to a straight line pointing to
the origin. More specifically,
= e?tf(z).
3
Proof. Computing df(?(t)) using the chain rule gives:
df(? (t))			_____
dt
dt
which is a differential equation in t for the function p(t) = f(?z(t)). Given
the initial condition, p(O) = f(z), it has the unique solution p(t) = e?tf(z),
?e. f(?z(t)) = e--Htf(z). ?
Using the properties of ? one can show the following
Lemma 2 For every z E C--H (Vf uv,1), ? is defined on a maximum interval
containing 0, which is of the following type, for some a, b E ?:
1. (--Hoo, +oo), and the flow comes in from infinity and goes to a root of
2. (--Hoo, a) and the flow comes in from infinity and goes to a root off';
3. (a, b) and the flow comes in from a root of f' and goes to another root
off';
?. (a, +00) and the flow comes in from a root of fI and goes to a root of
Proof. NJ is a C1 function W C where W = --H (Vf U Vff), and by
Theorem in 8.5 in [3], all flow must leave any compact set of W. By Lemma
1, the (maximum) interval of ? is unbounded upwards iff the flow goes to
a root of f. The same argument shows that the interval is not downward
bounded iff the flow comes in from 00. Since the flow leaves any compact set
of W the only other limit points are in Vft.
Definition 3 The Newtonian Graph of f ? C[z] is the plane graph G =
(V, E) with vertices V = V1 u V1 and directed edges being the curves of flow
between vertices, where these exist. (Zi
4
Figure 2: The graph of the field in Figure 1. The solid dots are the roots of
the polynomial, the hollow ones are the roots of the derivative.
We note that the graph is not just a combinatorial structure, as the edges
come with an embedding defined by ?.
Figure 2 shows the Newtonian graph of the field shown in Figure 1. An-
other example in Figure 3 shows that there can be connections between two
roots of f'.
We observe that under f, every edge maps onto a straight line segment
pointing to 0 in ? with at least one endpoint in ff(c) : = 01. This is an
immediate consequence of Lemma 1 and the fact that edges are curves of flow.
We can conversely look at the pre-images (under f) of such line segments, and
we get finitely many curves (at most n(n --H 1), since f is an n to 1 mapping,
and f' has at most n --H 1 roots). We will use this observation later, that the
graph is contained in the pre-image f-1(fmf(c) : f'(c) = 0 0 ? m < 11).
Thus the graph has finitely many edges. Furthermore [7] show that the graph
is connected and go on to classify the possible types of graphs that can arise.
A basin of attraction is a connected region where the flow comes in from
infinity and goes to one particular root of f. A basin boundary is the bound-
ary of two basins. It is not hard to show that there must be a root of f' on
every basin boundary, because it requires a discontinuity of Nf for the flow
to "split" into two directions, and these are only at the roots of f'. Also the
basin boundaries must be curves of flow themselves, so we conclude that every
basin boundary is flow into a root of f'. In particular this means that basin
boundaries are contained in the pre-image f--H'(fmf(c) : J'(c) = 0,1 < ml).
5
Figure 3: The field and the graph of a degree 3 (real) polynomial where the
two derivative roots are linked.
3 Computing Basins and Graph Edges
compute the basin boundaries and the edges
first we need a few preliminaries on cylindric
We will give an algorithm to
of the Newtonian graph. But
algebraic decomposition.
3.1 Cell Decomposition
We describe cylindric algebraic cell decomposition briefly. For more detailed
description, see [2] or [1].
Definition 4 A decomposition of Rk is a finite partition fCt?iEI such that
each Ci is connected, Ci fl Cj = ? if i # j and UjEl C --H Th?k For k = 1 such a
decomposition is cylindrn'c if the each Cj is either a point or an interval. For
k > 1, the decomposition is cylindric if for all r 1 < r < k, ??ri,..,r(Ci) : i C
I? is a cylindric decomposition of ?r ?
Definition 5 Given polynomial equations, fi(x1, ..., x?) = 0, i = 1,..., n,
with f? E ?ki, ..., x?], a Cylindric Algebraic Decomposition (CAD) of ik? is
6
a data structure D with the following properties.
o+ D contains a graph, where the nodes correspond to subsets (cells) of
??, each cell being homeomorphic to ?d for some d, and the cells are
a decomposition of ?m
o+ For all i = 1,.., n, sign(f?) is constant on every cell. Each cell is labelled
with the signs that the f? take on that cell.
o+ Every node contains an oracle such that given any c ? ?rn? the oracle
can answer if c is in the subset corresponding to the node.
o+ Every node contains dimension information, corresponding to the di-
mension of the cell.
o+ The edges of the graph correspond to adjacency of the cells in ?m? ?.e.
there is an edge (v, v) if the subsets that v and v represent are adjacent.
o+ The decomposition is cylindric.
Algorithms have been developed to compute (parts of) such a cell decom-
position dating back to Tarski in 1948 [8]. Collins [2] has a double exponential
algorithm, although it lacks some of the adjacency information. Ben-Or et
al. [1] developed a parallel algorithm giving the same kind of decomposition
(the BKR algorithm), and Kozen and Yap [4] extended that algorithm to ob-
tain the adjacency information as well (here after named the extended BKR
algorithm).
We note that due to the cylindric condition and adjacency information,
an algorithm computing such a decomposition can be used on a set of poly-
nomials with quantifiers, projecting down the result. If the input is a system
of polynomials of the form
??1, ..Yk :			f1(x1, ..., Xm, Yi, ., Yk) = 0
fn(xi, ..., Xm, Yi, ., ?) = 0
then using CAD on we can project the solution down to ??, by treating
the partitions according to Yi, .., Yk as insignificant. The resulting structure
7
can be used for answering questions of the form: Given c e ??, does there
exist ?1, ., Yk E ?k such that Yi, .., Yk, c is a solution to the system?
We note that the order of variables is important with respect to the
cylindric condition.
3.2 The Algorithm
Recall that every basin boundary and every edge is mapped by f onto a
straight line. Also, the basin boundaries and edges have a root of /` as a
limit. Thus, all these "interesting" curves of flow satisfy, for every z on the
curve,
C,m ? : f(z) = mf(c)
= 0
(1)
Any point z on a basin boundary or an edge must satisfy these two conditions.
We note that the converse is not true; z ? C can be a solution to (1) without
being on an edge or a basin boundary.
We proceed in two steps. First we find a decomposition of C describing
where we have solutions z to (1). Then we prune that output, because we
may get solution curves which do not correspond to basin boundaries or
edges.
To find the solutions to (1), we can feed the equations
f(z) = mf(c)			(2)
= 0
into our favorite cylindric algebraic decomposition algorithm. The resulting
structure would be a decomposition of If? x C x ? describing regions where
such m, c, z exist, along with the dimension of each region and adjacency
information. Projecting m and c, we get curves in C for which there exists a
solution to (1).
First let us note that algorithms such as Collins' and the extended BKR do
decomposition over the reals. But we can split the equations into a real and
imaginary parts, and get a decomposition of Ik? ? x ? x ? corresponding
to the equations
fR(x,y) = mf?(ci,c2)
8
f1(x,y) = mfi(c1,c2)
f?'(ci,c?) =			0
=			0
where f(x+iy) = f?(x,y) +ifi(x,y) with fR, fi C If?[x,y]. We get a decom-
position on ?? which corresponds to a decomposition on I? x ? x C.
We then project the dimensions corresponding to c = c1 t ic2 and again
project m, obtaining a decomposition of C corresponding to z for which there
exist m and c satisfying equation (2).
This decomposition will contain the basin boundaries and graph edges.
These may be partitioned into segments (bounded i-cells) and 0-cells between
such segments. There may be other i-cells present which are not solutions
to the system (introduced by the CAD algorithm to get a finer partition).
However, we can always identify the segments which are a solution to this
system, because all of them are labelled with the signs of the input polyno-
mials. The curves which are actual solutions to the system (2) will all show
f(z) = mf(c). Hence a solution curve to the system can be reconstructed
by linking such adjacent cells.
But not all solution curves are edges or basin boundaries. The following
lemma classifies the types:
Lemma 6 The output from the process above contains at most 0(n2) i-cells
which are solutions curves for the input system. They are of the following
types:
1. Adjacent to two 0-cells, one of which describes a root of f' and one
which describes a root of either f or
2. Adjacent only to one 0-cell which describes a root of f';
8. Adjacent only to one 0-cell which describes a root of f.
Cells of type 1 are edges of the Newtonian graph, a cell of type 2 is a basin
boundary and cells of type 8 are extraneous solutions to the system.
Proof. The only i-dimensional cells that can be solutions to the system
correspond to curves of flow. Then the classification is obvious from the
definition of edges and the properties of basin boundaries.
9
There are at most 0(n2) solution curves for f(z) = mf(c) with c a root
of f'(c) and m ? ?, because there are at most n --H 1 roots of f' and f is an
n to 1 mapping. ?
The cells of type 1 and 2 are the ones we are interested in and we must
tell these apart from the extraneous cells of type 3.
Since f is a part of the input, the signs of f are given on every cell. In
particular this allows us to verify if a curve ends at a root of f.
Depending on which algorithm we use, we may or may not have all the
information needed. The extended BKR guarantees that if f is a part of the
input, then the signs of f' will be provided on each cell. If we don't have this
guarantee, we can always add f'(z) = 0 to our input equations and get the
same information that way.
At this point we can determine the types of the solution curves. Now
it is easy to implement the promised "pruning" step. We simply eliminate
all cells of type 3. More precisely, we mark them as parts of the adjacent
2-dimensional cells (which are the basin that this cell lies in).
Now the structure can be used in answering queries. Two points are in
the same basin if they are in the same 2-cell or if they are separated only by
"fake" i-cells (of type 3).
4			Improvements
Recall we did cylindric decomposition on If?5 ?? ? ?2 of the equations
f(z) = mf(c)
= 0
and projected the solution onto ? This can be simplified by defining
g(m, z) = Res?(f(z) --H mf(c), f'(c)),
where Res? denotes the univariate resultant of two inputs, considered as poly-
nomials in c. (Here we view f(z) --H mf(c) and f'(c) as univariate polynomials
in C[z,m]kJ).
Then g has the property that g(m, z) = 0 iff ?c E f(z) --H mf(c) =
0 = f'(c). Hence, a decomposition of I? x C with respect to g is the same as
10
the projection of the decomposition of ? x ? x  with respect to the original
two equations.
The only thing we must be aware of is how to obtain the necessary signs
of f and f' on cells, in order to identify and link up solution curves and prune
off the redundant ones. One way would be to add the equation f(z) = 0 (and
= 0, if we are not using the extended BKR), and do a decomposition
with respect to f (f') and g. This is already an improvement in terms of
dimensions, since we are only working with 3 real variables (x = R?e(z),y =
Im(z) and m) instead of 5 before.
The asymptotic complexity remains the same, but the constants are
clearly much better. The extended BKR gives an NC circuit of depth
2o(d2) ?0gO(d) n where d is the number of variables and n is the maximum
of either the number of polynomials or their degrees. In our case the circuit
will be of depth o(logO(l) n) where n is the degree of the input polynomial
5 Applications
The Newtonian graph is of its own interest, as it describes the arrangement
of the roots of f and f'. Our computation gives a complete topological
information of both the graph and the basins of the Newtonian field.
The relation to Newton's method gives an interesting method of approx-
imating all the roots of f simultaneously, guaranteeing convergence. If we
start with a point z0 in a basin, we can apply modified Newton's method,
where the iteration Zk+1 & Zk t Nf(zk) is replaced by:
a 1,
repeat
Zk+1 & Zk t QNf(zk)
a a/2
until (Zk+1 is in the same basin as zk)
Le. we scale down the step in order to ensure that we stay within the same
basin. Here we use our pre-computed structure of basins to determine if two
points are in the same basin. If we furthermore require that a progress is
made at each step, (i.e. that If(zk+1)I < if(zk)I, then we are guaranteed that
the sequence fzkl will eventually converge to the root in the basin of z0.
11
We remark that this does not guarantee quadratic convergence every-
where, it only ensures that the method will converge.
6 More General Newtonian Graphs
In this extended abstract, we will briefly describe how the Newtonian graph
and its computation extend to more general vector fields on
The definition of a Newtonian graph (Definition 3) only uses the fact the
the function f : W ? is C2 on an open subset W c ? which makes
N1(z) = --Hf(z)/f'(z) a C1 vector field on W. Lemma 1 still holds, but
Lemma 6 now allows curves parameterized (--Hoo, +00) coming in from either
infinity or a pole of f and going to either infinity or a root of f.
Figure 4: Flow in the Newtonian field of a rational function with three roots
and four poles. Curves of fixed color indicate curves of flow.
12
In particular, most of the same observations hold for rational functions.
The function f maps the curves of flow onto straight line segments pointing
to the origin. The natural definition of a Newtonian graph for the rational
function is the graph whose edges are the directed curves from poles or infinity
to roots of f'; between the roots of f' and from roots of f' to roots of f or
to infinity. This differs only slightly from Smale's definition [6] in that for
polynomials we now count the basin boundaries as a part of the graph. A
nice property of the graph is that it is symmetric in poles and roots, i.e. the
graphs of Np/q and Nq?p are identical, except the directions of the edges are
reversed.
Wflte f(z) = p(z)/q(z) with p, q ? C[z] with gcd(p, q) = 1. If c is a root
of f'(c) = 0 then f maps any edge into c onto a ray ?mf(c) : 0 < m <
Again we can consider a system of equations
C,m ? fF?? : f(z) = mf(c)
= 0
which now is equivalent to the system of polynomials
? C, m ? : p(z)q(c) = mq(z)p(c)
p'(c)q(c) --H p(c)q'(c) = 0.
This we can solve with algebraic decomposition as before and determine the
actual solution curves which are edges. As before we can reduce the number
of variables by using resultants. Let
g(m, z) = Res?(p(z)q(c)mq(z)p(c), p'(c)q(c) --H
Then the previous system is equivalent to am g(m, z) = 0. As before we
have devised an NC algorithm to compute the Newtonian graph.
The key property used in the computation is that the flow satisfies f(?(t)) =
e--Htf(z). In the case of a rational function f(z) = p(z)/q(z) this translated
into the polynomial equation
--H			= 0.
This equation also implicitly defines the flow ?? satisfying ?z(0) =
13
In more generality, consider any flow (t)z defined on W
is almost all of ? Assume for some polynomial g with
cients and for all z ? W, ?z(O) = z and g(??(t),et,z) = 0.
Dg(??(t),et,z) = 0, ze.
D1g(??(t), et, Z)?'z(t) t D29(?zU), et, z) = 0
which is equivalent to
= --HD2g(??(t), et, z)/Dig(?z(?), et, z).
C C where W
complex coeffi-
Then we have
(3)
If a basin of attraction is the area where all the flow goes from one pole to
one root then as before there is a discontinuity in the field somewhere along
every basin boundary. In particular, ? will be undefined at such points.
Equation (3) shows that &(?) extends continuously to all of C except the
points where D1g(??(t),et,z) = 0 and D2g(?z(t),et,z) # 0. These points
can be computed as being the x ? C for which there exist m, m' ? ? and
z,x E?wfth
g(x,m,z) =			0
g(x',m',z) =			0
D1g(x',m',z) =			0.
The first two conditions force x and x' to be on the same curve whereas the
third condition places x' at a discontinuity of the field. It is easy to verify
that for rational functions we get the same equations as before.
Again we can do CAD on ?? and project down for the x variable to obtain
the solution curves, which then are the edges of the Newtonian graph.
7 Acknowledgments
This work was supported in part by the National Science Foundation and
in part by the United States Army Research Office through the Army Cen-
ter of Excellence for Symbolic Methods in Algorithmic Mathematics (AC-
SyAM), Mathematical Sciences Institute of Cornell University under contract
DAALO3-91-C-0027.
14
References
[1] Michael Ben-Or, Dexter Kozen and John Reif The Complexity of Elemen-
tary Algebra and Geometry Journal of Computer and System Sciences, Vol.
32 No. 2, April 1986
[2] G. E. Collins, Quantifler Elimination for Real Closed Fields by Cylindrical
Algebraic Decomposition, Lecture Notes in Computer Science, vol. 33, pp.
134-183, Springer-Verlag, Berlin Heidelberg New York 1975.
[3] Morris W. Hirsch and Stephen Smale Differential Equations, Dynamical
Systems, and Linear Algebra, Academic Press, 1974
[4] Dexter Kozen, Chee-Keng Yap Algebraic Cell Decomposition in NC, Pro-
ceedings of the 26th Annual Symposium on Foundations of Computer Sci-
ence (FOCS) 1985, pp. 515-521
[5] Morris Marden The geometry of the Zeroes of a Polynomial in a Complex
Variable, American Mathematical Society, New York, 1966
[6] Stephen Smale On the Efficiency of Algorithms of Analysis, Bulletin of
the American Mathematical Society, Vol. 13 No. 2, pp. 87-122, October
1985
[7] Michael Shub, David Tischler and Robert F. William The Newtonian
Graph of a Complex Polynomial, SlAM J. Math. Anal., Vol. 19, No. 1,
January 1988
[8] Alfred Tarski, A Decision Method for Elementary Algebra and Geometry,
Univ. of Calif. Press Berkeley, 1948; 2nd ed. 1951.
15
