Extracting a Gröbner Basis from Inexact Input ∗Kosaku Nagasaka
Kobe University †E-mail:nagasaka@main.h.kobe-u.ac.jp

Abstract
Computing a Gröbner basis with inexact input is one of the challenging problems
in symbolic-numeric computations, and has been studied by many researchers from thefollowing points of view: 1) numerical stability and 2) mathematical correctness. However,
it is not yet studied from the data mining point of view: how to extract a meaningful
result from the given inexact input when the amount of noise is not small or we do nothave enough information about the input. In this paper, we introduce a stabilization
method for extracting a numerical Gröbner basis from the given inexact input, by basic
numerical linear algebra, low rank approximations and known numerical Gröbner basis
algorithms.
Key words: Approximate Gröbner Bases, Symbolic-Numeric Computations

1

IntroductionGröbner basis is one of the most powerful tools in mathematical computations: solving algebraic
equations, simplifying algebraic relations, automatic theorem proving and so on. Recently,
computing a Gröbner basis for polynomials with inexact coeﬃcients has been studied by several
researchers ([17, 18, 10, 19, 20, 23, 25, 27]) since in practice we have to operate with empirical
data. In Sasaki and Kako [17], this problem is classiﬁed into the ﬁrst and second kinds ofproblem. The ﬁrst kind is computing a Gröbner basis for the ideal generated by the given
polynomials with exact coeﬃcients by numerical arithmetic (e.g. ﬂoating-point arithmetic).
The second kind is for the given polynomials with inexact coeﬃcients having a priori errors.
In this case, we have to operate with a priori errors whether we compute a basis by exact
arithmetic or not. In this paper, we focus on the second kind of problem.
The research aspects for the second kind are mainly numerical stability and mathematicalmeaning, of the results. For the former aspect, several numerical analyses are given by Sasaki
and Kako [18, 17] and an elaboration of the idea of a stabilized Gröbner basis computation
for zero-dimensional cases is studied by Kondratyev and others [10]. For the latter aspect,∗†

This work was supported in part by MEXT KAKENHI (22700011).
3-11 Tsurukabuto, Nada-ku, Kobe 657-8501 JAPAN.

1

Weispfenning introduced a solution by comprehensive Gröbner bases [27, 28, 26]. However,in general, a comprehensive Gröbner system has a huge number of segments and its computation time is quite slow (see [14] for example). Though Weispfenning [27] tried to decrease
the time-complexity by using only a single parameter to represent the inexact parts, whosebounding error mechanism is very similar to interval arithmetic and Traverso and Zanoni [25]pointed out that an interval easily becomes too large when we compute a Gröbner basis byinterval arithmetic. Moreover, Sasaki pointed out that it is very diﬃcult to represent errors as
parameters or intervals in practical situations, in personal conversations.
In this paper, from the data mining point of view which is diﬀerent from these known results,
we introduce a new approach for the second kind of problem by combining the following known
methods.
• Reduced row echelon form ([11, 4, 6, 10, 15, 24]) by which we can compute a Gröbner
basis by the well-known linear algebraic ways.
• Structured low rank approximations including Structured total least squares([13, 12, 9, 3] and citations therein) by which we can compute a rank deﬁcient structured
matrix satisfying some useful properties.
After the introductory subsection below, we review some relevant results of computing aGröbner basis by reduced row echelon form in Section 2. Our deﬁnition for the second kindof numerical Gröbner basis and our new method to compute it are given in Sections 3 and 4,
respectively. In Section 5, we give some useful remarks.

1.1

The Problem to be Solved

We show some examples of numerical diﬃculties to be solved in this paper. For example,
suppose that we compute a Gröbner basis w.r.t. the graded lexicographic order for the ideal
generated by the following polynomials.
Ff p = { 1.01x2 −2.09y 2 +0.002, 4.03x2 y+3.06xy, 2.04x2 y+0.504x2 +1.504xy−1.02y 2 }. (1.1)
In this case, by the algorithm (appGröbner, [17]), we get the following Gröbner bases Gaz,15 andGaz,2 that are computed with initial precision εInit = 10−16 and εInit = 10−3 and approximatezero threshold εZ = 10−15 and εZ = 10−2 , respectively. Please note that 1) we use our implementation1 with Mathematica’s error tracking system instead of eﬀective ﬂoating-point numbers
used in the original paper, 2) we show only limited leading ﬁgures (rounded) in the all examples
in this paper, and 3) the resulting Gröbner basis would be similar to the following ( though it
is depending on the speciﬁed tolerance) if we use other algorithms with some approximate zerotest (e.g. the present author [15]).
Gaz,15 = {1.0} ,Gaz,2 = {1.0x2 − 2.07y 2 + 0.00198, 1.0y 3 + 0.36xy} .
However, using any approximate zero test may make a diﬀerence between the ideal generated
by the resulting Gröbner basis and the given ideal, and we do not have any method to conﬁrm
1

We use the Buchberger algorithm with the Gebauer and Möller criteria and the sugar strategy.

2whether the resulting basis is suitable or not. This may cause a big problem. In fact, we haveto choose a suitable result between Gaz,15 and Gaz,2 . Moreover, there are so many possiblecombinations of initial precision and approximate-zero threshold though we use only two ofthem in this example.
On the other hand, we can use the stabilization techniques for algebraic algorithms ([19, 20,21]) if we assume that the given polynomial set is exact up to their input precisions. Though
we may be able to get some theoretical results, they just converge to Gst = {1} which is justa result by ordinary exact arithmetic. Please note again that the stabilization techniques are
introduced for exact inputs and not designed for inexact inputs.
Moreover, we can compute a comprehensive Gröbner system for Ff p after rationalizing it
and adding parameters on all the coeﬃcients as follows.
101
1
403
Fε = { ( 100+ ε1 )x2 + (− 209+ ε2 )y 2 + ( 500+ ε3 ), ( 100+ ε4 )x2 y + ( 153+ ε5 )xy,
100
50
51
63
188
51
2
2
( 25 + ε6 )x y + ( 125 + ε7 )x + ( 125 + ε8 )xy + (− 50 + ε9 )y 2 }.

In this case, it is very diﬃcult to get reasonable result and we have so many segments inresulting comprehensive Gröbner system (the number of segments is 263 even for very similar
but easier problem). In addition to this problem, even if we can deal with the large numberof segments, we are not able to determine which segment is preferable for this input since ingeneral we do not have any information on a priori errors (that may be continuous as well) .
Therefore, roughly speaking, any known method can not extract a meaningful result fromthe given inexact input in both mathematical and practical aspects at once. To solve thisproblem, from the data mining point of view, we go back to the starting point: “When do wecompute a Gröbner basis with inexact inputs?”. The present author thinks that we compute
it when it seems that there are some algebraic structures on inexact data. That is, we shouldﬁnd a Gröbner basis with lower entropy (unfortunately enlarged by some errors), which may
be hidden by a priori errors. In this paper, we introduce an algorithm from this point of view.
For example, by our algorithm, we have the following numerical Gröbner basis for Ff p .
G = {1.00000x2 − 2.06419y 2 + 1.32909 × 10−14 , 1.00000y 3 + 0.362958xy − 6.43879 × 10−15 y}.
A notable diﬀerence between G and Gaz,2 is the constant term of the ﬁrst basis polynomial,and our algorithm suggests that 0.00198 in Gaz,2 is generated by a priori errors and should bezero or tiny.

2Gröbner basis by RREF

Some researchers studied computing a Gröbner basis by using the reduced row echelon form
(RREF for short, [11, 4]). However, this is not eﬃcient since we have to operate with large
matrices. Using matrix operations partially like F4 and F5 ([6], [7], [10]) may be the best choice
if we want to decrease the computation time. However, it may be useful since we can use somany results from numerical linear algebra for the situation where we must inevitably operatewith a priori errors.
We assume that we compute a Gröbner basis or its variants for the ideal I ⊆ C[~x] generated
by a polynomial set F = {f1 , . . . , fk } ⊂ C[~x] where C[~x] is the polynomial ring in variables3

~x = x1 , . . . , x` over the complex number ﬁeld C. We note that we use the following deﬁnition
though there are several equivalents (see [1] or other text books).
Definition 1 (Gröbner Basis) G = {g1 , . . . , grG } ⊆ I \{0} is a Gröbner basis for I w.r.t.
a fixed term order  if for any f ∈ I \ {0}, there exists gi ∈ G such that ht(gi )|ht(f ) whereht(p) denotes the head term of p(~x) ∈ C[~x] w.r.t. .
/→
We consider the linear map φT : C[~x]T → C1×mT such that φT (ti ) = −ei where C[~x]T isthe submodule of C[~x] generated by an ordered set (the left-most element is the highest) of→terms T = {t1 , . . . , tmT } and −ei (i = 1, . . . , mT ) denotes the canonical basis of C1×mT . The→→−→
coeﬃcient vector −p of p(~x) ∈ C[~x] is deﬁned to be satisfying −p = φT (p) and p(~x) = φ−1
T ( p ).
With a ﬁxed T , we consider the following subset FT of I.
}
{ k
∑
FT =
si (~x)fi (~x) si (~x)fi (~x) ∈ C[~x]T , si (~x) ∈ C[~x] .
i=1

The Buchberger algorithm guarantees that G ⊆ FT if T has a large enough number of elements(however, note that T must include some required elements depending on the term order).
To compute a Gröbner basis for I, we construct the matrix MT (F ) whose each row vector→−p = φT (p) is corresponding to each element of p(~x) ∈ PT (f ) for each f (~x) ∈ F where−→PT (f ) = {ti × f (~x) ∈ C[~x]T | ti = φ−1
T ( ei ), i = 1, . . . , mT }.
By this deﬁnition, FT and the linear space VT generated by the row vectors of MT (F ) are
isomorphic. We note that a matrix is said to be in the reduced row echelon form if it satisﬁes
the following four conditions.
1. All nonzero rows appear above zero rows.
2. Each leading element of a row is in a column to the right of the leading element of therow above it.
3. The leading element in any nonzero row is 1.
4. Every leading element is the only nonzero element in its column.
For the sake of completeness we present the following two lemmas with proofs introducedin [15] by the present author.
Lemma 1 Let MT (F ) be the reduced row echelon form of MT (F ). If gi (~x) ∈ FT for a−→→/
fixed i ∈ {1, . . . , rG }, MT (F ) has a row vector −p satisfying ht(gi ) = ht(φ−1T ( p )).
Proof Since the linear map φT is deﬁned by the ordered set T , each leading element of→→a row vector −p of MT (F ) is corresponding to ht(φT−1 (−p )). The lemma follows from the facts
that FT and VT are isomorphic and all the leading entries of nonzero rows are disjoint sinceMT (F ) is in the reduced row echelon form.

4Lemma 2 Let MT (F ) be the reduced row echelon form of MT (F ). If T has a large enoughnumber of elements, the following GT is a Gröbner basis for I.
{
}
−
→
−
→
GT = φ−1
(
p
)
p
is
a
row
vector
of
M
(F
)
.
T
T
/
Proof The Buchberger algorithm guarantees that G ⊆ FT if T has a large enough numberof elements. Therefore, GT satisﬁes the condition of Deﬁnition 1 since we have gi (~x) ∈ GT ,i = {1, . . . , rG } by Lemma 1.
Example 1 We compute the reduced Gröbner basis w.r.t. the graded lexicographic orderfor the ideal generated by the following polynomials.
13
Fex = {x2 − 2y 2 , 4x2 y + 3xy, 2x2 y + x2 + xy − y 2 }.
2
2

(2.1)

We construct the following matrix MT (Fex ) with T = {x4 , x3 y, x2 y 2 , xy 3 , y 4 , x3 , x2 y, xy 2 , y 3 ,x2 , xy, y 2 } (this set is large enough and not the minimum set for Fex ) and compute MT (Fex ),the reduced row echelon form of MT (Fex ).
(Fex ) =MT (Fex ) = MT 

0 0
1 0 −2 0 0 0 0 0 0 0 0 0
1 0 0 0 0 0 0 32 0 0
 0 1 0 −2 0 0 0 0 0 0 0 0 
0 1 0 0 0 0 0 0 0 0 − 9 0
16




 0 0 1 0 −2 0 0 0 0 0 0 0 
0 0 1 0 0 0 0 3 0 0
0
0
4




 0 0 0 0 0 1 0 −2 0 0 0 0 
0 0 0 1 0 0 0 0 0 0 − 9 0
32




 0 0 0 0 0 0 1 0 −2 0 0 0 
0 0 0 0 1 0 0 3 0 0
0
0
8




 0 0 0 0 0 0 0 0 0 1 0 −2 
 0 0 0 0 0 1 0 −2 0 0
0 0

,
.

3
0 4 0 0 0 0 3 0 0 0 0 0
0 0 0 0 0 0 1 0 0 0
0
4




3
0 0 4 0 0 0 0 3 0 0 0 0
0 0 0 0 0 0 0 0 1 0
0
8




0 0 0 0 0 0 4 0 0 0 3 0
0 0 0 0 0 0 0 0 0 1
0 −2 
 0 2 0 0 0 1 3 −1 0 0 0 0 
0 0 0 0 0 0 0 0 0 0
0
0




2 2
1
3
0 0 2 0 0 0
0 0 0 0 0 0 0 0 0 0−1 0 0 0 
0 0
2
2
0 0 0 0 0 0 2 0 0 12 32 −1
0 0 0 0 0 0 0 0 0 0
0 0
99xy, x2 y 2 + 43 xy 2 , xy 3 − 32
xy, y 4 + 83 xy 2 , x3 − 2xy 2 , x2 y + 34 xy,
We have GT = {x4 + 32 xy 2 , x3 y − 163223
y + 8 xy, x − 2y } hence we obtain the following reduced Gröbner basis if we delete all the
redundant elements.
{}
G = x2 − 2y 2 , y 3 + 38 xy .

We note that in general we must check whether G is actually a Gröbner basis for I, or not,since Lemma 2 is valid only when T has a large enough number of elements.
/

5

2.1

Relation between Rank and Priori Errors

It seems that Ff p in (1.1) at the beginning of Section 1.1 is one of perturbed sets of Fex in (2.1)
at Example 1 by some errors. We note that we can not say this is correct in practice. In general,
no one knows the actual polynomial set (with no errors) of Ff p . This means that Fex is just oneof (inﬁnitely many) possible cases of Ff p . The aim of this paper is adding mathematical andpractical meanings to the situation that we compute a Gröbner basis for Ff p with this property.
To make the situation much simpler, we consider the following polynomial set F̃f p .
1
F̃f p = {x2 − 2y 2 , 4x2 y + 3xy, 2x2 y + x2 + 1.500001xy − y 2 }.
2

(2.2)

The diﬀerence between polynomial sets Fex and F̃f p is just 0.000001xy, and we have the following MT (F̃f p ) which is very similar to MT (Fex ). However, ranks of MT (Fex ) and MT (F̃f p )are not the same: rank(MT (Fex )) = 9 and rank(MT (F̃f p )) = 11.

1 0 −2 0 0 0
0
0
0 0
0
0
 0 1 0 −2 0 0
0
0
0 0
0
0 


 0 0 1 0 −2 0

0
0
0
0
0
0


0 0 0 0 0 1
0
−2
0 0
0
0 


0 0 0 0 0 0

1
0
−2
0
0
0


0 0 0 0 0 0

0
0
0
1
0
−2
.
MT (F̃f p ) = 
0 4 0 0 0 0
3
0
0 0
0
0 


0 0 4 0 0 0

0
3
0
0
0
0


0 0 0 0 0 0
4
0
0 0
3
0 


 0 2 0 0 0 1 1.500001

−1
0
0
0
0


2
1
0 0 2 0 0 0
1.500001
−1
0
0
0
2
1
0 0 0 0 0 02
0
0 2 1.500001 −1
Moreover, the singular values of MT (F̃f p ) are followings: {6.90335, 6.06098, 4.16844, 2.52968,
2.44446, 2.24685, 1.9732, 1.86836, 0.96183, 7.04525 × 10−7 , 2.89525 × 10−7 }.
By Lemma 2 and the rank of MT (F̃f p ), if T also has a large enough number of elements
for F̃f p , a Gröbner basis for the ideal generated by F̃f p has two extra polynomials against aGröbner basis GT for Fex . These two polynomials must be generated by the tiny diﬀerence
between the matrices and can be thought to be aftereﬀects of noise. This means that thediﬀerence ( 0.000001xy ) hides the well structured result GT and G for Fex . Therefore, if we’dlike to unveil the well structured result from the given inexact input, ﬁnding the rank deﬁcientmatrix MT (Fex ) from MT (F̃f p ) is the problem to be solved.
Moreover, it easily happens that the rank of a matrix will increase by some perturbations(errors) though it also happens but rarely that the rank of a matrix may decrease, since a rankincrease or decrease corresponds to an increase or decrease in entropy, respectively. Hence,the diﬀerence of ranks can be thought to be aftereﬀects of a priori errors. This situation isvery similar to several algorithms for approximate polynomial GCD, factorization and so on[9, 3, 29, 8]. Therefore, computing a Gröbner basis with inexact inputs also may be done bycomputing a rank deﬁcient (structured) matrix and some known algorithms using ﬂoating point
arithmetic, based on the conventional Buchberger algorithm.
6Remark 1 One may think that there is a big diﬀerence between exact and numerical ranks:a small element in absolute value may be ignored even if its relative error is much smaller thanothers. Hence a resulting Gröbner basis can be far from the original. This possibly happens
but inevitable for the second kind of problem since any relativeness of errors is easily lost evenif we just do only a single addition of polynomials (see the discussion in [15]). Again, please
note that the aim of this paper is ﬁnding a meaningful result hidden by a priori errors fromthe given inexact input, and in general it is diﬃcult as shown in Section 1.1 that we get some
practical results by known methods (e.g. CGS).
/

3

Structured Gröbner basis

For approximate GCD and factorization, we usually assume that inexact input polynomials(with a priori errors) and unknown exact polynomials (with no error) have the same degree at
most. This means that input and hidden desired polynomials are supposed to have the samestructure. In order to state this kind of constraint in its full generality, we deﬁne a structured
polynomial set. Consider a mapping Si (i = 1, . . . , k) from a parameter space Cni to a set ofpolynomials C[~x]. A polynomial set F = {f1 , . . . , fk } is called S-structured if each elementfi (~x) of the set is in the image of Si , i.e., if there exists a parameter p~i ∈ Cni , such thatfi (~x) = Si (~pi ).
In this paper, we try to compute a numerical Gröbner basis by ﬁnding a rank deﬁcientmatrix as noted in the previous section. From this point of view, we introduce the following
structured Gröbner basis (SGB for short).
Definition 2 (Structured Gröbner basis) We say G is a S-structured Gröbner basis
for F with tolerance ε ∈ R≥0 , rank deficiency d ∈ Z≥0 and set of terms T if they satisfy thefollowing conditions:
1. G is a Gröbner basis for the ideal generated by the following Fst = {fst,1 , . . . , fst,k } ∈ C[~x].
2. F and Fst are S-structured polynomial sets, i.e., there exists parameters p~i , p~st i ∈ Cni ,such that fi (~x) = Si (~pi ) and fst,i (~x) = Si (p~st i ).
3. k(~p1 . . . p~k ) − (p~st 1 . . . p~st k )k= ε where k·k denotes a suitable vector norm.
4. rank(MT (Fst )) = rank(MT (F )) − d.
/By this deﬁnition, the aim of this paper becomes to ﬁnd a method to compute a S-structuredGröbner basis for the given polynomial set F with the tolerance ε, the rank deﬁciency d, theset of terms T and S-structured polynomial set Fst .
Example 2 We show a S-structured Gröbner basis for the following F̃f p introduced in(2.2) at Section 2.1.
1
F̃f p = {x2 − 2y 2 , 4x2 y + 3xy, 2x2 y + x2 + 1.500001xy − y 2 }.
2
7

We take the following structure speciﬁcation S for example.
S1 :(p1 p2 ) →7
p1 x 2 + p2 y 2 ,S2 :(p3 p4 ) →7
p3 x2 y + p4 xy,
S3 : (p5 p6 p7 p8 ) →7
p5 x2 y + p6 x2 + p7 xy + p8 y 2 .
In this setting, the following G is a S-structured Gröbner basis for F̃f p w.r.t. the gradedlexicographic order with tolerance ε = 0.000001 in 2-norm, rank deﬁciency d = 2 and set ofterms T = {x4 , x3 y, x2 y 2 , xy 3 , y 4 , x3 , x2 y, xy 2 , y 3 , x2 , xy, y 2 }.
{}G = x2 − 2y 2 , y 3 + 83 xy .
Moreover, for this idealized case, we have
1
3
Fst = {x2 − 2y 2 , 4x2 y + 3xy, 2x2 y + x2 + xy − y 2 }.
2
2Note that S-structured Gröbner basis is not unique even for the ﬁxed S, ε, d and T . The
above G is just one of them.
/4

Algorithm for computing SGB

In this section, we introduce an algorithm for computing S-structured numerical Gröbner
basis. The algorithm is not for exact S-structured Gröbner basis but for numerical ones sincein general it is diﬃcult to compute a nearby exact rank deﬁcient structured matrix. We givesome remarks about this diﬃculty in Section 5.
To ﬁnd a S-structured polynomial set Fst we compute a rank deﬁcient structured matrix.
This can be done by several known methods for the following problem (SLRA: Structured Low
Rank Approximation).
Definition 3 (SLRA) Given a structure specification S : Rnα → Rm×n , a parameter
vector p~ ∈ Rnα , a vector norm k·k, and an integer r, 0 < r < min{m, n}, find a vector p~∗ suchthat
min k~p − p~∗ k and rank(S(p~∗ )) ≤ r.
p~∗/
The SLRA problem can be solved by the lift-and-project method or solvers for the STLS
problem (Structured Total Least Squares) under some convergent conditions (see Chapter 5
in [2]) and easily extended to complex numbers (see Chapter 2 in [12]). We note that thediﬀerence between the SLRA and STLS problems is the objective rank deﬁciency. The STLS
problem is a special case of the SLRA problem with r = min{m, n} − 1, however in some cases
we can convert the SLRA problem to the STLS problem though this conversion (see [16, 2]) is
not easy task for our problem according to our experiments.
By the above discussions, we have the following algorithm for computing a S-structured
numerical Gröbner basis.
8

Algorithm 1 (S-structured numerical Gröbner basis)
Input: F = {f1 (~x), . . . , fk (~x)} ⊂ C[~x], a term order  and a structure speciﬁcation S.
Output: a S-structured Gröbner basis G for F with the tolerance ε, the rank deﬁciency d, the setof terms T and S-structured polynomial set Fst , or failed.
1. Compute a numerical Gröbner basis G̃ for the ideal generated by F , by some known
algorithms and determine a suitable set of terms T based on the result.
2. Construct MT (F ) ∈ Cm×n , compute its non-zero singular values σi (i = 1, . . . , rorg ) anddetermine a suitable rank deﬁciency d (Take the largest d such that σrorg −d+1 /σrorg −d <10−2 for example). If there is no such d found, output failed.
3. By a solver for the SLRA problem, ﬁnd a S-structured polynomial set Fst with the rank
deﬁciency d satisfying rank(MT (Fst )) = rank(MT (F )) − d and compute the tolerance ε.
If there is no such Fst found, output failed.
4. Compute a numerical Gröbner basis G for the ideal generated by Fst , by some known
algorithms and output {G, ε, d, T , Fst }.
/The above algorithm is very simple however it works as shown in the following examples andfor Ff p in (1.1) at the beginning of Section 1.1.
Example 3 Suppose that we compute a S-structured numerical Gröbner basis w.r.t. thelexicographic order for the following polynomial set F .
F = {x3 + x2 y 2 , x2 y 2 − y 3 , −x2 y + 1.000001x2 + xy 2 + 0.999999y 2 }.
This is the following polynomial set Fδ with δ = 10−6 , which is introduced by Sasaki and Kako
(Example 5,[17]).
Fδ = {x3 + x2 y 2 , x2 y 2 − y 3 , −x2 y + (1 + δ)x2 + xy 2 + (1 − δ)y 2 }.
In this case, by the algorithm (appGröbner, [17]) with initial precision2 εInit = 10−16 andapproximate-zero threshold εZ = 10−15 , we get the following Gröbner basis and the set ofterms T .
Gaz = {1.0x2 − 2.000001y 3 + 0.999998y 2 , 1.0xy 2 + 1.000001y 3 , 1.0y 4 − 1.000001y 3 },T = {x4 y 4 , x4 y 3 , x4 y 2 , x4 y, x3 y 5 , x3 y 4 , x3 y 3 , x3 y 2 , x3 y, x3 , x2 y 5 , x2 y 4 , x2 y 3 ,x2 y 2 , x2 y, x2 , xy 6 , xy 5 , xy 4 , xy 3 , xy 2 , y 6 , y 5 , y 4 , y 3 , y 2 }.
We construct MT (F ) ∈ C30×26 and determine the rank deﬁciency d = 4 since we get thefollowing singular values and σ26−4+1 /σ26−4 = 1.47591 × 10−6 .
{2.12971, 1.84588, . . . (18 elements snipped) . . . ,0.282599, 0.238185, 3.51538 × 10−7 , 3.04624 × 10−7 , 1.66608 × 10−16 , 1.66608 × 10−16 }.
2

In this paper, we assume that we do not have any information on a priori errors hence we suppose that all
the input coeﬃcients uniformly have the machine-precision.

9

We take the following structure speciﬁcation S.
S1 :(p1 p2 ) →7
p 1 x3 + p 2 x2 y 2 ,S2 :(p3 p4 ) →7
p 3 x2 y 2 + p 4 y 3 ,S3 : (p5 p6 p7 p8 ) 7→ p5 x2 y + p6 x2 + p7 xy 2 + p8 y 2 .
By the lift-and-project method (a solver for the SLRA problem), we get the following Sstructured polynomial set Fst with tolerance ε = 7.2326×10−7 in 2-norm, such that rank(MT (Fst ))= rank(MT (F )) − d.
Fst = {1.00000023x3 + 0.9999998x2 y 2 , 1.00000018x2 y 2 − 0.9999998y 3 ,− 0.99999968x2 y + 1.00000073x2 + 1.00000027xy 2 + 0.99999932y 2 }.
Finally, by the algorithm (appGröbner, [17]), we get the following S-structured numericalGröbner basis for F .
G = {1.0x2 − 0.99999802y 5 − 0.99999908y 4 + 0.99999859y 2 ,1.0xy 2 + 0.99999954y 4 , 1.0y 6 − 1.0000006y 3 }.
We note that this result is compatible to the following comprehensive Gröbner system for F if
we think that δ represents a priori errors.
 2{x − y 5 − y 4 + y 2 , xy 2 + y 4 , y 6 − y 3 }(δ = 0),
{2x2 + xy 2 − y 3 , xy 3 + y 3 , y 4 − y 3 }(δ = 1),{x3 + y 3 , x2 y + y 3 − 2y 2 , xy 2 + y 3 , y 4 − y 3 }(δ = −1),{(1 + δ)2 x2 + 2xy 2 − 2δy 3 + (1 − δ 2 )y 2 , −(1 − δ)(xy 2 + y 3 ), y 4 − y 3 } (δ 3 6= δ).
/Example 4 Suppose that we compute a S-structured numerical Gröbner basis w.r.t. thegraded lexicographic order for the following polynomial set F .
F = {1.99974x2 yz + 0.999742x2 , 3.00263xyz 2 − 2.00013yz,
2.99826x2 z − 2.00053x − 2.00006yz − 0.99943}.
This is generated by adding small random perturbations to the following polynomial set Fex
introduced by Sasaki and Kako (Example 6,[17]).
Fex = {2x2 yz + x2 , 3xyz 2 − 2yz, 3x2 z − 2x − 2yz − 1}.
In this case, by the algorithm (appGröbner, [17]) with initial precision εInit = 10−16 andapproximate-zero threshold εZ = 10−15 , we get the following Gröbner basis and the set ofterms T .
Gaz = {1.0x − 0.142550, 1.0y + 0.106986, 1.0z − 4.67292},
T = {x5 yz 3 , x4 y 2 z 3 , x4 yz 4 , x3 y 2 z 4 , x2 y 2 z 5 , x5 yz 2 , x4 y 2 z 2 ,. . . (42 elements snipped) . . . ,2x z, xyz, xz 2 , y 2 z, yz 2 , x2 , xy, xz, yz, z 2 , x, y, z, 1}.
10We construct MT (F ) ∈ C75×63 and determine the rank deﬁciency d = 7 since we get thefollowing singular values and σ63−7+1 /σ63−7 = 0.00724885.
{2.14386, 2.08212, . . . (52 elements snipped) . . . , 0.04441, 0.0384606, 0.000278795,
0.000234535, 0.000181527, 0.000128085, 0.0000768688, 1.78766 × 10−16 , 1.78766 × 10−16 }.
We take the following structure speciﬁcation S.
S1 :(p1 p2 ) 7→ p1 x2 yz + p2 x2 ,S2 :(p3 p4 ) 7→ p3 xyz 2 + p4 yz,
S3 : (p5 p6 p7 p8 ) 7→ p5 x2 z + p6 x + p7 yz + p8 .
We get the following S-structured polynomial set Fst with tolerance ε = 0.00200393 in 2-norm,such that rank(MT (Fst )) = rank(MT (F )) − d, by the lift-and-project method (a solver for the
SLRA problem).
Fst = {1.99978x2 yz + 0.99966x2 , 3.00199xyz 2 − 2.00109yz,
2.99915x2 z − 1.9992x − 1.99991yz − 0.999725}.
Finally, by the algorithm (appGröbner, [17]), we get the following S-structured numericalGröbner basis for F .
G = {1.0yz + 0.499885, 1.0x + 1.33349y}.
We note that the result is very similar to the following exact S-structured Gröbner basis whichis the reduced Gröbner basis for the original polynomial set without perturbations though as
noted in former sections we can not say that this is the unique correct answer. However, theresult shows that we extracted a meaningful numerical Gröbner basis from the given inexact
input.
1
4
Gex = {yz + , x + y}.
2
3/
As shown in the above examples, Algorithm 1 can compute a S-structured numericalGröbner basis. However, we note that the algorithm works only as a preconditioning methodand does not improve the numerical stability. It is the diﬀerent problem. For example, thealgorithm outputs failed and does not give us any new useful information on the following
polynomial set used by Kondratyev and others [10].
{h1 = 1.027748y 2 − 0.467871xy + 2.972252x2 + 0.662026y + 0.0785252x − 3.888889,
h2 = 3.958378y 2 + 0.701807xy + 1.041622x2 − 0.0785252y + 0.662026x − 3.888889.
The reason is very simple. This polynomial set does not have any other S-structured polynomial
set with any small tolerance ε and rank deﬁciency d > 0 if we take T = {y 6 , xy 5 , y 5 , x2 y 4 , xy 4 ,
y 4 , x3 y 3 , x2 y 3 , xy 3 , y 3 , x4 y 2 , x3 y 2 , x2 y 2 , xy 2 , y 2 , x5 y, x4 y, x3 y, x2 y, xy, y, x6 , x5 , x4 , x3 , x2 , x,1} for example. In fact, h1 and h2 are introduced by Stetter [22], by decimal approximations
hence they do not have a priori errors from the data mining point of view in this paper.

115Remarks

In general the SLRA problem is NP-hard except for a few special cases (see [13] for details). Ouralgorithm fully depends on the SLRA problem hence Algorithm 1 is not eﬃcient in computation
time. However, the lift-and-project method can ﬁnd a local optimum in terms of the diﬀerence
ε between F and Fst if each iteration constructs a distinct projection (see [5] for details). Ouralgorithm is to ﬁnd a rank deﬁcient matrix from the matrix having a large gap among its
singular values hence this helps the lift-and-project method to start with a matrix near its local
optimum.
One may be interested in the eﬀectiveness of Algorithm 1 and think that the lift-and-project
method converges only for toy examples. To see this, we generated 3 cases (A,B and C) of 100
polynomial sets in two variables and apply the algorithm with the graded lexicographic order.
Each polynomial set is randomly generated as follows: 1) generate a Gröbner basis with 3
basis polynomials of total degree tb at most having 2,3 or 4 terms with integer coeﬃcients∈ [−10, 10], 2) generate a polynomial set with 4 syzygies of total degree ts at most having 2,3
or 4 terms with integer coeﬃcients ∈ [−10, 10], and 3) perturb the all coeﬃcients by relative
noises ∈ [−10−5 , 10−5 ].
In Table 1, ’#~p’ denotes the average of the number of parameters of structure speciﬁcations
used in the algorithm, ’w/o’ denotes the number of good results by our appGröbner implementation, ’with’ denotes the number of good results by Algorithm 1, ’ε’ denotes the average oftolerances of good results in 2-norm, and ’diﬀ’ denotes the average of diﬀerences between the
resulting Gröbner basis by Algorithm 1 and the Gröbner basis for the exact polynomial set in
2-norm, where “good result” means that its resulting ’diﬀ’ is enough small.
It shows that our algorithm can extract a meaningful basis for the inexact input which the
known algorithm can not work for. Note that our appGröbner implementation can work more if
we have enough information on a priori errors (e.g. suitable initial precision and approximatezero threshold). Moreover, for most of the failure polynomial sets, Algorithm 1 can compute
meaningful results if we specify suitable rank deﬁciencies though this is diﬃcult since in suchcases the gaps among singular values are not large.
case
A
B
C

ts
1
1
2

tb
2
3
2

#~p w/o with
33.9
0
96
40.7
3
45
42.9 1589

ε
diﬀ
0.000434 0.000043
0.000458 0.000034
0.000345 0.000040

Table 1: Experimental Results
For polynomials in more than two variables or polynomials of higher degrees, it becomesdiﬃcult to determine a suitable set of terms T since most of generated syzygies during computations in Step 1 becomes unreliable. This forces us to use much larger T hence the lift-and-project
method has to operate with a huge matrix and its convergence speed becomes very slow. A
simple workaround is replacing Step 4 with the following.
4. Compute a numerical Gröbner basis G for the ideal generated by Fst , by some known
algorithms. Output {G, ε, d, T , Fst } if the suitable set of terms T is unchanged,
otherwise go to Step 2 with a new suitable set of terms T .
12However, this step tends to cause an inﬁnite loop and not so eﬃcient according to our experiments hence this does not resolve the problem.
Therefore, the followings are open questions for further work.
• Finding a method to determine a suitable set of terms T in terms of S-structured Gröbner
basis for the input.
• Finding a better method for our SLRA problem.
• Finding a method to compute an optimal S-structured numerical Gröbner basis w.r.t.
the tolerance ε.

Acknowledgements
The author wishes to thank Prof. Nabeshima and Mr. Kurata for their contributions incomputing the large comprehensive Gröbner system (e.g. 263 segments). The information isvery helpful to note the introductory section.

References[1] T. Becker and V. Weispfenning. Gröbner bases, volume 141 of Graduate Texts in Mathematics. Springer-Verlag, New York, 1993. A computational approach to commutative
algebra, In cooperation with Heinz Kredel.
[2] P. Boito. Structured Matrix Based Methods for Approximate GCD. Ph.D. Thesis. Department of Mathematics, University of Pisa, Italia, 2007.
[3] B. Botting. Structured total least squares for approximate polynomial operations. Master’s
Thesis. School of Computer Science, University of Waterloo, 2004.
[4] M. Byröd, K. Josephson, and K. Åström. Fast optimal three view triangulation. In Y. Yagi,
I. S. Kweon, S. B. Kang, and H. Zha, editors, Asian Conference on Computer Vision, 2007.
[5] M. T. Chu, R. E. Funderlic, and R. J. Plemmons. Structured low rank approximation.
Linear Algebra Appl., 366:157–172, 2003. Special issue on structured matrices: analysis,
algorithms and applications (Cortona, 2000).
[6] J.-C. Faugére. A new eﬃcient algorithm for computing Gröbner bases (F4 ). J. Pure Appl.
Algebra, 139(1-3):61–88, 1999.
[7] J.-C. Faugère. A new eﬃcient algorithm for computing Gröbner bases without reduction tozero (F5 ). In Proceedings of the 2002 International Symposium on Symbolic and Algebraic
Computation, pages 75–83 (electronic), New York, 2002. ACM.
[8] S. Gao, E. Kaltofen, J. May, Z. Yang, and L. Zhi. Approximate factorization of multivariate
polynomials via diﬀerential equations. In ISSAC 2004, pages 167–174. ACM, New York,
2004.
13

[9] E. Kaltofen, Z. Yang, and L. Zhi. Approximate greatest common divisors of several polynomials with linearly constrained coeﬃcients and singular polynomials. In ISSAC ’06:
Proceedings of the 2006 international symposium on Symbolic and algebraic computation,
pages 169–176, New York, NY, USA, 2006. ACM Press.
[10] A. Kondratyev, H. J. Stetter, and S. Winkler. Numerical computation of gröbner bases.
In Proceedings of CASC2004 (Computer Algebra in Scientific Computing), pages 295–306,
2004.
[11] D. Lazard. Gröbner bases, Gaussian elimination and resolution of systems of algebraic
equations. In Computer algebra (London, 1983), volume 162 of Lecture Notes in Comput.
Sci., pages 146–156. Springer, Berlin, 1983.
[12] P. Lemmerling. Structured total least squares: analysis, algorithms and applications. Ph.D.
Thesis. Faculty of Applied Sciences, K.U. Leuven, Belgium, 1999.
[13] I. Markovsky. Structured low-rank approximation and its applications. Automatica J.
IFAC, 44(4):891–909, 2008.
[14] K. Nabeshima. A speed-up of the algorithm for computing comprehensive gröbner systems. In ISSAC 2007: Proceedings of the 2007 international symposium on Symbolic andalgebraic computation, pages 299–306, New York, NY, USA, 2007. ACM.
[15] K. Nagasaka. A study on gröbner basis with inexact input. In Proceedings of CASC 2009,
volume 5743 of Lecture Notes in Comput. Sci., pages 247–258. Springer, Berlin, 2009.
[16] H. Park, L. Zhang, and J. B. Rosen. Low rank approximation of a Hankel matrix bystructured total least norm. BIT, 39(4):757–779, 1999.
[17] T. Sasaki and F. Kako. Computing ﬂoating-point gröbner bases stably. In Proceedings of
SNC 2007, pages 180–189. ACM, New York, 2007.
[18] T. Sasaki and F. Kako. Floating-point gröbner basis computation with ill-conditionedness
estimation. In Proceedings of ASCM 2007, volume 5081 of Lecture Notes in Comput. Sci.,
pages 278–292. Springer, Berlin, 2008.
[19] K. Shirayanagi. An algorithm to compute ﬂoating point gröbner bases. In Proceedings ofthe Maple summer workshop and symposium on Mathematical computation with Maple V: ideas and applications, pages 95–106, Cambridge, MA, USA, 1993. Birkhauser Boston
Inc.
[20] K. Shirayanagi. Floating point gröbner bases. In Selected papers presented at the international IMACS symposium on Symbolic computation, new trends and developments,
pages 509–528, Amsterdam, The Netherlands, The Netherlands, 1996. Elsevier Science
Publishers B. V.
[21] K. Shirayanagi and M. Sweedler. A theory of stabilizing algebraic algorithms. Technical
Report 95-28, pages 1–92, 1995. http://www.ss.u-tokai.ac.jp/s̃hirayan/msitr95-28.pdf.
14[22] H. J. Stetter. Stabilization of polynomial systems solving with Groebner bases. In Proceedings of the 1997 International Symposium on Symbolic and Algebraic Computation (Kihei,HI), pages 117–124 (electronic), New York, 1997. ACM.
[23] H. J. Stetter. Approximate gröbner bases – an impossible concept? In Proceedings of SNC
2005 (Symbolic-Numeric Computation), pages 235–236, 2005.
[24] A. Suzuki. Computing gröbner bases within linear algebra. In Proceedings of CASC 2009,
volume 5743 of Lecture Notes in Comput. Sci., pages 310–321. Springer, Berlin, 2009.
[25] C. Traverso and A. Zanoni. Numerical stability and stabilization of groebner basis computation. In ISSAC 2002: Proceedings of the 2002 international symposium on Symbolicand algebraic computation, pages 262–269, New York, NY, USA, 2002. ACM.
[26] V. Weispfenning. Comprehensive Gröbner bases. J. Symbolic Comput., 14(1):1–29, 1992.
[27] V. Weispfenning. Gröbner bases for inexact input data. In Proceedings of CASC 2003
(Computer Algebra in Scientific Computing), pages 403–411, 2002.
[28] V. Weispfenning. Canonical comprehensive Gröbner bases. J. Symbolic Comput., 36(34):669–683, 2003. International Symposium on Symbolic and Algebraic Computation (ISSAC 2002) (Lille).
[29] L. Zhi. Displacement structure in computing approximate GCD of univariate polynomials.
In Computer mathematics, volume 10 of Lecture Notes Ser. Comput., pages 288–298. World
Sci. Publ., River Edge, NJ, 2003.
The preliminary implementation of Algorithm 1 on Mathematica 7.0, used in this paper canbe found at “ http://wwwmain.h.kobe-u.ac.jp/~nagasaka/research/snap/sgb.nb ”.

15