NUMERICAL SCHUBERT CALCULUS

BIRKETT HUBER, FRANK SOTTILE, AND BERND STURMFELS

Abstract. We develop numerical homotopy algorithms for solving systems of polynomial equations arising from the classical Schubert calculus. These homotopies are optimal in that generically no paths diverge. For problems defined by hypersurface Schubert conditions we give two algorithms based on extrinsic deformations of the Grassmannian: one is derived from a Gröbner basis for the Plücker ideal of the Grassmannian and the other from a SAGBI basis for its projective coordinate ring. The more general case of special Schubert conditions is solved by delicate intrinsic deformations, called Pieri homotopies, which first arose in the study of enumerative geometry over the real numbers. Computational results are presented and applications to control theory are discussed.

Contents

1. Introduction 1 2. Linear equations in Plücker coordinates 3 2.1. An example 4 2.2. Gröbner homotopy 5 2.3. SAGBI homotopy 7 3. Special Schubert conditions 8 3.1. Basics on Schubert varieties 10 3.2. Pieri homotopy algorithm 12 3.3. Definition of the moving cycles Z_{R,k}(t) 13 3.4. An example 14 3.5. Proof of correctness 16 4. Homotopy continuation of overdetermined systems 17 5. Applications 18 6. Computational results 20 References 20

Date: 9 June 1997. 1991 Mathematics Subject Classification. 65H10, 14N10, 14M15, 14Q99, 05E10. Key words and phrases. Grassmannian, Homotopy Continuation, Gröbner Basis, Overdetermined

System.

Huber was supported in part by NSF grant DMS-9508742. Sottile was supported in part by NSF grant DMS-9022140. Sturmfels was supported in part by a David and Lucille Packard Fellowship and an NSF National Young Investigator Fellowship.

1. Introduction

Suppose we are given linear subspaces K_{1}, …, K_{n} of C^{m+p} with dim K_{i} = m+1−k_{i} and k_{1} + ⋯ + k_{n} = mp. Our problem is to find all p-dimensional linear subspaces of C^{m+p} which meet each K_{i} nontrivially. When the given linear subspaces are in general position, the condition k_{1} + ⋯ + k_{n} = mp guarantees that there is a finite number d = d(m, p, k_{1}, …, k_{n}) of such p-planes. The classical Schubert calculus [20] gives the following recipe for computing the number d. Let h_{1}, …, h_{m} be indeterminates with degree(h_{i}) = i. For each integer sequence λ_{1} ≥ ⋯ ≥ λ_{p+1} we define the following polynomial:

S_{λ} := det(h_{λ_{i}+j−i})_{1≤i,j≤p+1}. (1)

Here h_{0} := 1 and h_{i} := 0 if i < 0 or i > m. Let I be the ideal in Q[h_{1}, …, h_{m}] generated by those S_{λ} with m ≥ λ_{1} and λ_{p+1} ≥ 1. The quotient ring A_{m,p} := Q[h_{1}, …, h_{m}]/I is the cohomology ring of the Grassmannian of p-planes in C^{m+p}. It is Artinian with one-dimensional socle in degree mp. In the socle we have the relation

d · (h_{m})^{p} − h_{k_{1}}h_{k_{2}} ⋯ h_{k_{n}} ∈ I. (2)

Thus we can compute the number d by normal form reduction modulo any Gröbner basis for I. More efficient methods for computing in the ring A_{m,p} are implemented in the Maple package SF [37].

In the important special case k_{1} = ⋯ = k_{n} = 1 there is an explicit formula for d:

d = \frac{1! 2! 3! ⋯ (p−2)! (p−1)! · (mp)!}{m! (m+1)! (m+ 2)! ⋯ (m+p−1)!}. (3)

The integer on the right hand side is the degree of the Grassmannian in its Plücker embedding. This formula is due to [31]; see also [16, XIV.7.8] and Section 2.3 below.

The objective of this paper is to present semi-numerical algorithms for computing all d solution planes from the input data K_{1}, …, K_{n}. This amounts to solving certain systems of polynomial equations. Our algorithms are based on the paradigm of numerical homotopy methods [26, 1, 2].

Homotopy methods have been developed for the following classes of polynomial systems:

1. complete intersections in affine or projective spaces [11, 14],

2. complete intersections in products of projective spaces [27],

3. complete intersections in toric varieties [41, 17].

In these cases the number of paths to be traced is optimal and equal to the standard combinatorial bounds:

1. the Bézout number (= the product of the degrees of the equations)

2. the generalized Bézout number for multihomogeneous systems

3. the BKK bound [4, 21, 18] (= mixed volume of the Newton polytopes)

None of these known homotopy methods is applicable to our problem, as the following simple example shows: Take m = 3, p = 2, and k_{1} = ⋯ = k_{6} = 1, that is, we seek the 2-planes in C^{5} which meet six general 3-planes nontrivially. By formula (1.3) there are d = 5 solutions. Formulating this in Plücker coordinates gives 11 homogeneous equations in ten variables, the five quadrics in display (6) below and six linear equations (7). A formulation in local coordinates (9) has 6 quadratic equations in 6 unknowns, giving a Bézout bound of 64. These 6 equations all have the same Newton polytope, which has normalized volume 17, giving a BKK bound of 17.

In Section 2 we give two homotopy algorithms which solve our problem in the special case k_{1} = ⋯ = k_{n} = 1, when the number of solutions equals (1.3). The first algorithm is derived from a Gröbner basis for the Plücker ideal of a Grassmannian and the second from a SAGBI basis for its projective coordinate ring. (See [10] or [39, Ch. 11] for an introduction to SAGBI bases). Both the Gröbner homotopy and SAGBI homotopy are techniques for finding linear sections of Grassmannians in their Plücker embedding.

In Section 3 we address the general case of our problem, that is, we describe a numerical method for solving the polynomial equations defined by special Schubert conditions. This is accomplished by applying a sequence of delicate intrinsic deformations, called Pieri homotopies, which were introduced in [35]. Pieri homotopies first arose in the study of enumerative geometry over the real numbers [34, 33]. For the experts we remark that it is an open problem to find Littlewood-Richardson homotopies, which would be relevant for solving polynomial equations defined by general Schubert conditions.

A main challenge in designing homotopies for the Schubert calculus is that one is not dealing with complete intersections: there are generally more equations than variables. In Section 4 we discuss some of the numerical issues arising from this challenge, and how we propose to resolve them. In Section 5 we discuss applications of these algorithms to control theory and real enumerative geometry. Finally, in Section 6 we present computational results.

In closing the introduction let us emphasize that all homotopies described in this paper are optimal in the sense the number of paths to be traced equals the number d. This means that for generic input data K_{1}, …, K_{n} no paths diverge.

2. Linear equations in Plücker coordinates

The set of p-planes in C^{m+p}, Grass(p, m + p), is called the Grassmannian of pplanes in C^{m+p}. This complex manifold of dimension mp is naturally a subvariety of the complex projective space P^{\binom{m+p}{p}−1}. To see this, represent a p-plane in C^{m+p} as the column space of an (m + p) × p-matrix X = (x_{ij}). The Plücker coordinates of that p-plane are the maximal minors of X, indexed by the set \binom{[m+p]}{p} of sequences α : 1 ≤ α_{1} < α_{2} < ⋯ < α_{p} ≤ m + p:

[α] ⟶ det [ \begin{matrix}x_{α_{1} 1} & ⋯ & x_{α_{1} p} \\ ⋮ & ... & ⋮ \\ x_{α_{p} 1} & ⋯ & x_{α_{p} p} \end{matrix} ] . (4)

This section deals with the “k_{i} = 1” case of the problem stated in the Introduction. Given mp general m-planes K_{1}, …, K_{mp}, we wish to find all p-planes X which meet K_{1}, …, K_{mp} nontrivially. This geometric condition translates into linear equations in the Plücker coordinates: Represent X as an (m+p)×p-matrix as above, represent K_{i} as an (m + p) × m-matrix, and form the (m + p) × (m + p)-matrix [ X | K_{i}]. Then

X ∩ K_{i} ≠ {0} if and only if det [ X | K_{i}] = 0.

Laplace expansion with respect to the first p columns gives

det [ X | K_{i}] = ∑_{α∈\binom{[m+p]}{p}} C_{α}^{i} · [α], (5)

where C_{α}^{i} is the correctly signed maximal minor of K_{i} complementary to α. Hence our problem is to solve mp linear equations (5) on the Grassmannian. The number of solutions is the degree of the Grassmannian in its Plücker embedding, which is given in (1.3).

The Grassmannian is represented either implicitly, as the zero set of polynomials in the Plücker coordinates, or parametrically, as the image of the polynomial map (4). These two representations lead to two different numerical homotopies. The implicit representation gives the Gröbner homotopy in Section 2.2 while the parametric representation gives the SAGBI homotopy in Section 2.3. The first is conceptually simpler but the second is more efficient. In both methods the number of paths to be traced equals the optimal number in (3).

2.1. An example. We describe the two approaches for the case (m, p) = (3, 2). The Grassmannian of 2-planes in C^{5} has dimension 6 and is embedded into P^{9}. Its degree (3) is five. The Gröbner homotopy works directly in the ten Plücker coordinates:

[12], [13], [14], [15], [23], [24], [25], [34], [35], [45].

The ideal I_{3,2} of the Grassmannian in the Plücker embedding is generated by five quadrics:

\begin{matrix}\underline{[14][23]} − [13][24] + [12][34], \\ \underline{[15][23]} − [13][25] + [12][35], \\ \underline{[15][24]} − [14][25] + [12][45], \\ \underline{[15][34]} − [14][35] + [13][45], \\ \underline{[25][34]} − [24][35] + [23][45]. \end{matrix} (6)

This set is the reduced Gröbner basis for I_{3,2} with respect to any term order which selects the underlined terms as leading terms (see Proposition 2.1 below).

Our problem is to compute all 2-planes which meet six sufficiently general 3-planes K_{1}, …, K_{6} nontrivially. This amounts to solving (6) together with six linear equations

\begin{matrix}C_{12}^{i} · [12] + C_{13}^{i} · [13] + C_{14}^{i} · [14] + C_{15}^{i} · [15] + C_{23}^{i} · [23] \\ + C_{24}^{i} · [24] + C_{25}^{i} · [25] + C_{34}^{i} · [34] + C_{35}^{i} · [35] + C_{45}^{i} · [45] = 0, \end{matrix} (7)

for i = 1, …, 6. This is an overdetermined system of 11 equations in 10 homogeneous variables. To solve it we introduce a parameter t into (6) as follows:

\begin{matrix}[14][23] − t · [13][24] + t^{2} · [12][34] = 0, \\ [15][23] − t^{2} · [13][25] + t^{4} · [12][35] = 0, \\ [15][24] − t · [14][25] + t^{5} · [12][45] = 0, \\ [15][34] − t^{2} · [14][35] + t^{4} · [13][45] = 0, \\ [25][34] − t · [24][35] + t^{2} · [23][45] = 0. \end{matrix}  (6′)

We call (6′) the Gröbner homotopy because this is an instance of the flat deformation which exists for any Gröbner basis; see [12, Theorem 15.17]. The flatness of this deformation ensures that, for almost every complex number t, the combined system (6′)&(7) has five roots.

For t=0 the equations (6′) are square-free monomials. We decompose their ideal:

(8) \begin{matrix}⟨ [14][23], [15][23], [15][24], [15][34], [25][34] ⟩ \\ = ⟨ [23], [24], [34] ⟩ ∩ ⟨ [15], [23], [34] ⟩ ∩ ⟨ [15], [23], [25] ⟩ \\ ∩ ⟨ [14], [15], [34] ⟩ ∩ ⟨ [14], [15], [25] ⟩. \end{matrix}

The five distinct solutions for t = 0 are computed by setting each listed triple of variables to zero and then solving the six linear equations (7) in the remaining seven variables. Thereafter we trace the five solutions from t = 0 to t = 1 by numerical path continuation. At t = 1 we get the five solutions to our original problem.

We next describe the SAGBI homotopy. For this we choose the local coordinates

X = [ \begin{matrix}1 & 0 \\ x_{21} & x_{22} \\ x_{31} & x_{32} \\ x_{41} & x_{42} \\ 0 & 1 \end{matrix} ]

on the Grassmannian. Substituting X into (7) we get six polynomials in six unknowns:

\begin{matrix}C_{23}^{i} · (x_{21}x_{32} − x_{22}x_{31}) + C_{24}^{i} · (x_{21}x_{42} − x_{22}x_{41}) + C_{34}^{i} · (x_{31}x_{42} − x_{32}x_{41}) \\ + C_{25}^{i} · x_{21} + C_{35}^{i} · x_{31} + C_{45}^{i} · x_{41} + C_{12}^{i} · x_{22} + C_{13}^{i} · x_{32} + C_{14}^{i} · x_{42} + C_{15}^{i}. \end{matrix} (9)

To solve these six equations we introduce a parameter t as follows:

\begin{matrix}C_{23}^{i}(x_{21}x_{32} − t · x_{22}x_{31}) + C_{24}^{i}(x_{21}x_{42} − t^{2} · x_{22}x_{41}) + C_{34}^{i}(x_{31}x_{42} − t · x_{32}x_{41}) \\ +C_{25}^{i} · x_{21} + C_{35}^{i} · x_{31} + C_{45}^{i} · x_{41} + C_{12}^{i} · x_{22} + C_{13}^{i} · x_{32} + C_{14}^{i} · x_{42} + C_{15}^{i}. \end{matrix} (9′)

The system (9′) has five complex roots for almost all t ∈ C. For t = 0 we get a generic unmixed sparse system (in the sense of [17]) with support

A = { 1, x_{21}, x_{22} , x_{31} , x_{32} , x_{41} , x_{42} , x_{21}x_{32} , x_{21}x_{42} , x_{31}x_{42} }.

We identify A with a set of ten points in Z^{6}. Their convex hull conv(A) is a 6dimensional polytope with normalized volume five. We can therefore solve (9′) for t = 0 using the homotopy method in [17] or [41], provided the input data K_{1}, …, K_{6} are sufficiently generic. Tracing the five roots from t = 0 to t = 1 by numerical path continuation, we obtain the five solutions to our original problem.

2.2. Gröbner homotopy. We next describe a quadratic Gröbner basis for the defining ideal of Grass(p, m+p). Let S be the polynomial ring over C in the variables [α] where α ∈ \binom{[m+p]}{p}. We define a partial order on these variables as follows: [α] ≤ [β] if and only if α_{i} ≤ β_{i} for i = 1, …, p. This partially ordered set is called Young’s poset. Figure 1 shows Young’s poset for (m, p) = (3, 2).

Figure 1. Young’s poset for (m, p) = (3, 2).

Fix any linear ordering on the variables in S which refines the ordering in Young’s poset, and let ≺ denote the induced degree reverse lexicographic term order on S. Let I_{m,p} be the ideal of polynomials in S which vanish on the Grassmannian, that is, I_{m,p} is the ideal of algebraic relations among the maximal minors of a generic (m + p) × p-matrix X. The Gröbner homotopy is based on the following well-known result.

Proposition 2.1. The initial ideal in_{≺}(I_{m,p}) is generated by all quadratic monomials [ α][ β] where α_{i} < β_{i} and α_{j} > β_{j} for some i, j ∈ {1, …, p}.

In other words, in_{≺}(I_{m,p}) is generated by products of incomparable pairs in Young’s poset. Let C_{m,p} denote the set of all maximal chains in Young’s poset. For example,

C_{3,2} = { {[12], [13], [14], [15], [25], [35], [45]} , {[12], [13], [14], [24], [25], [35], [45]} ,

{[12], [13], [14], [24], [34], [35], [45]} , {[12], [13], [23], [24], [34], [35], [45]} ,

{[12], [13], [23], [24], [25], [35], [45]} }

A standard result in combinatorics [36] states that the cardinality of C_{m,p} equals the number (3). From Proposition 2.1 we read off the following prime decomposition which generalizes (8):

in_{≺}(I_{m,p}) = ⋂_{C∈C_{m,p}}⟨ [α] : [α] ∉ C ⟩. (10)

For a proof of Proposition 2.1 see [16, §XIV.9] or [6, Theorem (4.3)] or [38, §3.1]. In these references one finds an explicit minimal Gröbner basis for I_{m,p}, which is classically called the set of straightening syzygies. In the special case p = 2 the straightening syzygies coincide with the reduced Gröbner basis:

Proposition 2.2. If p = 2 then the reduced Gröbner basis of I_{m,p} consists of the three-term Plücker relations \underline{[il][kj]} − [ik][jl] + [ij][kl] where 1 ≤ i < j < k < l ≤ m + p = m + 2.

For p ≥ 3 the straightening syzygies and the reduced Gröbner basis do not coincide, and they are complicated to describe. For our purposes the following coarse description suffices. Let Std be the set of all quadratic monomials in S which do not lie in in_{≺}(I_{m,p}). The reduced Gröbner basis consists of elements of the form

[α][β] − ∑_{[γ][δ]∈Std} E_{γ,δ}^{α,β} · [γ][δ] (11)

where [α][β] runs over all generators of in_{≺}(I_{m,p}). The constants E_{γ,δ}^{α,β} are integers which can be computed by substituting (4) into (11) and solving linear equations.

The term order ≺ can be realized for the ideal I_{m,p} by the following choices of weights. We define the weight of the variable [α] = [α_{1}α_{2} ⋯ α_{p}] to be

v_{α} := −\frac{1}{2} ∑_{1≤i<j≤p}(α_{j} − α_{i} − 1)^{2}. (12)

If we replace each variable [α] in (11) by [α]·t^{v_{α}} and clear t-denominators afterwards, then we get the Gröbner homotopy:

[α][β] − ∑_{[γ][δ]∈Std} E_{γ,δ}^{α,β} · [γ][δ] · t^{v_{γ}+v_{δ}−v_{α}−v_{β}}  (11′)

It can be checked that all exponents v_{γ} + v_{δ} − v_{α} − v_{β} appearing here are positive integers. The special case (m, p) = (3, 2) is presented in (6′).

In the Gröbner homotopy algorithm, we first solve systems of mp linear equations, one for each chain C ∈ C_{m,p}. These systems consist of the mp equations (5), one for each K_{i}, and the \binom{m+p}{p} − mp − 1 equations

[α] for [α] ∉ C,

suggested by the prime decomposition (10). Once this is accomplished, we trace each of these d solutions from t = 0 to t = 1 in the Gröbner homotopy (6′).

Clearly, the weights v_{α} of (12) are not best possible for any specific value of m and p. Smaller weights can be found using Linear Programming, as explained e.g. in the proof of [39, Proposition 1.11]. Another method would be to adapt the “dynamic” approach in [40] to our situation. This is possible since the Gröbner basis in (11) is reverse lexicographic: first deform the lowest variable to zero, then deform the second lowest variable to zero, then the third lowest variable, and so on.

2.3. SAGBI homotopy. Let X = (x_{ij}) be an (m+p)×p-matrix of indeterminates. We identify the coordinate ring of the Grassmannian with the C-algebra generated by the p×p-minors of X. Call this algebra R and write [α](x_{ij}) for the minor indexed by α. Reinterpreting classical results in [16, §XIV.9], it was shown in [38, Theorem 3.2.9] that these generators form a SAGBI basis with respect to the lexicographic term order induced from x_{11} > x_{12} > ⋯ > x_{1p} > x_{21} > ⋯ > x_{m+p,p}. This means that the initial algebra C[ in_{>}(f) : f ∈ R] is generated by the main diagonal terms of the p × p-minors,

in_{>}([α](x_{ij})) = x_{α_{1},1} x_{α_{2},2} x_{α_{3},3} · · · · · · x_{α_{m},m}. (13)

The resulting flat deformation can be realized by replacing x_{ij} with x_{ij}t^{(i−1)(p−j)} for t → 0 in the matrix X. If we expand [α](x_{ij}t^{(i−1)(p−j)}) as a polynomial in t, then the lowest term equals t^{w_{a}} times the main diagonal monomial (13), where

w_{α} := ∑_{j=1}^{p} (α_{j} − 1)(p − j).

In what follows we multiply that polynomial by t^{−w_{α}}. For any t ∈ C consider the algebra

R_{t} := C [ t^{−w_{α}} · [α](x_{ij}t^{(i−1)(p−j)}) : α ∈ \binom{[m+p]}{p} ] .

Then R_{1} is the coordinate ring of the Grassmannian, and R_{0} is the algebra generated by the monomials (13). This is a flat deformation of C-algebras; see [10] and [39, §11].

The SAGBI homotopy is the following system of mp equations:

∑_{α∈\binom{[m+p]}{p}} C_{α}^{i} · t^{−w_{α}} · [α](x_{ij}t^{(i−1)(p−j)}) = 0 (i = 1, …, mp). (14)

We reduce the number of variables to mp by introducing local coordinates as follows: x_{ii} = 1 for i = 1, …, p and x_{ij} = 0 for i < j or i > m + j. Our original problem is to solve the system (14) for t = 1.

The flatness of the family of algebras R_{t} guarantees that the system (14) has the same finite number of complex solutions (counting multiplicities) for almost every t ∈ C. For t = 0 we get a system of linear equations in R_{0}:

∑_{α∈\binom{[m+p]}{p}} C_{α}^{i} · x_{α_{1},1}x_{α_{2},2}x_{α_{3},3} ⋯ x_{α_{p},p} (i = 1, …, mp). (15)

In order to solve these equations we apply the symbolic-numeric algorithm in [17], while taking advantage of the following combinatorial structures described in [39, Remark 11.11]. The common Newton polytope of the equations (15) equals the order polytope of the product of an m-chain and a p-chain (Sturmfels, Remark 11.11). We have the following combinatorial result.

Proposition 2.3. The following five numbers coincide:

• the right hand side of (3),

• the number of linear extensions of the product of a m-chain and a p-chain,

• the number of maximal chains in Young’s poset,

• the normalized volume of the order polytope, and

• the number of roots in (C^{∗})^{mp} of a generic system (15).

The equality of (4) and (5) is a special case of Kouchnirenko’s Theorem [21], as all the equations have the same Newton polytope. The order polytope has a distinguished unimodular regular triangulation with simplices indexed by the chains in Young’s poset. This regular triangulation is induced by the system of weights given in (12). We may use these weights to define a numerical homotopy for finding all isolated solutions of (15). Once this is accomplished, we trace these roots from t = 0 to t = 1 in the homotopy (14).

3. Special Schubert conditions

We first describe a purely combinatorial method for computing the number d of solution planes. Instead of the algebraic relation (1.2) we shall make use of Young’s poset which was introduced in Section 2.2. A cover [α] ⋖ [β] in Young’s poset determines a unique index j = j(α, β) for which

α_{j} + 1 = β_{j} and α_{i} = β_{i} for i ≠ j.

A chain [α^{0}] ⋖ [α^{1}] ⋖ ⋯ ⋖ [α^{l}] in Young’s poset is increasing at i if either i = 1, or else i > 1 and j(α^{i−2}, α^{i−1}) ≤ j(α^{i−1}, α^{i}). For instance, [12] ⋖ [13] ⋖ [14] ⋖ [24] is increasing at 1 and 2, but decreasing at 3.

Given positive integers r_{0}, …, r_{a}, the Pieri tree T (r_{0}, …, r_{a}) consists of all chains of length r_{0} +⋯+r_{a} in Young’s poset which start at the bottom element [1, 2, …, p] and which increase everywhere, except possibly at r_{0} + 1, r_{0} + r_{1} + 1, …, r_{0} + ⋯ + r_{a−1} + 1. Here we include all initial segments of such chains and we order the chains by inclusion. Label a node in the Pieri tree by the endpoint of the chain which that node represents. Then the sequence of labels from the root to that node is the chain which that node represents. For example, here is T (2, 2) when m = 5, p = 2:

To compute the number d, partition the integer sequence k_{1}, …, k_{n} into three parts r_{0}, …, r_{a}, r′_{0}, …, r′_{a′}, and q. Then d is the number of pairs (R, S) where R is a leaf of T (r_{0}, …, r_{a}), S is a leaf of T (r′_{0}, …, r′_{a′}), and the endpoints [α] of R and [α′] of S satisfy Pieri’s condition:

α′_{1} ≤ m + p + 1 − α_{p} < α′_{2} ≤ ⋯ < α′_{p} ≤ m + p + 1 − α_{1} (16)

Call this set of pairs Sols = Sols(r_{0}, …, r_{a}; r′_{0}, …, r′_{a′}).

For instance, d = 6 for the sequence 2, 2, 2, 2, 2 with (m, p) = (5, 2): Of the 9 pairs (R, S) of leaves of T (2, 2), only the following 6 satisfy (16): (Here we represent a leaf by its label.)

{([25], [34]), ([34], [25]), ([25], [25]), ([25], [16]), ([16], [25]), ([16], [16])} . (17)

This combinatorial rule for the number d gives the same answer as the algebraic rule (1.2) because the Pieri tree and Pieri’s condition (16) encode the structure of the cohomology ring A_{m,p} with respect to its Schur basis {S_{λ} | p ≥ λ_{1} ≥ ⋯ ≥ λ_{p+1} = 0}; see [25, §I] or [13, §9.4]. Specifically,

∏_{i=0}^{a} h_{r_{i}} = ∑_{λ_{1}+···+λ_{p}=r_{0}+···+r_{a}} K_{λ,(r_{0},... ,r_{a})} · S_{λ},

where K_{λ,(r_{0},... ,r_{a})} is the number of leaves in T (r_{0}, …, r_{a}) with label

[α(λ)] := [λ_{p} + 1, λ_{p−1} + 2, …, λ_{1} + p].

The numbers K_{λ,(r_{0},... ,r_{a})} are called Kostka numbers. In A_{m,p} we calculate

∏_{i=1}^{n} h_{k_{i}} = (∏_{i=0}^{a} h_{r_{i}}) · (∏_{j=0}^{a′} h_{r′_{j}}) · h_{q}
= ∑_{λ,μ} K_{λ,(r_{0},... ,r_{a})}K_{μ,(r′_{0},... ,r′_{a′})} · S_{λ} · S_{μ} · h_{q}.

We evaluate this expression with Pieri’s formula (Proposition 3.1 below): If λ_{1} +⋯+ λ_{p} + μ_{1} + ⋯ + μ_{p} + q = mp, then S_{λ} · S_{μ} · h_{q} is either (h_{m})^{p} or 0 depending upon whether or not

λ_{p} ≤ m − μ_{1} ≤ λ_{2} ≤ ⋯ ≤ λ_{1} ≤ m − μ_{p}.  (16′)

The condition (16′) is equivalent to (16) under the transformation λ ↔ [α(λ)].

These methods correctly enumerate the p-planes which meet each K_{1}, …, K_{n} nontrivially because, under the isomorphism between A_{m,p} and the cohomology ring of Grass(p, m + p), the indeterminate h_{k_{i}} corresponds to the cohomology class Poincaré dual to Ω_{K_{i}}, the set of p-planes which meet K_{i} nontrivially. Moreover, (h_{m})^{p} represents the class dual to a point. The Pieri tree models certain intrinsic deformations (described in §3.3 and §3.5) of the Grassmannian which establish this isomorphism, and which we shall use for computing the p-planes which meet each K_{1}, …, K_{n} nontrivially.

3.1. Basics on Schubert varieties. For vectors f_{1}, …, f_{j} in C^{m+p}, let ⟨f_{1}, …, f_{j}⟩ be their linear span. Fix the columns e_{1}, …, e_{m+p} of the identity matrix as a standard basis for C^{m+p}. For α ∈ \binom{[m+p]}{p}, define α^{∨} ∈ \binom{[m+p]}{p} by α_{j}^{∨} := m + p + 1 − α_{p+1−j}.

A sequence α ∈ \binom{[m+p]}{p} determines a Schubert variety

Ω_{α} := {X ∈ Grass(p, m + p) | dim X ∩ ⟨e_{1}, …, e_{α_{j}^{∨}}⟩ ≥ j for 1 ≤ j ≤ p}.

This variety has complex codimension |α| := α_{1} − 1 + α_{2} − 2 + ⋯ + α_{p} − p. Similarly define

Ω′_{α} := {X ∈ Grass(p, m+p) | dim X∩⟨e_{α_{j}}, …, e_{m+p}⟩ ≥ p+1−j for 1 ≤ j ≤ p}.

For a linear subspace N of C^{m+p} of dimension m + 1 − q, define the special Schubert variety

Ω_{N} := {X ∈ Grass(p, m + p) | dim X ∩ N ≥ 1}.

This has codimension q. If N = ⟨e_{1}, …, e_{m+1−q}⟩, then Ω_{N} = Ω_{[1,... ,p−1,p+q]}.

The special Schubert variety Ω_{N} is cut out by the system of \binom{m+p}{q−1} polynomial equations:

(18) X ∈ Ω_{N} ⟺ all maximal minors of [X | N] are zero,

where X ∈ Grass(p, m + p) is represented by a (m + p) × p-matrix. The Laplace expansion of these equations in terms of the Plücker coordinates of X define Ω_{N} as a subscheme of Grass(p, m + p). These equations are redundant: select m+1−q rows of [X | N] such that the corresponding maximal minor of N is invertible. Consider the set of maximal minors of [X | N] which cover all the rows selected. This gives \binom{p}{q−1} polynomial equations which generate the same ideal as all \binom{m+p}{q−1} minors of [X | N]. For a purely set-theoretic (but not scheme-theoretic) representation of Ω_{N} a further substantial reduction in the number of equations is possible using the results of [5].

An intersection Y ∩ Z of subvarieties is generically transverse if every component of Y ∩ Z has an open subset along which Y and Z meet transversally. In this case the following identity in the cohomology ring holds:

[Y ∩ Z] = [Y] · [Z],

where [W] denotes the cycle class of a subvariety W. By Kleiman’s Transversality Theorem [19], subvarieties of Grass(p, m + p) in general position meet generically transversally. Transversality and generic transversality coincide when Y ∩ Z is finite.

Proposition 3.1 (Hodge and Pedoe, 1952, Theorem III in §XIV.4). Let α, α′ ∈ \binom{[m+p]}{p} with |α|+|α′|+q = mp and let N be a linear subspace of C^{m+p} with dimension m+1−q none of whose Plücker coordinates vanish. Then the intersection

Ω_{α} ∩ Ω′_{α′} ∩ Ω_{N} (19)

either is transverse consisting of a single p-plane or is empty, depending upon whether or not (16) holds.

Proof and Algorithm: The intersection Ω_{α} ∩ Ω′_{α′} is nonempty if and only if α′_{j} ≤ α_{j}^{∨} for j = 1, …, p. These are the weak inequalities in (16). We shall assume that they hold in what follows. The p-planes in Ω_{α} ∩ Ω′_{α′} are represented by (m + p) × p-matrices X = (x_{ij}) such that

x_{i,j} = 0 for i < α′_{j} or α_{j}^{∨} < i. (20)

Consider the nonzero coordinate subspaces C_{j} := ⟨e_{α′_{j}}, …, e_{α_{j}^{∨}}⟩, set C := C_{1} + ⋯ + C_{p}, and note that p + q = ∑_{j} dim(C_{j}) ≥ dim(C). From (20) we see that X ∈ Ω_{α} ∩ Ω′_{α′} implies X ⊆ C and hence N ∩ X ⊆ N ∩ C. Therefore the triple intersection (19) is nonempty only if the following equivalent conditions hold:

dim(C ∩ N) ≥ 1 ⟺ dim(C) = p + q ⟺ the sum C = C_{1} + ⋯ + C_{p} is direct

⟺ α_{j}^{∨} < α′_{j+1} for j = 1, …, p − 1 ⟺ (3.1) holds.

In this case we determine C ∩ N by computing vectors g_{j} ∈ C_{j} such that C ∩ N = ⟨g_{1} ⊕ g_{2} ⊕ ⋯ ⊕ g_{p}⟩. (This computation is the “algorithm” part in this proof.) The desired p-plane X satisfies C ∩ N = X ∩ N, and, in view of (20), this implies X = ⟨g_{1}, …, g_{p}⟩. Transversality of (19) is verified in local coordinates for Ω_{α} ∩ Ω′_{α′} by considering p + q − 1 independent linear forms which vanish on N.

For α ∈ \binom{[m+p]}{m} define λ(α) by λ(α)_{j} := α_{j} −j for 1 ≤ j ≤ p. Then S_{λ(α)} represents the cycle class of Ω_{α} (equivalently, of Ω′_{α}). If dim N = m + 1 − q, then h_{q} is the cycle class of Ω_{N}. Suppose |α| + |α′| + q = mp. Then Proposition 3.1 implies the following identity in A_{m,p}:

S_{λ(α)} · S_{λ(α′)} · h_{q} = \{ \begin{matrix}(h_{m})^{p} & if (16) holds \\ 0 & otherwise \end{matrix} .

This identity implies (via Poincaré duality) that

S_{λ(α)} · h_{q} = ∑ S_{λ(β)},

the sum over all β with |β| = |α| + q for which α, β^{∨} satisfy Pieri’s condition (16). Call this set α ∗ q, which is also the set of endpoints of increasing chains of length q in Young’s poset that begin at α.

This last form has geometric content. In [35], explicit deformations were given that transform the irreducible intersection Ω_{α} ∩ Ω_{N} into the cycle ∑_{β∈α∗q} Ω_{β}. Moreover, the branching of the components of the cycles in these deformations reflects the branching among these increasing chains above α. This process may be iterated to transform an intersection of several special Schubert varieties into a sum of triple intersections of the form (19), indexed by pairs (R, S) ∈ Sols. From this sum, we obtain a set of start solutions indexed by Sols. Also, every intermediate cycle in these deformations consists of the same number (counting multiplicities) of p-planes. The Pieri homotopy begins with one of the start solutions and uses numerical path continuation to trace the sequence of curves defined by these deformations which connect that start solution to a solution of the original problem.

3.2. Pieri homotopy algorithm. Given linear subspaces K_{1}, …, K_{n} in general position with dim K_{i} = m + 1 − k_{i} and k_{1} + … + k_{n} = mp, first partition K_{1}, …, K_{n} into three lists:

L_{0}, …, L_{a},  L′_{0}, …, L′_{a′},  N

where dim L_{i} = m + 1 − r_{i}, dim L′_{i} = m + 1 − r′_{i}, and dim N = m + 1 − q. Construct the Pieri trees T (r_{0}, …, r_{a}) and T (r′_{0}, …, r′_{a′}), and form the set Sols. Change coordinates so that L_{0} = ⟨e_{1}, …, e_{m+1−r_{0}}⟩ and L′_{0} = ⟨e_{p+r′_{0}}, …, e_{m+p}⟩. Set τ := max{r_{1} + ⋯ + r_{a}, r′_{1} + ⋯ + r′_{a′}}.

Given a chain R in the Pieri tree and a positive integer k, let R(k) be the kth element in that chain, or, if k exceeds the length of R, then let R(k) be the endpoint of R. For each (R, S) ∈ Sols and k from τ to 0 we shall construct (in Definition 3.2 below) one-parameter families Z_{R,k}(t) and Z′_{S,k}(t) of pure-dimensional subvarieties of Grass(p, m + p) with the following properties:

1. Z_{R,k}(t) ⊂ Ω_{R(r_{0}+k)} and Z′_{S,k}(t) ⊂ Ω′_{S(r′_{0}+k)}.

2. For t = 0 or 1 and each k, Z_{R,k}(t)∩Z′_{S,k}(t)∩Ω_{N} is transverse and 0-dimensional.

3. Z_{R,τ}(t) = Ω_{R(r_{0}+τ)} and Z′_{S,τ}(t) = Ω′_{S(r_{0}+τ)} .

4. Z_{R,k+1}(1) is a component of Z_{R,k}(0). Likewise, Z′_{S,k+1}(1) is a component of Z′_{S,k}(0).

5. Z_{R,0}(1) = Ω_{L_{0}} ∩ ⋯ ∩ Ω_{L_{a}} and Z′_{S,0}(1) = Ω_{L′_{0}} ∩ ⋯ ∩ Ω_{L′_{a′}}.

Property 4 is a consequence of Proposition 3.4, the others follow from the assumption of genericity and the definition (Definition 3.2) of the families Z_{R,k}(t) and Z′_{S,k}(t).

By 2, the family W_{(R,S),k}(t) over C whose fibre at general t (including t = 0 and t = 1) is

W_{(R,S),k}(t) := Z_{R,k}(t) ∩ Z′_{S,k}(t) ∩ Ω_{N}

consists of a finite number of curves. In fact, for general t (including t = 0 and t = 1), W_{(R,S),k}(t) has the following general form (see Definition 3.2 for the precise form):

W_{(R,S),k}(t) = Ω_{α} ∩ Ω′_{α′} ∩ Ω_{M_{1}} ∩ ⋯ ∩ Ω_{M_{s}},

where M_{1}, …, M_{s} are linear subspaces with M_{s} = N, which depend upon R, S, k, the subspaces L_{1}, …, L_{a}, L′_{1}, …, L′_{a′}, and at most two of the M_{i} depend upon t. Also, α and α′ depend upon R, S, and k with the typical case being α = R(r_{0} + k) and α′ = S(r′_{0} + k).

The numerical homotopy defined by the curves W_{(R,S),k}(t) may be expressed in a parameterization X = (x_{i,j}) of an open subset of Ω_{α} ∩ Ω′_{α′}:

x_{i,j} = 0 if i < α′_{j} or α_{j}^{∨} < i  and  x_{δ_{j}, j} = 1, (21)

where δ := S(r′_{0} + τ). The equations for W_{(R,S),k}(t) are then

maximal minors [X | M_{i}] = 0 i = 1, …, s.

The curves of W_{(R,S),k}(t) define the sequences of homotopies in the Pieri homotopy algorithm as follows: For (R, S) ∈ Sols, let X_{(R,S),τ} be the (unique by Proposition 3.1) p-plane in Ω_{R(r_{0}+τ)} ∩Ω′_{S(r′_{0}+τ)} ∩Ω_{N} = W_{(R,S),τ}(1). By 3 and 4, X_{(R,S),τ} ∈ W_{(R,S),τ−1}(0) and hence lies on a unique curve in W_{(R,S),τ−1}(t). Use numerical path continuation to trace this curve from t = 0 to t = 1 to obtain X_{(R,S),τ−1}, which is a point of W_{(R,S),τ−2}(0), by 4. Then X_{(R,S),τ−1} lies on a unique curve in W_{(R,S),τ−2}(t), which we trace to find X_{(R,S),τ−2} ∈ W_{(R,S),τ−2}(1). Continuing this process, after tracing τ curves, we obtain X_{(R,S),0} ∈ W_{(R,S),0}(1), which is a solution to the original system, by 5. We show shall prove in Theorem 3.5 that {X_{(R,S),0} | (R, S) ∈ Sols} consists of all the solutions to the original system.

3.3. Definition of the moving cycles Z_{R,k}(t). The cycle Z_{R,k}(t) will depend upon the choice of a general upper triangular (m + p) × (m + p)-matrix F with 1’s on its anti-diagonal,

( \begin{matrix}∗ & ∗ & 1 \\ ∗ & ... \\ 1 & & 0 \end{matrix} ) ,

the kth link in the chain R, and the data L_{1}, …, L_{a}. The key ingredient of this definition of Z_{R,k}(t) is the construction of a one-parameter family of linear subspaces Λ_{i}(t) in Definition 3.3, which depends upon F. The matrix F is fixed throughout the algorithm, its purpose is that ⟨e_{m+p−j}, …, e_{m+p}⟩ equals the span of the first j columns of F, and these columns are in general position with e_{1}, …, e_{m+p}. The subtle linear degeneracies of Λ_{i}(t) as t → 0 are at the heart of this homotopy algorithm, as well as the explicit proof of Pieri’s formula [35, Theorem 3.6], which we state below (Proposition 3.4).

Definition 3.2.

1. If r_{1} + ⋯ + r_{a} ≤ k, then set Z_{R,k}(t) := Ω_{R(r_{0}+τ)}.

2. Otherwise, define c by r_{1} + ⋯ + r_{c−1} ≤ k < r_{1} + ⋯ + r_{c}, and set i := k − r_{1} − ⋯ − r_{c−1}, α := R(r_{0} + r_{1} + ⋯ + r_{c−1}), and β := R(r_{0} + k).

(a) If i > 0 and β_{p} > α_{p}, then β + (0, …, 0, r_{c} − i) = R(r_{0} + ⋯ + r_{c}), and we

set

Z_{R,k}(t) := Ω_{R(r_{0}+···+r_{c})} ∩ Ω_{L_{c+1}} ∩ ⋯ ∩ Ω_{L_{a}}.

(b) Otherwise, let Λ_{i}(t) be the 1-parameter family of linear subspaces given by Definition 3.3, where we let L := L_{c} and r := r_{c}. If i = 0, then Λ_{0}(1) = L_{c}, α = β = R(r_{0} + k), and we set

Z_{R,k}(t) := Ω_{R(r_{0}+k)} ∩ Ω_{Λ_{0}(t)} ∩ Ω_{L_{c+1}} ∩ ⋯ ∩ Ω_{L_{a}}.

If i > 0, let j be maximal such that β_{j} > α_{j}. Then j < p as j = p is case 2(a). Set

Z_{R,k}(t) := Ω_{R(r_{0}+k)} ∩ Ω_{Λ_{i}(t)∩⟨e_{1},... ,e_{β{∨ \atop p+1−j}} ⟩} ∩ Ω_{L_{c+1}} ∩ ⋯ ∩ Ω_{L_{a}}.

We define Z′_{S,k}(t) similarly, but with the matrix F replaced by a lower triangular matrix with 1’s on its diagonal, and ⟨e_{1}, …, e_{β_{p+1−j}^{∨}}⟩ replaced by ⟨e_{β_{j}}, …, e_{m+p}⟩.

Definition 3.3. Let F be an upper triangular matrix with 1’s on the anti-diagonal and L be a general (m + 1 − r)-plane, represented as a (m + p) × (m + 1 − r)-matrix with columns l_{1}, …, l_{m+1−r}. Construct a (m + p) × p-matrix U = (u_{1}, …, u_{p}) as follows: Reverse the last m + p − α_{p} columns of F, then remove the columns indexed by α_{1}^{∨}, …, α_{p}^{∨}.

For each 0 ≤ i < r, define a one-parameter family of (m+p)×(m+1−r)-matrices Λ_{i}(t) for t ∈ C as follows:

1. If i = 0, then the bth column of Λ_{0}(t) is t · l_{b} + (1 − t) · u_{b}.

2. For 0 < i < r, the bth column of Λ_{i}(t) is

\begin{matrix}t · u_{b+i−1} + (1 − t) · u_{b+i} & b + i − 1 < α_{p} − p \\ t · u_{b+i−1} + (1 − t) · u_{p+1+i−r} & b + i − 1 = α_{p} − p \\ u_{b+i−1} & b + i − 1 > α_{p} − p \end{matrix}

3.4. An example. We give an example illustrating these definitions and the Pieri homotopy algorithm. Let L_{0}, L_{1}, L′_{0}, L′_{1}, and N be general 4-planes in C^{7}. We give a sequence of homotopies W_{(R,S),k}(t) for k = 2, 1, 0 for finding one of the six 2-planes which meet each of the five given 4-planes nontrivially.

Here, (m, p) = (5, 2) and k_{1} = ⋯ = k_{5} = 2 so that τ = 2. Construct the set Sols as in (17). Let (R, S) ∈ Sols be the following two sequences:

R := [12] ⋖ [13] ⋖ [14] ⋖ [24] ⋖ [25], S := [12] ⋖ [13] ⋖ [14] ⋖ [15] ⋖ [16].

Let e_{1}, …, e_{7} be the columns of a 7 × 7-identity matrix, a basis for C^{7}. Suppose that L_{0} = ⟨e_{4}, e_{5}, e_{6}, e_{7}⟩ and L′_{0} = ⟨e_{1}, e_{2}, e_{3}, e_{4}⟩ and represent L_{1} and L′_{1} as 7×4-matrices. Then Ω_{[14]} ∩ Ω_{L} ∩ Ω′_{[14]} ∩ Ω_{L′} ∩ Ω_{N} is the set of 2-planes which meet all five linear subspaces nontrivially.

We first find the plane X_{(R,S),2} ∈ Ω_{[25]} ∩Ω′_{[16]} ∩Ω_{N}, using the algorithm in the proof of Proposition 3.1. Suppose that N has the form

[ \frac{n}{I} ] ,

where I is the 4 × 4 identity matrix, and n is a 3 × 4-matrix. In this case, C_{1} = ⟨e_{1}, e_{2}, e_{3}⟩ and C_{2} = ⟨e_{6}⟩, hence C = ⟨e_{1}, e_{2}, e_{3}, e_{6}⟩. Thus the intersection C ∩ N is generated by the third column of N, and so X_{(R,S),2} is represented by the matrix:

X_{(R,S),2} = [ \begin{matrix}n_{13} & 0 \\ n_{23} & 0 \\ n_{33} & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 1 \\ 0 & 0 \end{matrix} ] .

Following Definition 3.3, we have:

Λ_{0}(t) = ⟨tl_{1} + (1 − t)u_{1}, tl_{2} + (1 − t)u_{2}, tl_{3} + (1 − t)u_{3}, tl_{4} + (1 − t)u_{4}⟩,
Λ_{1}(t) = ⟨tu_{1} + (1 − t)u_{2}, tu_{2} + (1 − t)u_{5}, u_{3}, u_{4}⟩.

Λ′_{i}(t) is defined similarly.

We describe the families W_{(R,S),k}(t) for k = 2, 1, 0 in local coordinates for Ω_{[14]}∩Ω′_{[14]} determined by the sequence [16]:

X = [ \begin{matrix}1 & 0 \\ x_{21} & 0 \\ x_{31} & 0 \\ x_{41} & x_{42} \\ 0 & x_{52} \\ 0 & 1 \\ 0 & x_{72} \end{matrix} ] .

The family W_{(R,S),2}(t) is the constant family {X_{(R,S),2}} = Ω_{25} ∩ Ω′_{16} ∩ Ω_{N}. Assuming that n_{13}, which is the [1457]th Plücker coordinate of N, is non-zero, then X_{(R,S),2} may be expressed in these local coordinates:

X_{(R,S),2} = [ \begin{matrix}1 & 0 \\ n_{23}/n_{13} & 0 \\ n_{33}/n_{13} & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 1 \\ 0 & 0 \end{matrix} ] .

When k = 1, first consider the definition of Z_{R,1}(t). Here we are in case 2(b) with β = [24] and i > 0, so that β^{∨} = [46]. Since Λ_{1}(t) ⊂ ⟨e_{1}, …, e_{6}⟩, we have Z_{R,1}(t) = Ω_{[24]} ∩ Ω_{Λ_{1}(t)}. For the definition of Z′_{S,1}(t), we are in case 2(a), so that Z′_{S,1}(t) = Ω′_{[16]}. Hence

W_{(R,S),1}(t) = Ω_{[24]} ∩ Ω_{Λ_{1}(t)} ∩ Ω′_{[16]} ∩ Ω_{N}.

This has 3 linear equations x_{42} = x_{52} = x_{72} = 0, which describe Ω_{[24]} ∩ Ω′_{[16]}, and 7 non-trivial equations, the vanishing of the maximal minors of [X | Λ_{1}(t)] and [X | N], which describe Ω_{Λ_{1}(t)} ∩ Ω_{N}.

For k = 0, i = 0 and we are in case 2(b) for both Z_{R,0}(t) and Z′_{S,0}(t) so that

W_{(R,S),0}(t) = Ω_{[14]} ∩ Ω_{Λ_{0}(t)} ∩ Ω′_{[14]} ∩ Ω_{Λ′_{0}(t)} ∩ Ω_{N}.

This has 21 non-trivial equations, the vanishing of the maximal minors of [X | Λ_{0}(t)], [X | Λ′_{0}(t)], and [X | N].

3.5. Proof of correctness. We describe the Pieri deformations linking the families Z_{R,k}(t) for k from r_{1} + ⋯ + r_{c−1} to r_{1} + ⋯ + r_{c}, which establishes Property 4 of Z_{R,k}(t) in §3.2. We also show that the set {X_{(R,S),0} | (R, S) ∈ Sols} consists of all the solutions to the original system.

Consider the dynamic part of Z_{R,k}(t), namely whichever of

Ω_{R(r_{0}+···+r_{c})}, Ω_{R(r_{0}+k)} ∩ Ω_{Λ_{0}(t)}, or Ω_{R(r_{0}+k)} ∩ Ω_{Λ_{i}(t)∩⟨e_{1},... ,e_{β{∨ \atop p+1−j}} ⟩},

appeared in the definition of Z_{R,k}(t). We call this cycle Y_{α,β,L}(t), where L := L_{c} and α = R(r_{0} + ⋯ + r_{c−1}), and β = R(r_{0} + k).

For β ∈ α ∗ i and γ ∈ α ∗ (i + 1) write β ≺_{α} γ if γ covers β and j(β, γ) ≥ j(α, β) := max{j | β_{j} > α_{j}}. This partitions α ∗ (i + 1) into sets {γ | β ≺_{α} γ} for β ∈ α ∗ i.

Proposition 3.4. [35, Theorem 3.6] Let α, β, i, r, L, Λ_{i}(t), and Y_{α,β,L}(t) be as above. Then

1. For all t, Ω_{α} ∩ Ω_{Λ_{0}(t)} is generically transverse.

2. Y_{α,β,L}(t) is free of multiplicities for all t and irreducible for t ≠ 0.

3. If i ≠ r − 1, then Y_{α,β,L}(0) = ∑_{β≺_{α}γ} Y_{α,γ,L}(1).

4. If β ∈ α ∗ (r − 1), then Y_{α,β,L}(0) = ∑_{β≺_{α}γ} Ω_{γ}.

By 3, the cycle class of ∑_{β∈α∗i} Y_{α,β,L}(t) is independent of i and t, and it equals the cycle class of Ω_{α} ∩ Ω_{Λ_{0}(1)} = Ω_{α} ∩ Ω_{L}. By 4, we see that the cycle classes of Ω_{α} ∩ Ω_{L} and ∑_{β∈α∗r} Ω_{β} coincide, furnishing another proof of Pieri’s formula. Property 4 of Z_{R,k}(t) follows from assertion 3.

Theorem 3.5. When K_{1}, …, K_{n} are generic, the Pieri homotopy algorithm finds all p-planes which meet each K_{1}, …, K_{n} nontrivially. That is,

{X_{(R,S),0} | (R, S) ∈ Sols} = Ω_{K_{1}} ∩ ⋯ ∩ Ω_{K_{n}}.

Proof. Note that for any R, S ∈ Sols, the families Z_{R,k}(t), Z′_{S,k}(t), and W_{(R,S),k}(t) depend only upon the initial segments R(0), …, R(r_{0} + k) and S(0), …, S(r′_{0} + k) of R and S.

By construction, the original system Ω_{K_{1}} ∩ ⋯ ∩ Ω_{K_{l}} coincides with W_{(R,S),0}(1), for any (R, S) ∈ Sols. We inductively construct chains R ∈ T (r_{0}, …, r_{a}) and S ∈ T (r′_{0}, …, r′_{a′}), and p-planes X_{k} for 0 ≤ k ≤ τ such that

X_{k+1} ∈ W_{(R,S),k}(0) ∩ W_{(R,S),k+1}(1)

and X_{k}, X_{k+1} lie on the same curve of W_{(R,S),k}(t). Then X_{τ} is the start solution X_{(R,S),τ}, which shows that X_{0} ∈ {X_{(R,S),0} | (R, S) ∈ Sols}.

First set R(0), …, R(r_{0}) to be the unique chain from [1, …, p] to [1, …, p−1, p+ r_{0}], and similarly for S(0), …, S(r′_{0}). Then X_{0} ∈ W_{(R,S),0}(1) and hence lies on a unique curve in W_{(R,S),0}(t). Let X_{1} be the point on that curve with t = 0. By Proposition 3.4 (4),

X_{1} ∈ Y_{(R(r_{0}),R(r_{0}),L_{1})}(0) = ∑_{β∈R(r_{0})∗1} Y_{(R(r_{0}),β,L_{1})}(1).

Let R(r_{0} +1) be the index β such that X_{1} ∈ Y_{(R(r_{0}),β,L_{1})}(1). Define S(r′_{0} +1) similarly. Then X_{1} ∈ W_{(R,S),1}(1).

In general, suppose that we have constructed R(0), …, R(r_{0} +k), S(0), …, S(r′_{0} + k), and X_{k} ∈ W_{(R,S),k}(1). Then X_{k} lies on a unique curve in W_{(R,S),k}(t). Let X_{k+1} be the point on that curve at t = 0. Let c be minimal subject to k < r_{1} + ⋯ + r_{c} and set α = R(r_{0} + ⋯ + r_{c−1}) and β = R(r_{0} + k). If k + 1 < r_{1} + ⋯ + r_{c}, then by Proposition 3.4 (3), there is a unique index γ ∈ β ∗ 1 such that X_{k+1} ∈ Y_{α,γ,L_{c}}(1). If k +1 = r_{1} +⋯+r_{c} then by Proposition 3.4 (4), there is a unique index γ ∈ β ∗1 such that X_{k+1} ∈ Ω_{γ}. Set R(r_{0} + k + 1) = γ and likewise define S(r′_{0} + k + 1). Continuing in this fashion, we construct the chains R and S, and X_{j} for 0 ≤ j ≤ τ.

We show that R increases everywhere, except possibly at 1, r_{0} + 1, …, r_{0} + ⋯ + r_{a−1}+1, and hence R ∈ T (r_{0}, …, r_{a}). Similar arguments show that S ∈ T (r′_{0}, …, r′_{a′}), which will complete the proof. Suppose k +1 ∉ {1, r_{1} +1, …, r_{1} +⋯+r_{a−1} +1}. Let c be minimal subject to k < r_{1} +⋯+r_{c} and let α = R(r_{0} +⋯+r_{c−1}), β = R(r_{0} +k), and γ = R(r_{0} + k + 1). Then by Proposition 3.4 (3) and (4), β ≺_{α} γ. The condition j(β, γ) ≥ j(α, β) in the definition of β ≺_{α} γ ensures that R increases at k + 1.

4. Homotopy continuation of overdetermined systems

Numerical homotopy continuation is a method for finding the isolated solutions of a system

F (X) = 0 (22)

where F = (f_{1}, …, f_{n}) are polynomials in the variables X = (x_{1}, …, x_{N}). First, a homotopy H(X, t) is found with the following properties:

1. H(X, 1) = F (X).

2. The isolated solutions of H(X, 0) = 0 are known.

3. The system H(X, t) = 0 defines finitely many (rational) curves σ_{i}(t), and each isolated solution of (22) is connected to an isolated solution σ_{i}(0) of H(X, 0) = 0 by one of these curves.

Given such a homotopy, numerical path continuation is used to trace these curves from solutions of H(X, 0) = 0 to solutions of the original system (22). When there are fewer solutions to F (X) = 0 than to H(X, 0) = 0, some curves will diverge or become singular as t → 1, and it is expensive to trace such a curve.

When N = n, the system (22) is square and the homotopy

H(X, t) := tF (X) + (1 − t)G(X), (23)

where G(X) = (x_{1}^{d_{1}} − a_{1}, …, x_{N}^{d_{N}} − a_{N}) with d_{i} := deg(f_{i}) and a_{i} ≠ 0, gives ∏ d_{i} curves. This is the Bézout bound for a generic dense system F.

In practice, F (X) = 0 may have fewer than ∏ d_{i} solutions and we desire a homotopy with no divergent curves. Methods for such deficient systems which reduce the number of divergent curves are developed in [23, 22, 24]. When the polynomials f_{1}, …, f_{n} have special forms [27, 28], then such homotopies (23) are constructed where G(X) shares this special form. When the polynomials f_{1}, …, f_{n} are sparse, polyhedral methods [41, 17] give a homotopy. The SAGBI homotopy algorithm (§2.3) is in the same spirit. We exploit a special feature of the coordinate ring of the Grassmannian to obtain a homotopy between the system (2.2) we wish to solve and one (2.12) whose solution may be obtained using polyhedral methods. Moreover, there are generically no divergent curves to be followed.

The overdetermined situation of n > N is more delicate. One difficulty is finding a homotopy H(X, t) for an overdetermined system as generic perturbations of F have no solutions. In [32, §2], this difficulty is avoided as follows: The system F = (f_{1}, …, f_{n}) is replaced by N random linear combinations of the f_{1}, …, f_{n} yielding a square system whose isolated solutions include all isolated solutions of F (X) = 0, but typically many more. They then find all isolated solutions of this random square subsystem.

For the Gröbner and Pieri homotopy algorithms, we gave (in §§2.2 and 3.2–3) homotopies H(X, t) = (h_{1}(X, t), …, h_{n}(X, t)) and solutions σ_{i}(0) at t = 0 as above. For these, there are generically no divergent curves. To efficiently follow the curves σ_{i}(t), we select a square subsystem MH(X, t) of H(X, t) (M is an (N × n)-matrix). If the Jacobian of MH at each σ_{i}(0) has the same rank (N) as does the Jacobian of H(X, t), then the curves σ_{i}(t) remain components of the algebraic set defined by the equations

MH(X, t) = 0.

Moreover, other components of this set meet the curves σ_{i}(t) in at most finitely many points t in C − {0}. Thus, we may use the square subsystem MH(X, t) to trace the curves σ_{i}(t) along some path from 0 to 1 in the complex plane. We remark that in practice, M may be chosen at random.

5. Applications

The algorithms of Sections 2 and 3 are useful for studying both the pole assignment problem in systems theory [8] and real enumerative geometry [33].

We describe the connection to the control of linear systems following [8]. Suppose we have a system (for example, a mechanical linkage) with inputs u ∈ R^{m} and outputs y ∈ R^{p} for which there are internal states x ∈ R^{n} such that the evolution of the system is governed by the first order linear differential equation

\begin{aligned}ẋ & = & Ax + Bu, \\ y & = & Cx. \end{aligned} (24)

If the input is controlled by constant output feedback, u = F y, then we obtain

ẋ = (A + BF C)x.

The natural frequencies of the controlled system are the roots s_{1}, …, s_{n} of

φ(s) := det(sI − A − BF C). (25)

The pole assignment problem asks, given a system (24) and a polynomial φ(s) of degree n, which feedback laws F satisfy (25)?

A standard transformation (cf. [8, §2]) transforms the input data A, B, C into matrices N(s), D(s) of polynomials with det(D(s)) = det(sI −A) and N(s)D(s)^{−1} = C(sI − A)^{−1}B such that

φ(s) = det [ \begin{matrix}F & D(s) \\ I & N(s) \end{matrix} ] . (26)

Here I is the p × p-identity matrix and the feedback law F is an m × p-matrix. If we let

X := [ \begin{matrix}F \\ I \end{matrix} ] and K(s) := [ \begin{matrix}D(s) \\ N(s) \end{matrix} ] ,

then F gives local coordinates on Grass(p, m + p) and (26) is equivalent to

X ∩ K(s_{i}) ≠ {0} for i = 1, …, n.

These conditions are independent for generic A, B, C and distinct s_{i}, hence n ≤ mp is necessary for there to be any feedback laws F. The critical case of n = mp is an instance of the situation in §2.

In [9] homotopy continuation was used to solve a specific feedback problem when (m, p) = (3, 2). From this result, they deduced that the pole assignment problem is not in general solvable by radicals. Despite this success, only few non-trivial examples have been computed in the control theory literature [30].

An important question is whether a given system may be controlled by real output feedback [42, 7, 29]. That is, if all roots of φ(s) are real, are there real feedback laws F satisfying (26)? Real enumerative geometry [33] asks a similar question: Are there real linear subspaces K_{1}, …, K_{n} in general position with dim K_{i} = m + 1 − k_{i} and k_{1} + ⋯ + k_{n} = mp such that all p-planes meeting each K_{i} nontrivially are real? When either m or p is 2 [34], n ≤ 5 [35], or when m = p = 3 and the k_{i} = 1 [33], the answer is yes. In fact, the Pieri homotopies arose from these investigations.

B. Shapiro and M. Shapiro give a precise conjecture relating both applications. Suppose

K_{i}(s) := [γ(s), γ′(s), γ′′(s), …, γ^{(m+1−k_{i})}(s)],

where γ(s) is a parameterization of a rational normal (non-degenerate) curve in R^{m+p} of degree m + p − 1. One such choice is

(27) γ(s) = transpose[1, s, s^{2}, …, s^{m+p−1}].

Geometrically, K_{i}(s) is the (m + 1 − k_{i})-plane which osculates the curve γ(s) at s. Such osculating m-planes have been used to prove non-degeneracy results in control theory.

Conjecture 5.1 (B. Shapiro and M. Shapiro). Let s_{1}, …, s_{n} be distinct real numbers and suppose K_{i}(s_{i}) osculates γ at s_{i} and k_{1} + ⋯ + k_{n} = mp. Then each of the finitely many p-planes X ⊂ C^{m+p} satisfying X ∩ K_{i}(s_{i}) ≠ {0} for i = 1, …, n is defined over the reals.

6. Computational results

Our algorithms have been tested successfully in MATLAB, finding all 14 solutions in the case (m, p) = (4, 2) for both the SABGI and Gröbner homotopy algorithm, and all 15 solutions when (m, p) = (6, 2) and k_{1} = ⋯ = k_{6} = 2 for the Pieri homotopy algorithm.

At present, the SAGBI and Gröbner homotopy algorithms have been fully implemented. Some timings from trial runs of these algorithms on a Sparc 20 are displayed in Table 1. The input for these were mp random complex (m + p) × m-matrices.

We provide a comparison to methods based upon Gröbner bases. Table 1 also gives the time on the Sparc 20 for the system Singular [15] to compute a degree reverse lexicographic Gröbner basis for the polynomial systems:

det [X | K_{i}] = 0, for i = 3, …, K_{mp}.

Here X is expressed in local coordinates for Ω_{[13]} ∩ Ω′_{[13]} and the K_{3}, …, K_{mp} are (m + p) × m-matrices with random integral entries between −4 and 4. A degree reverse lexicographic Gröbner basis is the input for some alternative numerical polynomial systems solvers (e.g. eigenvalue methods [3]). We note that the Gröbner basis calculation did not terminate within one week in the case (m, p) = (6, 2).

m p d(m, p) SABGI homotopy Gröbner homotopy Gröbner basis
3 2  5  <1  < 0.5  < 0.5
4 2 14  47  6  19
5 2 42  373  408  149, 897
6 2 132  3, 364  8, 626  ∞

Table 1. Time (in seconds)

The final version of this paper will include data from implementations of the Pieri homotopy algorithm, as well.

References

[1] E. Allgower and K. Georg, Numerical Continuation Methods, An Introduction, no. 13 in Computational Mathematics, Springer-Verlag, 1990.

[2] E. Allgower and K. Georg, Numerical path following. To appear in the Handbook of Numerical Analysis, edited by P. G. Ciarles and J.L. Lions, North Holland, 1997.

[3] W. Auzinger and H. Stetter, An elimination algorithm for the computation of the zeros of a system of multivariate polynomial equations, in Numerical Mathematics, Singapore, 1988, R. Agarwal, Y. Chow, and S. Wilson, eds., vol. 86 of Internat. Schriftenreihe Numer. Math., Birkhäuser, 1988, pp. 11–30.

[4] D. N. Bernstein, The number of roots of a system of equations, Funct. Anal. Appl., 9 (1975), pp. 183–185.

[5] W. Bruns and R. Schwänzl, The number of equations defining a determinantal variety, Bull. London Math. Soc., 22 (1990), pp. 439–445.

[6] W. Bruns and U. Vetter, Determinantal Rings, Lecture Notes in Math., Vol. 1327, Springer-Verlag, 1988.

[7] C. I. Byrnes, Control theory, inverse spectral problems, and real algebraic geometry, in Diffferential Geometric Methods in Control Theory, R. W. Brockett, R. S. Millman, and H. J. Sussmann, eds., Birkhäuser, Boston, 1982, pp. 192–208.

[8] C. I. Byrnes, Pole assignment by output feedback, in Three Decades of Mathematical Systems Theory, H. Nijmeijer and J. M. Schumacher, eds., vol. 135 of Lecture Notes in Control and Inform. Sci., Springer-Verlag, Berlin, 1989, pp. 31–78.

[9] C. I. Byrnes and P. K. Stevens, Global properties of the root-locus map, in Feedback Control of Linear and Non-Linear Systems, D. Hinrichsen and A. Isidori, eds., vol. 39 of Lecture Notes in Control and Inform. Sci., Springer-Verlag, Berlin, 1982.

[10] A. Conca, J. Herzog, and G. Valla, Sagbi bases with applications to blow-up algebras, J. reine angew. Math., 474 (1996), pp. 113–138.

[11] F.-J. Drexler, Eine Methode zur Berechnung sämtlicher Lösungen von Polynomgleichungssystemen, Numer. Math., 29 (1977), pp. 45–58.

[12] D. Eisenbud, Commutative Algebra With a View Towards Algebraic Geometry, no. 150 in GTM, Springer-Verlag, 1995.

[13] W. Fulton, Young Tableaux, Cambridge University Press, 1996.

[14] W. I. García, C. B.; Zangwill, Finding all solutions to polynomial systems and other systems of equations, Math. Programming, 16 (1979), pp. 159–176.

[15] G.-M. Greuel, G. Pfister, and H. Schönemann, Singular: A system for computation in algebraic geometry and singularity theory. Available via anonymous ftp from

helios.mathematik.uni-kl.de, 1996.

[16] W. V. D. Hodge and D. Pedoe, Methods of Algebraic Geometry, vol. II, Cambridge University Press, 1952.

[17] B. Huber and B. Sturmfels, A polyhedral method for solving sparse polynomial systems, Math. Comp., 64 (1995), pp. 1541–1555.

[18] A. Khovanskii, Newton polyhedra and the genus of complete intersections, Funct. Anal. Appl., 12 (1978), pp. 38–46.

[19] S. Kleiman, The transversality of a general translate, Comp. Math., 28 (1974), pp. 287–297.

[20] S. Kleiman and D. Laksov, Schubert calculus, Amer. Math. Monthly, 79 (1972), pp. 1061– 1082.

[21] A. G. Kouchnirenko, A Newton polyhedron and the number of solutions of a system of k equations in k unknowns, Usp. Math. Nauk., 30 (1975), pp. 266–267.

[22] T. Li and T. Sauer, A simple homotopy for solving deficient polynomial systems, Japan J. Appl. Math., 6 (1989), pp. 409–419.

[23] T. Li, T. Sauer, and J. Yorke, The random product homotopy and deficient polynomial systems, Numer. Math., 51 (1987), pp. 481–500.

[24] T. Li and X. Wang, Nonlinear homotopies for solving deficient polynomial systems with parameters, SIAM J. Numer. Anal., 29 (1992), pp. 1104–1118.

[25] I. G. Macdonald, Symmetric Functions and Hall Polynomials, Oxford University Press, 1995. second edition.

[26] A. Morgan, Solving polynomial systems using continuation for engineering and scientific problems, Prentice-Hall, Inc., Englewood Cliffs, NJ, 1987.

[27] A. Morgan and A. Sommese, A homotopy method for solving general polynomial systems that respects m-homogeneous structures, Appl. Math. and Comp., 24 (1987), pp. 101–113.

[28] A. Morgan, A. Sommese, and C. Wampler, A product-decomposition bound for Bezout numbers, SIAM J. Num. Anal., 32 (1995), pp. 1308–11325.

[29] J. Rosenthal, J. M. Schumacher, and J. C. Willems, Generic eigenvalue assignment by memoryless real output feedback, Systems & Control Lett., 26 (1995), pp. 253–260.

[30] J. Rosenthal and F. Sottile, Some remarks on real and complex output feedback. MSRI preprint # 1997 - 012, 1997. http://www.nd.edu/˜ rosen/pole.

[31] H. Schubert, Beziehungen zwischen den linearen Räumen auferlegbarennn charakteristischen Bedingungen, Math. Ann., 38 (1891), pp. 588–602.

[32] A. Sommese and C. Wampler, Numerical algebraic geometry, in The Mathematics of Numerical Analysis (Park City, UT, 1995), no. 32 in Lectures in Appl. Math., Amer. Math. Soc., Providence, RI, 1996, pp. 749–763.

[33] F. Sottile, Enumerative geometry for real varieties, in Algebraic Geometry, Santa Cruz, 1995, J. Kollár, ed., vol. 56, Part 2 of Proc. Sympos. Pure Math., Amer. Math. Soc., 1997.

[34] , Enumerative geometry for the real Grassmannian of lines in projective space, Duke Math. J., 87 (1997), pp. 59–85.

[35] , Pieri’s formula via explicit rational equivalence. Can. J. Math., to appear, 1997.

[36] R. Stanley, Some conbinatorial aspects of the Schubert calculus, in Combinatoire et Représentation du Groupe Symétrique, D. Foata, ed., vol. 579 of Lecture Notes in Math., Springer-Verlag, 1977, pp. 217–251.

[37] J. Stembridge, A maple package for symmetric functions, J. Symb. Comp., 20 (1995), pp. 755– 768. ftp://ftp.math.umich.edu/pub/jrs/software.

[38] B. Sturmfels, Algorithms in Invariant Theory, Texts and Monographs in Symbolic Computation, Springer-Verlag, 1993.

[39] , Gröbner Bases and Convex Polytopes, vol. 8 of University Lecture Series, American Math. Soc., 1996.

[40] J. Verschelde, K. Gatermann, and R. Cools, Mixed volume computation by dynamic lifting applied to polynomial system solving, Discrete Comput. Geometry, 16 (1996), pp. 69–112.

[41] J. Verschelde, P. Verlinden, and R. Cools, Homotopies exploitating Newton polytopes for solving sparse polynomial systems, SIAM J. Num. Anal., 31 (1994), pp. 915–930.

[42] J. C. Willems and W. H. Hesselink, Generic properties of the pole placement problem, in Proc. of the 7th IFAC Congress, 1978, pp. 1725–1729.

Birkett Huber, Department of Mathematics, Texas A&M University, College Station, Texas 77843

  E-mail address: birk@math.cornell.edu

Frank Sottile, Department of Mathematics, University of Toronto, 100 St. George Street, Toronto. Ont., M5S 3G3, CANADA

  E-mail address: sottile@math.toronto.edu

Bernd Sturmfels, Department of Mathematics, University of California, Berkeley CA 94720, USA

  E-mail address: bernd@math.berkeley.edu