HOMOTOPIES FOR INTERSECTING SOLUTION COMPONENTS OF POLYNOMIAL SYSTEMS^{∗}

ANDREW J. SOMMESE^{†}, JAN VERSCHELDE^{‡}, AND CHARLES W. WAMPLER^{§}

Abstract. We show how to use numerical continuation to compute the intersection C = A∩B of two algebraic sets A and B, where A, B, and C are numerically represented by witness sets. Enroute to this result, we first show how to find the irreducible decomposition of a system of polynomials restricted to an algebraic set. The intersection of components A and B then follows by considering the decomposition of the diagonal system of equations u − v = 0 restricted to {u, v} ∈ A × B. One offshoot of this new approach is that one can solve a large system of equations by finding the solution components of its subsystems and then intersecting these. It also allows one to find the intersection of two components of the two polynomial systems, which is not possible with any previous numerical continuation approach. 2000 Mathematics Subject Classification. Primary 65H10; Secondary 13P05, 14Q99, 68W30. Key words and phrases. Components of solutions, embedding, generic points, homotopy continuation, irreducible components, numerical algebraic geometry, polynomial system.

1. Introduction. In a series of papers [12, 13, 14, 15, 17], we have proposed numerical continuation algorithms that use witness sets as the basic construct for representing solution components of a system of polynomial equations on ℂ^{N}. Witness sets are the central concept of a young subject that we call numerical algebraic geometry, which uses numerical continuation [1, 2] and generalizes earlier work in computing isolated solutions of polynomial systems [8, 9]. The main concern of this paper is to provide an algorithm for computing the intersection of two solution components A, B from two possibly identical polynomial systems f, g, whose witness sets have been given. It is important to realize that naively combining f, g into one system h = {f, g} is not sufficient, even if we were willing to put aside the potentially prohibitive size of the combined system. For example, suppose A is the line x_{2} = 0 as a solution component of f(x) = x_{1}x_{2} and B is the line x_{1} − x_{2} = 0 as a solution component of g(x) = x_{1}(x_{1} − x_{2}). Then, A ∩ B, which is the isolated point (0, 0), does not appear as an irreducible component of the system h = {f, g}.

Questions involving intersection of components arise naturally in applications. Just as a single polynomial in one variable has multiple roots, a system of polynomial equations in several variables can have multiple solution components; these components can even appear at different dimensions (points, curves, surfaces, etc.) from the same set of equations. We may wish to find the intersection of just one of those components with another algebraic set. In our new approach, only the degrees of the components being intersected come into play in the determination of the number of paths followed by the homotopies that we use. This is important since the degree of

∗The authors acknowledge the support of the Volkswagen-Stiftung (RiP-program at Oberwolfach). †Department of Mathematics, University of Notre Dame, Notre Dame, IN 46556-4618, USA

Email: sommese@nd.edu URL: http://www.nd.edu/˜ sommese. This material is based upon work supported by the National Science Foundation under Grant No. 0105653; and the Duncan Chair of the University of Notre Dame.

‡Department of Mathematics, Statistics, and Computer Science, University of Illinois at Chicago, 851 South Morgan (M/C 249), Chicago, IL 60607-7045, USA Email: jan@math.uic.edu or jan.verschelde@na-net.ornl.gov URL: http://www.math.uic.edu/˜ jan. This material is based upon work supported by the National Science Foundation under Grant No. 0105739 and Grant No. 0134611.

§General Motors Research and Development, Mail Code 480-106-359, 30500 Mound Road, Warren, MI 48090-9055, U.S.A. Email: Charles.W.Wampler@gm.com.

a component of a given system of polynomials is typically much less than the number of paths required to find even all isolated solutions of the given system.

Viewed another way, the intersection operation is required for a Boolean algebra of constructible algebraic sets; a complete Boolean algebra also requires the operations of union and complement. Suppose W is a witness set for a component X. There are several probability-one algorithms for deciding if a point x ∈ ℂ^{N} is a member of X, using numerical continuation and the data in W. (We review witness sets and membership tests in §2.) The complement operation is just the logical inversion of a membership test, and the union operation is just a union of witness sets, utilizing membership tests to eliminate duplications. However, the operation of intersection is more difficult.

In our previous work, we have shown how to find the solution set of a system of polynomial equations as a union of witness sets, and further, we have shown how to decompose these into witness sets for the irreducible components. Said another way, this solves the problem of intersecting a collection of hypersurfaces defined by polynomial equations. But this does not give us an effective means of computing the intersection of two components represented by witness sets.

Our first step in creating an algorithm for the intersection of components is to generalize an earlier algorithm for generating the witness sets for the solution set of a system of polynomial equations on ℂ^{N}. The generalization instead considers the polynomial equations restricted to an algebraic set. The intersection of components A and B then follows by considering the decomposition of the diagonal system of equations u − v = 0 restricted to {u, v} ∈ A × B. Hence, we call the intersection algorithm the diagonal homotopy.

This paper is organized as follows. First, in §2, we review the definition of a witness set and its role in finding the numerical irreducible decomposition of the solution set of a system of polynomial equations. In §3 we introduce a slight generalization of the randomization procedure of [17], and in §4 we give a general construction of homotopies. These sections give the basic definitions and results that will be needed later in the article.

The original algorithm for constructing witness supersets was given in [17]. A much more efficient algorithm for constructing witness supersets was given in [12] by means of an embedding theorem. In §5, we show how to carry out the generalization of [12] to the case of a system of polynomials on a pure N-dimensional algebraic set X ⊂ ℂ^{m}, i.e., an algebraic subset of ℂ^{m} all of whose irreducible components are N-dimensional. We call this the “abstract embedding theorem” because it does not rely on any specific numerical description of X. In this generality we lose some control of multiplicities. However, since our main objective is to find the underlying reduced algebraic solution components, this loss of multiplicity information is of minor importance.

In §6 we show how to implement the abstract embedding theorem numerically. We need only the information about X that would be produced by the algorithm for the numerical irreducible decomposition of a polynomial system f, for which X is an irreducible component of the solution set of f.

In §7 we specialize to the situation where we have two polynomial systems f and g on ℂ^{N} and we wish to describe the irreducible decompositions of A ∩ B where A is an irreducible component of V (f) and B is an irreducible component of V (g). Computational experiments are discussed in §8.

In Appendix §A, we give some further discussion of the method of constructing homotopies described in §4.

In Appendix §B, we give the proof of Theorem 5.1 from §5. We would like to thank the referees for their helpful suggestions.

2. Witness Sets. We begin by reviewing the basics of numerical algebraic geometry, wherein the most fundamental concept is a witness set. We refer the reader to [13, 15, 17] for more details on irreducible components, the irreducible decomposition, and reduced algebraic sets.

Given a system of polynomials on ℂ^{N}

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

we denote the underlying point set {x ∈ ℂ^{N} | f(x) = 0} by V (f), i.e., the algebraic set f^{−1}(0) (with all multiplicity information that comes with f^{−1}(0) ignored). A pure i-dimensional algebraic subset X of V (f) is a subset X ⊂ V (f) equal to the closure of a union of i-dimensional connected components of the smooth points of V (f). We emphasize that X is reduced, i.e., that we are ignoring the multiplicity of X within f^{−1}(0). We represent X numerically by a witness set, defined as follows.

Definition 2.1. A witness set for a pure i-dimensional algebraic set X ⊂ V (f) ∈ ℂ^{N} consists of:

1. the dimension, i, of X;

2. the polynomial system f(x);

3. a general (N − i)-dimensional affine linear subspace L_{N−i} ⊂ ℂ^{N}; and

4. the set of deg X distinct points X = L_{N−i} ∩ X.

In other words, a witness set W for X is the ordered set W = {i, f, L_{N−i}, X }. We use the notation V (W) to denote the component represented by W; in the current context V (W) = X.

This definition is useful because it allows us to numerically represent and manipulate the irreducible decomposition of the solution set of a polynomial system. Let us quickly review that concept before describing our new results. Everything we say is over the complex numbers, e.g., even if the polynomials have real coefficients, we always deal with the sets of solutions on complex Euclidean space.

We start with a system of polynomials f on ℂ^{N} as in (2.1) above. Let V (f) denote the set of solutions of f on ℂ^{N}, i.e., the set of points x ∈ ℂ^{N} such that f(x) = 0. The set Z := V (f) is an affine algebraic set and decomposes into a union of distinct irreducible components. Recall that an algebraic set X is irreducible if and only if the Zariski open dense subset of manifold points on X is connected. We have the decomposition

Z = ⋃_{i=1}^{dim Z} Z_{i} = ⋃_{i=1}^{dim Z} ( ⋃_{j∈I_{i}} Z_{i,j}) , (2.2)

where

1. For each i, Z_{i} := (∪_{j∈I_{i}}Z_{i,j});

2. the sets I_{i} are finite and each Z_{i,j} is irreducible of dimension i;

3. Z_{i,j} is not contained in a union of a collection of the Z_{a,b} unless Z_{i,j} occurs in the collection.

Any collection of irreducible components Z_{i,j} having the same dimension, i, can be numerically represented by a witness set. A “numerical irreducible decomposition” is a list having one witness set W_{i,j} for the reduction of each irreducible component Z_{i,j}.

In a series of papers [12, 13, 14, 15, 17], we showed how to compute a numerical irreducible decomposition of Z := V (f). The approach is to intertwine two numerical algorithms: a witness generating algorithm, which finds a superset of witness points for each pure-dimensional algebraic set Z_{i}, and a decomposition algorithm, which eliminates spurious points from the superset and breaks it into irreducible components. To be more precise, at each dimension i = 0, …, dim Z, the witness generating algorithm gives a finite set of points Ŵ_{i} satisfying L_{N−i} ∩Z_{i} ⊂ Ŵ_{i} ⊂ L_{N−i} ∩(∪_{j≥i}Z_{j}), where L_{N−i} ⊂ ℂ^{N} is a general (N −i)-dimensional affine linear subspace. The second algorithm decomposes the Ŵ_{i}. Precisely:

1. Ŵ_{i} decomposes into the disjoint union

J_{i} ∪ ( ⋃_{j∈I_{i}} Z_{i,j}) , (2.3)

where J_{i} ⊂ ∪_{k>i}Z_{k} and Z_{i,j} consists of the deg Z_{i,j} points of L_{N−i} ∩ Z_{i,j}; and

2. J_{dim Z} = ∅. The points Z_{i,j} along with the dimension i, the system of equations f, and the linear subspace L_{N−i}, form a witness set for the irreducible component Z_{i,j}.

The key theoretical advance of this paper is to observe that the previous algorithms for numerical irreducible decomposition still work with restrictions of a polynomial system to a pure-dimensional algebraic set. Only the first algorithm constructing the witness point supersets Ŵ_{i} needs to be generalized. The decomposition algorithms starting with the witness point supersets Ŵ_{i} are proved in the papers [13, 14, 15] in sufficient generality to cover the present situation.

The above implicitly assumed that the components Z_{i,j} are reduced, i.e., of multiplicity one in f^{−1}(0). The algorithms in [12, 17] in fact produce sets Z_{i,j} consisting of deg Z_{i,j} distinct points each repeated μ_{i,j} times, where μ_{i,j} is greater than or equal to the multiplicity of Z_{i,j} in f^{−1}(0). Moreover the multiplicity of Z_{i,j} is one if and only if μ_{i,j} = 1. Unfortunately, in the algorithm in this article we can only assert that μ_{i,j} > 0 for any irreducible component Z_{i,j}.

As we mentioned in the introduction, an important aspect of a witness set X is that we can use it to test a point for membership in the algebraic set X = V (X) that X represents. This stems from the fact that we can sample X by continuously perturbing the linear slice and numerically tracking its intersection with X, starting from the witness points in X. Several different membership tests can be employed. At one expensive extreme, by sampling and fitting, we might compute a set of polynomials, whose set of common zeroes is exactly X. A much more efficient “probability-one” test for a point x ∈ ℂ^{N} to be in X is whether the pullback from ℂ^{dim X+1} of a deg X defining polynomial for π(X) ⊂ ℂ^{dim X+1} is zero on x, where π is a general linear projection from ℂ^{N} to ℂ^{dim X+1}. Finally, a very different sort, and quite efficient test, for x ∈ ℂ^{N} to be in X is to see whether x is one of the images of the set X under the homotopy taking L_{N−i} to a general (N − i)-dimensional linear subspace of ℂ^{N} that contains x. This test depends on the real-one-dimensional path between the general linear subspaces to remain general, which occurs with probability one.

3. Randomizing Systems. Randomization is a key element of our approach. This section introduces some notation for randomized systems and gives a lemma describing their most important properties. Given a system of n equations defined on ℂ^{N}, as in Eq.(2.1), and a positive integer k, we define a randomization operation

ℜ(f(x); k) := Λf(x), Λ ∈ ℂ^{k×n}, (3.1)

where Λ is chosen generically from ℂ^{k×n}. Note that k does not have to equal n. Gaussian elimination does not change the ideal generated by Λf(x). Therefore, when k ≥ n, Λf(x) is equivalent to the system consisting of f(x) plus k − n identically zero equations. For the same reason, when k ≤ n, Λf(x) is equivalent to the system

[ I_{k} R ] f(x), R ∈ ℂ^{k×(n−k)}, (3.2)

where I_{k} is the k × k identity matrix. Consequently, we may without loss of generality assume that ℜ(f(x); k) is of the form of Eq.(3.2) with a generic choice of R ∈ ℂ^{k×(n−k)}. This form allows us to take some advantage of the original equations. For example, if k = N and the original equations had total degrees d_{1} ≥ d_{2} ≥ ⋯ ≥ d_{n}, then the total degree of the original form of ℜ(f(x); N) is d_{1}^{N}, but the total degree of the modified form is d_{1}d_{2} ⋯ d_{N}.

The following lemma gives the main properties of randomization.

Lemma 3.1. Let

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

be a system of restrictions of n polynomials on ℂ^{m} to a pure N-dimensional affine algebraic set X ⊂ ℂ^{m}. Assume that k ≤ min{n, N}. Assume that f does not vanish on any component of X. The following conclusions follow.

1. The dimension of any component of V (ℜ(f(x); k)) is ≥ N − k.

2. The irreducible components of V (ℜ(f(x); k)) and V (f) of dimension greater than N −k are the same, and the irreducible (N −k)-dimensional components of V (f) are components of V (ℜ(f(x); k)).

Proof. This variant of Bertini’s Theorem follows by the same type of reasoning as the analogous result in [12, 17] for systems of N polynomials on ℂ^{N}. For the convenience of the reader, we give a proof.

The first conclusion is simply [11, Corollary 3.14].

Since V (f) ⊂ V (ℜ(f(x); k)), the second conclusion will follow if we show that all irreducible components V (ℜ(f(x); k)) ∩ (X \ V (f)) have dimension N − k. Thus it suffices to show that if V (f) is empty, then it follows that all irreducible components V (ℜ(f(x); k)) have dimension N − k. This is immediate from Theorem B.1.

4. Construction of Homotopies. Our algorithm for intersecting algebraic varieties is based on constructing homotopies to solve a system of polynomial equations restricted to an algebraic set. This is a generalization of existing homotopies, which have until now always worked on complex Euclidean space, ℂ^{m}. Accordingly, in this section we give a very general construction for homotopies on varieties.

Let X ⊂ ℂ^{m} be an irreducible N-dimensional affine algebraic set and let Y be an irreducible r-dimensional smooth algebraic set with r ≥ 1. Let

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

be a system of N algebraic functions on X × Y. In practice, Y is a parameter space defining a family of systems of interest, and for any one member of the family, we wish to find its solution points in X.

More precisely, suppose we have some parameter value y^{∗} ∈ Y for which we want to find a finite set F^{∗} of solutions of the system f(x, y^{∗}) = 0, such that all the isolated solutions of f(x, y^{∗}) = 0 are contained in F^{∗}. A procedure to do this proceeds in a number of steps in the same manner as if Y is ℂ^{N}.

1. Choose a point y′ ∈ Y for which we can find the isolated solutions F′ of f(x, y′) = 0, and the number of isolated solutions is the maximum number D for any system f(x, y) = 0 as a system in the x variables. We assume here that y′ ≠ y^{∗}, since otherwise we are done.

2. Construct a smooth connected algebraic curve B ⊂ Y which contains y^{∗} and y′. (Typically Y is a Euclidean space and we choose B equal to the complex line joining the points y^{∗} and y′.)

3. Construct a differentiable mapping c : [0, 1] × Γ → B where Γ is an interval or the unit circle, c(0, Γ) = y^{∗}, c(1, Γ) = y′, and where there is a positive integer K such that given any point y′′ ∈ c([0, 1] × Γ) not equal to y^{∗} or y′, it follows that c^{−1}(y′′) has at most K inverse images.

4. Choose a random γ ∈ Γ and starting with the isolated solutions F′ of f(x, y′) = 0 use “homotopy continuation” of the system f(x, c(t, γ)) = 0 to continue from the solutions F′ at t = 1 to solutions F^{∗} at t = 0.

Let us show that if we can make the choices specified by this procedure, we will find a finite set F^{∗} of solutions of the system f(x, y^{∗}) = 0, such that all the isolated solutions of f(x, y^{∗}) = 0 are contained in F^{∗}. In Appendix A we show how to relax item (2) so that the procedure can be carried in all situations where Y is irreducible.

It may well happen that the solution sets of f(x, y) = 0 for some or all the y ∈ Y also contain positive dimensional solution components. Nevertheless, the number D in item (1) exists and is finite by general results, e.g., [10]. Now choose B as in item (2) above. Lemma A.1 guarantees that for all but a finite number of points ŷ ∈ B, f(x, ŷ) = 0 has D isolated solutions, and that the closure of the union of the isolated solutions of f(x, ŷ) = 0 as ŷ runs over B is an algebraic curve B which surjects generically D-to-one onto B. Since the set of points in B over which this mapping is not a covering is an algebraic set and hence finite, the procedure is seen to work.

Remark 4.1. Algebraic functions on an affine algebraic set X ⊂ ℂ^{m} are the restrictions of polynomials from ℂ^{m}. If by deg f_{i} we denote a degree of a polynomial on ℂ^{m} restricting to f_{i}, it follows that the number D above is at most deg X × Π_{i=1}^{N} deg f_{i}. From this it further follows that if we can find a y′ such that f(x, y′) = 0 has deg X × Π_{i=1}^{N} deg f_{i} nonsingular isolated solutions, we can use y′.

Remark 4.2. Lemma A.1, which justifies the above procedure, is strong enough to yield the algorithms we need to construct witness points. However, the lemma is too weak to relate the multiplicity of the points as they appear in these algorithms to the multiplicity of the components that they represent. See App. A for more details.

5. An Abstract Embedding Theorem. The object of this section is to present Theorem 5.1, a generalization of the main theorem of [12]. We are aiming for the same results as in that article except that ℂ^{N} is replaced by a pure N-dimensional affine algebraic set X. We call the generalization in this section “abstract,” because we do not specify an explicit description of X; a numerical version is the topic of the next section. Since the proof of Theorem 5.1 follows the same line of reasoning of [12], we only state and discuss the parts of that article that need changes. Before we can state the theorem, we need some notation.

5.1. Definitions. Let X ⊂ ℂ^{m} be a reduced pure N-dimensional affine algebraic set, i.e., an affine algebraic subset of ℂ^{m}, all of whose irreducible components are of multiplicity one and dimension N. We assume that we have a system of restrictions of polynomials on ℂ^{m} to X

f(x) := [ f_{1}(x) ; ⋮ ; f_{N}(x) ] (5.1)

We assume that f does not vanish identically on any irreducible component of X. We will occasionally abuse notation and use the same notation f_{i} to denote the polynomial on ℂ^{m} and its restriction to X. In line with this abuse, we let

x := [ x_{1} ; ⋮ ; x_{m} ] (5.2)

denote both the coordinates on ℂ^{m} and the restrictions of the coordinates to X. Let Y denote the matrix space ℂ^{N×(1+m+N)}, with submatrices denoted as

[ A_{0} A_{1} A_{2} ] , (5.3)

where A_{0} ∈ ℂ^{N×1}, A_{1} ∈ ℂ^{N×m}, and A_{2} ∈ ℂ^{N×N}. We have the stratification of Y

Y_{0} ⊂ Y_{1} ⊂ ⋯ ⊂ Y_{N}, (5.4)

where Y_{i} is the subspace of Y obtained by setting the last N − i rows of Y equal to 0, and we define π_{i} : Y → Y_{i} as the corresponding projection. Note in particular that Y_{N} is Y, whereas Y_{0} is a N × (1 + m + N) matrix of zeroes. Defining e_{i} as the N × N matrix of all zeros except a 1 in the i-th diagonal element and letting P_{i} = ∑_{j=1}^{i} e_{j}, we can explicitly write π_{i}(Y) = P_{i}Y. We define P_{0} to be the N × N matrix with all entries zero. This notation will be useful in defining a homotopy below.

We let

z := [ z_{1} ; ⋮ ; z_{N} ] (5.5)

denote coordinates on ℂ^{N}.

5.2. Embedding and Cascade. For Y = (A_{0}, A_{1}, A_{2}) ∈ Y, we define the system

E(f)(x, z, Y) := [ \begin{matrix}f(x) + A_{2}^{T} z \\ z − A_{0} − A_{1}x \end{matrix} ] , (5.6)

which admits the embeddings

E_{i}(f)(x, z, Y) = E(x, z, π_{i}(Y)) = E(x, z, P_{i}Y). (5.7)

We often refer to E_{i}(f)(x, z, Y) by E_{i} or E_{i}(f). We regard E_{i} as a family of systems of equations on X × ℂ^{N} parameterized by Y_{i}. Note that

1. the E_{i} are the restrictions of systems E_{N} to Y_{i}; and

2. the E_{i} can be identified with systems on X × ℂ^{i} with coordinates z_{1}, …, z_{i} on ℂ^{i}. (This is because z_{j} = 0 for j > i.) Thus, E_{0} is naturally identified with the system f.

For i from 1 to N and γ_{i} ∈ {γ ∈ ℂ | |γ| = 1}, we define a cascade of homotopies that connect the embedded systems:

H_{i}(x, z, t, Y, γ_{i}) := [ \begin{matrix}f(x) + A_{2}^{T} z \\ P_{i−1}(z − A_{0} − A_{1}x) \\ +e_{i}((1 − t)z + γ_{i}t(z − A_{0} − A_{1}x)) \\ +(I_{N} − P_{i})z \end{matrix} ] . (5.8)

The nonzero parts of the three terms in the lower block of this expression occupy separate rows, with only the i-th row depending on t. At t = 1, H_{i}(x, z, 1, Y, γ_{i}) is equivalent to E_{i}(x, z, Y) (they differ only in that the i-th row of the lower block has been scaled by γ_{i}), and at t = 0, H_{i}(x, z, 0, Y, γ_{i}) = E_{i−1}(x, z, Y). Homotopy H_{i} allows us to compute solutions to the embedded systems by continuation, as described in the next paragraph.

For i from 1 to N, F_{i} denotes the solutions to E_{i} = 0 with z ≠ 0. In the case of i = 0, we make the convention that F_{0} is the empty set. Of course, like E_{i}, F_{i} depends on Y ∈ Y. We do not emphasize the dependence since the thrust of the main result is that a generic choice of Y, which is done once and for all using a random number generator in implementations, has a number of nice properties:

1. the solutions F_{i} of E_{i} = 0 are nonsingular and isolated and equal to the set of solutions of E_{i} = 0 with z_{i} ≠ 0;

2. the solutions of E_{i} = 0 equal F_{i} for i > dim V (f); and

3. for all u ∈ F_{i} and but a finite number of γ_{i}, there is a unique continuous map s_{u}(t) : [0, 1] → X × ℂ^{i} such that

(a) s_{u}(1) = u;

(b) H_{i}(s_{u}(t), t, Y, γ_{i}) = 0; and

(c) the Jacobian of H_{i}(x, z, t, Y, γ_{i}) with respect to (x, z) is invertible at

(s_{u}(t), t) for t ∈ (0, 1].

4. The limits of the functions s_{u}(t) as t → 0, which exist by the last properties, consists of the set F_{i−1} plus a finite set Ŵ_{i−1}.

The collection of sets Ŵ_{i} for i = 1, …, N contains the witness points for the irreducible decomposition of f^{−1}(0). This is stated formally in the following theorem, a generalization of the main theorem of [12].

Theorem 5.1. Let f be the restriction of a system of N polynomials on ℂ^{m} to a pure N-dimensional affine algebraic set X ⊂ ℂ^{m}. Assume that f is not identically zero on any irreducible component of X and that (A_{0}, A_{1}, A_{2}) is chosen generically. If j is the largest integer with Ŵ_{j} nonempty, then the dimension of f^{−1}(0) is j. Moreover given any irreducible component W of f^{−1}(0) of dimension i ≤ j, then, the finite set, Ŵ_{i} contains deg(W_{red}) generic points of W_{red}, each counted ν_{W} times, where ν_{W} is a positive integer, and W_{red} is the reduction of W. The remaining points J_{i} ⊂ Ŵ_{i} lie on components of f^{−1}(0) of dimension > i.

Theorem 5.1 is a consequence of Lemma A.1 and Lemmas B.2 and B.4 in the appendices.

6. Numerical Embedding. In this section we show how to numerically implement the algorithm of §5. We assume that we have f and X ⊂ ℂ^{m} as in Theorem 5. We assume that we have a system of polynomials on ℂ^{m}

g(x) := [ g_{1}(x) ; ⋮ ; g_{n}(x) ] (6.1)

such that X is a union of dimension N irreducible components of V (g). Once and for all choose a randomized system of m − N polynomials G(x) := ℜ(g(x); m − N). By Lemma 3.1 we know that X is a union of dimension N irreducible components of V (G).

We further assume that we begin with a witness set for X; that is, we know its dimension N, and have found the deg X smooth isolated witness points W = L_{m−N} ∩ X for a general linear subspace L_{m−N} of dimension m − N. (This will be on hand after computing the numerical irreducible decomposition of g(x) = 0.)

To convert the “abstract” systems of the previous section to systems we can compute with, we append G. Thus regarding the f_{i} as polynomials on ℂ^{m}, we replace E_{i}(f) by

[ \begin{matrix}G(x) \\ E_{i}(f)(x, z) \end{matrix} ] , (6.2)

which by abuse we still call E_{i}(f). We let Ẽ_{i}(f) denote the original E_{i}(f) without the G(x): we only use this in Equations 6.6 and 6.7.

To start the algorithm we need to solve E_{N}(f) = 0. Assume the total degree of f_{i} as a polynomial on ℂ^{m} is d_{i} for each i.

Choosing d_{1} + ⋯ + d_{N} general linear forms

L_{1,1}(x), …, L_{1,d_{1}}(x), …, L_{N,1}(x), …, L_{N,d_{N}} (x) (6.3)

on ℂ^{m}, we want them to have the good property that for any choice of integers i_{j} in 1, …, d_{j} for each j in 1, …, N, the solution set S_{i_{1},...,i_{N}} of the system of restrictions to X of the linear equations

L_{1,i_{1}}(x) = 0 ; ⋮ ; L_{N,i_{N}} (x) = 0 (6.4)

consists of deg X nonsingular isolated solutions, and moreover S_{i_{1},...,i_{N}} ∩ S_{k_{1},...,k_{N}} = ∅ unless (i_{1}, …, i_{N}) = (k_{1}, …, k_{N}). Let π : ℂ^{m} → ℂ^{N} denote a general linear projection. As discussed in [14], π_{X} is proper and (deg X)-to-one. Let B be the proper algebraic subset such that π_{X} is an unramified cover when restricted to X\π^{−1}(B). By composing with π we have reduced to the straightforward observation that choosing d_{1} + ⋯ + d_{N} general linear functions L_{i} on ℂ^{N} for i from 1 to d_{1} + ⋯ + d_{N}, it follows that

1. the unique zero of the linear functions L_{i_{1}}, …, L_{i_{N}} for distinct i_{1}, …, i_{N} between 1 and d_{1} + ⋯ + d_{N} vanishes at a general point of ℂ^{N} \ B; and 2. given any N + 1 of the d_{1} + ⋯ + d_{N} linear functions, there are no solutions

on ℂ^{N}. The system

L̂(x) = [ L_{1,1}(x) ⋯ L_{1,d_{1}}(x) ; ⋮ ; L_{N,1}(x) ⋯ L_{N,d_{N}} (x) ] = 0 (6.5)

has d_{1} ⋯ d_{N} · deg X nonsingular isolated solutions w_{α} contained in X_{reg}, the Zariski open set of smooth points of X. By homotopy continuation tracking from L_{m−N} to each of the d_{1} ⋯ d_{N} linear systems that occur in the system L̂(x), we can compute all the solutions w_{α} of L̂(x) = 0.

Fix the homotopy

H(x, z, t) := [ \begin{matrix}G(x) \\ (1 − t)Ẽ_{N}(f)(x, z) + tγ [ \begin{matrix}L̂(x) \\ z \end{matrix} ] \end{matrix} ] = 0. (6.6)

where γ is any of all but a finite number of norm one complex numbers. The solutions of E_{N}(f) = 0 are the nonsingular limits as t → 0 on X of paths starting at t = 1 with the w_{α} and z_{i} = 0 for all i.

Remark 6.1. In actual practice we often have some estimate, say N − 1 of the largest dimension of any component of the solution set of f on X. This will happen, for example, in §7. In such a situation we need only start with E_{N}. In this case we can replace the homotopy 6.6 with

H(x, z, t) := [ \begin{matrix}G(x) \\ (1 − t)Ẽ_{N} (f)(x, z) + tγ [ \begin{matrix}L̂(x) \\ z \end{matrix} ] \end{matrix} ] = 0. (6.7)

Note that the smooth nonsingular solutions of E_{i}(f) on X are generic. Thus they miss E := \overline{(G^{−1}(0) \ X)} ∩ X except for a proper algebraic set of parameter values. Thus for a Zariski open dense set of the homotopy parameters the homotopies with G compute the abstract homotopies. Though E might well contain the limits of a homotopy, the value of the limit is not influenced by G.

It is important to realize that serious numerical difficulties can arise, even when we are dealing with a nice smooth reduced component C of the system f on X. These occur if C is contained in a component of V (g) other than those in X. If this happens path tracking to decompose the witness point superset containing generic points of C will be ‘singular,’ and require the path tracker used in [15].

7. Diagonal Homotopies. Assume that A is an irreducible component of the solution set of polynomial system f_{A}(u) = 0 in u ∈ ℂ^{k} of dimension a > 0, and that B is an irreducible component of the solution set of polynomial system of a polynomial system f_{B}(v) = 0 in v ∈ ℂ^{k} of dimension b > 0. An important special case of this is when f_{A} and f_{B} are the same system, and A and B are distinct irreducible components. After renaming if necessary, we assume a ≥ b. Moreover we assume that B is not contained in A, since we would check this at the start of the algorithm and terminate if B was contained in A. Thus all components of A∩B are of dimension at most b − 1.

We wish to compute the irreducible decomposition of A∩B. Note that the product X := A × B ⊂ ℂ^{2k} is irreducible of dimension a + b. The theory of the preceding sections applies with m = 2k and N = a + b. The intersection of A and B can be identified, e.g., [4, Ex. 13.15], with X ∩ Δ where Δ is the diagonal of ℂ^{2k} defined by the system on X

δ(u, v) := [ u_{1} − v_{1} ; ⋮ ; u_{k} − v_{k} ] = 0. (7.1)

Remark 7.1. Notice that δ(u, v) plays the role of f in (5.1) in §5.

If a + b ≥ k set D(u, v) equal to δ(u, v) with a + b − k identically zero equations adjoined. If k < a + b, fix a randomization D(u, v) := ℜ(δ(u, v); a + b) once and for all. Note that the smallest dimensional nonempty component of A ∩ B is of dimension at least max{0, a+b−k}. Thus by Lemma (3.1), we can find the irreducible decomposition of A ∩ B by finding the irreducible decomposition of D(u, v) = 0 on X.

Fix randomizations F_{A}(u) := ℜ(f_{A}(u); k − a) and F_{B}(v) := ℜ(f_{B}(v); k − b) once and for all. We assume that we have already processed f_{A} and f_{B} through our numerical irreducible decomposition. So our data for A consists of a generic system L_{A}(u) = 0 of a = dim A linear equations and the deg A solutions {α_{1}, …, α_{deg A}} ∈ ℂ^{k} of the system

[ \begin{matrix}F_{A}(u) \\ L_{A}(u) \end{matrix} ] = 0, (7.2)

and the data for B consists of a generic system L_{B}(v) = 0 of b = dim B linear equations and the deg B solutions {β_{1}, …, β_{deg B}} ∈ ℂ^{m} of the system

[ \begin{matrix}F_{B}(v) \\ L_{B}(v) \end{matrix} ] = 0. (7.3)

Remark 7.2. We are not assuming that A and B occur with multiplicity one. If the multiplicity is greater than one, we must use a singular path tracker [15].

Note that A × B is an irreducible component of the solution set of the system

F(u, v) := [ \begin{matrix}F_{A}(u) \\ F_{B}(v) \end{matrix} ] = 0. (7.4)

In the following paragraphs, we write z_{h:k} to mean the column vector of variables z_{h}, …, z_{k}.

Since we know that all components of A ∩ B are of dimension at most b − 1, the first system of the cascade of homotopies is

E_{b}(u, v, z_{1:b}) = [ \begin{matrix}F(u, v) \\ ℜ(D(u, v), z_{1}, …, z_{b}; a + b) \\ z_{1:b} − ℜ(1, u, v; b) \end{matrix} ] = 0. (7.5)

This system consists of k − a + k − b + a + b + b = 2k + b equations in 2k + b variables.

To start the cascade, we must find the solutions of Eq. (7.5). Recall a ≥ b. Specializing the system from the end of §6, we have the homotopy

[ \begin{matrix}F(u, v) \\ (1 − t) [ \begin{matrix}ℜ(D(u, v), z_{1}, …, z_{b}; a + b) \\ z_{1:b} − ℜ(1, u, v; b) \end{matrix} ] + tγ [ \begin{matrix}L_{A}(u) \\ L_{B}(v) \\ z_{1:b} \end{matrix} ] \end{matrix} ] = 0. (7.6)

At t = 1, solution paths start at the deg A × deg B nonsingular solutions

{(α_{1}, β_{1}), …, (α_{deg A}, β_{deg B})} ⊂ ℂ^{2k} (7.7)

obtained by combining the witness points for A and B. At t = 0, the solution paths terminate at the desired start solutions for Eq. (7.5).

Since A ∩ B ≠ ∅ implies that

dim A ∩ B ≥ a + b − k, (7.8)

we see that when a + b ≥ k, we do not have to continue the cascade beyond level a + b − k. We can codify this into the numerics by noting that the system E_{b} is, with probability one, the same as the system

Ê_{b}(u, v, z_{b−ā+1}, …, z_{b}) = [ \begin{matrix}F(u, v) \\ ℜ(δ(u, v), z_{b−ā+1}, …, z_{b}; k) \\ ℜ(1, u, v; b − ā) \\ z_{(b−ā+1):b} − ℜ(1, u, v; ā) \end{matrix} ] = 0, (7.9)

where ā = k − a. This system has (k − a) + (k − b) + k + (k − a) + (a + b − k) = 3k − a equations in 3k − a variables. Notice that a + b ≥ k implies 3k − a ≤ 2k + b. To appreciate this, consider the case when a and b are both k − 1 and f_{A} and f_{B} are each a single equation. In this case the first system of the cascade is

Ê_{1}(u, v, z_{1}) = [ \begin{matrix}f_{A}(u) \\ f_{B}(v) \\ u − v + R_{k×1}z_{1} \\ ℜ(1, u, v; k − 2) \\ z_{1} − ℜ(1, u, v; 1) \end{matrix} ] = 0, (7.10)

where R_{k×1} is a generic complex k-vector.

In the important case when a + b ≥ k, we want to compute the start solutions for Eq. (7.9). Then, letting ā = k − a, the homotopy (7.6) reduces with probability one to

[ \begin{matrix}F(u, v) \\ (1 − t) [ \begin{matrix}ℜ(δ(u, v), z_{b−ā+1}, …, z_{b}; k) \\ ℜ(1, u, v; b − ā) \\ z_{(b−ā+1):b} − ℜ(1, u, v; ā) \end{matrix} ] + tγ [ \begin{matrix}L_{A}(u) \\ L_{B}(v) \\ z_{(b−ā+1):b} \end{matrix} ] \end{matrix} ] = 0 (7.11)

8. Computational Experiments. The diagonal homotopies are implemented in the software package PHCpack [18], recently upgraded to deal with positive dimensional solution components.

To compute witness points on all positive dimensional components of the intersection, we distinguish three stages:

1. given witness points on the two components, construct the top dimensional system in the cascade and the start system to start the cascade;

2. use polynomial continuation to compute the solutions at the start of the cascade; and

3. follow all paths defined by the cascade, in b stages, until all slack variables in z_{1:b} are eliminated or until no more paths are left to trace. When a + b ≥ k, we need work only with z_{(b−ā+1):b}.

The complexity of this procedure thus depends on

1. the number of variables (and equations) in the top dimensional system in the cascade;

2. the number of paths it takes to compute the solutions at the start of the cascade; and

3. the number of paths defined by the cascade. Although we will mention timings of runs done on a 2.4 Ghz Linux machine, the numbers describing the complexity are less transient.

8.1. An illustrative example. Consider the following example:

f(x, y, z, w) = [ \begin{matrix}xz \\ xw \\ yz \\ yw \end{matrix} ] = 0. (8.1)

There are two solution components of dimension two, characterized by the equations {x = 0, y = 0} and {z = 0, w = 0}. Pretending we do not know the two components intersect in the origin, we will set up a cascade of homotopies to compute the intersection of the two components.

Since we start out with four variables (k = 4), and work with two dimensional components (a = b = 2), the total number of variables at the start of the cascade is 2k + b = 10. The components are characterized by one witness point each, so there is only one path to trace. Tracing one path to start the cascade only takes 80 milliseconds CPU time, and gives a point with z_{2} ≠ 0, z_{1} ≠ 0. In the first stage of the cascade, we take z_{2} to zero, but z_{1} remains nonzero, showing that there is not a 1-dimensional component. The second stage of the cascade takes z_{1} to zero and yields the origin as the point of intersection of the two components, as expected. The two stages of the cascade together take just 30 milliseconds.

8.2. Intersection of a cylinder with a sphere. In Figure 8.1 we see a sphere intersected by a cylinder. The curve C defined by this intersection is

C := { (x, y, z) | x^{2} + y^{2} − 1 = 0 ∩ (x + 0.5)^{2} + y^{2} + z^{2} − 1 = 0 }. (8.2)

The total user CPU time of all path tracking is about a tenth of a second. First we track two paths to find a witness set for the cylinder, which takes 20 milliseconds.

Fig. 8.1. Intersection of a sphere with a cylinder. At the right we see the curve of degree four defined by the intersection.

Then it takes also 20 milliseconds to compute a witness set for the sphere. We have a = b = 2 and k = 3, thus a + b > k and the diagonal homotopy requires 7 variables, as 7 = 3k − a. Tracking the 2 × 2 paths defined by the diagonal homotopy takes 70 milliseconds CPU user time. At the end of the paths we find four points in the witness set for the curve C.

We may now move the slicing plane of the witness set to find the intersection of C with any desired plane. For example, to find the points on C of the form (x, x, z), we move the slice in a continuous fashion to x − y = 0. Tracking the four solutions in the witness set to this special plane takes only 10 milliseconds CPU time and gives two real and two complex-conjugate solutions.

8.3. Adding an extra leg to a moving platform. In this section we give an application of the important case where one of the components is a hypersurface. We consider a special case of a Stewart-Gough platform proposed by Griffis and Duffy [6]. When further specialized to have equilateral upper and lower triangles connected by six legs in cyclic fashion from a vertex of one triangle to a midpoint of an edge of the other triangle, and vice versa, the platform permits motion. This property was first identified and analyzed by Husty and Karger [7] and subsequently re-examined by the authors of this paper in [15].

When the legs of the mechanism described above have general lengths, a formulation of the kinematic equations using Study coordinates has one curve of degree 28 and 12 lines [15]. The lines are mechanically irrelevant, so we ignore them. Suppose we form a tetrahedron by adding a fourth point in general position to the base triangle and similarly for the upper triangle and then add a seventh leg of known length connecting these two points. The condition for assembling the mechanism is equivalent to intersecting the motion curve of degree 28 for the first six legs with a quadratic hypersurface that equates the length of the seventh leg to the distance between its points of connection. This hypersurface is of the same form as the main equations in the system defining the curve. With the addition of the seventh leg, the platform will no longer move, but will have instead a finite number of fixed postures.

The number of variables and equations in the original system is eight (k = 8). We intersect a one dimensional component with a hypersurface, for k = 8, this hypersurface is of dimension seven. Since a ≥ b, we have a = 7 and b = 1. So the cascade starts with 17 variables, as 2k + b = 2 × 8 + 1 = 17. The hypersurface is represented by 2 witness points and the curve we intersect has 28 witness points. To start the cascade, we trace 2 × 28 = 56 paths in dimension 17, using 20.3 seconds user CPU time. The cascade just has to remove one hyperplane to arrive at the 40 intersection points (16 of the 56 paths diverge), which requires 14.4 seconds user CPU time. Interestingly, a general Stewart-Gough platform also has 40 solution points.

Finally, we point out that the CPU time spent on the diagonal homotopy is considerably less than solving the system directly. For the direct approach the input is a system in 9 equations and 8 variables. Before giving it to the blackbox solver of PHCpack, we add to every equation one monomial, which is a new slack variable multiplied with a random constant. The mixed volume of this new 9-dimensional system is 164. The computation of the mixed volume and tracking of all 164 paths takes 108.5 seconds (1.8 minutes) CPU time. At the end we find the same 40 intersection points, the other 124 paths diverged to infinity. Notice that in the diagonal homotopy, only 16 paths diverged.

9. Conclusions. In this paper, we extend the cascade of [12] to compute witness points on all components of the intersection of two irreducible varieties. This is done by computing the irreducible decomposition of the diagonal of the product of the two irreducible varieties, and so we call the new procedure a “diagonal homotopy.” The procedure is justified as a special case of a method, also described herein, for the irreducible decomposition of the solution set of any polynomial system restricted to an irreducible algebraic set.

The diagonal homotopy given here always has at least twice the number of variables as the ambient space of the varieties being intersected. In a sequel to this paper, we will describe a modification to the diagonal homotopy that avoids the explicit doubling of the system, which leads to more efficient computation.

REFERENCES

[1] E.L. Allgower and K. Georg. Numerical Continuation Methods, an Introduction, volume 13 of Springer Ser. in Comput. Math. Springer–Verlag, 1990. To appear in the SIAM Classics in Applied Mathematics Series.

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

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

[4] D. Eisenbud. Commutative Algebra with a View Toward Algebraic Geometry, volume 150 of Graduate Texts in Mathematics. Springer–Verlag, New York, 1995.

[5] W. Fulton. Intersection theory, Ergebnisse der Mathematik und ihrer Grenzgebiete 3, second edition. Springer–Verlag, New York, 1998.

[6] M. Griffis and J. Duffy. Method and apparatus for controlling geometrically simple parallel mechanisms with distinctive connections. US Patent 5,179,525, 1993.

[7] M.L. Husty and A. Karger. Self-motions of Griffis-Duffy type parallel manipulators Proc. 2000 IEEE Int. Conf. Robotics and Automation, CDROM, San Francisco, CA, April 24–28, 2000.

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

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

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

[11] D. Mumford. Algebraic Geometry I: Complex Projective Varieties, volume 221 of Grundlehren der mathematischen Wissenschaften 221. Springer–Verlag, Berlin, 1976.

[12] A.J. Sommese and J. Verschelde. Numerical homotopies to compute generic points on positive

dimensional algebraic sets. Journal of Complexity 16(3):572–602, 2000.

[13] A.J. Sommese, J. Verschelde and C.W. Wampler. Numerical decomposition of the solution sets of polynomial systems into irreducible components. SIAM J. Numer. Anal. 38(6):2022– 2046, 2001.

[14] A.J. Sommese, J. Verschelde, and C.W. Wampler. Using monodromy to decompose solution sets of polynomial systems into irreducible components. In Application of Algebraic Geometry to Coding Theory, Physics and Computation, edited by C. Ciliberto, F. Hirzebruch, R. Miranda, and M. Teicher. Proceedings of a NATO Conference, February 25 - March 1, 2001, Eilat, Israel. Pages 297–315, Kluwer Academic Publishers.

[15] A.J. Sommese, J. Verschelde and C.W. Wampler. Symmetric functions applied to decomposing solution sets of polynomial systems. SIAM J. Numer. Anal. 40(6):2026–2046, 2002.

[16] A.J. Sommese, J. Verschelde, and C.W. Wampler. A method for tracking singular paths with application to the numerical irreducible decomposition. In Algebraic Geometry, a Volume in Memory of Paolo Francia, edited by M.C. Beltrametti, F. Catanese, C. Ciliberto, A. Lanteri, C. Pedrini. W. de Gruyter, pages 329-345, W. de Gruyter, 2002.

[17] 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, 749–763, 1996. Proceedings of the AMS-SIAM Summer Seminar in Applied Mathematics, Park City, Utah, July 17-August 11, 1995, Park City, Utah.

[18] J. Verschelde. Algorithm 795: PHCpack: A general-purpose solver for polynomial systems by homotopy continuation. ACM Transactions on Mathematical Software 25(2): 251–276,

          1999. Software available at http://www.math.uic.edu/~jan.

Appendix A. Homotopy on an Algebraic Set.

In §4, we give a procedure for constructing homotopies to solve a system of parameterized polynomials restricted to an algebraic set. In that procedure, Y is the parameter space, B is a smooth curve in Y, and we compute solution paths along a one-real-dimensional curve in B. Both Y and B are irreducible algebraic sets.

While choosing a smooth B is a difficulty when Y is irreducible and singular, it is always easy to find an irreducible curve B that contains y′ and y^{∗} with y′ ∉ B_{Sing}. If y′ is not in the singular set Y_{Sing} of Y, then B is not contained in Y_{Sing}. This is more than enough for the procedure to find a finite set F^{∗} of solutions of the system f(x, y^{∗}) = 0, such that all the isolated solutions of f(x, y^{∗}) = 0 are contained in F^{∗}.

In fact, the procedure given in §4, works with item (2) relaxed to finding an irreducible curve B containing y′ and y^{∗} such that

1. B is not contained in the singular set Y_{Sing} of Y; and

2. y′ ∉ B_{Sing}.

That the procedure finds a finite set F^{∗} of solutions of the system f(x, y^{∗}) = 0, such that all the isolated solutions of f(x, y^{∗}) = 0 are contained in F^{∗} may be shown by reducing to the nonsingular case:

a) Desingularize Y, i.e., let π : \overline{Y} → Y denote a surjective birational morphism which gives an isomorphism from \overline{Y} \ π^{−1}(Y_{Sing}) to Y \ Y_{Sing}.

b) Note that since B ⊄ Y_{Sing}, there is an algebraic curve B′ ⊂ \overline{Y} that maps generically one-to-one and onto B. By using embedded resolution of \overline{Y}, it can be further assumed that B′ is smooth.

c) By composition with π we get algebraic functions f′_{i} on X × \overline{Y}.

d) Note that given y^{∗} ∈ B ⊂ Y, there is a point y^{′∗} ⊂ B′ that maps onto y^{∗}.

e) Note that the result shown holds for X, \overline{Y} , f′, B′, and that (with the obvious identifications) the system f(x, y^{∗}) = 0 on X ×{y^{∗}} is identical to the system f′(x, y^{′∗}) on X × {y^{′∗}}.

The lemma which justifies the homotopy is as follows.

Lemma A.1. Let X ⊂ ℂ^{m} be an irreducible N-dimensional affine algebraic set and let Y be an irreducible smooth algebraic set. Let

f(x, y) = [ f_{1}(x, y) ; ⋮ ; f_{N}(x, y) ] = 0 (A.1)

be a system of N algebraic functions on X × Y. Let x^{∗} be an isolated solution of 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. Then there exists a neighborhood V of y^{∗} ∈ Y such that for any y ∈ V there exists at least one isolated solution of x ∈ O of f(x, y) = 0.

Proof. This result is a special case of a basic general result from complex algebraic geometry, e.g., [11, (3.10)]. Any irreducible component of f(x, y) = 0 is of dimension ≥ dim Y. Choose such a component C through (x^{∗}, y^{∗}). Consider \overline{C} ⊂ \overline{X} × Y where we close up X within ℙ^{m}. Since the induced projection π : \overline{C} → Y is proper, and there is a Euclidean neighborhood O of x^{∗} as in the lemma with π_{O}^{−1}(y^{∗}) = x^{∗}, we conclude that there is a neighborhood V ⊂ Y of y^{∗} such that π_{C∩(O×V)} : C ∩ (O × V) → V is proper. By the proper mapping theorem the image is an algebraic subset, and by the upper semicontinuity of fiber dimension it must be surjective. This proves the lemma.

We conclude this appendix with a few remarks on multiplicity. Lemma (A.1) is strong enough yield the algorithms we need to construct witness points, but unfortunately too weak for us to relate the multiplicity of x^{∗} as a solution of f(x, y^{∗}) = 0 to the multiplicity of the projection map from C to Y at (x^{∗}, y^{∗}). If X was a local complete intersection, then it would follow that C was Cohen-Macaulay in a neighborhood of (x^{∗}, y^{∗}), and we could use the stronger result [12, Lemma 6], and conclude the two multiplicities are the same, and thus have the same multiplicity statements as in [12, Theorem 3].

In the situation when we apply Lemma A.1 we know a bit more information, i.e., that for a general point y′ near y^{∗}, the solutions of f(x, y′) = 0 near (x^{∗}, y′) are nonsingular. It is worth noting in this case, e.g., using [11, Appendix to Chapter 6] that when we choose a sufficiently generic smooth curve in Y through y^{∗}, e.g., a generic line through y^{∗} when Y is Euclidean space, the number of paths coming into (x^{∗}, y^{∗}) is the multiplicity of the local ring of X at x^{∗} with respect to the ideal generated by the functions f_{i}(x, y^{∗}). Unfortunately, this multiplicity is in general only bounded by the multiplicity of (x^{∗}, y^{∗}) as a solution of f(x, y^{∗}) = 0.

Appendix B. Proof of the Main Theorem. For algebraic sets, there is a very strong version of Sard’s Theorem, e.g., [11, Theorem 3.7]. This result has a large number of consequences, going under the name Bertini’s Theorem, asserting that the zero set of a suitably general function inherits properties of the set the function is defined on. For the convenience of the reader we collect in one place a Bertini Theorem of sufficient generality to cover the needs of this article. Given a complex vector space V we let V^{r} denote the Cartesian product of V with itself r times, i.e., the space of are r-tuples of elements of V. V^{r} has a natural vector space structure given by addition of r-tuples and multiplication of an r-tuple by a complex number being defined as the r-tuple obtained by componentwise multiplication of elements of the r-tuple by the complex number. V^{r} with this vector space structure is denoted by V^{⊕r}. We use notation close to that of Fulton [5, Lemma B.9.1] in the Theorem B.1.

Theorem B.1 (Bertini’s Theorem). Let X denote an irreducible algebraic subset of ℂ^{k}. Let Z_{1}, …, Z_{q} denote a finite number of irreducible algebraic subsets of X (with one of the X_{i} possibly equal to X). Let 𝕍 denote the restriction to X of a finite dimensional vector space of polynomial functions on ℂ^{k}. Assume that for each point of X at least one element of 𝕍 does not evaluate to zero. Then for any integer r > 0, there is a Zariski open dense set U ⊂ 𝕍^{r} such that for f := (f_{1}, …, f_{r}) ∈ U it follows for each Z_{i} that

1. if V (f) ∩ Z_{i} is nonempty, then dim V (f) ∩ Z_{i} = dim Z_{i} − r; and

2. letting Sing(Z_{i}) denote the singular set of Z_{i}, V (f)∩(Z_{i}\Sing(Z_{i})) is smooth. Proof. We first apply [5, Lemma B.9.1]. For the vector bundle E in [5, Lemma

B.9.1] take X × ℂ^{r}; take p = 1 with C_{1} = X × {0}; for Γ take 𝕍^{⊕r}, i.e., take 𝕍^{r}. The conclusion from [5, Lemma B.9.1] is the existence of a Zariski open dense set Γ^{o} of Γ such that for f := (f_{1}, …, f_{s}) ∈ Γ^{o}, it follows that if V (f) ∩ Z_{i} is nonempty then

dim(V (f) ∩ Z_{i}) ≤ dim Z_{i} − r.

The opposite inequality is a property of zero sets of functions, e.g., [11, Corollary 3.14].

Since the intersection of a finite number of Zariski open and dense sets is Zariski open and dense, it suffices to show that there is a Zariski open dense set U_{i} ⊂ 𝕍^{⊕r} such that for f ∈ U_{i}, V (f) ∩ (Z_{i} \ Sing(Z_{i})) is smooth. For this we use [3, Theorem 1.7.1.1]. Restricting to (Z_{i} \ Sing(Z_{i})), we conclude from [3, Theorem 1.7.1.1] that there is a Zariski open dense set O_{1} ⊂ 𝕍 such that for f_{1} ∈ O_{1} we have that V (f_{1}) ∩ (Z_{i} \ Sing(Z_{i})) is smooth and if nonempty of dimension dim Z_{i} − 1. Applying [3, Theorem 1.7.1.1] to the restriction of 𝕍 to V (f_{1}) ∩ (Z_{i} \ Sing(Z_{i})), we conclude that there is a Zariski open dense set O_{2} ⊂ 𝕍 such that for f_{2} ∈ O_{2} we have that V (f_{1}, f_{2})∩ (Z_{i} \ Sing(Z_{i})) is smooth and if nonempty of dimension dim Z_{i} − 2. Proceeding this way for j going to r, U_{i} := O_{1} × ⋯ × O_{r} ⊂ 𝕍^{⊕r} is the desired Zariski open dense set.

Lemma B.2. Let f and X be as in Theorem 5.1. Assume further that Z is an algebraic subset of X of dimension < N. Assume that f does not vanish on any component of X or of Z. There is a Zariski open and dense set U ⊂ Y = ℂ^{N×(1+m+N)} such that

1. the solutions F_{i} of the system E_{i}(f)(x, z, Y) for Y ∈ U with z ≠ 0 are isolated nonsingular solutions and lie in the set (X \ Z) × ℂ^{i};

2. U ∩ Y_{i} is Zariski open and dense for each i < N; and

3. the solutions of E_{i}(f)(x, z, Y) for Y ∈ U with z ≠ 0 are the same as those with z_{i} ≠ 0.

Proof. Since the following result follows almost verbatim from the reasoning in the first half of the proof of [12, Lemma 2], we give only a brief sketch of the proof. As discussed in §5, we regard E_{i} as a system on X × ℂ^{i}.

Consider the vector space V_{1} of functions on X × ℂ^{i} generated by

f_{1}, …, f_{N}, z_{1}, …, z_{i}.

The common zeroes of the functions in V_{1} are the points

V (V_{1}) := {(x, 0) ∈ X × ℂ^{i}|f(x) = 0} .

From this we conclude, using Theorem B.1, that for a choice of a system S in a nonempty Zariski open set of the vector space V_{1}^{⊕N}, it follows that the common zeroes Z_{S} of S on X ×ℂ^{i} \V (V_{1}) is pure i-dimensional with singular set of dimension ≤ i−1. Similarly Z_{S} meets Z ×ℂ^{i} \V (V_{1}) in a set of dimension at most dim Z +i−N ≤ i−1.

Now let V_{2} be the vector space of functions on X × ℂ^{i} generated by

1, x_{1}, …, x_{m}, z_{1}, …, z_{i}.

Since 1 ∈ V_{2}, there are no common zeroes of the functions in V_{2}. Using Theorem B.1 again, we conclude that for a generic choice of a system S′ in a nonempty Zariski open set the vector space V_{2}^{⊕i}, it follows that the common zeroes of S′ on Z_{S} with z ≠ 0 is a finite set of isolated smooth points not contained in Z × ℂ^{i}. The above system

[ \begin{matrix}S \\ S′ \end{matrix} ] = 0, (B.1)

of N + i equations is of the form

\begin{matrix}B [ f_{1}(x) ; ⋮ ; f_{N}(x) ] + C [ z_{1} ; ⋮ ; z_{i} ] = 0 \\ D + E [ x_{1} ; ⋮ ; x_{m} ] + F [ z_{1} ; ⋮ ; z_{i} ] = 0, \end{matrix} (B.2)

where B is an N × N complex matrix, C is an N × i complex matrix, D is an i × 1 complex matrix, E is i × m complex matrix, and F is an i × i complex matrix. The above Bertini type results show that the set of

(B, C, D, E, F) ∈ ℂ^{N×(N+i)+i×(1+m+i)}

giving rise to systems of the form (B.2) with only isolated nonsingular solutions on X × ℂ^{i} \ X × {0} is dense in ℂ^{N×(N+i)+i×(1+m+i)} with respect to the usual Euclidean topology. The set of such (B, C, D, E, F) such that the maximal number of isolated solutions of the associated system (B.2) on X × ℂ^{i} \ X × {0} occurs is a dense constructible set, and thus by Chevalley’s Theorem, e.g., [11, Proposition 2.31], contains a dense Zariski open set O. Moreover we know that the systems of the form (B.2) with only isolated solutions on X×ℂ^{i}\X×{0} form a constructible set C of (B, C, D, E, F). By the density of the systems (B.1) in the usual Euclidean topology, we conclude that C is a dense Zariski constructible set and thus contains a dense Zariski open set O′. The systems arising with parameters from the set U′_{i} = O ∩ O′ have the properties required for the first assertion of the lemma. The matrices (B, C, D, E, F) giving rise to systems with the desired properties are invariant under the action

G_{1} × G_{2} × (B, C, D, E, F) → (G_{1}^{−1}B, G_{1}^{−1}C, G_{2}^{−1}D, G_{2}^{−1}E, G_{2}^{−1}F),

where G_{1} is an invertible N × N complex matrix and G_{2} is an invertible i × i complex matrix. Thus we can assume that U′_{i} is invariant under this action. Since the matrices (B, C, D, E, F) with B and F invertible form a Zariski open dense set invariant under the same action, we can assume that the Zariski open set U′_{i} is chosen so that all (B, C, D, E, F) in the set have B and F invertible. Since (I_{N}, B^{−1}C, F^{−1}D, F^{−1}E, I_{i}) is in U′_{i}, we see that the set

U_{i} := U′_{i} ∩ {(B, C, D, E, F)|B and F invertible }

is the desired set for conclusion 1) of the lemma.

We have natural projections π_{i} : Y → Y_{i} obtained by setting the last N − i rows of an element of Y to 0. The set U := ∩_{i=1}^{N}π^{−1}(U_{i}) is Zariski open and dense. Noting that since the maps π_{i} are surjective, the images are Zariski dense constructible sets, we have, upon redefining U_{i} := U ∩ Y_{i}, the first two assertions of the lemma.

For the last assertion we can assume without loss of generality that i ≥ 2. The desired assertion will follow if we show that the condition that the set of Y ∈ U for which there are solutions of E_{i}(f)(x, z, Y) with z_{i} = 0, but z ≠ 0, is not Zariski dense. Assume it was Zariski dense. Then, for a general Y ∈ U and a general (a_{i}, a_{1,1}, …, a_{i,N}) ∈ ℂ^{N+1}, the system

[ \begin{matrix}E_{i−1}(f)(x, z, Y) \\ a_{i} + a_{i,1}x_{1} + ⋯ + a_{i,m}x_{m} \end{matrix} ] = 0

has a solution with (z_{1}, …, z_{i−1}) ≠ 0. This is absurd, since we have already shown that for a general Y ∈ U, there are only a finite number of solutions of E_{i−1}(f)(x, z, Y) with (z_{1}, …, z_{i−1}) ≠ 0.

Remark B.3. The condition in Lemma B.2 that the F_{i} lie in (X − Z) × ℂ^{i} is important because we will typically not have defining polynomials for X, but only know that X is an irreducible component of V (g) for a system of polynomials g. Taking Z equal to the union of the intersections of X with other components of V (g) guarantees with probability-one that g will be a set of defining equations for X on a Zariski open set large enough so that all the homotopy continuations that are given in this article will be well defined.

We need some information about the isolated solutions of E_{i}(f)(x, z, Y) with z = 0. This is the generalization of the last assertion of [12, Lemma 2].

Lemma B.4. There is a Zariski open and dense set U ⊂ Y = ℂ^{N×(1+m+N)} such that the solutions of the system E_{i}(f)(x, z, Y) for Y ∈ U with z = 0 consist of

1. positive dimensional components all contained in components of V (f) of dimension greater than i; plus

2. for each dimension i irreducible component W of f^{−1}(0), isolated solutions consisting of deg(W_{red}) generic points of W_{red}, the reduction of W, each occurring the same number of times.

Proof. When z = 0, the system E_{i}(f)(x, z, Y) reduces to

[ \begin{matrix}f(x) \\ A_{0} + A_{1} · x \end{matrix} ] . (B.3)

The assertion is contained in the discussion in [17].

The remaining result from [12] that needs modification is the “Local Extension Lemma” [12, Lemma 6]. We use Lemma A.1 in its place.