Contents lists available atScienceDirect
Linear
Algebra
and
its
Applications
www.elsevier.com/locate/laa
On
the
null
spaces
of
the
Macaulay
matrix
Kim Batselier∗,1, Philippe Dreesen2, Bart De Moor3DepartmentofElectricalEngineeringESAT-STADIUSCenterforDynamical Systems,SignalProcessingandDataAnalytics,KULeuven/IBBTFutureHealth Department,3001Leuven,Belgium
a r t i c l e i n f o a b s t r a c t
Articlehistory:
Received25October2013 Accepted24July2014
Availableonline14August2014 SubmittedbyE.Zerz MSC: 54B05 12E05 65F99 Keywords: Macaulaymatrix Multivariatepolynomials Rowspace Nullspace Syzygies Roots
In this article both the left and right null space of the Macaulaymatrixaredescribed.Theleftnullspaceisshown tobelinkedwiththeoccurrenceofsyzygiesinitsrowspace. Itis also demonstrated how the dimension of the left null space is described by a piecewise function of polynomials. Wepresenttwoalgorithmsthatdeterminethesepolynomials.
Furthermore weshow how the finiteness of the number of
basissyzygiesresultsinthenotionofthedegreeofregularity. Thisconceptplaysacrucialroleindescribingabasisforthe rightnullspaceoftheMacaulaymatrixintermsofdifferential functionals.WedefineacanonicalnullspacefortheMacaulay matrixintermsoftheprojectiverootsofapolynomialsystem and extend the multiplication property of this canonical basis to the projective case. This results in an algorithm todeterminetheuppertriangularcommutingmultiplication matrices.Finally,wediscusshowStetter’seigenvalueproblem to determinethe rootsof a multivariate polynomial system
* Correspondingauthor.
E-mailaddress:kim.batselier@esat.kuleuven.be(K. Batselier). 1
KimBatselierisaresearchassistantattheKULeuven,Belgium. 2
PhilippeDreesenis supportedbytheInstitute forthePromotionof Innovationthrough Scienceand TechnologyinFlanders(IWT-Vlaanderen).
3 BartDeMoorisafullprofessorattheKULeuven,Belgium.ResearchsupportedbyResearchCouncil KUL:GOA/10/09MaNet,CoEPFV/10/002(OPTEC);PhD/Postdocgrants;FlemishGovernment:IOF: IOF/KP/SCORES4CHEM,iMindsMedicalInformationTechnologiesSBO2014,BelgianFederalScience PolicyOffice:IUAPP7/19(DYSCO,Dynamicalsystems,controlandoptimization,2012–2017).
http://dx.doi.org/10.1016/j.laa.2014.07.035
canbeextendedtothecasewhereamultivariatepolynomial systemhasbothaffinerootsandrootsatinfinity.
© 2014ElsevierInc.All rights reserved.
1. Introduction
Manyproblems incomputerscience,physicsandengineeringrequireamathematical modelling step and multivariate polynomials are a naturalmodelling tool [1]. This re-sults inproblemsinwhichoneneedstocomputetherootsofamultivariatepolynomial system,dividemultivariatepolynomials,eliminatevariables, computegreatestcommon divisors, etc. The area of mathematics in which multivariate polynomials are studied is algebraic geometry and hasa rich historyspanning manycenturies [2]. Most meth-ods to solve these problems are symbolical and involve the computation of aGröbner basis inarbitrary high precision [3,4]. Estimating themodel parametersof polynomial blackboxmodelsfrommeasureddataleadstopolynomialsystemswherethecoefficients are directly related to themeasurements [5–7]. Hence, these coefficients are subjectto noise and known only with a limited accuracy. In this case, reporting high precision results obtainedfrom acomputeralgebrasystem isnot meaningfuland mighteven be misleading. This motivates the development of a framework of numerical methods to solveproblemssuchaspolynomialroot-finding,elimination,etc.Althoughsome numer-ical methods are available, there was no unifying framework. For example, the most important andknownnumericalmethodforsolvingmultivariatepolynomialsystemsis numerical polynomial homotopy continuation(NPHC) [8–11]. Homotopy methodscan also be used for elimination [12], but it is not possible to compute greatest common divisors or do polynomial divisions in this framework. It is possible to solve many of these problems withmultivariate polynomials inanumerical linearalgebraframework. An important milestone in this respect was the discovery of Stetter that finding the affine roots of a multivariate polynomial system is equivalent with solving an eigen-value problem [13,14]. In his approach however, itis still necessary to first compute a Gröbner basisbefore theeigenvalue problem canbe writtendown.We havedeveloped a numerical linear algebraframework where no symbolical computations are required. Problemssuchaseliminationofvariables[15],computinganapproximategreatest com-mon divisor [16], computing a Gröbner and border basis [17], finding the affine roots of a polynomial system are all solved numerically in this Polynomial Numerical Lin-ear Algebra (PNLA)framework. TheMacaulaymatrix plays acentralrole inallthese problems [18,19].
The Macaulaymatrixis definedforacertaindegreed and itis importantto under-standhowitssizeanddimensionsofitsfundamentalsubspaceschangeasafunctionofd.
Theseareinfactdescribedbydifferentpolynomialsind.Wewillexplainhowthiscomes aboutinthisarticlethroughtheanalysisoftheleftnullspaceofM (d).Thisanalysiswill leadus tothedefinitionofthedegreeofregularityofamultivariatepolynomialsystem
f1,. . . ,fs.Inaddition,twoalgorithmsthatdeterminethepolynomialexpressionl(d) for
thedimension of theleft nullspace of M (d) are presented. Once thepolynomials that describe the dimensions of the fundamentalsubspaces are understood, we move on to describetheright nullspace oftheMacaulaymatrix.Theaffineroots ofamultivariate polynomialsystemf1,. . . ,fsare usuallydescribedbythedual vectorspaceofthe quo-tient spaceCn/f
1,. . . ,fs. Inthis articlewewill describethe rightnullspace ofM (d) as theannihilatorof itsrow spaceand express this bymeansof afunctional basis.An importantobservationhereisthatalsorootsatinfinityaredescribed bythis functional basisand thatmultiplicitiesofroots resultinacertainmultiplicity structure.Stetter’s methodto find the affineroots of amultivariate polynomial systemis bymeans of an eigenvalueproblemthatdescribesamonomial multiplicationwithinthequotient space
Cn/f
1,. . . ,fs.Wewillextendthis monomialmultiplicationpropertyto theprojective caseand present analgorithm thatderives the corresponding projective multiplication matrices.
Before defining the Macaulay matrix, we first discuss the numerical linear algebra frameworkinwhichwe willdescribemultivariate polynomialsandintroducethe mono-mialordering thatwillbe used.Mostof thealgorithms describedinthisarticleare im-plementedinMATLAB[20]/Octave[21]andarefreelyavailableathttps://github.com/ kbatseli/PNLA_MATLAB_OCTAVE.
2. Macaulaymatrix
BeforedefiningtheMacaulay matrix,wefirst discusssomebasicdefinitionsand no-tation. The ring of multivariate polynomials in n variables with complex coefficients is denotedby Cn. It is easy to show thatthe subset of Cn, containing all multivariate polynomials of total degreesfrom 0up to d forms avector space. We will denote this vectorspacebyCdn.Weconsidermultivariatepolynomialsthatoccurincomputerscience andengineeringapplications andlimitourselvestherefore,withoutloss ofgenerality,to multivariatepolynomialswith onlyrealcoefficients.Throughoutthisarticlewewill use amonomial basisasabasisforCdn.Since thetotalnumberof monomialsinn variables
fromdegree0upto degreed isgivenby
q(d) =
d + n n
it follows that dimCn
d = q(d).The total degree of amonomial xa = x a1
1 . . . xann is de-fined as |a| = ni=1ai. The degree of a polynomial p, deg(p), then corresponds with the degreeof the monomial of p with highestdegree. It is possible to order the terms of multivariate polynomials in different ways and results such as acomputed Gröbner basiswilldependonwhichordering isused.Forexample,itiswell-knownthata Gröb-nerbasiswithrespecttothelexicographicmonomialorderingistypicallymorecomplex (moretermsandofhigherdegree)thenwithrespect tothereverselexicographic order-ing[4,p. 114].It isthereforeimportantto specifywhichordering isused. Foraformal
definition of monomial orderings together with adetailed descriptionof some relevant orderings in computational algebraic geometry see [3,4]. The monomial ordering used in this article is the graded xel ordering [15, p. 3], which is sometimes also called the degree negative lexicographic monomial ordering. This ordering is graded because it first comparesthedegreesof thetwomonomials a,b andapplies thexelorderingwhen there isatie.Byconventionacoefficientvectorwillalwaysbearowvector.Depending on the context we will use the label f for both a polynomial and its coefficient vec-tor. (.)T will denote thetranspose ofamatrixor vector. Pointsat infinitywill playan important role in this article and these are naturally connectedto homogeneous poly-nomials. A polynomial of degree d is homogeneous when every term is of degree d. A non-homogeneous polynomial can easily be made homogeneous by introducing an extravariable x0.
Definition 2.1. Letf ∈ Cn
d ofdegreed, then itshomogenizationfh ∈ Cdn+1 is the poly-nomial obtainedbymultiplyingeachterm off with apower ofx0 suchthatits degree
becomesd.
The vectorspace ofall homogeneouspolynomials inn+ 1 variables and ofdegree d
is denoted byPn
d. This vectorspace is spanned by allmonomials inn+ 1 variablesof degreed andhence
dimPdn= d + n n .
In order to describe solutionsets of systems of homogeneous polynomials, the pro-jective space needs to be introduced. First,an equivalence relation ∼ on the non-zero points ofCn+1 isdefinedbysetting
x0, . . . , xn∼ (x0, . . . , xn)
ifthere isanon-zeroλ∈ C suchthat(x0,. . . ,xn)= λ(x0,. . . ,xn).
Definition 2.2. (See [4, p. 368].) The n-dimensional projective space Pn is the set of equivalenceclassesof∼ onCn+1− {0}.Eachnon-zero(n+ 1)-tuple(x
0,. . . ,xn) defines apointp inPn,andwesaythat(x
0,. . . ,xn) arehomogeneouscoordinatesofp.
The origin (0,. . . ,0) ∈ Cn+1 is not a point in the projective space. Because of the equivalence relation∼,aninfinitenumberofprojectivepoints(x0,. . . ,xn) canbe
asso-ciated with 1affine point(x1,. . . ,xn). Theaffinespace Cn canbe retrievedas a‘slice’
of theprojective space:
Cn =(1, x
1, . . . , xn)∈ Pn
This meansthatgiven aprojectivepoint p= (x0,. . . ,xn) with x0= 0, itsaffine
coun-terpart is (1,x1
x0,. . . ,
xn
x0). The projective points for which x0 = 0 are called points at
infinity.
2.1. Definition
Wenow introducethemain object ofthis article, theMacaulaymatrix, and discuss howitsdimensionsgrowas afunctionofthedegreed.
Definition2.3.Givenasetofpolynomialsf1,. . . ,fs∈ Cn,eachofdegreedi(i= 1,. . . ,s), thentheMacaulaymatrixofdegreed isthematrixcontainingthecoefficientsof
M (d) = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ f1 x1f1 .. . xd−d1 n f1 f2 x1f2 .. . xd−ds n fs ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ (1)
where each polynomialfi is multiplied with allmonomials from degree 0up to d− di foralli= 1,. . . ,s.
When constructing the Macaulay matrix, it is more practical to start with the co-efficient vectors of the original polynomial system f1,. . . ,fs, after which all the rows
correspondingtomultipliedpolynomialsxaf
i uptoadegreemax(d1,. . . ,ds) areadded. Then one can add the coefficient vectors of all polynomials xaf
i of one degree higher and so forth until thedesired degree d is obtained.This is illustrated in thefollowing example.
Example2.1.ForthefollowingpolynomialsysteminC2 2
f
1: x1x2− 2x2= 0,
f2: x2− 3 = 0,
we have that max(d1,d2) = 2 and we want to construct M (3). The first 2 rows then
correspond withthecoefficientvectorsof f1,f2. Sincemax(d1,d2)= 2 and d2 = 1,the
next2 rows correspond to the coefficient vectors of x1f2 and x2f2 of degree 2. Notice
thatthese first4rowsmakeupM (2) whenthecolumnsarelimitedtoallmonomialsof degree0upto2.Thenextrowsthatareaddedarethecoefficientvectorsofx1f1,x2f1
andx21f2,x1x2f2,x22f2which areallpolynomialsof degree3.Thisway ofconstructing
M (3) = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 x1 x2 x21 x1x2 x22 x31 x21x2 x1x22 x32 f1 0 0 −2 0 1 0 0 0 0 0 f2 −3 0 1 0 0 0 0 0 0 0 x1f2 0 −3 0 0 1 0 0 0 0 0 x2f2 0 0 −3 0 0 1 0 0 0 0 x1f1 0 0 0 0 −2 0 0 1 0 0 x2f1 0 0 0 0 0 −2 0 0 1 0 x2 1f2 0 0 0 −3 0 0 0 1 0 0 x1x2f2 0 0 0 0 −3 0 0 0 1 0 x2 2f2 0 0 0 0 0 −3 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ .
Each row of the Macaulay matrix contains the coefficients of one of the fi’s. The multiplicationofthefi’swiththemonomialsxaresultsintheMacaulaymatrixhavinga quasi-Toeplitz structure,inthesenseofbeing almostor nearlyToeplitz. TheMacaulay matrix depends explicitly on the degree d for which it is defined, hence the notation
M (d). The reason (1) is called the Macaulay matrix is because it was Macaulay who introduced this matrix, drawing from earlier work by Sylvester [22], in his work on elimination theory,resultantsand solvingmultivariate polynomialsystems[23,24].Itis infactageneralizationoftheSylvestermatrixton variablesandanarbitrarydegreed.
The MATLAB/Octave routineinthe PNLA frameworkthatreturns M (d) for agiven polynomialsystemanddegreed isgetM.m.
Foragivendegreed,thenumberofrowsp(d) ofM (d) isgivenbythepolynomial
p(d) = s i=1 d− di+ n n = s n!d n+ Odn−1 (2)
and thenumberofcolumns q(d) by
q(d) = d + n n = 1 n!d n+ Odn−1. (3)
From these two expressions it is clear that the number of rows will grow faster than the numberofcolumns as soonas s> 1.Wedenote therankof M (d) byr(d) and the dimensionof itsleftandright nullspacebyl(d) and c(d) respectively.Therank-nullity theorems forM (d)T andM (d) arethenexpressedas
q(d) = r(d) + c(d), p(d) = r(d) + l(d).
This shows that r(d), l(d), c(d) are also polynomials over all positive integers d >
max(d1,. . . ,ds).ThispolynomialincreaseofthedimensionsofM (d) isduetothe combi-natorialexplosionofthenumberofmonomialsandisthemainbottleneckwhensolving
problems inpractice.The following exampleillustrates this polynomialnatureof r(d), l(d), c(d),together with theinterestingobservation thatthedegreeof c(d) is linkedto thedimensionoftheaffinesolutionset off1,. . . ,fs.
Example 2.2. Consider the Macaulay matrix M (d) of one multivariate polynomial
f ∈ Cn
d1.Thestructureofthematrixensuresthatitisalwaysoffullrowrank(l(d)= 0).
Hence r(d) = p(d) = d− d1+ n n =d n n! + n(n− 2d1+ 1) 2n! d n−1+ Odn−2, and c(d) = q(d)− r(d) = d n n! + n(n + 1) 2n! d n−1+ Odn−2−dn n! − n(n− 2d1+ 1) 2n! d n−1− Odn−2 = d1 (n− 1)!d n−1+ Odn−2.
Aninterestingobservation from Example 2.2is thatthe dimensionof theright null spaceisapolynomialofdegreen− 1.Thiscorrespondsintuitivelywiththedimensionof theaffinesolutionset.Forexample,thesurfaceofaballis2-dimensionalandisdescribed byonepolynomialin3variables(n= 3)ofdegree2.Theconnectionbetweenthedegree ofc(d) andthedimensionofthesolutionsetwillbe mademoreexplicitinSection4.
2.2. Rowspace
We will now present two interpretations of the row space of the Macaulay matrix. Bothinterpretations willbeimportantlateronwhenwediscuss thenullspaceofM (d)
andM (d)T.FirstwediscusstheaffineinterpretationoftherowspaceofM (d).Therow spaceofM (d),denotedbyMd,containsalln-variatepolynomials
Md= s i=1 hifi : hi∈ Cd−dn i (i = 1, . . . , s) . (4)
Apolynomialidealf1,. . . ,fs isdefinedastheset
f1, . . . , fs = s i=1 hifi : h1, . . . , hs∈ Cn .
Itisnowtempting tohavethefollowing interpretation
or in words: the row space of M (d) contains all polynomials of the ideal f1,. . . ,fs from degree0upto d.Thisisnotnecessarilyvalid.Md doesnotingeneralcontainall polynomials ofdegreed thatcanbewrittenas apolynomialcombination(4).
Example 2.3.Considerthefollowing polynomialsysteminC2 1 f1: x21+ 2x1+ 1 = 0, f2: x21+ x1+ 1 = 0. From (−1 − x1)f1+ (2 + x1)f2= 1 (5)
itfollowsthat1∈ f1,f2.However,1∈ M/ 2.Infact,(5)tellsusthat1∈ M3.Deciding
whether a given multivariate polynomial p lies in the polynomial ideal generated by
f1,. . . ,fs is called the ideal membership problem. A numerical algorithm that solves this problemisalsopresented in[17].
As Example 2.3 shows, thereason thatnotall polynomials ofdegree d lie inMd is thatit ispossible thatapolynomialcombination ofadegreehigherthand isrequired. There isadifferentinterpretationoftherowspaceofM (d) suchthatallpolynomialsof degree d arecontained init.This requires thenotionofhomogeneous polynomialsand will becrucialinSection4tounderstandthenullspaceoftheMacaulaymatrix.Itwill turn outthatthedimensionofthenullspaceofM (d) isrelatedto thetotalnumberof projectiverootsofthepolynomialsystem.Thisincludesrootsatinfinityandinthisway homogeneous polynomials are relevant. Given a set of non-homogeneous polynomials
f1,. . . ,fs wecanalsointerpret Mdas thevectorspace
Md= s i=1 hifih: hi ∈ Pdn−di (i = 1, . . . , s) , (6) wherethefh
i’sarehomogeneousversionsoff1,. . . ,fsandthehi’sarealsohomogeneous. Thecorrespondinghomogeneousidealisdenotedbyfh
1,. . . ,fsh.Thehomogeneity guar-anteesthatallhomogeneouspolynomials of degreed arecontainedinMd.Orinother words,
Md=
f1h, . . . , fshd,
where fh
1,. . . ,fshd are all homogeneouspolynomials of degreed contained inthe ho-mogeneousidealfh
1,. . . ,fsh.Animportantconsequenceisthenthat dimf1h, . . . , fshd = r(d).
The homogenization of f1,. . . ,fs typically introduces extra roots that satisfy x0 = 0
and at least one xi = 0 (i = 1,. . . ,s). These points are roots at infinity. We revisit
Example 2.1to illustratethispoint.
Example 2.4. The homogenization of the polynomial system in Example 2.1 is f1h =
x1x2−2x2x0= 0,f2h= x2−3x0= 0.Allhomogeneouspolynomials
2
i=1hifihofdegree 3belongto therowspaceof
⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ x30 x1x20 x2x20 x21x0 x1x2x0 x22x0 x31 x12x2 x1x22 x32 x0f1 0 0 −2 0 1 0 0 0 0 0 x2 0f2 −3 0 1 0 0 0 0 0 0 0 x0x1f2 0 −3 0 0 1 0 0 0 0 0 x0x2f2 0 0 −3 0 0 1 0 0 0 0 x1f1 0 0 0 0 −2 0 0 1 0 0 x2f1 0 0 0 0 0 −2 0 0 1 0 x2 1f2 0 0 0 −3 0 0 0 1 0 0 x1x2f2 0 0 0 0 −3 0 0 0 1 0 x2 2f2 0 0 0 0 0 −3 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ,
whichequals M (3) fromExample 2.1.Note thatthenon-homogeneous polynomial sys-temhad only1root ={(2,3)}.Afterhomogenization,theresultingpolynomialsystem
fh
1,f2hhas2nontrivialroots={(1,2,3),(0,1,0)}.
Thehomogeneousinterpretationis ineffectnothingbutarelabellingof thecolumns androwsofM (d).Bothoftheseinterpretationsareusedinthisarticle.Whendiscussing theleft nullspace ofM (d) we will employtheaffine interpretation,while forthe right nullspace thehomogeneousinterpretationisimportant.
3. Leftnullspace
InthissectionwepresentadetailedanalysisoftheleftnullspaceofM (d).Themain focus will be to derive the polynomial expression l(d) for a given polynomial system
f1,. . . ,fs.Thiswillnaturallyleadtothenotionofthedegreeofregularity,whichwillbe importantfordescribingtherightnullspace.TheleftnullspaceofM (d),null(M (d)T), isthevectorspace
The vectors h are not to be interpreted as polynomials but rather as s-tuples of mul-tivariate polynomials.Indeed, from theaffine row spaceinterpretation we seethat the expression hM (d)= 0 is equivalentwith
s i=1
hifi= 0. (7)
The vector h therefore contains the coefficients of all polynomials hi. A polynomial combination such as (7) is called asyzygy [25], from the Greek word συζυγια,which refers to an alignment of three celestial bodies. In our case, the polynomials hi are thought tobeinsyzygywiththepolynomials fi,hencetheirpolynomialcombinationis zero.Thisbringsus totheinterpretationofthedimensionoftheleftnullspace,l(d).It simplycounts thetotalnumberofsyzygiesthatoccurinMd. Itisthereforepossibleto identify witheach syzygyalinearly dependentrowof M (d). Thelineardependence of this particularrowisthenwithrespect totheremainingrowsofM (d).
3.1. Expressing l(d) intermsofbasissyzygies
Eachelement oftheleftnullspacecorrespondswithasyzygyofmultivariate polyno-mials and with alinearlydependent rowof theMacaulay matrixM (d). Algorithm 3.1
findsamaximalsetofsuchlinearlydependentrowsl foragivenMacaulaymatrixM (d),
startingfrom thetoprowr1 whereri standsfortheith rowoftheMacaulaymatrix.
Algorithm3.1. Findamaximalsetoflinearlydependentrows Input:MacaulaymatrixM (d)
Output:amaximalsetoflinearlydependentrowsl l← ∅
if r1= 0 then
l← [l,r1]
end if
for i= 2: 1: p(d) do
if ri linearlydependentwithrespect to{r1,. . . ,ri−1} then
l← [l,ri] end if end for
Once thelinearlydependentrowsof M (d) areidentifiedwith Algorithm 3.1, itthen becomes possible to write down a polynomialexpression for l(d). Indeed, suppose the firstelementl1ofl isfoundusingAlgorithm 3.1.Thisrowisthenlinearlydependentwith
respect to all rowsabove it. Infact, this linear dependence expresses acertain syzygy s
i=1hifi= 0.Therowl1thenalsocorrespondswithacertainmonomialmultiplexαifk sinceitisarowoftheMacaulaymatrix.Observenowthat
xβj
s i=1
hifi = 0,
whichmeansthatallrowscorrespondingwithxβjxα
ifkwillalsobelinearlydependent.We willcalll1 inthiscaseabasissyzygyanditsmonomialmultiplesxβjl1,derivedsyzygies.
Thedegreeofabasissyzygyistakentobethemaximaldegreeoveritsterms.Theabove observationcannowbe summarizedinthefollowinglemma.
Lemma3.1. If abasissyzygyl hasadegree dl,thenitintroduces aterm d− dl+ n n (8) tothepolynomiall(d).
Proof. Thisfollows fromxβjsi=1hifi= 0 andthefactthatthetotalnumberof mono-mialsxβj atadegreed≥ dl isgivenby(8). 2
Itcanbeshownthatthenumberofbasissyzygiesisfinite.Thisisinfactlinkedwith thefinitenessoftheGröbnerbasisforapolynomialideal[3,p. 223].Wewillnowassume thatallbasissyzygieswerefound,usingforexampleAlgorithm 3.1forasufficientlylarge degreed,andexplainhowallbasissyzygiescanbe usedto deriveanexpressionforthe polynomial l(d). As mentioned above, each linearly dependent row can be labelled as amonomial multiple of oneof the polynomials f1,. . . ,fs. The first stepof the syzygy
analysis is todivide thebasis syzygiesinto groupsaccordingto the polynomialthatis multiplied.Thefollowingexamplewillbeusedthroughoutthewholesubsectioninorder toderivethealgorithm todeterminel(d) bymeansofbasissyzygies.
Example3.1. ConsiderthefollowingpolynomialsysteminC3
⎧ ⎪ ⎨ ⎪ ⎩ f1: x2y2+ z = 0, f2: xy− 1 = 0, f3: x2+ z = 0,
where x1 = x, x2= y, x3 = z. Thefirst basissyzygy isfound, using Algorithm 3.1, in
M (4) andcorresponds withtherow
xyf3.
TheremainingbasissyzygiesareallfoundinM (6) andcorrespond withtherows
Wecannowdividethesebasissyzygiesintothefollowingtwo groups {xyf3} and x3yf2, x2y2f2, xy2zf2 ,
whichisonegroupforf3 andonegroupforf2.
Thekeyobservationhereisthateachofthesegroupscanbeanalysedseparatelysince nointerferencebetweenrowsofdifferentgroupsispossible(indeed,theyinvolvedifferent polynomials).WewillnowcontinueExample 3.1andshowhowallcontributionsofbasis syzygiestol(d) are describedbybinomialcoefficients.
Example 3.2. The first group {xyf3} has only oneelement and describes a syzygy of
degree4.Lemma 3.1tellsus thenthatthiswill introduceaterm
d− 4 + 3
3
to l(d).Wecanthereforewrite
l(d) = d− 4 + 3 3 =1 6d 3− d2+11 6d− 1 (d ≥ 4). (9)
The second group has three basis syzygies, {x3yf2,x2y2f2,xy2zf2}, of degree 6 and
therefore introduces3terms
d− 6 + 3
3
.
Wecanthereforeupdatel(d) to
l(d) = d− 4 + 3 3 + 3 d− 6 + 3 3 =2 3d 3− 7d2+76 3d− 31 (d ≥ 4). (10)
Expression (10) for l(d) is still valid for degrees d ≥ 4, since the 3 extra binomial coefficient terms correspond with polynomials that have roots at d ∈ {3,4,5}. The differencebetween(9)and(10)isonlyvisiblethereforeford≥ 6.Wehavenotyetfound the final expression for l(d) however. The3binomial coefficient termsat degree 6will count too many contributions.Take forexamplethe basissyzygiescorresponding with the rows x3yf
2 and x2y2f2. Their least commonmultiple is x3y2f2, whichmeans that
the linearlydependentrowx3y2f
2 will becounted twiceby (10). It shouldhoweverbe
Example 3.2 shows that the analysis of all syzygies is reduced to a combinatorial problem:withinagroupofbasissyzygies,oneneedstocountthetotalnumberoflinearly dependentrowsthese basissyzygies‘generate’.This combinatorialproblem issolvedby theInclusion–Exclusionprinciple.
Theorem 3.1 (Inclusion–Exclusion Principle). (See [4, p. 454].) Let A1,. . . ,An be a
collection offinite setswith |Ai| the cardinality ofAi.Then n i=1 Ai = n k=1 (−1)k+1 1≤i1<...<ik≤n |Ai1∩ . . . ∩ Aik| . (11)
TheInclusion–ExclusionPrincipleis thefinal componentthatallowsusto conclude theanalysis of allsyzygies.Wenow apply Theorem 3.1and derivethe final expression forl(d) inExample 3.2.
Example3.3.Ifwedenotethesetofallmonomialmultiplesofx3yf
2 byA1 andlikewise
A2,A3 forx2y2f2,xy2zf2respectively,thenapplyingTheorem 3.1 onthese setsresults
inthe final expression for l(d). Note thatall termsof (11) for k = 1 arethe binomial coefficientsofLemma 3.1,whichalreadyhavebeenaddedtol(d) in(10).Theremaining analysis is hence on all terms of (11) for k ≥ 2. The cardinality of the intersections
A1∩ A2, A1∩ A3,A2∩ A3 aredescribedbythebinomialcoefficients
d− 7 + 3 3 , d− 7 + 3 3 , d− 8 + 3 3 ,
whichwilleachcontributeto l(d) withaminussignsincek = 2.Thedegreesforwhich thesetermsareintroducedarethedegreesoftheleastcommonmultiplesbetweenx3yf
2
andx2y2f
2,between x3yf2 andxy2zf2 and betweenx2y2f2 andxy2zf2.These degrees
are 7, 8 and 7 respectively. The next intersection, A1∩ A2∩ A3, corresponds with a
binomial term introduced at the degree of the least common multiple of all 3 basis syzygiesinf2.Thisleastcommonmultipleisx3y2zf2withadegreeof8.Thisconcludes
theanalysisofalllinearlydependentrowsofM (d) andwecanthereforewrite
l(d) = d− 4 + 3 3 + 3 d− 6 + 3 3 − 2 d− 7 + 3 3 − d− 8 + 3 3 + d− 8 + 3 3 = d− 4 + 3 3 + 3 d− 6 + 3 3 − 2 d− 7 + 3 3 = 1 3d 3− 2d2+2 3d + 9 (d≥ 4).
Notethatthetermsatdegree8havecancelledoneanother.Sincethetermd−7+33 has rootsatd∈ {4,5,6},this expressionforl(d) isvalidforalld≥ 4.
Animportantobservationisthatoncethepolynomialexpressionl(d) isknown,then the rankof M (d) and the dimension of its nullspace are also fully determined for all
d≥ 4 by: r(d) = p(d)− l(d) = 1 6d 3+ d2+5 6d− 10, c(d) = q(d)− r(d) = d + 11.
Algorithm 3.2 summarizestheiterative syzygyanalysis outlined abovetodetermine thepolynomialexpressionforl(d).
Algorithm3.2. Findl(d)
Input:polynomialsystemf1,. . . ,fn ∈ Cn Output: l(d)
d← max(deg(f1),deg(f2),. . . ,deg(fs))
l← ∅ l(d)← ∅
while not allbasissyzygiesfound do
findnewbasissyzygiesl usingAlgorithm 3.1onM (d)
update l(d) usingLemma 3.1andtheInclusion–Exclusionprinciple
d← d+ 1 end while
A stop criterionis needed to be ableto decide whetherallbasis syzygieshavebeen found. Thisis intimatelylinked withtheoccurrenceofaGröbner basisinMd.Indeed, it is a well-knownresult thatall basis syzygiesof f1,. . . ,fs can be determined from a Gröbner basis [3, p. 223]. The reduction to zero of every S-polynomial of a pair of polynomials inaGröbner basisprovidesabasissyzygy.Thisimpliesthatitisrequired toconstructM (d) foradegreewhichcontainsalltheseS-polynomials.Thelinkbetween aGröbner basis andMd isdescribed indetailin[17]. InSection4, wewill exclusively discusstherightnullspaceofM (d) forpolynomialsystemswithafinitesetofprojective roots.LazardshowedthatinthiscasethemaximaldegreeofareducedGröbnerbasisis atmostd1+ . . . + dn+1− n+ 1 withdn+1= 1 ifs= n.Thisprovidesanupperboundon thedegreeatwhichallbasissyzygiescanbefound.Thefinitenessoftheamountofbasis syzygieshasaveryimportantconsequence.ItensuresthatAlgorithm 3.2stopsandthat the polynomialexpressionsforl(d), r(d) and c(d) donotchangeanymoreafter afinite amount ofsteps. As aconsequence,the domainof the final polynomialexpressionsfor
l(d),r(d) andc(d) hasaparticularlowerboundd.InExample 3.3,thedomainofthese polynomialswaslowerboundedby4.Wewillcallthislowerboundonthedomainofthe final polynomialexpressions forl(d),r(d) andc(d) thedegreeofregularity.
Definition 3.1. The minimal degree d ∈ N for which the output of Algorithm 3.2 de-scribesthedimensionoftheleftnullspaceofM (d) iscalledthedegreeofregularity.
This degreeof regularity isof vitalimportanceinthenextsection where wediscuss the right nullspace of M (d).From the analysis abovewe see thatinorder to find the
degree of regularity d, it is required to construct the Macaulay matrix for a degree largerthan d. Forexample, thedegree of regularity d = 4 of the polynomialsystem inExample 3.1wasfoundfrom M (6).Thefollowingexampleillustratesthatthedegree forwhichallbasissyzygiesarefoundcanbe quitehighandconsequentlythatthetotal numberofbinomialtermsofl(d) canbeverylarge.
Example3.4.ConsiderthefollowingpolynomialsysteminC4
⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ f1: x22x3+ 2 x1x2x4− 2 x1− x3= 0, f2: −x31x3+ 4 x1x22x3+ 4 x21x2x4+ 2 x32x4+ 4 x21− 10 x22 + 4 x1x3− 10 x2x4+ 2 = 0, f3: 2 x2x3x4+ x1x24− x1− 2 x3= 0, f4: −x1x33+ 4 x2x23x4+ 4 x1x3x24+ 2 x2x34+ 4 x1x3 + 4 x2 3− 10 x2x4− 10 x24+ 2 = 0,
with degreesd1 = d3 = 3,d2 = d4 = 4.Thefirst basissyzygy groupis foundinM (7)
andcorrespondswith therowx2
2x3f2.Thenextbasissyzygiesgrouparetherows
x22x3f3, x32x4f3, x1x22x24f3, x31x2x23f3, x21x2x34f3, x41x3x44f3, x21x2x43x24f3,
x22x74f3, x21x2x63x4f3, x1x2x84f3, x21x2x83f3, x31x3x84f3
andarefoundforthedegrees
{6, 7, 8, 9, 9, 12, 12, 12, 13, 13, 14, 15}.
Thelastbasissyzygygrouparetherows
x22x3f4, x2x3x4f4, x1x2x24f4, x32x4f4, x31x24f4, x21x34f4, x41x3x4f4,
x31x2x23f4, x31x32x4f4, x21x33x4f4, x1x52f4, x51x23f4, x41x33f4, x31x43f4
andarefoundforthedegrees
{7, 7, 8, 8, 9, 9, 10, 10, 10, 10, 10, 11, 11, 11}.
Application oftheInclusion–ExclusionPrinciplefor eachofthese groupsresultsinthe finalexpression l(d) = d− 6 + 4 4 + 4 d− 7 + 4 4 + d− 8 + 4 4 − 6 d− 11 + 4 4 + 4 d− 13 + 4 4 − d− 14 + 4 4 =1 8d 4−13 12d 3−5 8d 2+19 12d + 105 (d≥ 10).
Observe thatl(d) consists infact of 10 241binomial terms, of which10 225 cancel out until only16termsareleft.
3.2. Recursivealgorithm todeterminel(d)
As Example 3.4shows, using Algorithm 3.2tofind theexpression forl(d) results in alargenumberofbinomialterms.A largenumberofcomputations areactuallywasted since most of these binomialterms cancel outand therefore do notcontribute to l(d).
Thisapproachhasthedisadvantagethatthetotalnumberofbinomialtermsthatneeded to be computed grows combinatorially, while in fact the majority of them cancel one another. Algorithm 3.3 remedies this disadvantage. This iterative algorithm does not need tocheckeachrowoftheMacaulaymatrix,norhastousetheInclusion–Exclusion Principle to findthe final expressionof l(d) and corresponding degreeof regularity d. Instead of finding basis syzygies and calculating how these will propagate to higher degrees,wewillsimplyupdatel(d) iteratively.Thealgorithmispresentedinpseudo-code in Algorithm 3.3. The main idea is to compute the numerical value p(d)− r(d) and
compare it with the evaluation of l(d) for eachdegree d.Here p(d) doesnot represent the polynomial expression in d but rather the evaluation of this polynomial for the valueofd inthealgorithm.Thepolynomialexpressionforr(d) isnotknownandhence cannot beevaluatedbutthisevaluation isfoundbythedeterminationofthenumerical rank of M (d). If p(d)− r(d) > l(d), then the polynomial l(d) needs to count p(d)− r(d)− l(d) additional linearly dependent rows. Furthermore, each of these additional rowswillpropagatetohigherdegreesandgiverisetoextrabinomialterms.Similarly,if
p(d)−r(d)< l(d),thentoomanylinearlydependentrowswerecountedandl(d) needsto beadjustedwithl(d)−p(d)+r(d) negativebinomialcontributions.Thep(d)−r(d)< l(d)
degrees for which positive contributions to l(d) are made are stored in the vector d+
and likewise for the l(d)− p(d)+ r(d) negative contributions in d−. All information on how to express l(d) interms of binomial coefficients is hencestored in d+ and d−.
Iterations need to start from a degree d = min(deg(f1),deg(f2),. . . ,deg(fs)) to make sure thatthe updatingof l(d) reflects the correctoccurrence of newbasis and derived syzygies.Obviously,whentherearepolynomialsfioff1,. . . ,fswithadegreehigherthan min(deg(f1),deg(f2),. . . ,deg(fs)),thentheyarenotincludedinM (d) forthatparticular
degree. Since thealgorithm iterates over thedegrees, a Singular Value Decomposition (SVD)-based recursive orthogonalization algorithm [26] can be used to determine the numericalvaluer(d) foreachiterationwithouttheneedtoconstructthewholeMacaulay matrix.Thebinomialterminl(d) thatappears atthehighestdegreedmax= max(d+,d−)
determinesthedegreeofregularity d.Indeed,thepolynomial
d− dmax+ n
n
haszerosford={dmax−n,dmax−n+ 1,. . . ,dmax−1} andthereforethefinalexpression
theoryofresultants.Macaulayshowedin[23]thatitispossible todeterminewhethera homogeneouspolynomialsystemf1h,. . . ,fnh ofdegrees d1,. . . ,dn hasanontrivial com-monrootbycomputingthedeterminantofasubmatrixofM (d) ford=ni=1di−n+ 1. Thisessentiallymeansthatd≤n
i=1di− n+ 1,whichresultsinamaximaldegreeof 1+ni=1di inAlgorithm 3.3.
Algorithm3.3.Findl(d) anddegreeofregularityd
Input:polynomialsystemf1,. . . ,fn∈ Cn Output:l(d) anddegreeofregularity d
d← min(deg(f1),deg(f2),. . . ,deg(fs))
d+← ∅ d−← ∅ l(d)← 0 r(d)← rankM (d) while d≤ 1+sideg(fi) do if p(d)− r(d)> l(d) then addd p(d)− r(d)− l(d) timesto d+ else if p(d)− r(d)< l(d) then addd l(d)− p(d)+ r(d) timesto d− end if l(d)←|d+| i=1 d−d+(i)+n n −|d−| i=1 d−d−(i)+n n d← d+ 1 r(d)← rankM (d) end while d← max(d +,d−)− n
Symbolical methods compute a Gröbner basis G of f1,. . . ,fs in order to describe allbasis syzygiesand find the degreeof regularity. Algorithm 3.3 doesnot requirethe computationof aGröbner basis.Instead,oneneedstodeterminethenumericalrankof
M (d) for increasing degreesd. Algorithm 3.3 is implementedinthe MATLAB/Octave routinealn.m.
Example3.5.WeillustrateAlgorithm 3.3withthepolynomialsystemfromExample 3.4: ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ f1: x22x3+ 2 x1x2x4− 2 x1− x3= 0, f2: −x31x3+ 4 x1x22x3+ 4 x21x2x4+ 2 x32x4+ 4 x21− 10 x22 + 4 x1x3− 10 x2x4+ 2 = 0, f3: 2 x2x3x4+ x1x24− x1− 2 x3= 0, f4: −x1x33+ 4 x2x23x4+ 4 x1x3x24+ 2 x2x34+ 4 x1x3 + 4 x2 3− 10 x2x4− 10 x24+ 2 = 0.
Theexpression for l(d) is initialized to 0.The Macaulay matrixis of full rowrank for degrees3,4and5. Ford= 6,we havethatp(6)− r(6)= 100− 99= 1> l(6)= 0.We
therefore setd+= 6 and l(d) = d− 6 + 4 4 .
After incrementingthedegreewefindthatp(7)− r(7)= 210− 201= 9> l(7)= 5 and we thereforeupdated+ andl(d) tod+={6,7,7,7,7} and
l(d) = d− 6 + 4 4 + 4 d− 7 + 4 4
respectively.Thealgorithmfinishesatd= 14 with
d+={6, 7, 7, 7, 7, 8, 13, 13, 13, 13},
and
d−={11, 11, 11, 11, 11, 11, 14}, whichindeedcorrespondswiththefinal expressionforl(d)
l(d) = d− 6 + 4 4 + 4 d− 7 + 4 4 + d− 8 + 4 4 + 4 d− 13 + 4 4 − 6 d− 11 + 4 4 − d− 14 + 4 4 = 1 8d 4−13 12d 3−5 8d 2+19 12d + 105 (d≥ 10).
ObservethatAlgorithm 3.3findsthedesiredexpression forl(d) atd= 14. Incontrast, thebasissyzygiesanalysis usingAlgorithm 3.2 requiredtheconstructionofM (15) and
the computationof 10 241binomial terms.From Algorithm 3.3 wecan derivethat the dimension oftheleft nullspace ofM (d) isdescribed bythefollowing piecewise-defined function l(d) = ⎧ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎩ 0, if 3≤ d ≤ 5 d4 24− 7d3 12 + 71d2 24 − 77d 12 + 5, if d = 6 d4 4 − 9d3 2 + 121d2 4 − 90d + 100, if 7≤ d ≤ 10 1 8d 4−13 12d 3−5 8d 2+19 12d + 105, if d≥ 10.
The minimal degree for which the final polynomial expression for l(d) is valid is 10. Hencethedegreeofregularityisd= 10.
4. Rightnullspace
Theknowledgethatc(d) isa polynomialind togetherwith thehomogeneous inter-pretationofMdallowsustolinkthenullspaceoftheMacaulaymatrixwiththenumber of projective roots of f1,. . . ,fs. The notion of the dual of the row space will play an importantroleindescribingtheseroots.Afterhavingdescribedthisdualvectorspacewe willextendthemonomialmultiplicationpropertyofStetter’seigenvalueproblem[13,14]
totheprojectivecase.
4.1. Link withprojective roots
Itis aclassicresultthatforapolynomialsystemf1h,. . . ,fsh with afinite numberof projectiveroots,thequotientringPdn/f1h,. . . ,fshd isafinite-dimensionalvectorspace
[3,4]. Thedimension ofthis vectorspace equals thetotalnumberof projective rootsof
fh
1,. . . ,fsh, counting multiplicities for a large enough degree d. From the rank-nullity theoremofM (d) itthenfollowsthat
c(d) = q(d)− r(d) = dimPdn− dimf1h, . . . , fshd = dimPdn/ f1h, . . . , fsh d (12)
Thisleadstothefollowingtheorem.
Theorem4.1. Fora zero-dimensionalhomogeneous idealf1h,. . . ,fsh withm projective roots(counting multiplicities) anddegree ofregularity d wehave that
c(d) = m ∀d ≥ d.
Proof. Thisfollows from(12)andDefinition 3.1. 2
Furthermore,when s= n,then c(d)= m= d1· · · ds accordingtoBézout’sTheorem
[3,p. 97]. This effectivelylinks thedegrees ofthe polynomials f1,. . . ,fs to the nullity oftheMacaulaymatrix.Theaffinerootscanberetrievedfromageneralizedeigenvalue problemasdiscussedin[14,27].Inthisarticlewewillextendthisgeneralizedeigenvalue approach to the projective case. Another interesting result is that the degree of the polynomialc(d) isthedimensionofprojectivevarietyoff1h,. . . ,fsh.
Definition4.1.Thepolynomial
c(d) = dimPdn/f1h, . . . , fshd ∀d ≥ d
iscalledtheprojectiveHilbertPolynomial[4,p. 462].Thedegreeofthispolynomialc(d)
Example 4.1. Forthepolynomialsystemfrom Example 3.1we hadthatc(d)= d+ 11. Since this is a polynomial of degree 1, it follows that the projective solution set of
fh
1,. . . ,fsh is one-dimensional. The number of affine solutions of f1h,. . . ,fsh is finite:
{(1,1,1,−1),(1,−1,−1,−1)},which implies thatthe non-zero-dimensionalpart of the solutionset lies‘atinfinity’.
From here on, we will only consider polynomial systems with a finite number of projective roots.
4.2. Dualvectorspace
As soonas d ≥ d and the number of projective roots is finite, then a basis of the nullspacecanbeexplicitlywrittendownintermsoftheroots.Thisrequiresthenotion ofthedual vectorspace.Wedenote thedual vectorspaceofCn
d byCn
d ,thedualofMd byMd andtheannihilatorofMdbyMod.Bydefinition,theelements ofModmap each element of Md to zero. Therefore, dimMod = c(d), which implies Mod ∼= null(M (d)). A basis forMod isdescribedbydifferentialfunctionals.
Definition 4.2. (See[14, p. 8].)Let j ∈ Nn
0 and z ∈ Cn,then thedifferential functional
∂j|z∈ Cn d isdefinedby ∂j|z≡ 1 j1! . . . jn! ∂j1+...+jn ∂xj1 1 . . . ∂x jn n z where |z standsforevaluationinz = (x1,. . . ,xn).
Being elements of the dual vector space, these differential functionals ∂j|z can be representedas vectors.Thisisillustratedinthefollowing simpleexample.
Example 4.2. In C32 the functionals ∂00|z, ∂10|z, ∂01|z, ∂20|z, ∂11|z and ∂02|z have the following coefficientvectors
⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ∂00|z ∂10|z ∂01|z ∂20|z ∂11|z ∂02|z 1 0 0 0 0 0 x1 1 0 0 0 0 x2 0 1 0 0 0 x2 1 2x1 0 1 0 0 x1x2 x2 x1 0 1 0 x2 2 0 2x2 0 0 1 x3 1 3x21 0 3x1 0 0 x2 1x2 2x1x2 x21 x2 2x1 0 x1x22 x22 2x1x2 0 0 x1 x32 0 3x22 0 0 3x2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ , (13)
withz = (x1,x2)∈ C2.Thehomogeneousinterpretationof Md impliesthatthe differ-ential functionals also havea homogeneousinterpretation. For theexample above,the coefficientvectorsofthecorrespondingdifferentialfunctionals inP2
3 are ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ∂000|z ∂010|z ∂001|z ∂020|z ∂011|z ∂002|z x3 0 0 0 0 0 0 x2 0x1 x20 0 0 0 0 x20x2 0 x20 0 0 0 x0x21 2x0x1 0 x0 0 0 x0x1x2 x0x2 x0x1 0 x0 0 x0x22 0 2x0x2 0 0 x0 x3 1 3x21 0 3x1 0 0 x21x2 2x1x2 x21 x2 2x1 0 x1x22 x22 2x1x2 0 0 x1 x3 2 0 3x22 0 0 3x2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ , (14)
with z = (x0,x1,x2) ∈ P2. The matrix in(13) canbe retrieved from (14) by setting
x0= 1.
We will make no further distinction between the linear functionals ∂j|z and their coefficient vectors. Notice that these coefficient vectors are column vectors, since they arethedualelementsoftherowspaceofM (d).Thevectors∂0|zintheaffinecasecanbe seenasageneralizationoftheVandermondestructuretothemultivariatecase.Applying thedifferential functional ∂j|z to theelements of Md is thensimply thematrix vector multiplicationM (d)∂j|z.
Weknowthatwhenapolynomialsystemfh
1,. . . ,fshhasafinitenumberofm projec-tiveroots,thendimMo
d= c(d)= m foralld≥ d.Hence,abasisforMod willconsistof differential functionals,evaluated ineach projective root and taking multiplicitiesinto account.This bringsustothedefinitionofthecanonicalnullspaceof M (d).
Definition 4.3. Let f1,. . . ,fs ∈ Cn with azero-dimensional projective solutionset and let m1,. . . ,mt be the multiplicities of the t projective roots zi (1 ≤ i ≤ t) such that t
i=1mi = m. Then for all d≥ d there exists amatrix K ofm linearly independent columnssuchthat
M (d)K = 0.
Furthermore,K canbe partitionedinto
suchthateachKiconsistsofmilinearcombinationsofdifferentialfunctionals∂j|zi∈ P n d (1≤ i≤ t).Wewill callthismatrixK thecanonicalnullspaceofM (d).
Definition 4.3 explicitlydepends onthe homogeneousinterpretationof Md.Indeed, it is only for the homogeneous case that the projective roots come into the picture. Defining themultiplicity ofa zerousing thedual space goes back to Macaulay[24]. It is also reminiscent of the univariate case. Remember that for aunivariate polynomial
f (x)∈ C1
d,azeroz withmultiplicitym meansthat ⎛ ⎜ ⎜ ⎜ ⎜ ⎝ f f D1 .. . f Dm−1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎠∂0|z= 0, (15)
where Di istheith orderdifferentialoperator.Orinotherwords, f (z)= f(z)= f=
. . . = f(m−1)(z)= 0.Alternatively,(15)canbewrittenas
f (∂0|z ∂1|z ∂2|z . . . ∂m−1|z) = 0. (16) As already mentioned inDefinition 4.3, the multivariate casegeneralizes this principle byrequiringlinearcombinationsofdifferentialfunctionals.
Example 4.3. Consider thefollowing polynomial systeminC2
2 withthe affine root z =
(1,2,3)∈ P2ofmultiplicity 4andnorootsat infinity
(x2− 3)2= 0,
(x1+ 1− x2)2= 0.
Thedegreeofregularity d is2andthecanonicalnullspaceis
K = ( ∂000|(1,2,3) ∂100|(1,2,3) ∂010|(1,2,3) ∂110|(1,2,3)− 2∂020|(1,2,3)). (17)
ThedifferentlinearcombinationsoffunctionalsneededtoconstructKiarecalledthe multiplicity structure of the root zi. Observethat the multiplicity structure of aroot is notunique. Indeed,for anynonsingular mi× mi matrix T we havethatthecolumn spaceofKiequalsthecolumnspaceofKiT .Findingthemultiplicitystructureforagiven root ofapolynomialsystemisanactiveareaofresearch[28,29].Iterativealgorithmsto compute the multiplicity structure ofa rootz such as in[29,30] exploittheclosedness propertyofthedifferentialfunctionals∂j|z[14,p. 330]toreducethesizeofthematrices ineveryiteration.Wewillnotfurtherdiscussthese algorithmshere.
4.3. ExtendingStetter’s eigenvalue problemtotheprojectivecase
Stetter’s approach to find the affine roots of multivariate polynomial systems is to phraseitasaneigenvalueproblem[13,14].Thiseigenvalueproblemexpressesthe multi-plicationofmonomialswithinthequotientspaceCn/f
1,. . . ,fs.Thestandardprocedure tofindallaffinerootsis
1. tocomputeaGröbnerbasisG forf1,. . . ,fs, 2. deriveamonomialbasisforCn/f
1,. . . ,fs from G,
3. solvetheeigenvalueproblem thatexpressesthemonomialmultiplicationwithinthe quotientspace Cn/f
1,. . . ,fs,
4. readofftheaffinesolutionsfrom theeigenvectors.
We now show how onecan write down Stetter’s eigenvalueproblem withoutthe com-putationofaGröbnerbasis.Instead, oneneedsto writedownthemultiplicationofthe differential functional ∂j|z with a monomial. Suppose that the polynomialsystem has only affine roots with no multiplicities. No homogeneous interpretation of the canon-ical null space K is required then. For a functional ∂0|z the following relationship holds: ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 x1 x2 .. . xd−1 n ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ x1= ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ x1 x2 1 x1x2 .. . x1xdn−1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ . (18)
Orinotherwords,theoperationofmultiplying∂0|zwithx1correspondswithaparticular
rowselectionofthesamefunctional.Inthiscase,thefirstrowbecomesthesecond,the second is mapped to rown+ 2, and so forth. If we wantto express the multiplication of functionals in Cn
d with monomials of degree 1, then only the rows corresponding withmonomials upto degreed− 1 areallowed to be multiplied. Indeed,monomials of degree d would be ‘shifted’ out of the coefficient vector. Hence (18) can be rewritten as
S10∂0|zx1= S01∂0|z, (19) where S10 selects at most all
d−1+n n
rows corresponding with monomials from de-gree 0 up to d− 1 and S01 selects the corresponding rows after multiplication with
Example 4.4.Writingdown(19)forfunctionalsinC22 resultsin ⎛ ⎝10 01 00 00 00 00 0 0 1 0 0 0 ⎞ ⎠ ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 x1 x2 x21 x1x2 x2 2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ x1= ⎛ ⎝00 10 00 01 00 00 0 0 0 0 1 0 ⎞ ⎠ ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 x1 x2 x21 x1x2 x2 2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ .
The selectionmatrixS10 selectsinthis caseallrowscorresponding withallmonomials
of degree0upto 1.
IfS10inExample 4.4 wouldhaveselectedanyoftherowscorresponding with
mono-mials of degree 2, then no corresponding S01 could have been constructed since the
functionals do not contain any monomials of degree 3. Observe that relations similar to (19)can be written downfor multiplication with any variable xi. Indeed, for every variable xi acorresponding rowselectionmatrixcanbe derived.Undertheassumption thatnoneofthem affinerootshasmultiplicities,(18)canbeextendedtoallfunctionals of affineroots K = (∂0|z1 . . . ∂0|zm) andany multiplication variablexi so thatwe can write
Si0KDi= S0iK, (20)
whereDi isasquarediagonalmatrixcontainingxi’s.Themeaningofthe0indexinthe rowselectionmatriceswillbecomeclearwhenwediscussthehomogeneouscase.Now,it will beshownhow(20)canbewrittenasastandardorgeneralizedeigenvalueproblem. Theq(d)× m matrixK cannotbedirectlycomputedfrom M (d).Itispossiblehowever, to computeanumericalbasisN forthenullspaceofM (d),usingforexampletheSVD. Since both N and K are bases for the null space, they are related by a nonsingular matrixV , orinotherwords,K = N V . WecanthereforereplaceK byN V in(20)and obtain
Si0N V Di= S0iN V,
BV Di= AV, (21)
where wehaveset B = Si0N and A= S0iN .Since A andB are overdetermined
matri-ces, (21)is not aneigenvalue problem yet. Oneway to transform itinto ageneralized eigenvalue problemwouldbe to chooseSi0 and S0i suchthatA andB become square.
Inaddition, B hasto be regularsinceweknow thatthediagonal ofDi containsthexi components of the affine roots. The second way is to transform (21)into an ordinary eigenvalueproblem bywritingitas
where B† istheMoore–Penrose pseudoinverseof B. Oncethe eigenvectorsV are com-puted, the canonical null space K can be reconstructed as N V . The affine roots are thensimplyread offfrom K.Theproceduretonumericallycompute theaffinerootsof
f1,. . . ,fswithoutaGröbner basisishence:
1. computeanumericalbasisforthenullspaceofM (d) ford≥ d, 2. chooseSij,Sji andformB,A,
3. solvetheeigenvalueproblemBV Di= AV , 4. reconstructK andreadofftheaffinesolutions.
More details on numerical affine root-finding using these two approaches, even when therearerootsat infinity,canbefoundin[19].
In order to extend the monomial multiplication property (18) to the homogeneous case,oneadjustmentneedsto bemade:amonomial multiplicationneedstobeinserted to its right-hand side. Indeed, since each component of the differential functional ∂j|z hasthesametotaldegree,∂j|zxi will notcorrespond withaparticularrow selectionof
∂j|z.Thefollowing exampleillustrates theextension of(18)tothehomogeneouscase. Example 4.5. The homogeneous extension of the monomial multiplication property in
Example 4.4is ⎛ ⎝10 01 00 00 00 00 0 0 1 0 0 0 ⎞ ⎠ ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ x2 0 x0x1 x0x2 x2 1 x1x2 x2 2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ x1= ⎛ ⎝00 10 00 01 00 00 0 0 0 0 1 0 ⎞ ⎠ ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ x2 0 x0x1 x0x2 x2 1 x1x2 x2 2 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ x0.
Theselectionmatrixontheleft-handsideisagaindenotedS10.Heretheleftmostindex
means thatinthis side of theequation wemultiply with x1 and with x0 on theother
side.LikewisefortherowselectionmatrixS01 ontheright-handside,herewemultiply
withx0andwithx1ontheotherside.Observethattheseselectionmatricesareidentical
tothoseof Example 4.4.
Consequently,thehomogeneousmonomial multiplicationpropertycanbewrittenas
SijKDi= SjiKDj. (22)
Choosingthetwo monomialsxi,xj with whichboth sidesof theequation willbe mul-tipliedrespectively completely determines theselection matrices Sij, Sji. Since (22) is validforthehomogeneouscase,thisrelationalsoholdsforrootsat infinity.
4.4. Derivingprojectivemultiplication matrices
Stetter’s Central Theorem of multivariate polynomialroot finding says thatthe Di matrices inthe eigenvalueproblems derived from (21) will ingeneral notbe aJordan normalform whentheaffinerootshavemultiplicities[14,p. 52]. Theyare howeverstill upper triangular. The sameis true for the homogeneouscase. We will now derive our algorithm tocomputetheseuppertriangularmultiplication matricesDiinthe homoge-neouscasebymeansofanexample.
Example 4.6.Supposewehavethefollowing polynomialsysteminC42
(x1− 2)x22= 0,
(x2− 3)2= 0.
It canbe easilyshownthatthere isanaffine rootz1 = (1,2,3) with multiplicity2and
arootat infinityz2= (0,1,0) withmultiplicity 4.Thedegreeofregularity is4and the
canonicalnullspaceK is
(K1 K2) = (∂000|z1 ∂100|z1+ 2∂010|z1 ∂000|z2 ∂100|z2 ∂010|z2 2∂200|z2+ 3∂101|z2).
Westartwiththeaffinerootz1and determineitsentriesofthemultiplicationmatrices
D0, D1. The first step is to write down the homogeneous multiplication property for
∂000|z1
S01∂000|z1x0= S10∂000|z1x1. (23)
Both ∂100|z1 and∂010|z1 areneededtodescribe thesecond columnofthecanonicalnull
space. Bytaking thepartial derivativeof(23)withrespect tox0we obtain
S01(∂000|z1+ ∂100|z1x0) = S10∂100|z1x1. (24)
Likewise, takingthe partial derivativeof (23)with respect to x1 andmultiplying both
sideswith 2resultsin
S012∂010|z1x0= S10(2∂000|z1+ 2∂010|z1x1). (25)
Combining (23), (24)and(25)andevaluatingx0,x1 resultsin
S01( ∂000|z1 ∂100|z1+ 2∂010|z1) 1 1 0 1 = S10( ∂000|z1 ∂100|z1+ 2∂010|z1) 2 2 0 2 . (26)
Observerhow thereisa2onthesuperdiagonalontheright-handside.This already in-dicatesthatthemultiplicationmatrixwillnotbeinJordannormalform.Theprocedure todeterminethepartofthemultiplicationmatricesfortherootatinfinityiscompletely analogoustotheanalysisabove.Startingoffwithwritingdownthemultiplication prop-ertyfor∂000|z2
S01∂000|z2x0= S10∂000|z2x1,
andtakingpartial derivativeswith respecttox0,
S01(∂000|z2+ ∂100|z2x0) = S10∂100|z2x1, (27)
andwithrespect tox1
S01∂010|z2x0= S10(∂000|z2+ ∂010|z2x1).
Differentialfunctionalsofseconddegreearealsoneeded.Takingthepartialderivativeof
∂100|z2 withrespect tox0 tocompute∂200|z2 resultsin
∂ ∂x0
(∂100|z2) = 2 ∂200|z2
dueto Definition 4.2.Anadditionalpartial derivativeof(27)withrespectto x0 results
inthedesiredequation
S01(∂100|z2+ ∂100|z2+ 2∂200|z2x0) = S102∂200|z2x1,
S01(2 ∂100|z2+ 2 ∂200|z2x0) = S102∂200|z2x1,
andlikewiseforthe∂101|z2 functional
S01(3∂001|z2+ 3∂101|z2x0) = S103∂101|z2x1.
Combiningtheresultsandevaluatingthecomponents,wecanwrite
S01( ∂000|z2 ∂100|z2 ∂010|z2) ⎛ ⎜ ⎜ ⎝ 0 1 0 0 0 0 0 2 0 0 0 3 0 0 0 0 ⎞ ⎟ ⎟ ⎠ = S10( ∂000|z2 ∂100|z2 ∂010|z2) ⎛ ⎜ ⎜ ⎝ 1 0 1 0 0 1 0 0 0 0 1 0 0 0 0 1 ⎞ ⎟ ⎟ ⎠ . (28)
Finally,combining(26)and(28)resultsinthehomogeneousmultiplicationproperty for thewholecanonicalnullspaceK
S01K ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 1 1 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 2 0 0 0 0 0 3 0 0 0 0 0 0 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ = S10K ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 2 2 0 0 0 0 0 2 0 0 0 0 0 0 1 0 1 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ .
Byabuseofnotation,wewillalsodenotetheuppertriangularmultiplicationmatrices byDi,Dj.Doingso,(22)alsoholdsforthecaseofrootswithmultiplicities.Algorithm 4.1 summarizesthewholeproceduretoderivethemultiplicationmatricesfortheprojective case.
Algorithm4.1.
Input:canonicalnullspaceK,multiplicationmonomialsxi,xj Output:projectivemultiplicationmatricesDi,Dj
determine rowselectionmatricesfrom xi, xj for each projectiverootzi do
write homogeneousmultiplicationpropertyforzi
apply partialderivativestoobtainhigherorderfunctionals make linearcombinations toreconstructthemultiplicity structure write outtheobtainedmultiplicationrelationinmatrixform end for
combinemultiplication matricesforeachroot intoDi,Dj
Since themultiplicationofmonomialsiscommutative,thisimpliesthatthe multipli-cation matrices Di, Dj will also commute. The occurrence of multiple roots, whether theyareaffineoratinfinity,posesaproblemtodeterminethemviaaneigenvalue com-putation. Indeed,inorderto beabletowrite
SijKDi= SjiKDj
as aneigenvalueproblem,eitherDi orDj hastobeinvertible.Foraffinerootsthiswill alwaysbethecaseforD0.Rootsatinfinitywillneedamorecarefulchoiceofxi,xj.Let us supposethatDj isinvertible,wecanthenwrite
SijKDiDj−1= SjiK.
Again,substitutingK byN V andsettingSijN = B,SjiN = A weget
Let
Jij = T−1DiD−1j T
betheJordannormalform ofDiD−1j .Thisallowsusto rewrite(29)as
BV T Jij = AV T. (30)
ThiscanagainbeconvertedintoageneralizedeigenvalueproblembymakingA,B square
and B regularorinto an ordinaryeigenvalueproblem bycomputing thepseudoinverse ofB.TheretrievedeigenvectorsareinthiscaseV T ,fromwhichwecannotreconstruct the canonicalnull space K = N V since T is unknown. Anadditional difficultylies in thenumericallystablecomputationoftheJordannormalformJ .Alternatively,onecan computethenumericallystable Schurdecomposition
B†A = QU Q−1,
whereQ isunitaryand U uppertriangular.Theeigenvaluesxi/xj canthenbe readoff fromthediagonalofU .Intheprojectivecasethere isalwaysat leastoneDj invertible. This meansthatthe n differentxi/xj componentsof theprojective roots canbe com-putedfromn Schur factorizations.Afterwardsthesecomponentsneedtobematched to formtheprojectiveroots
x0 xj ,x1 xj , . . . , 1, . . . ,xn xj ,
wherethe1appearsinthejth position.
5. Conclusions
In this article we provided a detailed analysis of the left and right null spaces of the Macaulay matrix.The left null space was shownto be linked with the occurrence of syzygies in the row space. It was also demonstrated how the dimension of the left null space is described by a piecewise function of polynomials. Two algorithms were presented that determine these polynomial expressions for l(d). The finiteness of the numberofbasissyzygiesresultedinthenotionofthedegreeofregularity,whichplayed acrucialroleindescribingabasisfortherightnullspaceofM (d) intermsofdifferential functionals. The canonical null space K of M (d) was defined and the multiplication property of this canonical basis was extended to the projective case. This resulted in uppertriangularcommutingmultiplicationmatrices.Finally,wediscussedhowStetter’s eigenvalue problem canbe extended to the case where both affine roots and roots at infinityarepresent.
