Certifying isolated singular points and their multiplicity structure
aa r X i v : . [ m a t h . AG ] S e p Certifying isolated singular points and their multiplicitystructure
Jonathan D. Hauenstein ∗ Department of Applied andComputational Mathematicsand Statistics,University of Notre DameNotre Dame, IN, 46556, USA [email protected]
Bernard MourrainInria Sophia Antipolis2004 route des Lucioles,B.P. 93,06902 Sophia Antipolis,Cedex France [email protected]
Agnes Szanto † Department of Mathematics,North Carolina State UniversityCampus Box 8205,Raleigh, NC, 27965, USA. [email protected]
10 July 2015
Abstract
This paper presents two new constructions related to singular solutions of polynomialsystems. The first is a new deflation method for an isolated singular root. This constructionuses a single linear differential form defined from the Jacobian matrix of the input, and definesthe deflated system by applying this differential form to the original system. The advantagesof this new deflation is that it does not introduce new variables and the increase in thenumber of equations is linear instead of the quadratic increase of previous methods. Thesecond construction gives the coefficients of the so-called inverse system or dual basis, whichdefines the multiplicity structure at the singular root. We present a system of equations in theoriginal variables plus a relatively small number of new variables. We show that the roots ofthis new system include the original singular root but now with multiplicity one, and the newvariables uniquely determine the multiplicity structure. Both constructions are “exact” in thatthey permit one to treat all conjugate roots simultaneously and can be used in certificationprocedures for singular roots and their multiplicity structure with respect to an exact rationalpolynomial system.
Our motivation for this work is twofold. On one hand, in a recent paper [1], two of the co-authorsof the present paper studied a certification method for approximate roots of exact overdeterminedand singular polynomial systems, and wanted to extend the method to certify the multiplicitystructure at the root as well. Since all these problems are ill-posed, in [1] a hybrid symbolic-numeric approach was proposed, that included the exact computation of a square polynomialsystem that had the original root with multiplicity one. In certifying singular roots, this exactsquare system was obtain from a deflation technique that added subdeterminants of the Jacobianmatrix to the system iteratively. However, the multiplicity structure is destroyed by this deflationtechnique, that is why it remained an open question how to certify the multiplicity structure ofsingular roots of exact polynomial systems.Our second motivation was to find a method that simultaneously refines the accuracy of a singularroot and the parameters describing the multiplicity structure at the root. In all previous numerical ∗ Research partly supported by DARPA YFA, NSF grant ACI-1460032, and Sloan Research Fellowship. † Research partly supported by NSF grant CCF-1217557.
Related work.
The treatment of singular roots is a critical issue for numerical analysis and there is a huge literatureon methods which transform the problem into a new one for which Newton-type methods convergequadratically to the root.Deflation techniques which add new equations in order to reduce the multiplicity have alreadybeen considered in [26], [25]: By triangulating the Jacobian matrix at the (approximate) root, newminors of the polynomial Jacobian matrix are added to the initial system in order to reduce themultiplicity of the singular solution.A similar approach is used in [10] and [8], where a maximal invertible block of the Jacobian matrixat the (approximate) root is computed and minors of the polynomial Jacobian matrix are addedto the initial system. In [8], an additional step is considered where the first derivatives of the inputpolynomials are added when the Jacobian matrix at the root vanishes.These constructions are repeated until a system with a simple root is obtained.In these methods, at each step, the number of added equations is ( n − r ) × ( m − r ), where n isnumber of variables, m is the number of equations and r is the rank of the Jacobian at the root.In [12], a triangular presentation of the ideal in a good position and derivations with respect to theleading variables are used to iteratively reduce the multiplicity. This process is applied for p-adiclifting with exact computation.In other approaches, new variables and new equations are introduced simultaneously. In [31],new variables are introduced to describe some perturbations of the initial equations and somedifferentials which vanish at the singular points. This approach is also used in [18], where it isshown that this iterated deflation process yields a system with a simple root.In [20], perturbation variables are also introduced in relation with the inverse system of the singularpoint to obtain directly a deflated system with a simple root. The perturbation is constructed froma monomial basis of the local algebra at the multiple root.In [13, 14], only variables for the differentials of the initial system are introduced. The analysisof this deflation is improved in [5], where it is shown that the number of steps is bounded by theorder of the inverse system. This type of deflation is also used in [17], for the special case wherethe Jacobian matrix at the multiple root has rank n − Contributions.
We propose a new deflation method for polynomial systems with isolated singular points, whichdoes introduce new parameters. At each step, a single differential of the system is considered basedon the analysis of the Jacobian at the singular point. A linear number of new equations is addedinstead of the quadratic increases of the previous deflations. The deflated system does not involvedany approximate coefficients and can therefore be used in certification methods as in [1].To approximate efficiently both the singular point and its inverse system, we propose a new defla-tion, which involves a small number of new variables compared to other approaches which rely onMacaulay matrices. It is based on a new characterization of the isolated singular point togetherwith its multiplicity structure. The deflated polynomial system exploits the nilpotent and commu-tation properties of the multiplication matrices in the local algebra of the singular point. We provethat it has a simple root which yields the root and the coefficients of the inverse system at thissingular point. Due to the upper triangular form of the multiplication matrices in a convenientbasis of the local algebra, the number of new parameters introduced in this deflation is less than n ( δ − δ where n is the number of variables and δ the multiplicity of the singular point. Theparameters involved in the deflated system are determined from the analysis of an approximation ofthe singular point. Nevertheless, the deflated system does not involve any approximate coefficientsand thus it can also be used in certification techniques as [1].In this paper we present two new constructions. The first one is a new deflation method for asystem of polynomials with an isolated singular root. The new construction uses a single lineardifferential form defined from the Jacobian matrix of the input, and defines the deflated systemby applying this differential form to the original system. We prove that the resulting deflatedsystem has strictly lower multiplicity and depth at the singular point than the original one. Theadvantage of this new deflation is that it does not introduce new variables, and the increase inthe number of equations is linear, instead of the quadratic increase of previous deflation methods.The second construction gives the coefficients of the so called inverse system or dual basis, whichdefines the multiplicity structure at the singular root. The novelty of our construction is thatwe show that the nilpotent and commutation properties of the multiplication matrices definesmoothly the singular points and its inverse system. We give a system of equations in the originalvariables plus a relatively small number of new variables, and prove that the roots of this newsystem correspond to the original multiple roots but now with multiplicity one, and they uniquelydetermine the multiplicity structure. The number of unknowns used to describe the multiplicitystructure is significantly smaller, compared to the direct computation of the dual bases from theso called Macaulay matrices. Both constructions are “exact” in the sense that approximations ofthe coordinates of the singular point are only used to detect numerically non-singular submatrices,3nd not in the rest of the construction. Thus these constructions would allow to treat all conjugateroots simultaneously, as well as to apply these constructions in the certification of the singularroots and the multiplicity structure of an exact rational polynomial system. Let f := ( f , . . . , f N ) ∈ K [ x ] N with x = ( x , . . . , x n ) for some K ⊂ C field. Let ξ = ( ξ , . . . , ξ n ) ∈ C n be an isolated multiple root of f . Let I = h f , . . . , f N i , m ξ be the maximal ideal at ξ and Q bethe primary component of I at ξ so that √ Q = m ξ .Consider the ring of power series C [[ ∂ ξ ]] := C [[ ∂ ,ξ , . . . , ∂ n,ξ ]] and we use the notation for β =( β , . . . , β n ) ∈ N n ∂ βξ := ∂ β ,ξ · · · ∂ β n n,ξ . We identify C [[ ∂ ξ ]] with the dual space C [ x ] ∗ by considering ∂ βξ as derivations and evaluations at ξ , defined by ∂ βξ ( p ) := ∂ β ( p ) | ξ := d | β | pdx β · · · dx β n n ( ξ ) for p ∈ C [ x ] . (1)Hereafter, the derivations “at x ” will be denoted ∂ β instead of ∂ β x . The derivation with respect tothe variable ∂ i is denoted d ∂ i ( i = 1 , . . . , n ). Note that1 β ! ∂ βξ (( x − ξ ) α ) = ( α = β , where we use the notation β ! = β ! ··· β n ! .For p ∈ C [ x ] and Λ ∈ C [[ ∂ ξ ]] = C [ x ] ∗ , let p · Λ : q Λ( p q ) . We check that p = ( x i − ξ i ) acts as a derivation on C [[ ∂ ξ ]]:( x i − ξ i ) · ∂ βξ = d ∂ i,ξ ( ∂ βξ )For an ideal I ⊂ C [ x ], let I ⊥ = { Λ ∈ C [[ ∂ ξ ]] | ∀ p ∈ I, Λ( p ) = 0 } . The vector space I ⊥ is naturallyidentified with the dual space of C [ x ] /I . We check that I ⊥ is a vector subspace of C [[ ∂ ξ ]], whichis stable by the derivations d ∂ i,ξ . Lemma 2.1. If Q is a m ξ -primary component of I , then Q ⊥ = I ⊥ ∩ C [ ∂ ξ ] . This lemma shows that to compute Q ⊥ , it suffices to compute all polynomials of C [ ∂ ξ ] which arein I ⊥ . Let us denote this set D = I ⊥ ∩ C [ ∂ ξ ]. It is a vector space stable under the derivations d ∂ i,ξ . Its dimension is the dimension of Q ⊥ or C [ x ] /Q , that is the multiplicity of ξ , denote it by δ ξ ( I ), or simply by δ if ξ and I is clear from the context.For an element Λ( ∂ ξ ) ∈ C [ ∂ ξ ] we define the order ord(Λ) to be the maximal | β | such that ∂ βξ appears in Λ( ∂ ξ ) with non-zero coefficient.For t ∈ N , let D t be the elements of D of order ≤ t . As D is of dimension d , there exists a smallest t ≥ D t +1 = D t . Let us call this smallest t , the nil-index of D and denote it by o ξ ( I ),or simply by o . As D is stable by the derivations d ∂ i,ξ , we easily check that for t ≥ o ξ ( I ), D t = D and that o ξ ( I ) is the maximal degree of the elements in D .4 Deflation using first differentials
To improve the numerical approximation of a root, one usually applies a Newton-type method toconverge quadratically from a nearby solution to the root of the system, provided it is simple. Inthe case of multiple roots, deflation techniques are employed to transform the system into anotherone which has an equivalent root with a smaller multiplicity or even with multiplicity one.We describe here a construction, using differentials of order one, which leads to a system with asimple root. This construction improves the constructions in [13, 5] since no new variables areadded. It also improves the constructions presented in [10] and the “kerneling” method of [8] byadding a smaller number of equations at each deflation step. Note that, in [8], there are smartpreprocessing and postprocessing steps which could be utilized in combination with our method.In the preprocessor, one adds directly partial derivatives of polynomials which are zero at the root.The postprocessor extracts a square subsystem of the completely deflated system for which theJacobian has full rank at the root.Consider the Jacobian matrix J f ( x ) = [ ∂ j f i ( x )] of the initial system f . By reordering properly therows and columns (i.e., polynomials and variables), it can be put in the form J f ( x ) := (cid:20) A ( x ) B ( x ) C ( x ) D ( x ) (cid:21) (2)where A ( x ) is an r × r matrix with r = rank J f ( ξ ) = rank A ( ξ ).Suppose that B ( x ) is an r × c matrix. The c columnsdet( A ( x )) (cid:20) − A − ( x ) B ( x )Id (cid:21) (for r = 0 this is the identity matrix) yield the c elementsΛ x = n X i =1 λ ,j ( x ) ∂ j , . . . , Λ x c = n X i =1 λ c,j ( x ) ∂ j . Their coefficients λ i,j ( x ) ∈ K [ x ] are polynomial in the variables x . Evaluated at x = ξ , theygenerate the kernel of J f ( ξ ) and form a basis of D . Definition 3.1.
The family D x = { Λ x , . . . , Λ x c } is the formal inverse system of order 1 at ξ . For i = { i , . . . , i k } ⊂ { , . . . , c } with | i | 6 = 0, the i - deflated system of order 1 of f is { f , Λ x i ( f ) , . . . , Λ x i k ( f ) } . By construction, for i = 1 , . . . , c ,Λ x i ( f ) = n X j =1 ∂ j ( f ) λ i,j ( x ) = det( A ( x )) J f ( x )[ λ i,j ( x )]has n − c zero entries. Thus, the number of non-trivial new equations added in the i -deflatedsystem is | i | · ( N − n + c ). The construction depends on the choice of the invertible block A ( ξ ) in J f ( ξ ). By a linear invertible transformation of the initial system and by computing a i -deflatedsystem, one obtains a deflated system constructed from any | i | linearly independent elements ofthe kernel of J f ( ξ ). Example 3.2.
Consider the multiplicity 2 root ξ = (0 ,
0) for the system f ( x ) = x + x and f ( x ) = x + x . Then, J f ( x ) = (cid:20) A ( x ) B ( x ) C ( x ) D ( x ) (cid:21) = (cid:20) x x x (cid:21) . − x T yields the elementΛ x = − x ∂ + ∂ . Since Λ x ( f ) = 0, the { } -deflated system of order 1 of f is (cid:8) x + x , x + x , − x x + 2 x (cid:9) which has a multiplicity 1 root at ξ .We use the following to analyze this deflation procedure. Lemma 3.3 (Leibniz rule) . For a, b ∈ K [ x ] , ∂ α ( a b ) = X β ∈ N n β ! ∂ β ( a ) d β ∂ ( ∂ α )( b ) . Proposition 3.4.
Let r be the rank of J f ( ξ ) . Assume that r < n . Let i ⊂ { , . . . , n } with < | i | ≤ n − r and f (1) be the i - deflated system of order of f . Then, δ ξ ( f (1) ) ≥ and o ξ ( f (1) ) < o ξ ( f ) , which also implies that δ ξ ( f (1) ) < δ ξ ( f ) .Proof. By construction, for i ∈ i , the polynomials Λ x i ( f ) vanish at ξ , so that δ ξ ( f (1) ) ≥
1. Byhypothesis, the Jacobian of f is not injective yielding o ξ ( f ) >
0. Let D (1) be the inverse systemof f (1) at ξ . Since ( f (1) ) ⊃ ( f ), we have D (1) ⊂ D . In particular, for any non-zero elementΛ ∈ D (1) ⊂ K [ ∂ ξ ] and i ∈ i , Λ( f ) = 0 and Λ(Λ x i ( f )) = 0.Using Leibniz rule, for any p ∈ K [ x ], we have Λ(Λ x i ( p )) = Λ n X j =1 λ i,j ( x ) ∂ j ( p ) = X β ∈ N n n X j =1 β ! ∂ β ξ ( λ i,j ( x )) d β ∂ξ (Λ) ∂ j,ξ ( p )= X β ∈ N n n X j =1 β ! ∂ β ξ ( λ i,j ( x )) ∂ j,ξ d β ∂ξ (Λ)( p )= X β ∈ N n ∆ i,β d β ∂ξ (Λ)( p ) where ∆ i, β = n X j =1 λ i,j, β ∂ j,ξ ∈ K [ ∂ ξ ] and λ i,j, β = 1 β ! ∂ β ξ ( λ i,j ( x )) ∈ K . The term ∆ i, is P nj =1 λ i,j ( ξ ) ∂ j,ξ which has degree 1 in ∂ ξ since [ λ i,j ( ξ )] is a non-zero element ofker J f ( ξ ). For simplicity, let φ i (Λ) := P β ∈ N n ∆ i, β d β ∂ (Λ).For any Λ ∈ C [ ∂ ξ ], we have d ∂ j,ξ ( φ i (Λ)) = X β ∈ N n λ i,j,β d β ∂ (Λ) + ∆ i,β d β ∂ ( d ∂ j,ξ (Λ))= X β ∈ N n λ i,j,β d β ∂ (Λ) + φ i ( d ∂ j,ξ (Λ)) . Moreover, if Λ ∈ D (1) , then by definition φ i (Λ)( f ) = 0. Since D and D (1) are both stable byderivation, it follows that ∀ Λ ∈ D (1) , d ∂ j,ξ ( φ i (Λ)) ∈ D (1) + φ i ( D (1) ). Since D (1) ⊂ D , we know D + φ i ( D (1) ) is stable by derivation. For any element Λ of D + φ i ( D (1) ), Λ( f ) = 0. We deducethat D + φ i ( D (1) ) = D . Consequently, the order of the elements in φ i ( D (1) ) is at most o ξ ( f ). Thestatement follows since φ i increases the order by 1, therefore o ξ ( f (1) ) < o ξ ( f ).6e consider now a sequence of deflations of the system f . Let f (1) be the i -deflated system of f . Weconstruct inductively f ( k +1) as the i k +1 -deflated system of f ( k ) for some choices of i j ⊂ { , . . . , n } . Proposition 3.5.
There exists k ≤ o ξ ( f ) such that ξ is a simple root of f ( k ) .Proof. By Proposition 3.4, δ ξ ( f ( k ) ) ≥ o ξ ( f ( k ) ) is strictly decreasing with k until it reachesthe value 0. Therefore, there exists k ≤ o ξ ( I ) such that o ξ ( f ( k ) ) = 0 and δ ξ ( f ( k ) ) ≥
1. This impliesthat ξ is a simple root of f ( k ) .To minimize the number of equations added at each deflation step, we take | i | = 1. Then, thenumber of non-trivial new equations added at each step is at most N − n + c .We described this approach using first order differentials arising from the Jacobian, but this canbe easily extended to use higher order differentials. Before describing our results, we start this section by recalling the definition of orthogonal primal-dual pairs of bases for the space C [ x ] /Q and its dual. The following is a definition/lemma: Lemma 4.1 (Orthogonal primal-dual basis pair) . Let f , ξ , Q , D , δ = δ ξ ( f ) and o = o ξ ( f ) beas above. Then there exists a primal-dual basis pair of the local ring C [ x ] /Q with the followingproperties: • The primal basis of the local ring C [ x ] /Q has the form B := { ( x − ξ ) α , ( x − ξ ) α , . . . , ( x − ξ ) α δ } . (3) We can assume that α = 0 and that the monomials in B are connected to 1 (c.f. [24]).Define the set of exponents in B E := { α , . . . , α δ } . • There is a unique dual basis Λ ⊂ D orthogonal to B , i.e. the elements of Λ are given in thefollowing form: Λ = ∂ α ξ = 1 ξ Λ = 1 α ! ∂ α ξ + X | β |≤ oβ E ν α ,β ∂ βξ ... (4)Λ δ − = 1 α δ ! ∂ α δ ξ + X | β |≤ oβ E ν α δ ,β ∂ βξ , • We have ) ≤ · · · ≤ ord(Λ δ − ) , and for all ≤ t ≤ o we have D t = span { Λ j : ord(Λ j ) ≤ t } , where D t denotes the elements of D of order ≤ t , as above.Proof. Let ≻ be the graded reverse lexicographic ordering in C [ ∂ ] such that ∂ ≺ ∂ ≺ · · · ≺ ∂ n . Weconsider the initial In( D ) = { In(Λ) | Λ ∈ D } of D for the monomial ordering ≻ . It is a finite set ofincreasing monomials D := { ∂ α , ∂ α , . . . , ∂ α δ − } , which are the leading monomials of the elementsof Λ = { Λ , Λ , . . . , Λ δ − } ⊂ D . As 1 ∈ D and is the lowest monomial ≻ , we have Λ = 1. As ≻
7s refining the total degree in C [ ∂ ], we have ord(Λ i ) = | α i | and 0 = ord(Λ ) ≤ · · · ≤ ord(Λ δ − ).Moreover, every element in D t reduces to 0 by the elements in Λ . As only the elements Λ i oforder ≤ t are involved in this reduction, we deduce that D t is spanned by the elements Λ i withord(Λ i ) ≤ t .Let E = { α , . . . , α δ − } . The elements Λ i are of the formΛ i = 1 α i ! ∂ α i ξ + X | β |≺| α i | ν α i ,β ∂ βξ . By auto-reduction of the elements Λ i , we can even suppose that β E in the summation above,so that they are of the form (4).Let B ( ξ ) = { ( x − ξ ) α , . . . , ( x − ξ ) α δ − } ⊂ C [ x ]. As (Λ i (( x − ξ ) α j )) ≤ i,j ≤ δ − is the identity matrix,we deduce that B is a basis of C [ x ] /Q , which is dual to Λ .As D is stable by derivation, the leading term of dd∂ i (Λ j ) is in D . If dd∂ i ( ∂ α j ξ ) is not zero, then itis the leading term of dd∂ i (Λ j ), since the monomial ordering is compatible with the multiplicationby a variable. This shows that D is stable by division by the variable ∂ i and that B is connectedto 1. This ends the proof of the lemma.Such a basis of D can be obtained from any other basis of D by transforming first the coefficientmatrix of the given dual basis into row echelon form and then reducing the elements above thepivot coefficients. The integration method described in [20] computes a primal-dual pair, such thatthe coefficient matrix has a block row-echelon form, each block being associated to an order. Thecomputation of a basis as in Lemma 4.1 can be then performed order by order. Example 4.2.
Let f = x − x + x , f = x − x + x , which has a multiplicity 3 root at ξ = (0 , B = { , x , x } , ˜ Λ = (cid:26) , ∂ + ∂ , ∂ + 12 ∂ + ∂ ∂ + 12 ∂ (cid:27) . This primal dual pair does not form an orthogonal pair, since ( ∂ + ∂ )( x ) = 0. However, usinglet say the degree lexicographic ordering such that x > x , we easily deduce the primal-dual pairof Lemma 4.1: B = (cid:8) , x , x (cid:9) , Λ = ˜ Λ . Throughout this section we assume that we are given a fixed primal basis B for C [ x ] /Q suchthat a dual basis Λ of D satisfying the properties of Lemma 4.1 exists. Note that such a primalbasis B can be computed numerically from an approximation of ξ and using a modification of theintegration method of [20].Given the primal basis B , a dual basis can be computed by Macaulay’s dialytic method which canbe used to deflate the root ξ as in [14]. This method would introduce n + ( δ − (cid:0)(cid:0) n + on (cid:1) − δ (cid:1) newvariables, which is not polynomial in o . Below, we give a construction of a polynomial system thatonly depends on at most n + nδ ( δ − / multiplication matrices that we define next. Let M i : C [ x ] /Q → C [ x ] /Qp ( x i − ξ i ) p be the multiplication operator by x i − ξ i in C [ x ] /Q . Its transpose operator is M ti : D → D Λ Λ ◦ M i = ( x i − ξ i ) · Λ = dd∂ i,ξ (Λ) = d ∂ i,ξ (Λ)8here D = Q ⊥ ⊂ C [ ∂ ]. The matrix of M i in the basis B of C [ x ] /Q is denoted M i .As B is a basis of C [ x ] /Q , we can identify the elements of C [ x ] /Q with the elements of the vectorspace span C ( B ). We define the normal form N ( p ) of a polynomial p in C [ x ] as the unique element b of span C ( B ) such that p − b ∈ Q . Hereafter, we are going to identify the elements of C [ x ] /Q withtheir normal form in span C ( B ).For any polynomial p ( x , . . . , x n ) ∈ C [ x ], let p ( M ) be the operator of C [ x ] /Q obtained by replacing x i − ξ i by M i . Lemma 4.3.
For any p ∈ C [ x ] , the normal form of p is N ( p ) = p ( M )(1) and we have p ( M )(1) = Λ ( p ) 1 + Λ ( p ) ( x − ξ ) α + · · · + Λ δ − ( p ) ( x − ξ ) α δ − . This shows that the coefficient vector [ p ] of N ( p ) in the basis B of is [ p ] = (Λ i ( p )) ≤ i ≤ δ − .The following lemma is also well known, but we include it here with proof: Lemma 4.4.
Let B as in (3) and denote the exponents in B by E := { α , . . . , α δ } as above. Let E + := n [ i =1 ( E + e i ) with E + e i = { ( γ , . . . , γ i + 1 , . . . , γ n ) : γ ∈ E } and we denote ∂ ( E ) = E + \ E . The values ofthe coefficients ν α,β for ( α, β ) ∈ E × ∂ ( E ) appearing in the dual basis (4) uniquely determine thesystem of pairwise commuting multiplication matrices M i , namely, for i = 1 , . . . , n M ti = 0 ν α , e i ν α , e i · · · ν α d − , e i ν α ,α + e i · · · ν α d − ,α + e i ... ... ... · · · ν α d − ,α d − + e i · · · Moreover, ν α i ,α k + e j = ( if α i = α k + e j if α k + e j ∈ E, α i = α k + e j . Proof. As M ti acts as a derivation on D and D is closed under derivation, so the third propertyin Lemma 4.1 implies that the matrix of M ti in the basis of Λ of D has an upper triangular formwith zero (blocks) on the diagonal.For an element Λ j of order k , its image by M ti is M ti (Λ j ) = ( x i − ξ i ) · Λ j = X | α l |
Note that we can reduce the number of free parameters in the parametric multipli-cation matrices by exploiting the commutation rules of the multiplication matrices correspondingto a given primal basis B . For example, consider the breadth one case, where we can assume that E = { , e , e , . . . , ( δ − e } . In this case the only free parameters appear in the first columns of M ( µ ) , . . . , M n ( µ ), the other columns are shifts of these. Thus, it is enough to introduce ( n − δ − n − δ −
1) free parameters.10 efinition 4.7 (Parametric normal form) . Let K ⊂ C be a field. We define N z ,µ : K [ x ] → K [ z , µ ] δ p
7→ N z ,µ ( p ) := X γ ∈ N n γ ! ∂ γ z ( p ) M ( µ ) γ [1] . where [1] = [1 , , . . . ,
0] is the coefficient vector of 1 in the basis B . This sum is finite since for | γ | ≥ δ , M ( µ ) γ = 0, so the entries of N z ,µ ( p ) are polynomials in µ and z .Notice that this notation is not ambiguous, assuming that the matrices M i ( µ ) ( i = 1 , . . . , n ) arecommuting. The specialization at ( x , µ ) = ( ξ, ν ) is the vector N ξ,ν ( p ) = [Λ ( p ) , . . . , Λ δ − ( p )] t ∈ C δ . We can now characterize the multiplicity structure by polynomial equations.
Theorem 4.8.
Let K ⊂ C be any field, f ∈ K [ x ] N and let ξ ∈ C n be an isolated solution of f . Let M i ( µ ) for i = 1 , . . . n be the parametric multiplication matrices as in (8) and N ξ,µ be the parametricnormal form as in Defn. 4.7 at z = ξ . Then the ideal J ξ of C [ µ ] generated by the polynomialsystem ( N ξ,µ ( f k ) for k = 1 , . . . , N, M i ( µ ) · M j ( µ ) − M i ( µ ) · M i ( µ ) for i, j = 1 , . . . , n (9) is the maximal ideal m ν = ( µ α,β − ν α,β , ( α, β ) ∈ E × ∂ ( E )) where ν α,β are the coefficients of the dual basis defined in (4).Proof. As before, the system (9) has a solution µ α,β = ν α,β for ( α, β ) ∈ E × ∂ ( E ). Thus J ξ ⊂ m ν .Conversely, let C = C [ µ ] /J ξ and consider the mapΦ : C [ x ] → C δ , p
7→ N ξ,µ ( p ) . Let K be its kernel. Since the matrices M i ( µ ) are commuting modulo J ξ , we can see that K is anideal. As f k ∈ K , we have I := ( f k ) ⊂ K .Next we show that Q ⊂ K . By construction, for any α ∈ N n we have modulo J ξ N ξ,µ (( x − ξ ) α ) = X γ ∈ N n γ ! ∂ γξ (( x − ξ ) α ) M ( µ ) γ [1] = M ( µ ) α [1] . Using the previous relation, we check that ∀ p, q ∈ C [ x ],Φ( pq ) = p ( ξ + M ( µ ))Φ( q ) (10)where p ( ξ + M ( µ )) is obtained by replacing x i − ξ i by M i ( µ ). Let q ∈ Q . As Q is the m ξ -primarycomponent of I , there exists p ∈ C [ x ] such that p ( ξ ) = 0 and p q ∈ I . By (10), we haveΦ( p q ) = p ( ξ + M ( µ ))Φ( q ) = 0 . Since p ( ξ ) = 0 and p ( ξ + M ( µ )) = p ( ξ ) Id + N with N lower triangular and nilpotent, p ( ξ + M ( µ ))is invertible. We deduce that Φ( q ) = p ( ξ + M ( µ )) − Φ( pq ) = 0 and q ∈ K .Let us show now that Φ is surjective and more precisely, that Φ(( x − ξ ) α k ) = e k (abusing thenotation as here e k has length δ not n ). Since B is connected to 1, either α k = 0 or there exists11 j ∈ E such that α k = α j + e i for some i ∈ { , . . . , n } . Thus the j th column of M i ( µ ) is e k by(7). As { M i ( µ ) : i = 1 , . . . , n } are pairwise commuting, we have M ( µ ) α k = M j ( µ ) M ( µ ) α j , and if weassume by induction on | α j | that the first column of M ( µ ) α j is e j , we obtain M ( µ ) α k [1] = e k . Thus,for k = 1 , . . . , δ , Φ(( x − ξ ) α k ) = e k .We can now prove that m ν ⊂ J ξ . As M i ( ν ) is the multiplication by ( x i − ξ i ) in C [ x ] /Q , for any b ∈ B and i = 1 , . . . , n , we have ( x i − ξ i ) b = M i ( ν )( b ) + q with q ∈ Q ⊂ K . We deduce that for k = 0 , . . . , δ − Φ(( x i − ξ i ) ( x − ξ ) αk ) = M i ( µ )Φ(( x − ξ ) αk ) = M i ( µ )( e k ) = M i ( ν )( e k ). This shows that µ α,β − ν α,β ∈ J ξ for ( α, β ) ∈ E × ∂ ( E ) and that m ν = J ξ .In the proof of the next theorem we need to consider cases when the multiplication matrices donot commute. We introduce the following definition: Definition 4.9.
Let K ⊂ C be any field. Let C be the ideal of K [ z , µ ] generated by entries of thecommutation relations: M i ( µ ) · M j ( µ ) − M j ( µ ) · M i ( µ ) = 0, i, j = 1 , . . . , n . We call C the commutatorideal . Lemma 4.10.
For any field K ⊂ C , p ∈ K [ x ] , and i = 1 , . . . , n , we have N z ,µ ( x i p ) = x i N z ,µ ( p ) + M i ( µ ) N z ,µ ( p ) + O i,µ ( p ) . (11) where O i,µ : K [ x ] → K [ z , µ ] δ is linear with image in the commutator ideal C .Proof. N z ,µ ( x i p ) = P γ γ ! ∂ γ z ( x i p ) M ( µ ) γ [1]= x i X γ γ ! ∂ γ z ( p ) M ( µ ) γ [1] + X γ γ ! γ i ∂ γ − e i z ( p ) M ( µ ) γ [1]= x i X γ γ ! ∂ γ z ( p ) M ( µ ) γ [1] + X γ γ ! ∂ γ z ( p ) M ( µ ) γ + e i [1]= x i N z ,µ ( p ) + M i ( µ ) X γ γ ! ∂ γ z ( p ) M ( µ ) γ [1] ! + X γ γ ! ∂ γ z ( p ) O i,γ ( µ )[1]where O i,γ = M i ( µ ) M ( µ ) γ − M ( µ ) γ + e i is a δ × δ matrix with coefficients in C . Therefore, O i,µ : p P γ γ ! ∂ γ z ( p ) O i,γ ( µ )[1] is a linear functional of p with coefficients in C .The next theorem proves that the system defined as in (9) for general z has ( ξ, ν ) as a simple root. Theorem 4.11.
Let f ∈ K [ x ] N and ξ ∈ C n be as above. Let M i ( µ ) for i = 1 , . . . n be the parametricmultiplication matrices defined in (8) and N x ,µ be the parametric normal form as in Defn. 4.7.Then ( z , µ ) = ( ξ, ν ) is an isolated root with multiplicity one of the polynomial system in K [ z , µ ] : ( N z ,µ ( f k ) = 0 for k = 1 , . . . , N, M i ( µ ) · M j ( µ ) − M j ( µ ) · M i ( µ ) = 0 for i, j = 1 , . . . , n. (12) Proof.
For simplicity, let us denote the (non-zero) polynomials appearing in (12) by P , . . . , P M ∈ K [ z , µ ] , where M ≤ N δ + n ( n − δ − δ − /
4. To prove the theorem, it is sufficient to prove thatthe columns of the Jacobian matrix of the system [ P , . . . , P M ] at ( z , µ ) = ( ξ, ν ) are linearlyindependent. The columns of this Jacobian matrix correspond to the elements in C [ z , µ ] ∗ ∂ ,ξ , . . . , ∂ n,ξ , and ∂ µ α,β for ( α, β ) ∈ E × ∂ ( E ) , ∂ i,ξ defined in (1) for z replacing x , and ∂ µ α,β is defined by ∂ µ α,β ( q ) = dqdµ α,β (cid:12)(cid:12) ( z ,µ )=( ξ,ν ) for q ∈ C [ z , µ ] . Suppose there exist a , . . . , a n , and a α,β ∈ C for ( α, β ) ∈ E × ∂ ( E ) not all zero such that∆ := a ∂ ,ξ + · · · + a n ∂ n,ξ + X α,β a α,β ∂ µ α,β ∈ C [ z , µ ] ∗ vanishes on all polynomials P , . . . , P M in (12). In particular, for an element P i ( µ ) correspondingto the commutation relations and any polynomial Q ∈ C [ x , µ ], using the product rule for the lineardifferential operator ∆ we get∆( P i Q ) = ∆( P i ) Q ( ξ, ν ) + P i ( ν )∆( Q ) = 0since ∆( P i ) = 0 and P i ( ν ) = 0. By the linearity of ∆, for any polynomial C in the commutatorideal C , we have ∆( C ) = 0.Furthermore, since ∆( N z ,µ ( f k )) = 0 and N ξ,ν ( f k ) = [Λ ( f k ) , . . . , Λ δ − ( f k )] t , we get that ( a ∂ ,ξ + · · · + a n ∂ n,ξ ) · Λ δ − ( f k ) + X | γ |≤| α δ − | p γ ( ν ) ∂ γξ ( f k ) = 0 (13) where p γ ∈ C [ µ ] are some polynomials in the variables µ that do not depend on f k . If a , . . . , a n are not all zero, we have an element ˜Λ of C [ ∂ ξ ] of order strictly greater than ord(Λ δ − ) = o thatvanishes on f , . . . , f N .Let us prove that this higher order differential also vanishes on all multiples of f k for k = 1 , . . . , N .Let p ∈ C [ x ] such that N ξ,ν ( p ) = 0, ∆( N z ,µ ( p )) = 0. By (11), we have N ξ,ν (( x i − ξ i ) p )= ( x i − ξ i ) N ξ,ν ( p ) + M i ( ν ) N ξ,ν ( p ) + O i,ν ( p ) = 0and ∆( N z ,µ (( x i − ξ i ) p ))= ∆(( x i − ξ i ) N z ,µ ( p )) + ∆( M i ( µ ) N z ,µ ( p )) + ∆( O µ ( p ))= ∆( x i − ξ i ) N ξ,ν ( p ) + ( ξ i − ξ i )∆( N z ,µ ( p ))+ ∆( M i ( µ )) N ξ,µ ( p ) + M i ( ν )∆( N z ,µ ( p ))+ ∆( O i,µ ( p ))= 0 . As N ξ,ν ( f k ) = 0, ∆( N z ,µ ( f k )) = 0, i = 1 , . . . , N , we deduce by induction on the degree of themultipliers and by linearity that for any element f in the ideal I generated by f , . . . , f N , we have N ξ,ν ( f ) = 0 and ∆( N z ,µ ( f )) = 0 , which yields ˜Λ ∈ I ⊥ . Thus we have ˜Λ ∈ I ⊥ ∩ C [ ∂ ξ ] = Q ⊥ (by Lemma 2.1). As there is no elementof degree strictly bigger than o in Q ⊥ , this implies that a = · · · = a n = 0 . Then, by specialization at x = ξ , ∆ yields an element of the kernel of the Jacobian matrix of thesystem (9). By Theorem 4.8, this Jacobian has a zero-kernel, since it defines the simple point ν .We deduce that ∆ = 0 and ( ξ, ν ) is an isolated and simple root of the system (12).13he following corollary applies the polynomial system defined in (12) to refine the precision of anapproximate multiple root together with the coefficients of its Macaulay dual basis. The advantageof using this, as opposed to using the Macaulay multiplicity matrix, is that the number of variablesis much smaller, as was noted above. Corollary 4.12.
Let f ∈ K [ x ] N and ξ ∈ C n be as above, and let Λ ( ν ) , . . . , Λ δ − ( ν ) be its dualbasis as in (4). Let E ⊂ N n be as above. Assume that we are given approximates for the singularroots and its inverse system as in (4) ˜ ξ ∼ = ξ and ˜ ν α i ,β ∼ = ν α i ,β ∀ α i ∈ E, β E, | β | ≤ o. Consider the overdetermined system in K [ z , µ ] from (12). Then a random square subsystem of (12)will have a simple root at z = ξ , µ = ν with high probability. Thus, we can apply Newton’s methodfor this square subsystem to refine ˜ ξ and ˜ ν α i ,β for ( α i , β ) ∈ E × ∂ ( E ) . For ˜ ν α i ,γ with γ E + wecan use (6) for the update. Example 4.13.
Reconsider the setup from Ex. 3.2 with primal basis { , x } and E = { (0 , , (0 , } .We obtain M ( µ ) = (cid:20) µ (cid:21) and M ( µ ) = (cid:20) (cid:21) . The resulting deflated system in (12) is F ( z , z , µ ) = z + z µ + 2 z z + z µz + 2 z which has a nonsingular root at ( z , z , µ ) = (0 , ,
0) corresponding to the origin with multiplicitystructure { , ∂ } . Computations for the following examples, as well as several other systems, along with
Matlab code can be found at . We first consider a modification of [17, Example 3.1]. For any n ≥
2, the following system has n polynomials, each of degree at most 3, in n variables: x + x − x , x + x − x , . . . , x n − + x n − − x n , x n . The origin is a multiplicity δ := 2 n root having breadth 2 (i.e., the corank of Jacobian at the originis 2).We apply our parametric normal form method described in §
4. Similarly as in Remark 4.6, wecan reduce the number of free parameters to be at most ( n − δ −
1) using the structure of theprimal basis B = { x a x b : a < n − , b < } .The following table shows the multiplicity, number of variables and polynomials in the deflatedsystem, and the time (in seconds) it took to compute this system (on a iMac, 3.4 GHz Intel Core i7processor, 8GB 1600Mhz DDR3 memory). Note that when comparing our method to an approachusing the null spaces of Macaulay multiplicity matrices (see for example [6, 14]), we found that for n ≥ n − , so the size of the Macaulay multiplicity matrixis n · (cid:0) n − + n − n − (cid:1) × (cid:0) n − + nn (cid:1) . New approach Null space n mult vars poly time vars poly time2 4 5 9 1 .
476 8 17 2 . .
596 192 241 2084 16 49 100 19 .
698 7189 19804 > . N/A N/A N/A . N/A N/A N/A
We consider the Caprasse system [3, 29]: f ( x , x , x , x ) = x x − x x x − x x x − x x − x +10 x − x x + 10 x x − ,x x − x x x − x x x − x x − x x +10 x x − x + 10 x − ,x x + 2 x x x − x − x ,x x + 2 x x x − x − x at the multiplicity 4 root ξ = (2 , −√− , , √− × × § | i | = 1 step creates a deflated system consisting of 6 polynomials in4 variables. In fact, since the null space of the Jacobian at the root is 2 dimensional, adding twopolynomials is necessary and sufficient.Next, we consider the computation of both the point and multiplicity structure. Using an intrinsicnull space approach via a second order Macaulay matrix, the resulting system consists of 64polynomials in 37 variables. In comparison, using the primal basis { , x , x , x x } , the approachof § In our last set of examples, we consider simply deflating a root of the last three systems from [6, §
7] and a system from [12, § { x − x x x , x − x x x , x − x x x , x − x x x } at (0 , , ,
0) with δ = 131 and o = 10;2: { x , x y + y , z + z − x − x } at (0 , , −
1) with δ = 16 and o = 7;3: { x + 33 y − √ x + 4 xy + 4 y + 2) + √ x + 6 x y + 12 xy + 8 y , x − y − √ x − x y + 6 xy − y + 3 √ xy − x − y − } at Z ≈ (1 . , . δ = 5 and o = 4;4: { x + 2 x + 2 x + 2 x + x − , ( x + x − x − − x , (2 x + 5 x + 10 x + 5 x + 5) − x } at (0 , , −
1) with δ = 18 and o = 7. We compare using the following four methods: (A) intrinsic slicing version of [6, 13]; (B) isosingulardeflation [10] via a maximal rank submatrix; (C) “kerneling” method in [8]; (D) approach of § | i | = 1 step. We performed these methods without the use of preprocessing and15ostprocessing as mentioned in § References [1] T.A. Akoglu, J.D. Hauenstein, and A. Szanto. Certifying solutions to overdetermined andsingular polynomial systems over Q . Preprint, http://arxiv.org/abs/1408.2721 , 2014.[2] L. Alberti, B. Mourrain, and J. Wintz. Topology and arrangement computation of semi-algebraic planar curves. Computer Aided Geometric Design , 25(8):631–651, 2008.[3] H. Caprasse, J. Demaret, E. Schr¨ufer: Can EXCALC be used to Investigate High-DimensionalCosmological Models with Non-Linear Lagrangians? ISSAC 1988: 116-124.[4] R.M. Corless, P.M. Gianni, and B.M. Trager. A reordered Schur factorization method forzero-dimensional polynomial systems with multiple roots. In
ISSAC 1997 , ACM, New York,1997, pp. 133–140.[5] B.H. Dayton, T.-Y. Li, and Z. Zeng. Multiple zeros of nonlinear systems.
Math. Comput. ,80(276):2143–2168, 2011.[6] B.H. Dayton and Z. Zeng. Computing the multiplicity structure in solving polynomial systems.In
ISSAC’05 , ACM, New York, 2005, pp. 116–123.[7] M. Giusti, G. Lecerf, B. Salvy, and J.-C. Yakoubsohn. On location and approximation of clus-ters of zeros: Case of embedding dimension one.
Foundations of Computational Mathematics ,7:1–58, 2007.[8] M. Giusti and J.-C. Yakoubsohn. Multiplicity hunting and approximating multiple roots ofpolynomial systems.
Contemp. Math. , 604:105–128, 2013.[9] S. Graillat and P. Tr´ebuchet. A new algorithm for computing certified numerical approxima-tions of the roots of a zero-dimensional system. In
ISSAC2009 , pages 167–173, 2009.[10] J.D. Hauenstein and C. Wampler. Isosingular sets and deflation.
Found. Comput. Math. ,13(3):371–403, 2013.[11] W. Hao, A.J. Sommese, and Z. Zeng. Algorithm 931: an algorithm and software for computingmultiplicity structures at zeros of nonlinear systems.
ACM Trans. Math. Software , 40(1):Art.5, 2013.[12] G. Lecerf. Quadratic newton iterarion for systems with multiplicity.
Found. Comput. Math. ,(2):247–293, 2002. 1613] A. Leykin, J. Verschelde, and A. Zhao. Newton’s method with deflation for isolated singular-ities of polynomial systems.
Theor. Comput. Sci. , 359(1-3):111–122, 2006.[14] ——. Higher-order deflation for polynomial systems with isolated singular solutions. InA. Dickenstein, F.-O. Schreyer, and A. Sommese, editors,
Algorithms in Algebraic Geometry ,volume 146 of
The IMA Volumes in Mathematics and its Applications , Springer, New York,2008, pp. 79–97.[15] N. Li and L. Zhi. Computing isolated singular solutions of polynomial systems: case of breadthone.
SIAM J. Numer. Anal. , 50(1):354–372, 2012.[16] ——. Computing the multiplicity structure of an isolated singular solution: case of breadthone.
J. Symbolic Comput. , 47(6):700–710, 2012.[17] ——. Verified error bounds for isolated singular solutions of polynomial systems: Case ofbreadth one.
Theor. Comput. Sci. , 479:163–173, 2013.[18] ——. Verified error bounds for isolated singular solutions of polynomial systems.
SIAM J.Numer. Anal. , 52(4):1623–1640, 2014.[19] F.S. Macaulay.
The Algebraic Theory of Modular Systems . Cambridge Univ. Press, 1916.[20] A. Mantzaflaris and B. Mourrain. Deflation and certified isolation of singular zeros of poly-nomial systems. In
ISSAC 2011 , ACM, New York, 2011, pp. 249–256.[21] M.G. Marinari, T. Mora, and H. M¨oller. Gr¨obner duality and multiplicities in polynomialsystem solving. In
ISSAC ’95 , ACM, New York, 1995, pp. 167–179.[22] H.M. M¨oller and H.J. Stetter. Multivariate polynomial equations with multiple zeros solvedby matrix eigenproblems.
Numer. Math. , 70(3):311–329, 1995.[23] B. Mourrain. Isolated points, duality and residues.
J. Pure Appl. Alg. , 117-118:469–493, 1997.[24] ——. A new criterion for normal form algorithms.
LNCS , 1719:430–443, 1999.[25] T. Ojika. Modified deflation algorithm for the solution of singular problems. I. a system ofnonlinear algebraic equations.
J. Math. Anal. Appl. , 123(1):199–221, 1987.[26] T. Ojika, S. Watanabe, and T. Mitsui. Deflation algorithm for the multiple roots of a systemof nonlinear equations.
J. Math. Anal. Appl. , 96(2):463–479, 1983.[27] S. Pope and A. Szanto. Nearest multivariate system with given root multiplicities.
J. Symb.Comput. , 44(6):606–625, 2009.[28] H.J. Stetter. Analysis of zero clusters in multivariate polynomial systems. In
ISSAC ’96 ,ACM, New York, 1996, pp. 127–136.[29] C. Traverso. The posso test suite. with complement by D. Bini & B. Mourrain, 1993. .[30] X. Wu and L. Zhi. Determining singular solutions of polynomial systems via symbolic-numericreduction to geometric involutive form.
J. Symb. Comput. , 27:104–122, 2008.[31] N. Yamamoto. Regularization of solutions of nonlinear equations with singular jacobian ma-trices.
J. Inform. Proc. , 7(1):16–21, 1984.[32] Z. Zeng. Computing multiple roots of inexact polynomials.