# Using Gröbner bases for solving polynomial equations

In my attempts to understand just how computer algebra systems "do things", I tried to dig around a bit on Gröbner bases, which are described almost everywhere as "a generalization of the Euclidean algorithm and Gaussian elimination". I've tried to look for examples of Gröbner bases in action, but have been unable to find any (that can be easily understood).

I could go ahead and just ask for an explanation with examples from people, but I'll go one step further.

General plane conics can be represented by the Cartesian equation

$ax^2+2bxy+cy^2+dx+fy+g=0$

One common problem in dealing with conics is finding out if two conics intersect, and if so, where the intersection point(s) are. Usually one can do all this by eliminating variables accordingly.

Now I ask: how would, say, Buchberger's method, proceed on determining if two given conics intersect, and then find where they intersect?

-
I don't really understand how these computations work. But what I do understand, I learned from a paper of Mumford and Bayer, called "What can be computed in algebraic geometry?". It's really good! – Sam Lichtenstein Aug 28 '10 at 16:27
Since you ask for examples of applications of Gröbner bases, you might be interested to know that they can be used to give a systematic method for generating thermodynamical identities. – jbc Jun 21 '13 at 23:30
@jbc, would you consider writing an answer demonstrating this? – J. M. Jun 21 '13 at 23:56

In order to find the solutions of the system of conics you mention, it suffices to give a procedure to find the projection of the simultaneous vanishing set of the two conics onto some axis, for instance, the $y$-axis: the coordinates of the projection onto the $y$-axis are the $y$-coordinates of the intersection points. Knowing them allows you to substitute back into the equations the values of $y$ and solve a system of equations in a single variable $x$, thus making the problem simpler. In terms of ideals, suppose that we can find a non-zero element $r$ in the ideal generated by the two conics, that depends only on a single variable, say $y$. This means that every solution of the system has $y$-coordinate satisfying the polynomial $r$. Thus we have severely limited the choices for the $y$-coordinates of the intersection (namely, they all have to satisfy the polynomial $r$) and, provided we can actually solve the polynomial $r$ in one variable, we can then substitute the various values of $y$ we found back into the initial equations and solve for $x$, again using our algorithm for solving polynomials in a single variable. So far, so good, I hope! The question is how to produce the element $r$. This is where Gröbner bases come in.

In order to do the computation you ask, one would have to choose an appropriate monomial order. I will gloss over this, and simply "do the obvious". It is an exercise for you to figure out which order to use so that the computation I am going to make is actually a Gröbner bases computation. Scattered throughout the computation there will be also a couple of special cases that I will not deal with: again, you can treat those as exercises for you!

First, I am going to choose a generic basis, so that the first conic has the form $$x^2 + \alpha x + \beta ,$$ where $\alpha$ and $\beta$ are polynomials in $y$ only (I am trying to simplify the notation as much as possible; the only "assumption that I have made is that the coefficient of $x^2$ is non-zero, which can be arranged unless the "conic" is defined by a polynomial of degree at most one).

In this basis, the second equation can be assumed to have no $x^2$ term (indeed, eliminating the $x^2$ term using the first equation would be the first step in any reasonable Gröbner basis computation in which the term $x^2$ is the highest term in sight). Thus I am going to write the second conic as $$\gamma x + \delta$$ where, as before, $\gamma$ and $\delta$ are polynomials in $y$ only. Again, just to fix on a definite case, I am going to assume that the polynomial $\gamma$ is non-zero. (If $\gamma$ were zero, then we would have found a polynomial in the ideal generated by the two conics which only depends on $y$: this was our goal at the start anyway! Obviously, you should worry about the case $\gamma=\delta=0$, but I won't.)

To compute a Gröbner basis, you would now compute S-polynomials: let's do it here as well. I am going to try to eliminate the $x^2$ term from the first equation, by using the second equation. This is easy: multiply the first equation by $\gamma$ and the second one by $x$ and take the difference: we are left with the equation $$(\alpha \gamma - \delta) x + \beta \gamma .$$ Now we are going to eliminate the $x$ term from this last equation using again the second equation: multiply the second equation by $(\alpha \gamma - \delta)$, the last equation by $\gamma$ and subtract to obtain $$\delta^2 - \alpha \gamma \delta + \beta \gamma^2 .$$ We found an expression independent of $x$!! We are done... provided this expression is not identically zero. You can figure out what this would mean and what happens in this case. Note also that the final expression is what we would have obtained if, at the very beginning, we had "solved" $x=-\frac{\delta}{\gamma}$ using the second equation, substituted in the first equation and cleared the denominators.

I hope that this "hybrid" computation explains what is going on: Gröbner bases and Buchberger's algorithm are a systematic way of "solving" systems of equations. You do not have to do any thinking, once you set up the problem. But you need to set up the problem so that it computes what you want. In this case, you could have used several shortcuts to get to the answer, without following all the steps. In more complicated situations, Buchberger's algorithm might be the best way of keeping track of all the steps to be taken.

Let me also comment that, except in the case in which you have a computer doing the computations for you, it is highly unlikely that Gröbner bases will help you with a specific question, unless you could have also found simple tricks to solve it right away.

-
Very cogent. I'll have to try those out on paper myself (using a CAS isn't always helpful for gaining mathematical insight). Thanks! – J. M. Aug 28 '10 at 23:54
@J.M., The mathematical insight you get from learning the ideas behind the CAS: that is explained in many places, like in the book by Cox et al. I referred to. – Mariano Suárez-Alvarez Apr 18 '12 at 17:02
@Mariano: Yes (and I've managed to do it profitably in the months after I asked this question); I bumped up this thread since I wanted to point this out to somebody who was asking about Gröbner in mathematica.SE... (on another note, I think I've read Cox/Little/O'Shea three or four times already. :)) – J. M. Apr 18 '12 at 17:29

You have ro read through the book[David A. Cox, John B. Little and Don O'Shea, Ideals, Varieties, Algorithms].

In any case, if you start with a system of polynomial equations and compute a Groebner basis for the ideal they generate, you get a "maximally triangular" system of equations which is equivalent to the original one---that is why Groebner bases generalize Gaussian elimination.

Let me do a simple example using Macaulay 2. Consider the ring $\mathbb Q[x,y,z]$:

i1 : R = QQ[x, y, z, MonomialOrder => Lex]

o1 = R

o1 : PolynomialRing

the Fermat quintic surface $x^5+y^5+z^5=1$, whose ideal is

i2 : surface = ideal (x^5 + y^5 + z^5 - 1)

5    5    5
o2 = ideal(x  + y  + z  - 1)

o2 : Ideal of R

and the curve $x^2=y^2$, $y^3+1=z^3$:

i3 : curve = ideal (x^2 - y^2, z^3 - y^3 - 1)

2    2     3    3
o3 = ideal (x  - y , - y  + z  - 1)

o3 : Ideal of R

The ideal of the intersection of the surface and the curve is the sum of the ideals of the intersecands:

i4 : intersection = curve + surface

2    2     3    3       5    5    5
o4 = ideal (x  - y , - y  + z  - 1, x  + y  + z  - 1)

o4 : Ideal of R

The number of points on the intersection, counting multiplicities, is the degree of the ideal:

i5 : degree intersection

o5 = 30

if we now compute the degree of the radical of the ideal, we get a lower number, so not all the points are simple:

o6 = 25

Now look at the lexicographic Groebner basis of the ideal of the intersection:

i7 : ideal gens gb intersection

17      16      15      14      13      12       11       10      9      8       7       6      5      4      3
o7 = ideal (9z   + 27z   + 54z   + 50z   + 15z   - 63z   - 104z   - 108z   - 35z  + 35z  + 108z  + 104z  + 63z  - 15z  - 50z  -
----------------------------------------------------------------------------------------------------------------------------
2                 5            16      15      14      13      12       11      10      9       8      7       6       5
54z  - 27z - 9, 20y*z  - 20y + 27z   + 18z   + 45z   - 75z   - 35z   - 164z   + 79z   - 10z  + 230z  - 10z  + 119z  - 164z
----------------------------------------------------------------------------------------------------------------------------
4       3      2              3    3                                       16        15        14        13        12
- 35z  - 155z  + 45z  + 18z + 67, y  - z  + 1, 50x*z - 50x + 50y*z - 50y + 1242z   + 2718z   + 4770z   + 2220z   - 1235z   -
----------------------------------------------------------------------------------------------------------------------------
11        10        9        8        7        6        5        4        3        2                 2    2
7979z   - 7481z   - 6010z  + 1810z  + 5240z  + 9394z  + 5346z  + 1240z  - 4030z  - 4005z  - 2657z - 583, x  - y )

o7 : Ideal of R

The first generator depends only on $z$. The second one and the third, only on $z$ and $y$, and the fourth depends (linearly) on $x$ too. This is a "triangular" system of equations, from which you can compute the 30 points of intersection, assuming you get to solve the first equation for $z$, which is of degree 17...

-
I'll try finding it in a library when I go downtown. For now, I'm hoping somebody would demonstrate how one would proceed in the problem I posed using Gröbner bases. Still, thanks for the pointer! – J. M. Aug 28 '10 at 16:12
Mariano, you could choose a slightly more elaborate example: the "curve" you used is a union of six lines, and it will not be hard to solve the resulting equations! – damiano Aug 28 '10 at 17:26
@damiano, I added a fudge summand of $1$ :) – Mariano Suárez-Alvarez Aug 28 '10 at 17:35
May I suggest also adding a fudge factor to the first equation: that one also factors at the moment... ;) – damiano Aug 28 '10 at 17:40
Eugh, I was just asking for how things work on a quadratic, and I get a cubic and a quintic! :D In 3D! I'll have to pore through these; thanks again Mariano! – J. M. Aug 29 '10 at 0:05

The example that you propose is a bit too general to serve as an introductory pedagogical example (compare e.g. the famous Kahan SIGSAM problem 9 of determining conditions for an ellipse to lie inside the unit circle). Instead, begin with one of the earliest historical examples: Gauss's proof that every symmetric polynomial can be written uniquely as a polynomial in the elementary symmetric polynomials. This is one of the first and the simplest applications of rewriting reduction via a lexicographic order. For a nice presentation see Chapter 7 of Cox, Little, O'Shea: Ideals, Varieties and Algorithms. They also present generalizations to the ring of invariants of a finite matrix group $\rm G \subset GL(n,k)$. You can find online copies of said book via various ebook databases.

-

You might find Lecture 3 of James Carlson's CIMAT Lectures of value as well as other material on this page:

http://www.math.utah.edu/~carlson/cimat/

-
I'll look into it, thanks. – J. M. Aug 28 '10 at 23:54

This is a simplistic version of the Gröbner basis computation and complements Mariano's answer. The idea behind the Buchberger algorithm is to build a basis so that you are to be able to do "unique division" between nonlinear polynomials. This is done by computing (what is called the) S-polynomial between any 2 given polynomials until there are no more produced.

Decide some ordering (grlex, lex, revlex)
Initialize I = <f1, f2>
where we use the polynomials you gave above.
f_1=ax^2+2bxy+cy^2+dx+fy+g
f_2=ax^2+2bxy+cy^2+dx+fy+g
Set n=2

Repeat the following in a loop
Iteratively choose any two polynomials ( f_i, f_j ) from set I
n = n + 1
Compute f[n] = S-polynomial( f_i, f_j )
Until all S-polynomials return 0

At the end of the computation, you will have a set I = {f1,f2, ..., fk} such that one of the polynomials will be in only one term (say 'x'). Solve that equation to get the one coordinate of the intersection, and then work upwards through the rest of the polynomials.

The general question is related to Elimination theory, and AFAIK there are three basic methods ( I can't speak of their complexity )

1. Gröbner basis, as seen above
2. Use Resultants (example below)
3. Wu's method of characteristic sets, which I can't say much about.

If you wanted to use resultants to find the intersection of say

E={f1(x, y,z), f2(x, y, z), f3(x, y, z)},

you would first compute

E1 = {g1 = Resultant_in_x (f1, f2), g2 = Resultant_in_x(f2, f3)}

and then

E2 = { h1 = Resultant_in_y(g1, g2} }

The last polynomial h1 would be in one variable and amenable to solution.

-
Computing/using resultants is a much clearer concept to me, yes. At least "work upwards through the rest of the polynomials" tells why Gröbner generalizes Gaussian. I am supposing from your description that finding the order (e.g. lexicographic ordering) is a nontrivial matter much like pivoting is in Gaussian elimination? – J. M. Aug 29 '10 at 0:02
Also another question: how would Buchberger halt in the case of the two given conics not intersecting (or maybe more properly, having imaginary intersections)? – J. M. Aug 29 '10 at 2:01
The choice of Order is tied to running time/complexity. There are papers on this topic which I haven't read. On your second question, there is an underlying ideal-variety correspondence(i.e. comm. algebra to algebraic geometry). The Groebner basis Ideal I=<f1, ... fk> you end up with may not contain the initial two polynomials you started with (due to S-poly reductions), therefore Variety(Ideal) can be a null set. This gets decided by the Extension Theorem books.google.com/… – SandeepJ Aug 29 '10 at 22:30
– SandeepJ Aug 29 '10 at 22:40

At the OP's request here is a short description of how to use Groebner bases to generate thermodynamical identities. The basic setup for the thermodynamics of a Gibbsian substance is a pair of equations $T=f(p,V)$ and $S=g(p,V)$ which express the temperature and entropy as functions of $p$ and $V$ (the example that everybody knows is the ideal gas which, omitting physical constants, has $f(p,V) = pV$, $g(p,V)= \frac 1 {\gamma-1} (\ln p + \gamma \ln V)$). One can then express all thermodynamical quantities (e.g., the heat capacities--$C_p=\frac{fg_p}{f_p}$ and $C_V=\frac{fg_V}{f_V}$) in terms of the functions $f$, $g$ and their partials. This firstly provides a unified and systematic method to verify such identities but also one to generate them. The basic reason for this is that one can express a very large number of thermodynamic quantities (it runs into many tens of thousands) as simple algebraic expressions of a very few terms so that there are bound to be multiple relations between them (the law of small numbers). This is a situation which Groebner bases are tailor-made to handle. The details can be found in my article "A Systematic Approach to Thermodynamical Identities" (arXiv 1108.4760) which also contains a Mathematica programme which allows one to compute all thermodynamical quantities such as the the heat capacities at the click of a button (doing this by hand can be rather tedious for more elaborate models such as the van der Waals gas or the Feynman gas---it is usually only done in the secondary literature for very simple quantities and models---typically the ideal gas).

-
This is quite nice! I do some amount of work on equations of state (e.g. the Redlich-Kwong EOS and modifications thereof) so this is interesting to me. I recall tediously deriving $\Delta H$ and $\Delta G$ and inversion temperatures by hand... – J. M. Jun 22 '13 at 7:51

Here is an example taken from Golden Fields: A case for the heptagon by Peter Steinbach

By applying Ptolemy's thereom to the diagonals (a and b) of a regular heptagon one arrives at the equations:

{a^2-(1+b),a b-(a+b),b^2-(1+a b)} (={0,0,0})

Applying GroebnerBasis gives:

GroebnerBasis[{a^2-(1+b),a b-(a+b),b^2-(1+a b)},{b,a}]

{1-2 a-a^2+a^3,1-a^2+b}

(={0,0})

The first of which can be solved to find a:

{1.80194,-1.24698,0.445042}

Whence the second equation can be used to find b.

-

## protected by user26857Oct 31 at 9:02

Thank you for your interest in this question. Because it has attracted low-quality answers, posting an answer now requires 10 reputation on this site.