-- -*- Mode: M2; mode: auto-fill; coding: utf-8; fill-column: 100 -*-
Key
     "BeginningMacaulay2"
Headline
     Mathematicians' Introduction to  Macaulay2
Description
     Text
	  We assume you've installed {\em Macaulay2} and can type
     Code
          EXAMPLE { PRE ///M2/// }
     Text
	  on a command line to bring up the program. You should see something like:
     Code
     	  EXAMPLE {
	      PRE concatenate("Macaulay2, version ", version#"VERSION",
		  newline, "Type `help` to see useful commands")}
     Text
          We suggest you do that now, so that you can experiment while you read
	  this tutorial!


     Code
     	 SUBSECTION "Arithmetic with integers, rings and ideals"
     Text
	  You can immediately do arithmetic with integers:
     Example
	  2+2
	  107*431
	  25!
	  binomial(5,4)
	  factor 32004
     Text

	  Most {\em Macaulay2} applications involve polynomial rings over fields
	  and their quotient rings. Fields can be made in various ways:
     Example
	  ZZ/101
	  QQ
	  GF 2^5
	  k = toField (QQ[i]/(i^2+1))
     Text
	  After making {\tt k} we can compute in it:
     Example
	  1/i
     Text

	  Computation is often fastest and needs least
	  memory when performed over finite prime fields of the form 
	  $\ZZ/p$.
	  Fortunately, when the characteristic $p$ is not too small,
	  qualitative questions often have similar answers over
	  $\ZZ/p$ and over $\QQ$, so we mostly use the former.
	  In {\em Macaulay2} the prime $p$
	  can range up to 32749.

	  We make a polynomial ring in 5 variables over $\ZZ/101$:
     Example
	  kk=ZZ/101
	  S=kk[x_1..x_5]
     Text
	  Here is another way:
     Example
	  S=kk[a,b,c,d,e] 
     Text
	  One can do arithmetic on polynomials:
     Example
	  (3*a^2+1)^5
     Text

	  We make an ideal in $S$:
     Example
	  I=ideal(a^3-b^3, a+b+c+d+e)
     Text
	  Using this ideal, we can make a factor ring:
     Example
	  R=S/I
     Text
	  Another way to make an ideal, with more compact notation (familiar to anyone who used the
	  classic Macaulay) is:
     Example
	  use S
	  I=ideal"3(a+b)3, 4c5"
     Text
	  Note the command ``{\tt use S}'', which specifies 
	  that we want to work with the generators of the polynomial ring S again;
	  otherwise the variables {\tt a}, {\tt b}, and {\tt c}
	  would still have had values in $R$ instead of in $S$.

	  Algebraic operations on ideals are available:
     Example
	  I^2
	  I*I
	  I+ideal"a2"
     Text

	  In case you forget any of these things, @TO help@ is available! The most
	  useful way to get it is often to type something like:
     Code
     	  EXAMPLE { PRE ///viewHelp ideal/// }
     Text
     	  Then a browser window will pop up that contains documentation about the function @ TO
          ideal @ that we've been using; links on that page allow one to explore all of the {\em Macaulay2} documentation.

	  On the other hand, we might have wanted information about the @TO class@ of all ideals.
	  Not too surprisingly, this class is called @TO Ideal@. We could get information about
	  what functions create or use ideals by typing:
     Code                                                                                                                                                                                                
          EXAMPLE { PRE ///viewHelp Ideal/// }                                                                                                                                                           
     Text                                                                                                                                                                                                
          To see the names of classes, you can begin by looking at the output
	  of commands; the second line output (the one introduced by a colon) often contains the name of the 
	  class of the result.
	  
	  Here are some basic operations on matrices:
     Example
	  M= matrix{{a,b,c},{b,c,d},{c,d,e}}
	  M^2
	  determinant M
	  trace M
	  M-transpose M
     Text
	  The function @TO entries@ gives the entries of a matrix:
     Example
	  entries M
     Text
	  The result is a list of lists, one for each row of the matrix $M$.
	  The function @TO flatten@ can be used to merge the
	  lists into a single list:
     Example
	  flatten entries M
     Text
	  If you want a particular entry, say the one in the upper left corner, 
	  you can use the underscore operator @TO2 {(symbol _,Matrix,Sequence),"_"}@:
     Example
	  M_(0,0)
     Text
	  Here, as everywhere in {\em Macaulay2}, all indexing starts with 0.
	  For example:
     Example
	  I_0
     Text
	  is the first generator of I. You can list all the generators with:
     Example
	  I_*
     Text

	  A {\em module} can be defined as a cokernel, kernel, image, or even as a subquotient:
     Example
	  coker M     
	  image M
	  kernel matrix"a,b,0;0,a,b"
	  N = matrix{{a,b},{b,c},{c,d}}
	  (image M)/(image N)
	  subquotient(M,N)
     Text
	 Note that the matrix $N$ above was defined with an
	 alternate syntax, parallel to the alternate syntax for @TO ideal@.

	  Before going on, the reader might want to explore a bit. A good place to 
	  start is the top of the documentation tree, which can be reached, for
	  example, by typing:
     Code
     	  EXAMPLE { PRE "viewHelp Macaulay2Doc" }
     Text
          

     Code
     	  SUBSECTION "Properties of ideals and modules",
     Text
	  To compute the Gröbner basis of an ideal
	  $(x^2y,xy^2+x^3)$ in the polynomial ring in
	  four variables we proceed as follows.
	  First we make our favorite field:
     Example
	  kk = ZZ/32003
     Text
	  Then the polynomial ring:
     Example
	  R = kk[x,y,z,w]
     Text
	  And then the ideal:
     Example
	  I = ideal(x^2*y,x*y^2+x^3)
     Text
	  Now the punch line.  We compute the Gröbner basis with the @ TO groebnerBasis @ function:
     Example
	  J = groebnerBasis I
     Text

	  Gr\"obner bases are always computed with respect to a particular
	  monomial order on the ring. In fact, the ring we defined above has
	  a default monomial order, the graded reverse lex order. For many
	  other possibilities, see @ TO MonomialOrder @, or type:
     Code
     	  EXAMPLE { PRE "viewHelp MonomialOrder" }
     Text

	  The analogue of factorization in the theory of ideals
	  is primary decomposition.
	  For example, we can begin by intersecting three ideals:
     Example
	  I= intersect (ideal"x2,y3", ideal"y2,z3", (ideal"x,y,z")^4)
     Text
	  We can almost undo this operation by computing
	  a primary decomposition:
     Example
	  primaryDecomposition I
     Text
          Inspecting the output, we see that the first two ideals
	  are the same as the first two ideals we intersected, but the 
	  third one differs from the corresponding input ideal.
	  This is because only the primary components corresponding
	  to minimal primes (here, the first two) are unique. All three of the input ideals
	  are primary, so they constitute a primary decomposition of $I$
	  different from the one provided by {\em Macaulay2} on the output line.
	  
	  For larger examples, primary decomposition is computationally challenging!
	  Sometimes it is easier to compute just the minimal primes. To do
	  this we can use @ TO decompose @: 
     Example
	  decompose I
     Text

	  Using Gröbner bases we can compute 
	  codimensions, dimensions,
	  degrees, Hilbert
	  functions, and Hilbert polynomials.  
	  This will be more fun if we work with a
	  meaningful example.  We will use
	  the ideal defining the smooth
	  rational quartic curve in $\PP^3$ given
	  parametrically (in an affine representation)
	  by $$t \mapsto{} (t,t^3,t^4).$$
	  (The reader more interested in algebra than geometry
	  may simply treat the ideal given below as a 
	  gift from the gods.)
	  First we make the
	  polynomial ring in 4 variables, to serve as the
	  homogeneous coordinate ring of $\PP^3$:
     Example
	  R = kk[a..d]
     Text
	  We introduce the ring map $\phi: R \to kk[s,t]$ defined by 
     	  $(a,b,c,d) \mapsto{} (s^4, s^3 t, s t^3, t^4)$:	  
     Example
	  phi = map(kk[s,t],R,{s^4, s^3*t, s*t^3, t^4})
     Text
	  Here the syntax of the function @ TO2 {(map,Ring,Ring,List),"map"} @ has the target ring first and the source ring second:
	  maps in {\em Macaulay2} generally go from right to left!
	  The last input to the command is a 
	  list of the elements to which to send the variables of the source ring.
	  The ideal we want is the kernel of this map:
     Example
	  I = ker phi
     Text
	  Shortcut notation for this construction is provided by the function @ TO
	  monomialCurveIdeal @:
     Example
	  I = monomialCurveIdeal(R,{1,3,4})
     Text
	  We can compute the @ TO2 {(dim,Ideal),"dimension"}@, @ TO2{(codim,Ideal), "codimension"}@ (also called the
	  height) and @ TO2{(degree,Ideal),"degree"} @ of this ideal:
     Example
	  dim I
	  codim I
	  degree I
     Text
	  The Hilbert polynomial is obtained with the function @ TO hilbertPolynomial@:
     Example
	  hilbertPolynomial(R/I)
     Text
	  The output above may not be what the user expected:
	  the term ${\mathbf P}_m$ represents the Hilbert polynomial of
	  projective $m$-space.  Thus the output tells
	  us that the Hilbert polynomial of $M$ is
	  $i \mapsto{} -3*1+4*(i+1) = 4i + 1$.  Thus the degree
	  is four, the dimension of the projective variety
	  that is the support of $M$ is 1 (and so the affine
	  dimension is 2), and the (arithmetic) genus is 0 (obtained as 1 minus the
	  constant term of the polynomial.)

	  The more usual expression for the Hilbert polynomial can
	  be obtained as follows:
     Example
	  hilbertPolynomial(R/I, Projective => false)
     Text
	  The construction {\tt Projective => false} is our first example of
	  an {\em option} to a function: we specified that the option 
	  {\tt Projective} was to have the value {\tt false}.
	  The form we used first could also have been written this way:
     Example
	  hilbertPolynomial(R/I, Projective => true)
     Text
	  The Hilbert series of $M$ (the generating function
	  for the dimensions of the graded pieces of $M$) is
     	  obtained with:
     Example
	  hilbertSeries (R/I)
     Text
	  This generating function is expressed
	  as a rational function with denominator equal to (1-T)^n, where
	  n is the number of variables in R. 
	  Since R/I has dimension 2, it can also be written
	  with denominator (1-t)^2. To see it in this form, use @ TO reduceHilbert @:
     Example
	  reduceHilbert hilbertSeries (R/I)
     Text
	  It is possible to manipulate the numerator and denominator of this
	  expression. To learn how to do so, see @ TO hilbertSeries @ or type:
     Code
     	  EXAMPLE { PRE "viewHelp hilbertSeries" }
     Text

	  A great deal of subtle information about a module is visible using
	  free resolutions. For an example, we begin
	  by turning $R/I$ into a module. Here the code @ TT "R^1" @ produces the free module of
	  rank 1 over $R$, and @ TO2 {"OldChainComplexes :: resolution(Module)", "res"} @ computes a free resolution:
     Example
	  M=R^1/I
	  Mres = res M
     Text
	  To get more precise information about {\tt Mres},
	  we could compute its Betti table with @TO betti@:
     Example
	  betti Mres
     Text
	  The display is chosen for compactness. Each column of the
	  table corresponds
	  to a free module in the resolution. The column's heading
	  specifies the {\em homological degree} (the position of the free
	  module in the resolution).
	  The entry just below the homological degree
	  is the rank of the free module, also called the
	  {\em total betti 
	  number}. The remaining entries in the column
          tell us how many generators of each degree this free
	  module has: the number in the column labelled $j$ and in the row labelled $d$
	  tells how many generators of degree $j+d$ the $j$-th free module has.
	  Thus, in our case, the single
	  generator of the third (and last) free module in the
	  resolution has degree $3+2=5$.

	  Commonly computed homological invariants
	  such as @ TO2 { pdim, "projective dimension" } @
	  and @ TO2 { regularity, "regularity"} @ are (also) available directly:
     Example
	  pdim M
	  regularity M
     Text

     
     Code
     	  SUBSECTION "Division With Remainder"
     Text
	  A major application of Gröbner bases is
	  to give the normal form for an element modulo an
	  ideal, allowing one, for example, to decide whether
	  the element is in the ideal.
	  For example, we can decide which power of the trace
	  of a generic 3x3 matrix is expressible in terms of the entries of the 
	  cube of the matrix with the following code:
     Example
	  R = kk[a..i]
	  M = genericMatrix(R,a,3,3)
	  I = ideal M^3
     Text
	  This gives the ideal of entries of the matrix. In the expression
	  ``{\tt M = genericMatrix(R,a,3,3)}'' the arguments ``{\tt R,a,3,3}'' specify the
	  ring, the first variable to use, and the numbers of rows and columns
	  desired.
     Example
	  Tr = trace M 
	  for p from 1 to 10 do print (Tr^p % I)
     Text
	  The expression ``@ TT "Tr^p % I"@'' computes the normal form for the p-th power
	  of the trace {\tt Tr} with respect to the Gröbner basis of I.
	  The expression ``@ TT "for p from 1 to 10 do" @'' specifies a 
	  {\em for loop} that executes the following expression, ``@ TT "print (Tr^p % I)" @'',
	  with 10 consecutive values of {\tt p}. For more information on such loops see @ TO "for" @
	  or type:
     Code
     	  EXAMPLE { PRE ///viewHelp "for"/// }
     Text
	  Here we have put quotes around ``for'' because
	  ``for'' is a keyword in the {\em Macaulay2} language.  (In general, it's always safe to use
	  quotes with viewHelp.)

	  We see from the output of these commands that the 6-th power
	  of the trace is NOT in the ideal of entries of the cube of M,
	  but the 7-th power is. We can compute the coefficients in the expression for it 
	  using the division algorithm, denoted in this setting by @ TO2 {(symbol //, RingElement,
	  Matrix), "//"} @:
     Example
	  Tr^7//(gens I)
     Text
     
     
     Code
     	  SUBSECTION "Elimination Theory"
     Text
	  Consider the problem of projecting the
	  ``twisted cubic'', a curve in $\PP^3$ defined
	  by the three $2 \times{} 2$ minors of a certain
	  $2 \times{} 3$ matrix.  
	  We already have the simplest tools for solving
	  such a problem.
	  We first clear the earlier meaning of {\tt x} 
	  to allow it to be used as a subscripted variable:
     Example
	  x = symbol x
     Text
	  Since we are going to deal with a curve in $\PP^3$,
	  we begin with a polynomial ring in four variables:
     Example
	  R = kk[x_0..x_3] 
     Text
	  The ideal of the twisted cubic curve is generated by the $2 \times{} 2$
	  minors of a ``catalecticant" or ``Hankel" matrix, conveniently
	  defined as follows:
     Example
	  M = map(R^2, 3, (i,j)->x_(i+j))
	  I = minors(2,M)
     Text
	  As projection center we
	  take the point with homogeneous coordinates $(1,0,0,-1)$,
	  which is defined by the ideal:
     Example
	  pideal = ideal(x_0+x_3, x_1, x_2)
     Text
	  The ideal J of the image of the curve under the projection from this point
	  is the kernel of the ring map $S=kk[u,v,w] \to R/I$
	  sending the variables
	  of S to the generators of {\tt pIdeal},
	  regarded as elements of $R/I$.  This is the same as the more usual formulation:
	  $$J = I \cap{} kk[x_0+x_3, x_1, x_x]$$ 
	  To compute this we first substitute {\tt pIdeal} into $R/I$, and then form
	  the necessary ring map:
     Example
	  Rbar = R/I
	  pideal = substitute(pideal, Rbar)
	  S = kk[u,v,w]
	  J=kernel map (Rbar, S, gens pideal)
     Text
	  The ideal J defines a curve with one singular point.
	  We can compute the ideal of the singular locus with:
     Example
	  K = ideal singularLocus(J)
     Text
	  This doesn't look like the ideal of a reduced point! But
	  that's because it isn't yet saturated:
     Example
	  saturate K
     Text
	  We have just seen the @TO "Saturation :: saturate"@ function in its most
	  common use: to saturate with respect to the maximal ideal.
	  but we can also find the saturation of any ideal with
	  respect to another:
     Example
	  saturate (ideal"u3w,uv", ideal"u")
     Text
	  We can also take the ``ideal quotient'' I:J of an ideal I with
	  respect to another, J
	  defined as the set of elements f such that
	  f*J is contained in I:
     Example
	  ideal"u3w,uv":ideal"u"
     Text
     
     
     Code
     	  SUBSECTION "Defining functions and loading packages"
     Text
	  It is easy to define your own functions in {\em Macaulay2}, and this
	  can save a lot of typing. Functions are defined with the 
	  symbol ->. For example, the famous {\em Collatz Conjecture}
	  (also called the ``hailstone problem'') asks
	  about the following procedure: given an integer $n$, divide it
	  by 2 if possible, or else multiply by 3 and add 1. 
	  If we repeat this over and over,
	  does the process always reach 1?  Here is a function that 
	  performs the Hailstone procedure again and again,
	  producing a list of the intermediate results.
     Example
	  Collatz = n ->
	      while n != 1 list if n%2 == 0 then n=n//2 else n=3*n+1
     Text
	  For example:
     Example
	  Collatz 27
     Text
	  If you don't understand this code easily, see @ TO Function @ and @ TO "while" @, or try:
     Code
     	  EXAMPLE { 
	       PRE "viewHelp Function",
	       PRE ///viewHelp "while"///
	       }
     Text
	  In order to understand a process it is often useful to tabulate the 
	  results of applying it many times. One feature of the Collatz process
	  is how many steps it takes to get to 1. We can tabulate this statistic
	  for the first 25 values of n with the function @ TO tally @, as follows:
     Example
	  tally for n from 1 to 30 list length Collatz n
     Text
	  A line of the form
     Code
     	  EXAMPLE { 
	       PRE "            18 => 3"
	       }
     Text
	  in the result means that a Collatz sequence of length 18
	  was seen 3 times. 
	  To see the successive ``record-breakers'', 
	  that is, the numbers with longer Collatz sequences than any
	  number before them, we might try:
     Example
	  record = length Collatz 1
	  L = for n from 2 to 1000 list (
		  l := length Collatz n;
		  if l > record 
		    then (record = l; (n,l))
		    else continue)
     Text
	  If you want to see a list of just the successive records, 
	  you can apply the 
	  function @ TO last @ to each element of the list $L$. 
	  A convenient way to do this is with this syntax:
     Example
	  L/last
     Text

	  Note that in
	  writing functions of more than one expression (usually
	  there's one expression per line), the expressions must be
	  separated by semicolons. For example in the ``for'' loop
	  above, the first expression was ``{\tt l = length Collatz n}''. 
	  After the last expression of an input line or of a function body,
	  a semicolon suppresses output, useful when the output
	  would be large. 	  
	  
	  There are many packages of ready-made functions available for
	  your use, many written by other users (perhaps you'll contribute one
	  someday!) A list of ``installed'' packages can be found with:
     Code
     	  EXAMPLE { PRE ///viewHelp/// }
     Text
	  For example, there is a package called @TO "EdgeIdeals::EdgeIdeals"@. 
	  To load the package, use:
     Example
	  loadPackage "EdgeIdeals"
     Text
	  After loading it, you can view its documentation with 
     Code
     	  EXAMPLE { PRE "viewHelp EdgeIdeals" }
     Text
     	  or you can call its functions, 
	  such as @TO "EdgeIdeals::randomGraph"@ and @TO "EdgeIdeals::edgeIdeal"@:
     Example
	  R = kk[vars(0..10)]
	  G=randomGraph (R,20)
	  K=edgeIdeal G
	  hilbertSeries K
	  betti res K
     Text

	  When testing a conjecture  one sometimes wants to run a 
	  large number of randomly chosen
	  examples.
	  Here's some typical code that one might use to study
	  a random graph ideal.  First we use ``{\tt for ... list ...}'' to construct a list {\tt L}
	  and suppress its printing by ending the line that creates
	  it with a ``;''.  Each entry of {\tt L} is a triple consisting of the
	  codimension, degree, and Betti table of a random graph ideal
	  on 10 vertices having only 4 edges.
     Example
	  R = ZZ/2[vars(0..10)]
	  L=for j from 1 to 100 list(
	      I = edgeIdeal randomGraph (R,5);
	      (codim I, degree I, betti res I));
     Text
	  We can use @ TO tally @ to find out how many examples
	  were found with each combination of codimension and degree and Betti table.
     Example
	  tally L
     Text
	  We can determine how many distinct patterns were found:
     Example
	  #tally L
     Text
     
     
     Code
     	  SUBSECTION "Ext, Tor, and cohomology"
     Text
	  {\em Macaulay2} can compute the homology of complexes;
	  for example, let's compute the homology of a
	  Koszul complex that is not a resolution:
	  $$ {\mathbf K}(x^2, x y^2):\ \  0 \rightarrow{} S(-5) \rightarrow{} S(-2)\oplus S(-3) \rightarrow{} S \rightarrow 0 $$
	  The free module $S(-2) \oplus{} S(-3)$ can be defined with this
	  syntax:
     Example
	  S^{-2,-3} 
     Text
          Here is how we can define the maps in the Koszul complex:
     Example
	  S = kk[x,y]
	  phi1 = map(S^1, S^{-2,-3}, matrix"x2,xy2")
	  phi2 = map(S^{-2,-3}, S^{-5}, matrix"xy2;-x2")
     Text
	  Let's check that this is will really make a complex:
     Example
	  phi1*phi2
     Text
	  To get the homology we can, for example compute:
     Example
	  (ker phi1)/(image phi2)
     Text
	  We could also use the data type @TO "OldChainComplexes :: ChainComplex"@ 
	  and use a built-in facility to take homology (in our case $H_1$):
     Example
	  FF = chainComplex(phi1,phi2)
	  FF.dd
	  homology FF
	  presentation (homology FF)_1
     Text
	  Either way, the first homology is $((x^2):(xy^2)) / (x^2) \cong{} S/(x)$, in accord
	  with general theory.

     	  There are other ways to construct Koszul complexes.  One way is as the tensor product of
	  chain complexes of length 1:
     Example
     	  FF = chainComplex matrix {{x^2}} ** chainComplex matrix {{x*y^2}}
	  FF.dd
     Text
     	  Another way is by using the function @ TO koszul @, designed for that purpose:
     Example
     	  FF = koszul matrix {{x^2, x*y^2}}
	  FF.dd
     Text

	  Since {\em Macaulay2} can compute resolutions and homology, it can
	  compute things such as $Ext$, $Tor$ and sheaf cohomology, as in the 
	  following examples. The first uses Serre's formula to compute
	  the multiplicity with which a 2-plane meets the union
	  of two 2-planes in 4-space (this is the first case in which
	  the length of the intersection scheme is NOT the right answer.)
	  The notation ``{\tt M**N}'' denotes the tensor product of the modules $M$ and $N$.
	  We use the syntactical forms
	  ``{\tt for j from 0 to 4 list ...}'' to list some results and
	  ``{\tt sum(0..4, j -> ...)}'' to sum some results.
     Example
	  S=kk[a,b,c,d]
	  IX = intersect(ideal(a,b), ideal(c,d))
	  IY = ideal(a-c, b-d)
	  degree ((S^1/IX) ** (S^1/IY))
	  for j from 0 to 4 list degree Tor_j(S^1/IX, S^1/IY)
	  sum(0..4, j-> (-1)^j * degree Tor_j(S^1/IX, S^1/IY))
     Text
	  Similarly, we can compute Hom and Ext:
     Example
	  Hom(IX, S^1/IX)
	  Ext^1(IX, S^1/IX)
     Text
	  or the cohomology of the sheaf associated to a module. 

	  Here is how to compute
	  the first cohomology of the structure
	  sheaf twisted by $-2$ of the curve $Proj(S/IX)$, which
	  in this case is the disjoint union of two
	  lines in $\PP^3$:
     Example
	  HH^1 (sheaf (S^{-2}**(S^1/IX)))
-- Local Variables:
-- compile-command: "make -C /Users/de/src/M2/Macaulay2/packages PACKAGES=BeginningMacaulay2 RemakeAllDocumentation=true IgnoreExampleErrors=false"
-- End:
