Elementary uses of Groebner bases
 \def\P{\bf P}
 	In this tutorial we introduce a number
 of basic operations using Groebner bases, and
 at the same time become familiar with a range
 of useful Macaulay2 constructs. The sections are:
 A.  First steps; example with a monomial curve
 B.  Random regular sequences
 C.  Division with remainder
 D.  Elimination theory
 E.  Quotients and saturation
\beginsection{ A. First Steps; example with a monomial curve}\par
 To compute the Groebner basis of an ideal
 $(x^2y,xy^2+x^3)$ in the polynomial ring in
 four variables we proceed as follows:
 Our favorite field
  
    | i1 : KK = ZZ/31991
 o1 = KK
 
 o1 : QuotientRing
 | 
 The polynomial ring
  
    | i2 : R = KK[x,y,z,w]
 o2 = R
 
 o2 : PolynomialRing
 | 
 and the ideal
  
    | i3 : I = ideal(x^2*y,x*y^2+x^3)
 2    3      2
 o3 = ideal (x y, x  + x*y )
 
 o3 : Ideal of R
 | 
 now the punch line:
  
    | i4 : J = gens gb I
 o4 = {0} | x2y x3+xy2 xy3 |
 
 1       3
 o4 : Matrix R  <--- R
 | 
 From this we can for example compute the
 codimension, dimension,
 degree, and the whole Hilbert
 function and polynomial.  
 This will be more fun if we work with an
 example having some meaning.  We choose
 to work with the ideal defining the
 rational quartic curve in $\P^3$ given
 parametrically in an affine representation
 by 
        $$t \mapsto (t,t^3,t^4).$$
 (The reader who doesn't understand this
 terminology may ignore it for the moment,
 and treat the ideal given below as a 
 gift from the gods... .)
 We obtain the ideal by first making the 
 polynomial ring in 4 variables (the
 homogeneous coordinate ring of $\P^3$)
  
    | i5 : R = KK[a..d]
 o5 = R
 
 o5 : PolynomialRing
 | 
 and then using a function {\tt monomialCurve}, which we shall
 treat for now as a black box
  
    | i6 : I = monomialCurve(R,{1,3,4})
 3      2     2    2    3    2
 o6 = ideal (b*c - a*d, c  - b*d , a*c  - b d, b  - a c)
 
 o6 : Ideal of R
 | 
 From Macaulay's point of view, $I$ is an
 ideal, and the codimension of its support
 is 2, while its dimension is 2:
 This is the codimension of $R/I$ in $R$ 
 and the dimension of $R/I$.  We could work with
 the module $R/I$ as well.
 Precision requires writing $R^1$ instead
 of $R$ ($R$ is a ring, and $R^1$ is
 the free module of rank 1 over it)
  
    | i9 : codim (R^1/(I*R^1))
 o9 = 2
 | 
 We could also extract the generators of
 $I$ (as a matrix) and take the cokernel to
 get the same thing:
  
    | i10 : M = coker gens I
 o10 = cokernel {0} | bc-ad c3-bd2 ac2-b2d b3-a2c |
 
 1
 o10 : R - module, quotient of R
 | 
 And similarly for the degree:
 As one might expect, the degree of the quartic
 is 4 !
 
 The Hilbert polynomial is obtained by
  
    | i15 : hilbertPolynomial M
 o15 = - 3*P  + 4*P
 0      1
 
 o15 : ProjectiveHilbertPolynomial
 | 
 The term $\P^i$ represents the Hilbert polynomial of
 projective $i$-space.  This formula tells
 us that the Hilbert polynomial of $M$ is
 $H(m) = 4(m+1) - 3 = 4m + 1$.  Thus the degree
 is four, the dimension of the projective variety
 which is the support of $M$ is 1 (and so the affine
 dimension is 2),
 and that the (arithmetic) genus is 0 (1 minus the
 constant term, for curves).
 The Hilbert series of $M$ (the generating function
 for the dimensions of the graded pieces of $M$) is
  
    | i16 : hilbertSeries M
 3      2
 - $T  + 2$T  + 2$T + 1
 o16 = ----------------------
 2
 (- $T + 1)
 
 o16 : Divide
 | 
 The indeterminate in this expression is \$T.
 Another way to get information about
 the module $M$ is to see its free resolution
  
    | i17 : Mres = res M
 1      4      4      1
 o17 = R  <-- R  <-- R  <-- R
 
 0      1      2      3
 
 o17 : ChainComplex
 | 
 To get more precise information about Mres,
 we could do
  
    | i18 : betti Mres
 o18 = total: 1 4 4 1
 0: 1 . . .
 1: . 1 . .
 2: . 3 4 1
 
 o18 : Net
 | 
 The display is chosen for compactness:
 the first line gives the total betti 
 numbers, the same information given when
 we type the resolution.  The remaining
 lines express the degrees of each of the
 generators of the free modules in the
 resolution.  The $j$th column after the colons
 gives the degrees of generators of the
 $j$th module(counting from $0$); 
 an $n$ in the $j$th column in the
 row headed by ``$d$:'' means that the $j$th
 free module has $n$ generators of degree
 $n+j$.  Thus for example in our case, the
 generator of the third (last) free module in the
 resolution has degree $3+2=5$.
\beginsection{ B. Random regular sequences}\par
 An interesting and illustrative open problem
 is to understand the initial ideal (and 
 the Groebner basis) of a ``generic'' 
 regular sequence.  To study a very simple case
 we take a matrix of 2 random forms 
 in a polynomial ring in
 3 variables:
  
    | i19 : R = KK[x,y,z]
 o19 = R
 
 o19 : PolynomialRing
 | 
  
    | i20 : F = random(R^1, R^{-2,-3})
 o20 = {0} | -8514x2+5827xy-9613y2+7886xz+5702yz+11835z2 -15791x3+12712x2y-11866xy2-11819y3-2169x2z+1548xyz+10647y2z+11189xz2+6279yz2+8360z3 |
 
 1       2
 o20 : Matrix R  <--- R
 | 
 makes $F$ into a $1 \times 2$ matrix whose elements
 have degrees $2,3$ (that is, $F$ is a random map
 to the free module $R^1$, which has its one
 generator in the (default) degree, $0$, from
 the free module with generators in the listed
 degrees, $\{2,3\}$).  We now can compute
  
    | i21 : GB = gens gb F
 o21 = {0} | x2-14140xy+5100y2+12053xz+14443yz+7878z2 xy2-6589y3-2097xyz-8415y2z-5942xz2+4719yz2-9862z3 y4+2957y3z-7768xyz2+13545y2z2+7229xz3+5489yz3-4313z4 |
 
 1       3
 o21 : Matrix R  <--- R
 | 
  
    | i22 : LT = leadTerm gens gb F
 o22 = {0} | x2 xy2 y4 |
 
 1       3
 o22 : Matrix R  <--- R
 | 
  
    | i23 : betti LT
 o23 = total: 1 3
 0: 1 .
 1: . 1
 2: . 1
 3: . 1
 
 o23 : Net
 | 
 shows that there are Groebner basis elements
 of degrees 2,3, and 4.  This result is
 dependent on the monomial order in the ring $R$;
 for example we could take the lexicographic
 order
  
    | i24 : R = KK[x,y,z, MonomialOrder => Lex]
 o24 = R
 
 o24 : PolynomialRing
 | 
 (see {\tt help MonomialOrder} for other possibilities).
 We get 
  
    | i25 : F = random(R^1, R^{-2,-3})
 o25 = {0} | -9394x2+8312xy+12075xz-3147y2-11542yz+3907z2 -13293x3-1763x2y+10247x2z-8437xy2+2153xyz-4835xz2+2447y3-1786y2z+12876yz2+12592z3 |
 
 1       2
 o25 : Matrix R  <--- R
 | 
  
    | i26 : GB = gens gb F
 o26 = {0} | x2+2417xy+2168xz+13360y2-639yz+14449z2 xy2+198xyz-4686xz2+1906y3-585y2z+12211yz2+3780z3 xyz2-5641xz3+3811y4-12973y3z-2454y2z2+11675yz3+3823z4 xz4-3470y5-8776y4z-8700y3z2+1864y2z3-14123yz4-5748z5 y6-8183y5z+7264y4z2-9842y3z3-14142y2z4-12444yz5-15013z6 |
 
 1       5
 o26 : Matrix R  <--- R
 | 
  
    | i27 : LT = leadTerm gens gb F
 o27 = {0} | x2 xy2 xyz2 xz4 y6 |
 
 1       5
 o27 : Matrix R  <--- R
 | 
  
    | i28 : betti LT
 o28 = total: 1 5
 0: 1 .
 1: . 1
 2: . 1
 3: . 1
 4: . 1
 5: . 1
 
 o28 : Net
 | 
 and there are Groebner basis elements of degrees 
 $2,3,4,5,6.$
\beginsection{ C. Division With Remainder}\par
 A major application of Groebner bases is
 to decide whether an element is in a given
 ideal, and whether two elements reduce to
 the same thing modulo an ideal.  For
 example, everyone knows that the trace
 of a nilpotent matrix is 0. We can produce
 an ideal $I$ that defines the variety $X$ of 
 nilpotent $3 \times 3$ matrices by taking the cube
 of a generic matrix and setting the entries
 equal to zero.  Here's how:
  
    | i29 : R = KK[a..i]
 o29 = R
 
 o29 : PolynomialRing
 | 
  
    | i30 : M = genericMatrix(R,a,3,3)
 o30 = {0} | a d g |
 {0} | b e h |
 {0} | c f i |
 
 3       3
 o30 : Matrix R  <--- R
 | 
  
    | i31 : N = M^3
 o31 = {0} | a3+2abd+bde+2acg+bfg+cdh+cgi        a2d+bd2+ade+de2+cdg+afg+efg+dfh+fgi a2g+bdg+cg2+adh+deh+fgh+agi+dhi+gi2 |
 {0} | a2b+b2d+abe+be2+bcg+ach+ceh+bfh+chi abd+2bde+e3+bfg+cdh+2efh+fhi        abg+beg+bdh+e2h+cgh+fh2+bgi+ehi+hi2 |
 {0} | a2c+bcd+abf+bef+c2g+cfh+aci+bfi+ci2 acd+cde+bdf+e2f+cfg+f2h+cdi+efi+fi2 acg+bfg+cdh+efh+2cgi+2fhi+i3        |
 
 3       3
 o31 : Matrix R  <--- R
 | 
  
    | i32 : I = flatten N
 o32 = {0} | a3+2abd+bde+2acg+bfg+cdh+cgi a2b+b2d+abe+be2+bcg+ach+ceh+bfh+chi a2c+bcd+abf+bef+c2g+cfh+aci+bfi+ci2 a2d+bd2+ade+de2+cdg+afg+efg+dfh+fgi abd+2bde+e3+bfg+cdh+2efh+fhi acd+cde+bdf+e2f+cfg+f2h+cdi+efi+fi2 a2g+bdg+cg2+adh+deh+fgh+agi+dhi+gi2 abg+beg+bdh+e2h+cgh+fh2+bgi+ehi+hi2 acg+bfg+cdh+efh+2cgi+2fhi+i3 |
 
 1       9
 o32 : Matrix R  <--- R
 | 
 (actually this produces a 1 x 9 matrix of
 of forms, not the ideal: {\tt J = ideal I};
 the matrix will be more useful to us).
 But the trace is not in $I$!  This is obvious
 from the fact that the trace has degree $1$,
 but the polynomials in $I$ are of degree $3$.
 We could also check by division with
 remainder:
  
    | i33 : Tr = trace M 
 o33 = a + e + i
 
 o33 : R
 | 
  
    | i34 : Tr //I  -- the quotient, which is 0
 o34 = 0
 
 9       1
 o34 : Matrix R  <--- R
 | 
  
    | i35 : Tr % I  -- the remainder, which is Tr again
 o35 = a + e + i
 
 o35 : R
 | 
 (Here {\tt Tr} is an element of $R$, not a matrix.
 We could do the same thing with a $1 \times 1$ matrix
 with {\tt Tr} as its element.)
 This is of course because the entries of $I$ do
 NOT
 generate the ideal of all forms 
 vanishing on $X$ -- this we may find with
 {\tt J = radical ideal I},
 (but this takes a while: see the documentation for
 {\tt radical} on a faster way to find this)
 which shows that the radical is generated by
 the trace, the determinant, and the sum of 
 the principal $2 \times 2$ minors, that is, by the
 coefficients of the characteristic polynomial.
 In particular, we can try the powers of the
 radical:
  
    | i36 : Tr^2 % I
 2           2                  2
 o36 = a  + 2a*e + e  + 2a*i + 2e*i + i
 
 o36 : R
 | 
  
    | i37 : Tr^3 % I
 2                 2     3                                2               2                          2       2     3
 o37 = 3a e + 3b*d*e + 3a*e  + 3e  + 3b*f*g + 3c*d*h + 6e*f*h + 3a i + 6a*e*i + 3e i + 3c*g*i + 6f*h*i + 3a*i  + 3e*i  + 3i
 
 o37 : R
 | 
  
    | i38 : Tr^4 % I
 2 2         2       3     4                                                  2                   2 2      2                       2       3                                                               2 2          2     2 2         2          2       3        3     4
 o38 = 6a e  + 6b*d*e  + 6a*e  + 6e  + 6b*e*f*g + 6c*d*e*h - 6b*d*f*h + 6a*e*f*h + 12e f*h - 6c*f*g*h - 6f h  + 12a e*i + 12b*d*e*i + 12a*e i + 12e i + 12c*e*g*i + 6b*f*g*i + 6c*d*h*i + 6a*f*h*i + 36e*f*h*i + 6a i  + 12a*e*i  + 6e i  + 6c*g*i  + 12f*h*i  + 6a*i  + 12e*i  + 6i
 
 o38 : R
 | 
  
    | i39 : Tr^5 % I
 2 2           2         3       4         2         2                                       2                         2 2       2   2            2        2 2      3 2            2            2             2      2 3          3          3      2 3          3          3        4        4      5
 o39 = 30a e i + 30b*d*e i + 30a*e i + 30e i + 30c*e g*i + 30a f*h*i + 30b*d*f*h*i + 60a*e*f*h*i + 90e f*h*i + 30c*f*g*h*i + 30f h i + 30a e*i  + 30b*d*e*i  + 30a*e i  + 30e i  + 30c*e*g*i  + 60a*f*h*i  + 120e*f*h*i  + 30a i  + 30b*d*i  + 30a*e*i  + 30e i  + 30c*g*i  + 90f*h*i  + 30a*i  + 30e*i  + 30i
 
 o39 : R
 | 
  
    | i40 : Tr^6 % I
 2 2 2          2 2        3 2      4 2        2   2      2     2              2               2       2     2              2      2 2 2      2   3            3        2 3      3 3            3             3             3      2 4          4          4      2 4          4           4        5        5      6
 o40 = 90a e i  + 90b*d*e i  + 90a*e i  + 90e i  + 90c*e g*i  + 90a f*h*i  + 90b*d*f*h*i  + 180a*e*f*h*i  + 270e f*h*i  + 90c*f*g*h*i  + 90f h i  + 90a e*i  + 90b*d*e*i  + 90a*e i  + 90e i  + 90c*e*g*i  + 180a*f*h*i  + 360e*f*h*i  + 90a i  + 90b*d*i  + 90a*e*i  + 90e i  + 90c*g*i  + 270f*h*i  + 90a*i  + 90e*i  + 90i
 
 o40 : R
 | 
  
    | i41 : Tr^7 % I
 o41 = 0
 
 o41 : R
 | 
 The seventh power is the first one in the 
 ideal!  (Bernard Mourrain has worked out a
 formula for which power in general.)
 In this case
  
    | i42 : Tr^6 // I
 o42 = {3} | a3+6a2e+3bde+15ae2+22e3+3bfg+3cdh-12efh+6a2i+30aei+60e2i+3cgi-30fhi+15ai2+60ei2+4i3                                       |
 {3} | 18de2+18efg+18dei+18fgi                                                                                                   |
 {3} | 18deh+18egi+18dhi+18gi2                                                                                                   |
 {3} | -27abe-63be2-36ach+9ceh+18bfh-72abi-162bei+90chi                                                                          |
 {3} | -2a3+15a2e+21bde+6ae2+e3-6bfg-6cdh-36afh+15efh+60a2i+72bdi+30aei+6e2i+66cgi+141fhi-30ai2-75ei2-26i3                       |
 {3} | 63beg+36a2h+18bdh-18aeh+18e2h+9cgh+18fh2+117bgi-162ahi-27ehi+45hi2                                                        |
 {3} | -18a2c+18ace+18abf-18bef+36cfh+45aci+18cei-81bfi+153ci2                                                                   |
 {3} | -18cde-36a2f-36bdf+72aef+9e2f-45cfg-18f2h-108cdi+162afi+9efi+54fi2                                                        |
 {3} | 16a3+18abd-30a2e-42bde-30ae2-44e3+18acg-18ceg+3bfg+3cdh+18afh-57efh-75a2i-90bdi-60aei-75e2i-87cgi-165fhi-84ai2-84ei2-89i3 |
 
 9       1
 o42 : Matrix R  <--- R
 | 
 is not 0.  It is a matrix that makes the
 following true:
  
    | i43 : Tr^6 == I * (Tr^6 // I) + (Tr^6 % I)
 o43 = true
 | 
\beginsection{ D. Elimination Theory}\par
 Consider the problem of projecting the
 ``twisted cubic'', a curve in $\P^3$ defined
 by the three $2 \times 2$ minors of a certain
 $2 \times 3$ matrix into the plane.  
 Such problems can be solved in a very 
 simple and direct way using Groebner bases.
 The technique lends itself to many extensions,
 and in its developed form can be used to find
 the closure of the image of any map of 
 affine varieties.  
    In this section we shall first give a 
 simple direct treatment of the problem above,
 and then show how to use Macaulay2's 
 general tool to solve the problem.
 We first
 clear the earlier meaning of {\tt x} to make it
 into a subscripted variable
  
    | i44 : erase global x
 o44 = x
 
 o44 : Symbol
 | 
 and then set
  
    | i45 : R = KK[x_0..x_3] 
 o45 = R
 
 o45 : PolynomialRing
 | 
 the homogeneous coordinate ring of $\P^3$
 and
  
    | i46 : M = map(R^2, 3, (i,j)->x_(i+j))
 o46 = {0} | x_0 x_1 x_2 |
 {0} | x_1 x_2 x_3 |
 
 2       3
 o46 : Matrix R  <--- R
 | 
  
    | i47 : I = gens minors(2,M)
 o47 = {0} | -x_1^2+x_0x_2 -x_1x_2+x_0x_3 -x_2^2+x_1x_3 |
 
 1       3
 o47 : Matrix R  <--- R
 | 
 a matrix whose image is 
 the ideal of the twisted cubic.
 As projection center we
  take the point defined by
  
    | i48 : pideal = ideal(x_0+x_3, x_1, x_2)
 o48 = ideal (x  + x , x , x )
 0    3   1   2
 
 o48 : Ideal of R
 | 
 To find the image we must intersect the ideal
 $I$ with the subring generated by the 
 generators of {\tt pideal}.  We make a change of
 variable so that these generators become
 the last three variables in the ring; that
 is, we write the ring as $KK[y_0..y_3]$
 where 
  $$y_0 = x_0, y_1 = x_1, y_2 = x_2, y_3 = x_0+x_3$$
 and thus
 $x_3 = y_3-y_0$, etc. 
 We want the new ring to have an ``elimination
 order'' for the first variable.
  
    | i49 : erase global y
 o49 = y
 
 o49 : Symbol
 | 
  
    | i50 : S = KK[y_0..y_3,MonomialOrder=> Eliminate 1]
 o50 = S
 
 o50 : PolynomialRing
 | 
 Here is one way to make the substitution
  
    | i51 : I1 = substitute(I, matrix{{y_0,y_1,y_2,y_3-y_0}})
 o51 = {0} | y_0y_2-y_1^2 -y_0^2+y_0y_3-y_1y_2 -y_0y_1-y_2^2+y_1y_3 |
 
 1       3
 o51 : Matrix S  <--- S
 | 
 The elimination of 1 variable from the 
 matrix of Groebner basis elements proceeds
 as follows:
  
    | i52 : J = selectInSubring(1,gens gb I1)
 o52 = {0} | y_1^3+y_2^3-y_1y_2y_3 |
 
 1       1
 o52 : Matrix S  <--- S
 | 
 and gives (a matrix with element)
 the cubic equation of a rational
 curve with one double point in the plane.
 However, we are still in a ring with 4 
 variables, so if we really want a plane
 curve (and not the cone over one) we must
 move to yet another ring:
  
    | i53 : S1 = KK[y_1..y_3]
 o53 = S1
 
 o53 : PolynomialRing
 | 
  
    | i54 : J1 = substitute(J, S1)
 o54 = {0} | y_1^3+y_2^3-y_1y_2y_3 |
 
 1        1
 o54 : Matrix S1  <--- S1
 | 
 This time we didn't have to give so much
 detail to the {\tt substitute} command because of 
 the coincidence of the names of the variables.
 Having shown the primitive method, we
 now show a much more flexible and transparent
 one:  we set up a ring map from the polynomial
 ring in $3$ variables (representing the plane)
 to $R/I$, taking the variables $y$ to the three
 linear forms that define the projection 
 center.  Then we just take the kernel of
 this map!  (``Under the hood'', 
 Macaulay2 is doing a more refined version
 of the same computation as before.)
 Here is the ring map
  
    | i55 : Rbar = R/(ideal I)
 o55 = Rbar
 
 o55 : QuotientRing
 | 
  
    | i56 : f = map(Rbar, S1, matrix(Rbar,{{x_0+x_3, x_1,x_2}})
 )
 
 o56 = map(Rbar,S1,{x  + x , x , x })
 0    3   1   2
 
 o56 : RingMap Rbar <--- S1
 | 
 and the desired ideal
  
    | i57 : J1 = ker f
 3             3
 o57 = ideal(y  - y y y  + y )
 2    1 2 3    3
 
 o57 : Ideal of S1
 | 
\beginsection{ E. Quotients and saturation}\par
 Another typical application of 
 Groebner bases and syzygies is to the
 computation of ideal quotients and 
 saturations.  Again we give an easy example
 that we can treat directly, and then 
 introduce the tool used in Macaulay2 to 
 treat the general case.
 If $I$ and $J$ are ideals in a ring $R$, we define
 $(I:J)$, the ideal quotient, by 
   $$(I:J) = \{f \in R \mid fJ \subset I\}.$$
 In our first examples we consider 
 the case where $J$ is 
 generated by a single element $g$.
 This arises in practice, for example, in the
 problem of homogenizing an ideal.  Suppose
 we consider the affine space curve
 parametrized by
 $t \mapsto (t,t^2,t^3)$.  The ideal of polynomials
 vanishing on the curve is easily seen to
 be $(b-a^2, c-a^3)$ (where we have taken
 $a,b,c$ as the coordinates of affine space).
 To find the projective closure of the curve
 in $\P^3$, we must homogenize these equations 
 with respect to a new variable d, getting
 $dc-a^2, d^2c-a^3$.  But these forms do NOT
 define the projective closure! In general,
 homogenizing the generators of the ideal $I$ of
 an affine variety one gets an ideal $I_1$ that
 defines the projective closure UP TO
 a component supported on the hyperplane
 at infinity (the hyperplane $d=0$).  To see
 the ideal of the closure we must remove
 any such components, for example by
 replacing $I_1$ by the union $I_2$ of all the
 ideals $(I_1:d^n)$, where $n$ ranges over positive
 integers.  This is not so hard as it seems:
 First of all, we can successively compute
 the increasing sequence of ideals
 $(I_1:d), (I_1:d^2), \ldots $ until we get two 
 that are the same; all succeeding ones
 will be equal, so we have found the union.
 A second method involves a special property
 of the reverse lex order, and is much more
 efficient in this case.  We shall illustrate
 both. First we set up the example above:
  
    | i58 : R = KK[a,b,c,d]
 o58 = R
 
 o58 : PolynomialRing
 | 
  
    | i59 : I1 = ideal(d*c-a^2, d^2*c-a^3)
 2           3      2
 o59 = ideal (- a  + c*d, - a  + c*d )
 
 o59 : Ideal of R
 | 
 How to compute the ideal quotient:
 If $I$ is generated by $f_1,\ldots,f_n$, we see that
 $s \in (I:g)$ iff there are ring elements 
 $r_i$ such that 
        $$\sum_{i=1}^{n} r_i f_i + s g = 0. $$
 Thus it suffices to compute the kernel
 (syzygies) of the $1 \times (n+1)$ matrix
           $$(f_1, ... ,f_n, g)$$
 and collect the coefficients of $g$, that is,
 the entries of the last row of a matrix
 representing the kernel. 
 Thus in our case we may compute $(I_1:d)$
 by concatenating the matrix for $I_1$
 with the single variable $d$
  
    | i60 : I1aug = (gens I1) | matrix{{d}}
 o60 = {0} | -a2+cd -a3+cd2 d |
 
 1       3
 o60 : Matrix R  <--- R
 | 
  
    | i61 : augrelations = gens ker I1aug
 o61 = {2} | -a    d     |
 {3} | 1     0     |
 {1} | ac-cd a2-cd |
 
 3       2
 o61 : Matrix R  <--- R
 | 
 There are 3 rows (numbered 0,1,2 !) and
 2 columns, so to extract the last row we
 may do
  
    | i62 : I21 = submatrix(augrelations, {2}, {0,1})
 o62 = {1} | ac-cd a2-cd |
 
 1       2
 o62 : Matrix R  <--- R
 | 
 But this is not an ``ideal'' properly speaking:
 first of all, it is a matrix, not an ideal,
 and second of all its target is not $R^1$
 but $R(-1)$, the free module of rank 1 with
 generator in degree 1.  We can fix both
 of these problems by
  
    | i63 : I21 = ideal I21
 2
 o63 = ideal (a*c - c*d, a  - c*d)
 
 o63 : Ideal of R
 | 
 This is larger than the original ideal, having
 two quadratic generators instead of a 
 quadric and a cubic, so
 we continue.  Instead of doing the same
 computation again, we introduce the built-in
 command
  
    | i64 : I22 = I21 : d
 2                    2
 o64 = ideal (c  - c*d, a*c - c*d, a  - c*d)
 
 o64 : Ideal of R
 | 
 which is again larger than {\tt I21}, having
 three quadratic generators. Repeating,
  
    | i65 : I23 = I22 : d
 2                    2
 o65 = ideal (c  - c*d, a*c - c*d, a  - c*d)
 
 o65 : Ideal of R
 | 
 we get another ideal with three quadratic
 generators.  It must be the same as {\tt I21},
 but the generators are written differently
 because of the route taken to get it, so
 (being suspicious) we might check with
  
    | i66 : (gens I23) % (gens I22)
 o66 = 0
 
 1       3
 o66 : Matrix R  <--- R
 | 
 which returns 0, showing that {\tt I23} is 
 contained in (gives remainder 0 when divided
 by) {\tt I22}.  Thus the homogeneous ideal {\tt I2} of 
 the projective closure is equal to {\tt I23}
 (this is the homogeneous ideal of 
 the twisted cubic, already encountered above).
 A more perspicuous way of approaching the
 computation of the union of the $(I:d^n)$,
 which is called the saturation of $I$ with
 respect to $d$, and written $(I:d^\infty)$,
 is first to compute a reverse lex Groebner
 basis.
  
    | i67 : gens gb I1
 o67 = {0} | a2-cd acd-cd2 c2d2-cd3 |
 
 1       3
 o67 : Matrix R  <--- R
 | 
 This yields {\tt (a2-cd, acd-cd2, c2d2-cd3)},
 meaning 
  $$(a^2-cd, acd-cd^2, c^2d^2-cd^3).$$
 We see that the second generator is divisible
 by $d$, and the third is divisible by $d^2$.
 General theory says that we get the right
 answer simply by making these divisions,
 that is, the saturation is
  $$(a^2-cd, ac-cd, c^2-cd),$$
 as previously computed.  The same thing
 can be accomplished in one line by
  
    | i68 : I2 = divideByVariable(gens gb I1,d)
 o68 = {0} | a2-cd ac-cd c2-cd |
 
 1       3
 o68 : Matrix R  <--- R
 | 
 This saturation may be found directly in Macaulay2:
  
    | i69 : saturate(I1, d)
 2                    2
 o69 = ideal (c  - c*d, a*c - c*d, a  - c*d)
 
 o69 : Ideal of R
 | 




