NUMERICAL HOMOTOPIES TO COMPUTE GENERIC POINTS ON POSITIVE DIMENSIONAL ALGEBRAIC SETS

ANDREW J. SOMMESE AND JAN VERSCHELDE

Abstract. Many applications modeled by polynomial systems have positive dimensional solution components (e.g., the path synthesis problems for four-bar mechanisms) that are challenging to compute numerically by homotopy continuation methods. A procedure of A. Sommese and C. Wampler consists in slicing the components with linear subspaces in general position to obtain generic points of the components as the isolated solutions of an auxiliary system. Since this requires the solution of a number of larger overdetermined systems, the procedure is computationally expensive and also wasteful because many solution paths diverge. In this article an embedding of the original polynomial system is presented, which leads to a sequence of homotopies, with solution paths leading to generic points of all components as the isolated solutions of an auxiliary system. The new procedure significantly reduces the number of paths to solutions that need to be followed. This approach has been implemented and applied to various polynomial systems, such as the cyclic n-roots problem.

AMS Subject Classification. 65H10, 68Q40. Keywords. polynomial system, numerical homotopy continuation, components of solutions, numerical algebraic geometry, generic points, embedding.

Contents

1. Introduction 1 2. An Embedding of a Polynomial System 3 3. A Worked Out Example 8 4. Bertini’s Theorem and a Local Extension Theorem 10 5. Applications and Computational Experiences 12 5.1. A planar four-bar mechanism 13 5.2. Constructing Runge-Kutta formulas 13 5.3. Modular equations for special algebraic number fields 13 5.4. On Fourier transforms: the cyclic n-roots problem 14 6. Conclusions and Future Directions 15 References 16

Date: June 28, 1999. Research at MSRI is supported in part by NSF grant DMS-9701755.

1. Introduction

Let

f(x) := [ f_{1}(x_{1}, …, x_{n}) ; ⋮ ; f_{N}(x_{1}, …, x_{n}) ]  (1)

denote a system of N polynomials on ℂ^{n}. Positive dimensional components of the solution set of f(x) = 0 are a common occurrence, even when N = n. Sometimes they are an unpleasant side show that happens with a system generated using a model, for which only the isolated nonsingular solutions are of interest; and sometimes, the positive dimensional solution components are of primary interest. In either case, dealing with positive dimensional components, is usually computationally difficult. In [35], Sommese and Wampler presented a numerical algorithm, which uses auxiliary systems to numerically find sets of solutions of the original system. These sets of solutions, which let one numerically decide the dimension of the zero set of the original system, include at least the isolated solutions of the original system, plus “generic points” of each positive dimensional irreducible component of the solutions of the original solutions. Generic points are the basic numerical data which we are using to investigate the positive dimensional solution components.

The algorithm from [35], which is based on slicing with general linear spaces of different dimensions, leads to n auxiliary systems, which must be dealt with. In this paper we present an embedding of the system f(x) = 0 into a family of system of polynomials depending on 2n variables (x, z) ∈ ℂ^{2n}, and a large space of parameters. We then single out n + 1 of the systems, E_{i}(x, z) for i from n to 0 (here E_{0}(x, z) = 0 is equivalent to the system f(x) = 0) obtained by choosing particular values of the parameters, plus a homotopy H_{i} going from E_{i} to E_{i−1}. For simplicity we assume that f is not identically zero. The system E_{n} = 0 has isolated nonsingular solutions. We use polynomial continuation [1, 2], [23], [27] to implement the following algorithm

Algorithm 1. Cascade of homotopies between embedded systems.

Input: f, n.  system with solutions in ℂ^{n}
Output: (E_{i}, X_{i}, Z_{i})_{i=0}^{n}.  embeddings with solutions
E_{0} := f;  initialize embedding sequence
for i from 1 up to n do  slice and embed
E_{i} := Embed(E_{i−1}, z_{i});  z_{i} = new added variable
end for;  homotopy sequence starts
Z_{n} := Solve(E_{n});  all roots are isolated, nonsingular, with z_{n} ≠ 0
for i from n − 1 down to 0 do  countdown of dimensions
H_{i+1} := tE_{i+1} + (1 − t) ( \begin{matrix}E_{i} \\ z_{i+1} \end{matrix} );

homotopy continuation t : 1 → 0 to remove z_{i+1}

X_{i} := limits of solutions of H_{i+1}
as t → 0 with z_{i} = 0;  on component
Z_{i} := H_{i+1}(x, z_{i} ≠ 0, t = 0);  not on component: these solutions

are isolated and nonsingular

end for.

The routine ‘Embed’ will be defined in the next section. In section three we present a worked out example of the algorithm. Section four contains the mathematical background needed to prove our main results:

1. if i is the largest integer with X_{i} nonempty, then the dimension of f^{−1}(0) is i; and 2. given any irreducible component W of f^{−1}(0) of dimension i, then, counting multiplic-

ities, X_{i} contains deg(W) generic points of W.

The applications described in section five illustrate the performance of the new procedure. We end this paper with some directions for future research.

Acknowledgements. Both authors would like to thank MSRI for providing the stimulating environment where our collaboration began. We would also like to thank the Center for Applied Mathematics, and the Duncan Chair of the University of Notre Dame for their support of our work. The authors are grateful to Göran Björck for pointing at the results of Uffe Haagerup. The second author is supported by a post-doctoral fellowship at MSRI, and supported in part by the NSF under Grant DMS - 9804846 at the Department of Mathematics, Michigan State University.

2. An Embedding of a Polynomial System

Throughout this paper we work with algebraic functions. This allows the possibility of using rational functions, and not just polynomials. This extra flexibility will be needed in a sequel, where, even though we start with a system of polynomials on ℂ^{n}, it becomes necessary to work with rational functions on a Zariski open set of an associated Euclidean space. We assume in what follows that locally we have the same number of equations as unknowns. A procedure for reducing to this case is presented in [35].

We refer to [35] for a discussion of generic points.

Given a system of algebraic functions

f(x) := [ f_{1}(x) ; ⋮ ; f_{n}(x) ]  (2)

on a connected algebraic manifold X of dimension n embedded into a Euclidean space ℂ^{A} we have the following basic embedding into a family of systems. First we restrict n linear functions on ℂ^{A} to X. Thus for j from 1 to n we have

L_{j}(x) := a_{j} + a_{j,1}x_{1} + ⋯ + a_{j,A}x_{A},  (3)

where x_{i} is the restriction of the i-th coordinate function of ℂ^{A} to X. By abuse of notation we let L_{j} ∈ ℂ^{A+1} denote

( a_{j} a_{j,1} ⋯ a_{j,A} ) .  (4)

We fix linear coordinates z_{1}, …, z_{n} on a complex Euclidean space ℂ^{n}. Then we have the system of equations

E_{i}(f)(x, z, Λ_{1}, …, Λ_{i}, L_{1}, …, L_{i}) := [ f_{1}(x) + ∑_{j=1}^{i} λ_{1,j}z_{j} ; ⋮ ; f_{n}(x) + ∑_{j=1}^{i} λ_{n,j}z_{j} ; a_{1} + a_{1,1}x_{1} + ⋯ + a_{1,A}x_{A} + z_{1} ; ⋮ ; a_{i} + a_{i,1}x_{1} + ⋯ + a_{i,A}x_{A} + z_{i} ]  (5)

where we let

Λ_{i} := [ λ_{1,i} ; ⋮ ; λ_{n,i} ] .  (6)

We often refer to E(f)(x, z_{1}, …, z_{n}, Λ_{1}, …, Λ_{n}, L_{1}, …, L_{n}) by E_{n} or E_{n}(f). Further we let E_{i} or E_{i}(f) denote E(f)(x, z_{1}, …, z_{i}, Λ_{1}, …, Λ_{i}, L_{1}, …, L_{i}) on ℂ^{n+i}. Note that E_{0} is just f; and that the solutions (x, z_{1}, …, z_{i}) ∈ X × ℂ^{i} of

E_{i}(x, z_{1}, …, z_{i}, Λ_{1}, …, Λ_{i}, L_{1}, …, L_{i}) = 0  (7)

are naturally identified with the solutions (x, z_{1}, …, z_{i}, 0 …, 0) ∈ X × ℂ^{n} of the system

\{ \begin{matrix}E_{i}(x, z_{1}, …, z_{n}, Λ_{1}, …, Λ_{i}, 0, …, 0, L_{1}, …, L_{i}, 0, …, 0) & = & 0 \\ z_{i+1} & = & 0 \\ ⋮ & ⋮ \\ z_{n} & = & 0. \end{matrix}  (8)

We let Y denote the space ℂ^{n×(A+1)} × ℂ^{n×n} of parameters

[ \begin{matrix}a_{1} & a_{1,1} & ⋯ & a_{1,A} & λ_{1,1} & ⋯ & λ_{n,1} \\ ⋮ & ⋮ & ... & ⋮ & ⋮ & ... & ⋮ \\ a_{n} & a_{n,1} & ⋯ & a_{n,A} & λ_{1,n} & ⋯ & λ_{n,n} \end{matrix} ] ∈ ℂ^{n×(A+1)} × ℂ^{n×n}  (9)

for these systems. We have used the transpose of the λ_{i,j} for convenience in describing the stratification Y_{0} ⊂ Y_{1} ⊂ ⋯ ⊂ Y_{n} of the space Y given by defining Y_{n} := Y; and Y_{i}, for i = 0, …, n − 1, as the subset of Y with the coordinates

[ \begin{matrix}a_{i+1} & a_{i+1,1} & ⋯ & a_{i+1,A} & λ_{1,i+1} & ⋯ & λ_{n,i+1} \\ ⋮ & ⋮ & ... & ⋮ & ⋮ & ... & ⋮ \\ a_{n} & a_{n,1} & ⋯ & a_{n,A} & λ_{1,n} & ⋯ & λ_{n,n} \end{matrix} ] ∈ ℂ^{(n−i)×(A+1)} × ℂ^{(n−i)×n}  (10)

set equal to 0. Thus using the identification (8) above, we can regard Y_{i} as the parameter space of the systems of equations

E_{i}(x, z_{1}, …, z_{i}, Λ_{1}, …, Λ_{i}, L_{1}, …, L_{i}) = 0.  (11)

We consider homotopies H_{i}, for i from n to 1 and t from 1 to 0, defined by:

H_{i}(x, z_{1}, …, z_{i}, t) := [ f_{1}(x) + ∑_{j=1}^{i−1} λ_{1,j}z_{j} + tλ_{1,i}z_{i} ; ⋮ ; f_{n}(x) + ∑_{j=1}^{i−1} λ_{n,j}z_{j} + tλ_{n,i}z_{i} ; L_{1}(x) + z_{1} ; ⋮ ; L_{i−1}(x) + z_{i−1} ; tL_{i}(x) + z_{i} ] ,  (12)

with the convention that if i = 1, we mean

H_{1}(x, z_{1}, t) := [ f_{1}(x) + tλ_{1,1}z_{1} ; ⋮ ; f_{n}(x) + tλ_{n,1}z_{1} ; tL_{1}(x) + z_{1} ] .  (13)

Thus H_{i}(x, z_{1}, …, z_{i}, 1) = E_{i} and H_{i}(x, z_{1}, …, z_{i}, 0) = 0 is

\{ \begin{matrix}E_{i−1}(x, z_{1}, …, z_{i−1}) & = & 0 \\ z_{i} & = & 0 \end{matrix}  (14)

Note that using this convention, H_{i} can be rewritten as tE_{i} + (1 − t) ( \begin{matrix}E_{i−1} \\ z_{i} \end{matrix} ).

Lemma 2. There is a nonempty Zariski open set U of points

[ \begin{matrix}a_{1} & a_{1,1} & ⋯ & a_{1,A} & λ_{1,1} & ⋯ & λ_{n,1} \\ ⋮ & ⋮ & ... & ⋮ & ⋮ & ... & ⋮ \\ a_{n} & a_{n,1} & ⋯ & a_{n,A} & λ_{1,n} & ⋯ & λ_{n,n} \end{matrix} ] ∈ ℂ^{n×(A+1)} × ℂ^{n×n}  (15)

such that for each i = 1, …, n, 1. the solutions of E_{i}(x, z_{1}, …, z_{i}) = 0 with (z_{1}, …, z_{i}) ≠ 0 are isolated and nonsingular;

2. given any irreducible component W of f^{−1}(0) of dimension i, and counting multiplicities, the isolated solutions of E_{i}(x, z_{1}, …, z_{i}) = 0 with (z_{1}, …, z_{i}) = 0, contain deg(W) generic points of W; and

3. the solutions of E_{i}(x, z_{1}, …, z_{i}) = 0 with (z_{1}, …, z_{i}) ≠ 0 are the same as the solutions of E_{i}(x, z_{1}, …, z_{i}) = 0 with z_{i} ≠ 0.

Proof. For a fixed i ∈ {0, …, n} consider the system of equations

{ ∑_{j=1}^{n} α_{1,j}f_{j} + ∑_{j=1}^{i} β_{1,j}z_{j} = 0 ; ⋮ ; ∑_{j=1}^{n} α_{n,j}f_{j} + ∑_{j=1}^{i} β_{n,j}z_{j} = 0 ; γ_{1} + ∑_{j=1}^{A} γ_{1,j}x_{j} + ∑_{j=1}^{i} δ_{1,j}z_{j} = 0 ; ⋮ ; γ_{i} + ∑_{j=1}^{A} γ_{i,j}x_{j} + ∑_{j=1}^{i} δ_{i,j}z_{j} = 0.  (16)

By Bertini’s theorem (see section four for a convenient form of this result), there is a nonempty Zariski open set U of

(α, β, γ, δ) ∈ ℂ^{n×n} × ℂ^{n×i} × ℂ^{i×(A+1)} × ℂ^{i×i}  (17)

where

α := [ \begin{matrix}α_{1,1} & ⋯ & α_{1,n} \\ ⋮ & ... & ⋮ \\ α_{n,1} & ⋯ & α_{n,n} \end{matrix} ] ; β := [ \begin{matrix}β_{1,1} & ⋯ & β_{1,i} \\ ⋮ & ... & ⋮ \\ β_{n,1} & ⋯ & β_{n,i} \end{matrix} ]  (18)

and;

γ := [ \begin{matrix}γ_{1} & γ_{1,1} & ⋯ & γ_{1,A} \\ ⋮ & ⋮ & ... & ⋮ \\ γ_{i} & γ_{1,1} & ⋯ & γ_{i,A} \end{matrix} ] ; δ := [ \begin{matrix}δ_{1,1} & ⋯ & δ_{1,i} \\ ⋮ & ... & ⋮ \\ δ_{i,1} & ⋯ & δ_{i,i} \end{matrix} ] .  (19)

such that, if not empty, the zero set Z of the above system of equations on X × ℂ^{i} minus the set of common zeroes of f_{1}, …, f_{n}, z_{1}, …, z_{i} is smooth and of dimension 0. Left multiplying this n + i vector of equations with an invertible (n + i) × (n + i) matrix G of the form

[ \begin{matrix}A & 0 \\ 0 & B \end{matrix} ]  (20)

where A is an n × n matrix and B is a i × i matrix, results in the equivalent system of equations

\begin{aligned}A · α · f & + & A · β · [ z_{1} ; ⋮ ; z_{i} ] & = & 0 \\ B · γ · [ 1 ; x_{1} ; ⋮ ; x_{A} ] & + & B · δ · [ z_{1} ; ⋮ ; z_{i} ] & = & 0. \end{aligned}  (21)

Thus we can assume that the set V is invariant under this action by the matrices G. Thus we have a nonempty Zariski open set U of

[ \begin{matrix}a_{1} & a_{1,1} & ⋯ & a_{1,A} & λ_{1,1} & ⋯ & λ_{n,1} \\ ⋮ & ⋮ & ... & ⋮ & ⋮ & ... & ⋮ \\ a_{n} & a_{n,1} & ⋯ & a_{n,A} & λ_{1,n} & ⋯ & λ_{n,n} \end{matrix} ] ∈ ℂ^{n×(A+1)} × ℂ^{n×n}  (22)

such that for each i the equivalent system is of the form

\{ \begin{matrix}f_{1}(x) + ∑_{i=1}^{i} λ_{1,i}z_{i} & = & 0 \\ ⋮ & ⋮ \\ f_{n}(x) + ∑_{i=1}^{i} λ_{n,i}z_{i} & = & 0 \\ a_{1} + ∑_{j=1}^{A} a_{1,j}x_{j} + z_{1} & = & 0 \\ ⋮ & ⋮ \\ a_{i} + ∑_{j=1}^{A} a_{i,j}x_{j} + z_{i} & = & 0 \end{matrix} ,  (23)

and has smooth nonsingular zeroes when (z_{1}, …, z_{i}) ≠ 0. Thus we have the first assertion of the Lemma.

To see the second assertion of the Lemma, note that the solutions of E_{i} = 0 with (z_{1}, …, z_{i}) = 0 are naturally identified with the solutions of the system

\{ \begin{matrix}f_{1}(x) & = & 0 \\ ⋮ & & ⋮ \\ f_{n}(x) & = & 0 \\ a_{1} + ∑_{j=1}^{A} a_{1,j}x_{j} & = & 0 \\ ⋮ & & ⋮ \\ a_{i} + ∑_{j=1}^{A} a_{i,j}x_{j} & = & 0 \end{matrix}  (24)

The second assertion follows now from the Algorithm in [35, §3.1].

We prove the third assertion of the lemma by induction on i. If i = 1, then it is a tautology that the solutions of E_{i}(x, z_{1}, …, z_{i}) = 0 with (z_{1}, …, z_{i}) ≠ 0 are the same as the solutions of E_{i}(x, z_{1}, …, z_{i}) = 0 with z_{i} ≠ 0.

So we can assume that the result is true for k < i where i > 1. Note that a solution of E_{i} = 0 with z_{i} = 0 but (z_{1}, …, z_{i}) ≠ 0 is a solution of the system

{ E_{i−1}(x, z_{1}, …, z_{i−1}) = 0 ; L_{i}(x) = 0 ; z_{i} = 0 ; ⋮ ; z_{n} = 0.  (25)

Since the solutions of E_{i−1} = 0 with z_{i−1} ≠ 0 are isolated and nonsingular, a generic choice of L_{i}(x) will not be zero on any of the solutions. But this means that the solutions of E_{i} = 0 for generic L_{i}(x) will have no solutions with z_{i} = 0. □

Note that if we choose y generically in Y, then we have chosen the associated y_{i} generically in Y_{i} for each i from n to 1. Thus we can assume that with an initial generic choice of parameters y, the behavior for each of the systems E_{i} = 0 is the behavior we expect with a generic choice of parameters on Y_{i}. One minor point remains. Given a generic choice of parameters y ∈ Y, generic behavior might not occur for the homotopy H_{i} with t ∈ (0, 1]. This is easily dealt with by a trick of Morgan and Sommese [28, 29, 30].

Assume that we have chosen the parameters y for the systems in the nonempty Zariski open set U of

[ \begin{matrix}a_{1} & a_{1,1} & ⋯ & a_{1,A} & λ_{1,1} & ⋯ & λ_{n,1} \\ ⋮ & ⋮ & ... & ⋮ & ⋮ & ... & ⋮ \\ a_{n} & a_{n,1} & ⋯ & a_{n,A} & λ_{1,n} & ⋯ & λ_{n,n} \end{matrix} ] ∈ ℂ^{n×(A+1)} × ℂ^{n×n}.  (26)

We want our homotopies H_{i}(x, z_{1}, …, z_{i}, t) = 0 to define algebraic sets that are flat over a Zariski open set containing (0, 1]. Numerically this means that we want to have generic behavior for each t ∈ (0, 1], i.e., the same number of isolated points, no components that do not correspond to components in any other fiber. Exploiting generic flatness in this way is the underlying approach of the work of Morgan and Sommese [28, 29, 30], to which we refer for more details. Since for all but a finite number of values of t ∈ ℂ generic behavior occurs, we can conclude that given a set of parameters y ∈ U as above, the homotopy

H_{η,i}(x, z_{1}, …, z_{i}, t) := (x, z_{1}, …, z_{i}, ηt) = [ f_{1}(x) + ∑_{j=1}^{i−1} λ_{1,j}z_{j} + tηλ_{1,i}z_{i} ; ⋮ ; f_{n}(x) + ∑_{j=1}^{i−1} λ_{n,j}z_{j} + tηλ_{n,i}z_{i} ; L_{1}(x) + z_{1} ; ⋮ ; L_{i−1}(x) + z_{i−1} ; tηL_{i}(x) + z_{i} ] , (27)

defines an algebraic set flat over a Zariski open set of ℂ containing (0, 1], for all but a finite set of η ∈ ℂ with |η| = 1. Using openness of flatness, we can absorb the η into the parameters we use. We still have a dense open set of general parameters, but the set is only Zariski open in the underlying real algebraic structure.

By Lemma 2 if the L_{i}, Λ_{i} are chosen randomly, then with probability one, E_{n}(x, z) = 0 has a solution with z = 0 only if f(x) is identically zero on ℂ^{n}.

Recall from the introduction that for i from n to 1,

1. Z_{i} denotes the solutions to E_{i} = 0 with z_{i} ≠ 0; and

2. X_{i−1} denote the limits with z_{i−1} = 0 of the paths of the homotopy H_{i}(t), from t = 1 to t = 0, starting at points of Z_{i}. By convention, the condition z_{i−1} = 0 is empty when i = 1.

Theorem 3. Let f be as above. Assume that f is not identically zero and that the Λ_{i} are chosen generically. If i is the largest integer with X_{i} nonempty, then the dimension of f^{−1}(0) is i. Moreover given any irreducible component W of f^{−1}(0) of dimension i, then, counting multiplicities, X_{i} contains deg(W) generic points of W.

Thus our algorithm achieves the same numerical goal of the algorithm of [35], but much more efficiently. Proof of Theorem 3. We use the notation of the proof of Lemma 2. By Lemma 2, it follows that there is a nonempty Zariski open set of points y ∈ Y, such that for each i from n to 1, all elements of the set Z_{i} of solutions of E_{i}(x, z_{1}, …, z_{i}) = 0 with z_{i} ≠ 0 are nonsingular and isolated. Thus for a dense set (Zariski open in the underlying real algebraic structure of y ∈ Y) the paths over t ∈ (0, 1] are smooth if they start at nonsingular isolated solutions of E_{i} = 0. By Lemma 5 and our random choice of y ∈ Y, each nonsingular isolated solutions (x^{∗}, z_{1}^{∗}, …, z_{i−1}^{∗}, 0) of E_{i−1} = 0 must start from a nonsingular isolated solution (x′, z′_{1}, …, z′_{i}) of E_{i} = 0.

If z′_{i} = 0 for a solution, then (z′_{1}, …, z′_{i}) = 0 by Lemma 2. Thus f(x′) = 0 and hence H_{i}(x′, 0, …, 0, t) = 0 identically in t, and hence (x^{∗}, z_{1}^{∗}, …, z_{i−1}^{∗}, 0) = (x′, 0, …, 0). Thus we know from Lemma 2 that for a Zariski open set of Y_{i} the solutions of the system E_{i} = 0 with z_{i} = 0 (and hence (z_{1}, …, z_{i}) = 0) are on irreducible components of the algebraic set f(x) = 0 of dimension i. Thus (x^{∗}, z_{1}^{∗}, …, z_{i−1}^{∗}, 0) = (x′, 0, …, 0) is on an irreducible component of the algebraic set of points with E_{i−1} = 0, which has dimension ≥ 1. This contradicts (x^{∗}, z_{1}^{∗}, …, z_{i−1}^{∗}) being an isolated nonsingular solution of E_{i−1} = 0. Thus the isolated nonsingular solutions of E_{i−1} = 0 must be endpoints of paths of homotopy H_{i} starting at points of Z_{i}.

Using Lemma 5 the above argument can be modified to show that isolated, but possibly singular, solutions of (x^{∗}, z_{1}^{∗}, …, z_{i−1}^{∗}, 0) of E_{i−1} = 0 must start from a nonsingular isolated solution (x′, z′_{1}, …, z′_{i}) of E_{i} = 0 with (z′_{1}, …, z′_{i}) ≠ 0 (and hence by Lemma 2 with z′_{i} ≠ 0).

The property of a generic system of the form E_{n}(f) that it has isolated nonsingular solutions is useful. In this direction, see [24].

3. A Worked Out Example

Consider the polynomial system f(x) = 0 and start system g(x) = 0.

f(x) = \{ \begin{matrix}x_{1}^{2}x_{2} = 0 \\ x_{1}^{2}(x_{2}^{2} + x_{1}) = 0 \end{matrix}  g(x) = \{ \begin{matrix}x_{1}^{3} − c_{1} = 0 \\ x_{2}^{4} − c_{2} = 0 \end{matrix}  (28)

To solve f(x) = 0 we trace D = 3 × 4 = 12 solution paths starting at the solutions of g(x) = 0. One path diverges to infinity, three paths converge to (0, 0), and the remaining eight paths end at the solution component x_{1} = 0, for some x_{2} ≠ 0. The condition numbers at the end of the paths do not allow us to decide which solution is isolated.

To embed f(x) = 0 we take a random hyperplane L(x) = a_{0} + a_{1}x_{1} + a_{2}x_{2} = 0 and choose two random complex constants, λ_{1} and λ_{2}:

E_{1}(x, z) = \{ \begin{matrix}x_{1}^{2}x_{2} + λ_{1}z = 0 \\ x_{1}^{2}(x_{2}^{2} + x_{1}) + λ_{2}z = 0 \\ a_{0} + a_{1}x_{1} + a_{2}x_{2} + z = 0 \end{matrix}  (29)

We can solve this system by tracing D = 12 solution paths, using a standard linear homotopy with the equations g(x) = 0, z − 1 = 0 as start system, with g(x) = 0 as in (28). Five paths diverge to infinity. Two paths go to the same solution with z = 0, which reveals the degree of the solution component x_{1}^{2} = 0. Note that geometrically this component corresponds to a pair of lines. The five remaining paths go to regular solutions with z ≠ 0.

To compute the possible remaining isolated solutions, we trace five solution paths starting at the five regular solutions of the system (29). In going with t from 1 to 0, we use the homotopy

H_{1}(x, z, t) = t (\{ \begin{matrix}x_{1}^{2}x_{2} + λ_{1}z = 0 \\ x_{1}^{2}(x_{2}^{2} + x_{1}) + λ_{2}z = 0 \\ a_{0} + a_{1}x_{1} + a_{2}x_{2} + z = 0 \end{matrix} ) + (1 − t) (\{ \begin{matrix}x_{1}^{2}x_{2} = 0 \\ x_{1}^{2}(x_{2}^{2} + x_{1}) = 0 \\ z = 0 \end{matrix} ) (30)

Three paths converge to (0, 0), two paths go to solutions on the component with x_{1} = 0 and x_{2} ≠ 0. End games are still needed to decide whether the solutions are isolated.

With our new method, 17 solution paths instead of 24 were traced, as 24 would have been the number of paths with an iteration of the procedure in [35].

Polyhedral root counting methods provide a generically sharp root count for polynomial systems. In particular, the mixed volume of the Newton polytopes of the system equals the number of roots in (ℂ^{∗})^{n}, ℂ^{∗} := ℂ \ {0}, for a system with generic coefficients. When the system has only few monomials with nonzero coefficients, then the mixed volume provides a much lower root count than the Bézout bounds based on the degrees of the polynomials.

For our type of applications, the distinction between solutions with z = 0 and z ≠ 0 is instrumental in identifying components of solutions. This difference does depend on the values of the coefficients of the original system and is not neglected by the ordinary mixed volume. So we will not miss any solutions with z = 0, but we may miss solution components for which some x_{i} = 0. Fortunately, extensions of the polyhedral methods that allow to count and compute all affine roots (that is in ℂ^{n} instead of (ℂ^{∗})^{n}) are covered amply in the literature (see [15], [19], [22], [25], [32, 33] and [34]). The key idea [25] is to add a random constant to every equation to shift the roots with zero components away from the coordinate axes. In removing these constants by continuation, all affine roots lie at the end of some path that starts at a root in (ℂ^{∗})^{n}.

Consider again the polynomial system f(x) = 0 in (28). Because of the first equation we immediately see that there cannot be any solution with all components different from zero. The direct application of polyhedral methods does not yield anything, since the mixed volume for f(x) = 0 equals zero. With affine polyhedral methods, we consider the system

f^{(0)}(x) = \{ \begin{matrix}x_{1}^{2}x_{2} + γ_{1} = 0 \\ x_{1}^{2}(x_{2}^{2} + x_{1}) + γ_{2} = 0 \end{matrix}  (31)

where γ_{1} and γ_{2} are randomly chosen complex numbers. Here the mixed volume equals five. Note the difference with the total degree D = 12. Letting the γ’s go to zero, three of the five paths converge to the origin, and the other two paths go to other solutions on the component.

To obtain information about the components, with polyhedral methods we consider the embedding

E_{1}(x, z) = \{ \begin{matrix}x_{1}^{2}x_{2} + γ_{1} + λ_{1}z = 0 \\ x_{1}^{2}(x_{2}^{2} + x_{1}) + γ_{2} + λ_{2}z = 0 \\ a_{0} + a_{1}x_{1} + a_{2}x_{2} + z = 0 \end{matrix}  (32)

This system has mixed volume equal to seven. We use this as start system to solve it with γ_{1} = 0 and γ_{2} = 0, following the paths that start at the seven solution paths of E_{1}(x, z) = 0. From the seven paths, five paths go to solutions with z ≠ 0, and the other two paths go to the component ending with z = 0.

4. Bertini’s Theorem and a Local Extension Theorem

Here is a weak, but convenient form of Bertini’s theorem, e.g., Fulton [17, Example 12.1.11]. For a further discussion of Bertini theorems, see also [4, §1.7].

Theorem 4 (Bertini). Let X be an algebraic manifold of dimension n, e.g., complex Euclidean space or complex projective space. Let L_{1}, …, L_{n} be line bundles on X. For i from 1 to n, let {s_{i,j}|j = 1, …, r_{i}} be a set of sections of L_{i}. Let B_{i} denote the set of common zeroes of the sections {s_{i,j}|1 ≤ j ≤ r_{i}}; and let B := ∪_{i}B_{i}. Then given general real or complex numbers {λ_{i,j}|j = 1, …, r_{i}; i = 1, …, n}, all solutions on X − B, of the system of equations

{ ∑_{j=1}^{r_{1}} λ_{1,j}s_{1,j} = 0 ; ⋮ ; ∑_{j=1}^{r_{n}} λ_{n,j}s_{n,j} = 0  (33)

are isolated and nonsingular.

In the proof of Theorem 3, we need an extension theorem, giving conditions ensuring that if we have a system of equations f(x, y) = 0 depending on parameters y with a multiplicity μ isolated solution x^{∗} of the system for some point y^{∗} in the parameter space, then there is a neighborhood U of x^{∗} such that for points y near y^{∗} in the parameter space, there are μ solutions counting multiplicities of the system f(x, y) = 0. Special cases are well known in the context of all polynomial systems, but we do not know a general reference in the numerical analysis literature. This sort of result is standard for algebraic geometers in the algebraic context, or within the German school of several complex variables in the complex analytic context. So we are simply explaining why this sort of result follows immediately from standard results in these fields. We work locally. All open sets are in the usual Euclidean topology, i.e., not in the Zariski topology.

A possibly nonreduced complex analytic space, X, is said to be a local complete intersection if given any point x ∈ X, there is a set U ⊂ X, open in the usual complex topology, that contains x, and such that

1. there is an embedding ϕ : U → B of U into an open ball in B ⊂ ℂ^{N} for some N := n + m > 0;

2. the ideal of the complex space ϕ(U) is defined by holomorphic functions g_{1}, …, g_{n}; and 3. the dimension of the maximal dimensional irreducible component of X through x is m.

Local complete intersections are very well behaved. One elementary, but important fact about them is that all the irreducible component of X through a given x ∈ X are equal dimensional. A less elementary, but important result is that for surjective holomorphic maps (respectively algebraic morphisms) from possibly nonreduced complex analytic spaces (respectively possibly nonreduced algebraic varieties) to complex manifolds (respectively algebraic manifolds), flatness of the morphism is equivalent to the openness of the morphism (which is equivalent to all the irreducible components of all the fibers of the morphism restricted to a given connected component of X being equal dimensional). This last result is a standard result for complex analytic spaces X, whose structure sheaf consists of local rings which are Cohen-Macaulay, e.g., Fischer [16, page 158]. Local complete intersections are among the simplest examples after manifolds of complex analytic spaces with Cohen-Macaulay local rings.

Let X be a connected n-dimensional complex manifold. Let Y be a connected m-dimensional complex manifold. Let

f(x, y) = [ f_{1}(x, y) ; ⋮ ; f_{n}(x, y) ] = 0  (34)

be a system of n holomorphic functions. Let x^{∗} be an isolated solution of the system f(x, y^{∗}) = 0 for a fixed value y^{∗} ∈ Y, i.e., assume that there is an open set O ⊂ X containing x^{∗} with x^{∗} the only solution of f(x, y^{∗}) = 0 on O. Assume that the multiplicity of x^{∗} is μ. This number, which is 1 exactly when the Jacobian of f(x, y^{∗}) is invertible at x^{∗}, is equal to

dim_{ℂ} O_{X|x^{∗}}/J (f_{1}(x, y^{∗}), …, f_{n}(x, y^{∗}))  (35)

where O_{X|x^{∗}} is the local ring of convergent power series on X centered at the point x^{∗}, and J (f_{1}(x, y^{∗}), …, f_{n}(x, y^{∗})) is the ideal in O_{X|x^{∗}} generated by the functions f_{1}(x, y^{∗}), …, f_{n}(x, y^{∗}).

Lemma 5 (Local Extension Lemma). Let X, Y, f, x^{∗}, and y^{∗} be as above. There are open neighborhoods U of x^{∗} ∈ X and V of y^{∗} ∈ Y such that for any y ∈ V there exist μ isolated solutions (counting multiplicities) of f(x, y) = 0 on U.

Proof. Let dim_{{(x^{∗},y^{∗})}} Z denote the dimension of an analytic set Z ⊂ X × Y at the point (x^{∗}, y^{∗}). Let X′ denote the zero set of f_{1}, …, f_{n} on X × Y. Since there are n-functions, the dimension of each irreducible component of X′ is at least m, and in particular

dim_{(x^{∗},y^{∗})} X′ ≥ m.  (36)

Since x^{∗} is isolated we have that dim_{(x^{∗},y^{∗})}(X × {y^{∗}}) ∩ X′ = 0. Since X × Y is smooth, we have

dim_{(x^{∗},,y^{∗})}(X × {y^{∗}}) ∩ X′ ≥ dim_{(x^{∗},,y^{∗})} X × {y^{∗}}) + dim_{(x^{∗},,y^{∗})} X′ − dim_{(x^{∗},,y^{∗})} X × Y
= n + dim_{(x^{∗},,y^{∗})} X′ − n − m.

Thus we conclude that

dim_{(x^{∗},y^{∗})} X′ = m.  (37)

Thus the union X′′ of the components of X′ passing through (x^{∗}, y^{∗}) has pure dimension m. It follows, e.g., use Gunning [20, Theorem 16], that there are neighborhoods U of x^{∗} ∈ X and V of y^{∗} ∈ Y such that the projection π of X′′∩(U ×V) to V is proper and finite. By shrinking U and V we can assume that X′′∩(U ×V) = X′∩(U ×V) and that X′∩(U ×{y^{∗}}) = (x^{∗}, y^{∗}).

Now X′ ∩ (U × V) is a local complete intersection, and thus since the map π_{X′} : X′ → V is proper with finite fibers, we conclude by the discussion before the discussion of flatness, that this map is flat. Thus the direct image E := π_{∗}(O_{X′}) of the structure sheaf of X′ is locally free, e.g., Fischer [16, Proposition 3.13]. On U × V, X′ is defined by the functions f_{i}(x, y) = 0, for i = 1, …, n. By definition, at a point (x′, y′) ∈ (U × V) ∩ X′ we have O_{X′|{(x′,y′)}} is

O_{(U×V)|{(x′,y′)}}/J (f_{1}(x, y), …, f_{n}(x, y))  (38)

where O_{(U×V)|{(x′,y′)}} is the local ring of convergent power series on U ×V centered at {(x′, y′)} and J (f_{1}(x, y), …, f_{n}(x, y)) denotes the ideal in O_{(U×V)|{(x′,y′)}} generated by f_{1}, …, f_{n}.

The statement that E is locally free, is equivalent to ranks of E at different points of V being equal. The rank of E at a point y′ ∈ V is by definition

dim_{ℂ} E_{y′} := E/ (𝔪_{y′} ⊗_{O_{Y |y′}} E)  (39)

where 𝔪_{y′} is the maximal ideal generated O_{Y |y′} consisting of convergent power series vanishing at y′. Thus comparing to the definition of the multiplicity of the solution x^{∗} of f(x, y^{∗}) being μ, we see that the rank of E is μ. Unwinding the definition of π_{∗}, E_{y′} for any fixed y′ ∈ V is the direct sum of the modules

O_{X|x′}/J (f_{1}(x, y′), …, f_{n}(x, y′)),  (40)

with index set the set of distinct points x′ ∈ U with f(x′, y′) = 0. This proves the assertion of the Lemma. □

5. Applications and Computational Experiences

We have conducted a systematic set of experiments with PHCpack [36] on four case studies of familiar polynomial systems. Little symbolic manipulation of polynomials is needed to set up homotopies derived from the embedding presented in this paper.

To avoid wasting paper, we have omitted the algebraic formulations of the systems, which can either be found electronically in the database maintained on the web sites of the second author, or can be consulted in the cited literature. Although the systems are academic examples and only interesting for benchmarking purposes, we try to indicate the relevance to their application fields.

Unless stated otherwise, all reported timings concern a 450MHz Intel Pentium machine with 1Gb of main memory and 1Gb swap space, running Debian GNU/Linux. We observe that the initial stage is always the most expensive one. Dealing with components of solutions is a much harder problem than just approximating the isolated solutions.

With a faster computer one can attack larger problems, but also the quality of the software matters. Concerning this latter aspect we want to point out that no special-purpose software has been developed that exploits the structure of the embedded systems. The continuation could for instance go faster if one eliminated explicitly some variables using the linear hyperplanes that were introduced in the slicing.

5.1. A planar four-bar mechanism. Four-bar mechanisms are ubiquitous in mechanical design. The 4-variable polynomial system that was derived and solved in [31] has total degree 256 and lowest multi-homogeneous Bézout bound 96. The mixed volume equals 80, avoiding the calculation of zero component solutions. There is a solution component of dimension two, with sum of degrees equal to two. The lowest multi-homogeneous Bézout bound of E_{2} equals 240 whereas the mixed volume equals 96.

We summarize the computational experiments in Table 1, whose format goes as follows. The cascade of homotopies starts with solving the start system g by polyhedral methods since mixed volumes are sharper than the bounds based on the degrees. To solve the system E_{k} we must trace #paths solution curves. This number is partitioned into ones on the component (z = 0), regular finite solutions (z ≠ 0), and diverging ones (→ ∞). Observe the take over of the (z ≠ 0)-solutions to the next row. The number in the column with heading z ≠ 0 on the E_{0}-row is the number of isolated solutions of the original system. In the last column we list cpu times.

The information in Table 1 is useful to make some rough comparisons with other homotopy methods. For instance, if we are only interested in the isolated roots, the cost would be of the same magnitude as on the E_{0}-row. To compare with the procedure of Sommese and Wampler [35], we could add up the entries in the columns with headings z = 0 and → ∞ to the #paths in the next row, simulating the process of no solution recycling.

system #paths z = 0 z ≠ 0 → ∞  cpu time
g, start  96  0  96  0  10s 880ms
E_{2}  96  2  68  26 4m 47s 830ms
E_{1}  68  12  56  0  18s 370ms
E_{0}  56  20  36  0  11s 10ms
total  316  34 256  26 5m 28s 90ms

Table 1. The number of paths and timings for a planar four-bar mechanism.

We must note that there exists an isotropic formulation of this problem [37] for which the mixed volume is an exact root count for the 36 isolated solutions. The black-box solver of PHC needs only 11s 710ms to solve that system.

5.2. Constructing Runge-Kutta formulas. The following system is due to Butcher [10] and was used as test problem in [9]. The system has seven variables and a three-dimensional component. Its total degree is 4608 and mixed volume equals 24. After slicing and embedding thrice, the resulting 10-variable system has mixed volume 247. We summarize the running statistics in Table 2.

system #paths z = 0 z ≠ 0 → ∞  cpu time
g, start  247  0 247  0 17m 13s 770ms
E_{3}  247  3 193  51 9m 57s 60ms
E_{2}  193  15 161  17 1m 59s 470ms
E_{1}  161  11  72  78 4m 57s 600ms
E_{0}  72  0  4  68 2m 24s 650ms
total  920  29 677 214 34m 32s 550ms

Table 2. The number of paths and timings for Butcher’s system.

5.3. Modular equations for special algebraic number fields. The authors of [18] discuss the finding of the solution components of the system that was derived in [11]. That system has a component of dimension two, and in [11] the authors reasoned that any three of the four polynomials generate the same ideal. The challenge posed in [11] to computer algebra was to verify this claim.

Although we believe that nowadays, more than 10 years after [11], some computer algebra packages might be strong enough to compute a Gröbner basis for this problem, we would like to indicate how one could answer the challenge numerically. The approach we propose is to replace any of the four equations by the same random hyperplane. This operation cuts the dimension by one. So we have to embed each of the four systems with the same random constants. We solve one of them, this takes about 2 minutes, and check whether the generic points we find satisfy the other equations.

5.4. On Fourier transforms: the cyclic n-roots problem. In [7] an example was given which stems from the problem of finding all “bi-equimodular” vectors x ∈ ℂ^{n}, i.e.: all x with coordinates of constant absolute value such that the Fourier transform of x is a vector with coordinates of constant absolute value, see also [5] and [6]. This system was popularized in [12], and is by far the most notorious benchmark problem in polynomial system solving. Other references are [3], [8] and [14].

In [26], the following conjectures of Ralf Fröberg are mentioned. If n has a quadratic divisor, then there are infinitely many solutions. If the number of solutions is finite, then this number equals ( \begin{matrix}2n − 2 \\ n − 1 \end{matrix} ). Uffe Haagerup [21] has proven that for n prime, the number

of roots is always finite, and equals indeed \frac{(2n−2)!}{(n−1)!^{2}}.

In Table 3 we summarize the results of the embedding for the cyclic 8-roots problem, which has a solution component of dimension one. The original system has mixed volume equal to 2,560. The embedded system has mixed volume 4,176. See Table 3 for the summary of the runs. We remark that if one is only interested in the isolated roots, one only has to trace 2,560 solution paths.

system #paths z = 0 z ≠ 0 → ∞  cpu time
g, start 4,176  0 4,176  0 1h 30m 26s 800ms
E_{1}  4,176 144 3,975  57 1h 10m 3s 930ms
E_{0}  3,975 495 1,152 2,328 7h 7m 55s 580ms
total 12,327 639 9,303 2,385 9h 48m 26s 310ms

Table 3. The number of paths and timings for cyclic 8-roots.

Before developing the embedding method, we attacked the cyclic 8-roots problem first with the slicing method of Sommese and Wampler. The computation of generic points on the solution component required the tracing of 10,940 solution paths (which takes now 4,176 paths). For this task, the black-box solver of PHC needed 4 days and nights 15h 6m 49s 392ms cpu time on a SUN workstation. Although the Linux machine we used is somewhat faster than that SUN workstation, it certainly does not speed up things 50 times!

We end this section with another illustration on the importance of having good polynomial equations for the same problem. Thanks to a trick of John Canny, described in [13], there exists a reduced version of the cyclic n-roots problem. This version reduces the number of roots by eight and yields a significant saving in the path tracking. For instance, the mixed volume of the reduced cyclic 8-roots problem equals 320, instead of 2,560 for the original problem. Table 4 summarizes the computational experiments.

system #paths z = 0 z ≠ 0 → ∞  cpu time
g, start  775  0 775  0 10m 1s 860ms
E_{1}  775  18 744  13 11m 23s 680ms
E_{0}  744 445 144 155 14m 52s 340ms
total  2,294 463 1,663 168 36m 17s 860ms

Table 4. The number of paths and timings for reduced cyclic 8-roots.

Finally, this reduction also allows us to treat the cyclic 9-roots system that has a twodimensional component of solutions. The sum of the degrees of the reduced cyclic 9-roots problem equals two. The mixed volume of embedded system E_{2} equals 4,044, computed in cpu time 2h 48m 4s 320ms. The polyhedral continuation for the start system requires 1h 4m 5s 130ms and the path following to solve E_{2} takes 2h 15m 4s 540ms. In total, the black-box solver of PHCpack needs 6h 7m 45s 270ms to solve E_{2}.

6. Conclusions and Future Directions

In the paper we have presented an embedding of polynomial systems to compute generic points on components of solutions. The homotopies recycle solutions while moving down to the isolated solutions. The major advantage compared to the method in [35] is that fewer solution paths diverge. We provide practical evidence for the efficiency of the embedding technique.

Besides the reported progress on homotopy methods, we want to point out that our embedding technique might be of independent interest, in the sense that it allows us to focus on one component of a particular dimension. The main novelty of the embedding consists in the fact that the solutions of the embedded system that do not lie on the focused component provide meaningful information about the components of lower dimension. We hope that besides homotopy methods, other solvers for polynomial systems could benefit from our embedding technique.

It is sometimes said that every scientific problem solved generates three other questions. Here we briefly indicate three different directions for future research. First of all, the embedding leads to higher-dimensional systems with special structure. It would be worthwhile to exploit this structure to obtain a more efficient mixed-volume computation to alleviate the complexity of the first stage. Secondly, the computed generic points on the components are an interesting starting point for future computational explorations of the solution components. Note that continuation is a very natural tool to scan the component while wiggling the added hyperplanes. Lastly, in the final stage of the cascade of embedded homotopies, it remains difficult to separate solution paths converging to components from other possibly ill-conditioned isolated solutions. End games that produce certificates in the form of a random hyperplane through a nonisolated solution will have to be developed.

References

[1] E.L. Allgower and K. Georg. Numerical Continuation Methods, an Introduction, volume 13 of Springer Ser. in Comput. Math. Springer–Verlag, Berlin Heidelberg New York, 1990.

[2] E.L. Allgower and K. Georg. Numerical Path Following. In P.G. Ciarlet and J.L. Lions, editors, Techniques of Scientific Computing (Part 2), volume 5 of Handbook of Numerical Analysis, pages 3–203. North-Holland, 1997.

[3] J. Backelin and R. Fröberg. How we proved that there are exactly 924 cyclic 7-roots. In Proceedings of ISSAC-91, pages 101–111. ACM, 1991.

[4] M. Beltrametti and A.J. Sommese. The adjunction theory of complex projective varieties, volume 16 of Expositions in Mathematics. Walter De Gruyter, Berlin, 1995.

[5] G. Björck. Functions of modulus one on Z_{p} whose Fourier transforms have constant modulus. In Proceedings of the Alfred Haar Memorial Conference, Budapest, volume 49 of Colloquia Mathematica Societatis János Bolyai, pages 193–197. 1985.

[6] G. Björck. Functions of modulus one on Z_{n} whose Fourier transforms have constant modulus, and “cyclic n-roots”. In J.S. Byrnes and J.F. Byrnes, editors, Recent Advances in Fourier Analysis and its Applications, volume 315 of NATO Adv. Sci. Inst. Ser. C: Math. Phys. Sci., pages 131–140. Kluwer, [7] G. Björck and R. Fröberg. A faster way to count the solutions of inhomogeneous systems of algebraic

equations, with applications to cyclic n-roots. J. Symbolic Computation, 12(3):329–336, 1991.

[8] G. Björck and R. Fröberg. Methods to “divide out” certain solutions from systems of algebraic equations, applied to find all cyclic 8-roots. In M. Gyllenberg and L.E. Persson, editors, Analysis, Algebra and Computers in Math. research, volume 564 of Lecture Notes in Applied Mathematics, pages 57–70. Marcel Dekker, 1994.

[9] W. Boege, R. Gebauer, and H. Kredel. Some examples for solving systems of algebraic equations by calculating Groebner bases. J. Symbolic Computation, 2:83–98, 1986.

[10] C. Butcher. An application of the Runge-Kutta space. BIT, 24:425–440, 1984.

[11] H Cohn and J. Deutsch. An explicit modular equation in two variables for Q[sqrt(3)]. Math. Comp., 50:557–568, 1988.

[12] J Davenport. Looking at a set of equations. Technical report 87-06, Bath Computer Science, 1987.

[13] I.Z. Emiris. Sparse Elimination and Applications in Kinematics. PhD thesis, Computer Science Division, Dept. of Electrical Engineering and Computer Science, University of California, Berkeley, 1994.

[14] I.Z. Emiris and J.F. Canny. Efficient incremental algorithms for the sparse resultant and the mixed volume. J. Symbolic Computation, 20(2):117–149, 1995. Software available at

http://www.inria.fr/saga/emiris.

[15] I.Z. Emiris and J. Verschelde. How to count efficiently all affine roots of a polynomial system. Rapport de recherche no. 3212, INRIA, 1997. To appear in Discrete Applied Mathematics.

[16] G. Fischer. Complex Analytic Geometry, volume 538 of Lecture Notes in Mathematics. Springer–Verlag, New York, 1976.

[17] W. Fulton. Intersection Theory, volume (3) 2 of Ergeb. Math. Grenzgeb. Springer–Verlag, Berlin, 1984.

[18] A. Galligo and C. Traverso. Practical determination of the dimension of an algebraic variety. In E. Kaltofen and S.M. Watt, editors, Computers and Mathematics, pages 46–52. Springer–Verlag, Berlin, New York, 1989.

[19] T. Gao, T.Y. Li, and X. Wang. Finding isolated zeros of polynomial systems in C^{n} with stable mixed volumes. To appear in J. of Symbolic Computation.

[20] R.C. Gunning. Lectures on Complex Analytic Varieties: The Local Parameterization Theorem. Princeton University Press, Princeton, 1970.

[21] U. Haagerup. Orthogonal MASA’s in the n x n - matrices, complex Hadamard matrices and cyclic n-roots. Talk at the conference on “Operator Algebras and Quantum Field Theory”, Rome, July 1-6, 1996.

[22] B. Huber and B. Sturmfels. Bernstein’s theorem in affine space. Discrete Comput. Geom., 17(2):137–141, 1997.

[23] T.Y. Li. Numerical solution of multivariate polynomial systems by homotopy continuation methods. Acta Numerica, 6:399–436, 1997.

[24] T.Y. Li, T. Sauer, and J.A. Yorke. The cheater’s homotopy: An efficient procedure for solving systems of polynomial equations. SIAM J. of Numerical Analysis, 26:1241–1251, 1989.

[25] T.Y. Li and X. Wang. The BKK root count in C^{n}. Math. Comp., 65(216):1477–1484, 1996.

[26] H.M. Möller. Gröbner bases and numerical analysis. In B. Buchberger and F. Winkler, editors, Gröbner Bases and Applications, volume 251 of London Mathematical Lecture Note Series, pages 159–178. Cambridge University Press, 1998.

[27] A.P. Morgan. Solving polynomial systems using continuation for engineering and scientific problems. Prentice-Hall, Englewood Cliffs, N.J., 1987.

[28] A.P Morgan and A.J. Sommese. A homotopy for solving general polynomial systems that respects m-homogeneous structures. Appl. Math. Comput., 24:101–113, 1987.

[29] A.P Morgan and A.J. Sommese. Computation of all solutions to polynomial systems using homotopy continuation. Appl. Math. Comput., 24:115–138, 1987. Errata: Appl. Math. Comput. 51 (1992), p. 209.

[30] A.P. Morgan and A.J. Sommese. Coefficient-parameter polynomial continuation. Appl. Math. Comput., 29:123–160, 1989. Errata: Appl. Math. Comput. 51 (1992), p. 207.

[31] A.P. Morgan and C.W. Wampler. Solving a planar four-bar design problem using continuation. ASME J. of Mechanical Design, 112:544–550, 1990.

[32] J.M. Rojas. A convex geometric approach to counting the roots of a polynomial system. Theoret. Comput. Sci., 133(1):105–140, 1994. Extensions and corrections available at http://www.cityu.edu.hk/ma/staff/rojas/.

[33] J.M. Rojas. Toric intersection theory for affine root counting. Journal of Pure and Applied Algebra, 136(1):67–100, 1999.

[34] J.M. Rojas and X. Wang. Counting affine roots roots of polynomial systems via pointed Newton polytopes. J. Complexity, 12:116–133, 1996.

[35] A.J. Sommese and C.W. Wampler. Numerical algebraic geometry. In J. Renegar, M. Shub, and S. Smale, editors, The Mathematics of Numerical Analysis, volume 32 of Lectures in Applied Mathematics, pages 749–763, 1996. Proceedings of the AMS-SIAM Summer Seminar in Applied Mathematics, Park City, Utah, July 17-August 11, 1995, Park City, Utah.

[36] J. Verschelde. PHCpack: A general-purpose solver for polynomial systems by homotopy continuation. To appear in ACM Transactions on Mathematical Software. Paper and software available at

http://www.math.msu.edu/~jan.

[37] C.W. Wampler. Isotropic coordinates, circularity and Bezout numbers: planar kinematics from a new perspective. Proceedings of the 1996 ASME Design Engineering Technical Conference, Irvine, California August 18-22, 1996. CD-ROM edited by McCarthy, J.M., American society of mechanical engineers. Also available as GM Technical Report, Publication R&D-8188., 1996.

Department of Mathematics, University of Notre Dame, Notre Dame, IN 46556, U.S.A.

E-mail address: sommese@nd.edu and url:  http://www.nd.edu/∼sommese

Mathematical Sciences Research Institute, 1000 Centennial Drive, Berkeley, CA 94720-5070, U.S.A.

Department of Mathematics, Michigan State University, East Lansing, MI 48824-1027, U.S.A.

E-mail address: jan@msri.org, jan@math.msu.edu or jan.verschelde@na-net.ornl.gov, and urls:
http://www.msri.org/people/members/jan http://www.mth.msu.edu/∼jan