@comment -*-texinfo-*-
@comment this file contains the examples

@c The following directives are necessary for proper compilation
@c with emacs (C-c C-e C-r).  Please keep it as it is.  Since it
@c is wrapped in `@ignore' and `@end ignore' it does not harm `tex' or
@c `makeinfo' but is a great help in editing this file (emacs
@c ignores the `@ignore').
@ignore
%**start
\input texinfo.tex
@setfilename smallexamples.info
@node Top, Examples
@menu
* General concepts::
@end menu
@node Examples, Mathematical Background, Tricks and pitfalls, Top
@appendix Examples
%**end
@end ignore

@ifinfo
The following topics are treated:
@end ifinfo

@ifset singularmanual
@menu
* Programming::procedures and libraries, formatting output, etc.
* Computing Groebner and Standard Bases::GB conversion, slim GB
* Commutative Algebra::saturation, elimination, free resolution, factorization, primary decomposition, normalization, etc.
* Singularity Theory::singular and critical points, invariants of hypersurface singularities, classification, resolution of singularities, etc.
* Invariant Theory:: G_a-invariants, invariants of finite groups
* Geometric Invariant Theory:: Geometric Invariant Theory
* Non-commutative Algebra:: Groebner bases and applications in G-algebras
* Applications:: Solving systems of polynomial equations, AG codes
@end menu
@end ifset

@ifclear singularmanual
@menu
* Programming::
* Computing Groebner and Standard Bases::
* Commutative Algebra::
* Singularity Theory::
* Non-commutative Algebra::
* Applications::
@end menu
@end ifclear

@c ----------------------------------------------------------------------------
@c @node Start SINGULAR, Milnor and Tjurina,Examples, Examples
@c @section Start SINGULAR
@c @cindex Start SINGULAR

@c Call @sc{Singular} by typing @code{Singular} [return]

@c To use the online help type for instance:
@c    @code{help;} @code{help command;} @code{help General syntax;} @code{help ring;}...
@c Please note:  EVERY COMMAND MUST END WITH A SEMICOLON ";"

@c To leave @sc{Singular}, type one of the:
@c    @code{quit;} @code{exit;} @code{$}

@c The two characters @code{//} make the rest of the line a comment.


@c ----------------------------------------------------------------------------
@node Programming, Computing Groebner and Standard Bases, Examples, Examples
@section Programming
@cindex Programming

@ifset singularmanual
@menu
* Basic programming::
* Writing procedures and libraries::
* Rings associated to monomial orderings::
* Long coefficients::
* Parameters::
* Formatting output::
* Cyclic roots::
* Parallelization with ssi links::
* Dynamic modules::
@end menu
@end ifset

@ifclear singularmanual
@menu
* Basic programming::
* Writing procedures and libraries::
* Rings associated to monomial orderings::
* Parameters::
* Formatting output::
* Dynamic modules::
@end menu
@end ifclear

@c ----------------------------------------------------------------------------
@node Basic programming, Writing procedures and libraries, Programming, Programming
@subsection Basic programming
@cindex Basic programming
@cindex Basic programming

We show in the example below the following:
@itemize @bullet
@item
define the ring @code{R} of characteristic 32003, variables @code{x,y,z}, monomial
  ordering @code{dp}   (implementing F_32003[x,y,z])
@item
list the information about @code{R} by typing its name
@item
check the order of the variables
@item
define the integers @code{a,b,c,t}
@item
define a polynomial @code{f} (depending on @code{a,b,c,t}) and display it
@item
define the jacobian ideal @code{i} of @code{f}
@item
compute a Groebner basis of @code{i}
@item
compute the dimension of the algebraic set defined by @code{i} (requires
the computation of a Groebner basis)
@item
create and display a string in order to comment the result (text between
quotes "  "; is a 'string')
@item
load a library (see @ref{primdec_lib})
@item
compute a primary decomposition for @code{i} and assign the result to a
list @code{L} (which is a list of lists of ideals)
@item
display the number of primary components and the first primary and prime
components (entries of the list L[1])
@item
implement the localization of F_32003[x,y,z] at the homogeneous maximal
ideal (generated by x,y,z) by defining a ring with local monomial
ordering (@code{ds} in place of @code{dp})
@item
map i to this ring (see @ref{imap}) - we may use the same name @code{i},
since ideals are ring dependent data
@item
compute the local dimension of the algebraic set defined by @code{i} at
the origin (= dimension of the ideal generated by @code{i} in the localization)
@item
compute the local dimension of the algebraic set defined by @code{i} at
the point (-2000,-6961,-7944) (by applying a linear coordinate transformation)
@end itemize
For a more basic introduction to programming in @sc{Singular}, we refer
to @ref{Getting started}.

@smallexample
@c example
ring R = 32003,(x,y,z),dp;
R;
x > y;
y > z;
int a,b,c,t = 1,2,-1,4;
poly f = a*x3+b*xy3-c*xz3+t*xy2z2;
f;
ideal i = jacob(f);    // Jacobian Ideal of f
ideal si = std(i);     // compute Groebner basis
int dimi = dim(si);
string s = "The dimension of V(i) is "+string(dimi)+".";
s;
LIB "primdec.lib";     // load library primdec.lib
list L = primdecGTZ(i);
size(L);               // number of prime components
L[1][1];               // first primary component
L[1][2];               // corresponding prime component
ring Rloc = 32003,(x,y,z),ds;   // ds = local monomial ordering
ideal i = imap(R,i);
dim(std(i));
map phi = R, x-2000, y-6961, z-7944;
dim(std(phi(i)));
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Writing procedures and libraries, Rings associated to monomial orderings,  Basic programming, Programming
@end ifset
@ifclear singularmanual
@node Writing procedures and libraries, Rings associated to monomial orderings,  Basic programming, Programming
@end ifclear
@subsection Writing procedures and libraries
@cindex Procedures and libraries
@cindex Libraries

The user may add their own commands to the commands already available in
@sc{Singular} by writing @sc{Singular} procedures. There are basically
two kinds of procedures:
@itemize @bullet
@item procedures written in the @sc{Singular} programming language (which are
usually collected in @sc{Singular} libraries).
@item procedures written in C/C++ (collected in dynamic modules).
@end itemize

At this point, we restrict ourselves to describing the first kind of
(library) procedures, which are sufficient for most applications. The
syntax and general structure of a library (procedure) is described in
@ref{Procedures}, and @ref{Libraries}.

The probably most efficient way of
writing a new library is to use one of the official @sc{Singular}
libraries, say @code{ring.lib} as a sample. On a Unix-like operating
system, type @code{LIB "ring.lib";} to get information on where the
libraries are stored on your disk.

@sc{Singular} provides several commands and tools, which may be useful
when writing a procedure, for instance, to have a look at intermediate
results (see @ref{Debugging tools}).

If such a library should be contributed to @sc{Singular} some formal
requirements are needed:
@itemize
@item
the library header must explain the purpose of the library and
(for non-trivial algorithm) a pointer to the algorithm (text book, article, etc.)
@item
all global procedures must have a help string and an example which shows its usage.
@item
it is strongly recommend also to provide test scripts which test the functionality: one should test the essential functionality of the
  library/command in a relatively short time (say, in no more than 30s),
  other tests  should check the functionality of the
  library/command in detail so that, if possible, all relevant
  cases/results are tested. Nevertheless, such a test should not run
  longer than, say, 10 minutes.
@end itemize

We give short examples of procedures to demonstrate the following:

@itemize @bullet
@item
Write procedures which return an integer (ring independent), see also
@ref{Milnor and Tjurina number}. (Here we restrict ourselves to the main
body of the procedures).
  @itemize @minus
  @item
  The procedure @code{milnorNumber} must be called with one parameter, a
  polynomial.
  The name g is local to the procedure and is killed automatically when
  leaving the procedure.
  @code{milnorNumber} returns the Milnor number (and displays a comment).
  @item
  The procedure @code{tjurinaNumber} has no specified number of
  parameters. Here, the parameters are referred to by @code{#[1]} for
  the 1st, @code{#[2]} for the 2nd parameter, etc.
  @code{tjurinaNumber} returns the Tjurina number (and displays a comment).
  @item
  the procedure @code{milnor_tjurina} which returns a list consisting of two
  integers, the Milnor and the Tjurina number.
  @end itemize

@item
Write a procedure which creates a new ring and returns data dependent on
this new ring (two numbers) and an int. In this example, we also show
how to write a help text for the procedure (which is optional, but
recommended).
@end itemize


@smallexample
@c example
proc milnorNumber (poly g)
{
   "Milnor number:";
   return(vdim(std(jacob(g))));
}

proc tjurinaNumber
{
   "Tjurina number:";
   return(vdim(std(jacob(#[1])+#[1])));
}

proc milnor_tjurina (poly f)
{
   ideal j=jacob(f);
   list L=vdim(std(j)),vdim(std(j+f));
   return(L);
}

proc real_sols (number b, number c)
"USAGE: real_sols (b,c);  b,c number
ASSUME: active basering has characteristic 0
RETURN: list: first entry is an integer (the number of different real
        solutions). If this number is non-negative, the list has as second
        entry a ring in which the list SOL of real solutions of x^2+bx+c=0
        is stored (as floating point number, precision 30 digits).
NOTE:   This procedure calls laguerre_solve from solve.lib.
"
{
  def oldring = basering;  // assign name to the ring active when
                           // calling the procedure
  number disc = b^2-4*c;
  if (disc>0) { int n_of_sols = 2; }
  if (disc==0) { int n_of_sols = 1; }
  string s = nameof(var(1));  // name of first ring variable
  if (disc>=0) {
    execute("ring rinC =(complex,30),("+s+"),lp;");
    if (not(defined(laguerre_solve))) { LIB "solve.lib"; }
    poly f = x2+imap(oldring,b)*x+imap(oldring,c);
                        // f is a local ring-dependent variable
    list SOL = laguerre_solve(f,30);
    export SOL;         // make SOL a global ring-dependent variable
                        // such variables are still accessible when the
                        // ring is among the return values of the proc
    setring oldring;
    return(list(n_of_sols,rinC));
  }
  else {
    return(list(0));
  }
}

//
// We now apply the procedures which are defined by the
// lines of code above:
//
ring r = 0,(x,y),ds;
poly f = x7+y7+(x-y)^2*x2y2;

milnorNumber(f);
tjurinaNumber(f);
milnor_tjurina(f);     // a list containing Milnor and Tjurina number

def L=real_sols(2,1);
L[1];                  // number of real solutions of x^2+2x+1
def R1=L[2];
setring R1;
listvar(R1);           // only global ring-dependent objects are still alive
SOL;                   // the real solutions

setring r;
L=real_sols(1,1);
L[1];                  // number of reals solutions of x^2+x+1

setring r;
L=real_sols(1,-5);
L[1];                  // number of reals solutions of x^2+x-5
def R3=L[2];
setring R3; SOL;       // the real solutions
@c example
@end smallexample

Writing a dynamic module is not as simple as writing a library
procedure, since it does not only require some knowledge of C/C++, but
also about the way the @sc{Singular} kernel works.
See also @ref{Dynamic modules}.


@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Rings associated to monomial orderings, Long coefficients,  Writing procedures and libraries, Programming
@end ifset
@ifclear singularmanual
@node Rings associated to monomial orderings, Parameters,  Writing procedures and libraries, Programming
@end ifclear
@subsection Rings associated to monomial orderings
@cindex monomial orderings
@cindex localization
@cindex local rings, computing in

In @sc{Singular} we may implement localizations of the polynomial ring
by choosing an appropriate monomial ordering (when defining the ring by
the @code{ring} command). We refer to @ref{Monomial orderings} for a
thorough discussion of the monomial orderings available in
@sc{Singular}.

At this point, we restrict ourselves to describing the relation between
a monomial ordering and the ring (as mathematical object) which is
implemented by the ordering. This is most easily done by describing the
set of units: if > is a monomial ordering then precisely those elements
which have leading monomial 1 are considered as units (in all
computations performed with respect to this ordering).

In mathematical terms: choosing a monomial ordering @code{>} implements the
localization of the polynomial ring with respect to the multiplicatively
closed set of polynomials with leading monomial 1.

@tex
That is, choosing $>$ implements the ring
$$
  K[x]_> := \left\{{{f}\over{u}}\; \bigg|\; f, u \in K[x],
  \, LM(u) =
  1\right\}.
$$
@end tex

If > is global (that is, 1 is the smallest monomial), the implemented
ring is just the polynomial ring. If > is local (that is, if 1 is the
largest monomial), the implemented ring is the localization of the
polynomial ring w.r.t.@: the homogeneous maximal ideal. For a mixed
ordering, we obtain "something in between these two rings":

@smallexample
@c example
ring R = 0,(x,y,z),dp;     // polynomial ring (global ordering)
poly f = y4z3+2x2y2z2+4z4+5y2+1;
f;                         // display f in a degrevlex-ordered way
short=0;                   // avoid short notation
f;
short=1;
leadmonom(f);              // leading monomial

ring r = 0,(x,y,z),ds;     // local ring (local ordering)
poly f = fetch(R,f);
f;                         // terms of f sorted by degree
leadmonom(f);              // leading monomial

// Now we implement more "advanced" examples of rings:
//
// 1)   (K[y]_<y>)[x]
//
int n,m=2,3;
ring A1 = 0,(x(1..n),y(1..m)),(dp(n),ds(m));
poly f  = x(1)*x(2)^2+1+y(1)^10+x(1)*y(2)^5+y(3);
leadmonom(f);
leadmonom(1+y(1));         // unit
leadmonom(1+x(1));         // no unit

//
// 2)  some ring in between (K[x]_<x>)[y] and K[x,y]_<x>
//
ring A2 = 0,(x(1..n),y(1..m)),(ds(n),dp(m));
leadmonom(1+x(1));       // unit
leadmonom(1+x(1)*y(1));  // unit
leadmonom(1+y(1));       // no unit

//
// 3)  K[x,y]_<x>
//
ring A4 = (0,y(1..m)),(x(1..n)),ds;
leadmonom(1+y(1));       // in ground field
leadmonom(1+x(1)*y(1));  // unit
leadmonom(1+x(1));       // unit

@c example
@end smallexample

Note, that even if we implicitly compute over the localization of
the polynomial ring, most computations are explicitly performed with
polynomial data only.
In particular,  @code{1/(1-x);} does not return a power series
expansion or a fraction but 0 (division with remainder in polynomial ring).

See @ref{division} for division with remainder in the localization and
@ref{invunit} for a procedure returning a truncated power series expansion
of the inverse of a unit.

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Long coefficients, Parameters, Rings associated to monomial orderings, Programming
@subsection Long coefficients
@cindex Long coefficients
@cindex coefficients, long

The following innocent example produces in its standard basis
extremely long coefficients in char 0 for the lexicographical
ordering.
But a very small deformation does not (the undeformed
example is degenerated with respect to the Newton boundary).
This example demonstrates that it might be wise, for complicated
examples, to do the calculation first in positive char (e.g., 32003).
It has been shown, that in complicated examples, more than 95 percent of
the time needed for a standard basis computation is used in the
computation of the coefficients (in char 0).
The representation of long integers with real is demonstrated.

@smallexample
@c example
timer = 1;                              // activate the timer
ring R0 = 0,(x,y),lp;
poly f = x5+y11+xy9+x3y9;
ideal i = jacob(f);
ideal i1 = i,i[1]*i[2];                 // undeformed ideal
ideal i2 = i,i[1]*i[2]+1/1000000*x5y8;  // deformation of i1
i1; i2;
ideal j = std(i1);
j;
// Compute average coefficient length (=51) by
//   - converting j[2] to a string in order to compute the number
//   of characters
//   - divide this by the number of monomials:
size(string(j[2])) div size(j[2]);
vdim(j);
// For a better representation normalize the long coefficients
// of the polynomial j[2] and map it to real:
poly p=(1/12103947791971846719838321886393392913750065060875)*j[2];
ring R1=real,(x,y),lp;
short=0; // force the long output format
poly p=imap(R0,p);
p;
// Compute a standard basis for the deformed ideal:
setring R0; // return to the original ring R0
j = std(i2);
j;
vdim(j);
@c example
@end smallexample
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Parameters, Formatting output, Long coefficients, Programming
@end ifset
@ifclear singularmanual
@node Parameters, Formatting output, Rings associated to monomial orderings, Programming
@end ifclear
@subsection Parameters
@cindex Parameters

@ifset singularmanual
Let us deform the ideal in @ref{Long coefficients} by
introducing a parameter t and compute over the ground field Q(t).
We compute the dimension at the generic point,
@end ifset
@ifclear singularmanual
Let us now deform a given 0-dimensional ideal j by introducing a parameter t
and compute over the ground field Q(t).
We compute the dimension at the generic point,
@end ifclear
i.e.,
@tex
$dim_{Q(t)}Q(t)[x,y]/j$.
@end tex
@ifinfo
dim_Q(t) Q(t)[x,y]/j.
@end ifinfo
@ifset singularmanual
(This gives the
same result as for the deformed ideal above. Hence, the above small
deformation was "generic".)
@end ifset

For almost all
@tex
$a \in Q$
@end tex
@ifinfo
a in Q
@end ifinfo
this is the same as
@tex
$dim_Q Q[x,y]/j_0$,
@end tex
@ifinfo
dim_Q Q[x,y]/j0,
@end ifinfo
where
@tex
$j_0=j|_{t=a}$.
@end tex
@ifinfo
j_0=j_t=a
@end ifinfo

@smallexample
@c example
  ring Rt = (0,t),(x,y),lp;
  Rt;
  poly f = x5+y11+xy9+x3y9;
  ideal i = jacob(f);
  ideal j = i,i[1]*i[2]+t*x5y8;  // deformed ideal, parameter t
  vdim(std(j));
  ring R=0,(x,y),lp;
  ideal i=imap(Rt,i);
  int a=random(1,30000);
  ideal j=i,i[1]*i[2]+a*x5y8;  // deformed ideal, fixed integer a
  vdim(std(j));
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Formatting output, Cyclic roots, Parameters, Programming
@end ifset
@ifclear singularmanual
@node Formatting output, Dynamic modules, Parameters, Programming
@end ifclear
@subsection Formatting output
@cindex Formatting output
@cindex Output, formatting of

We show how to insert the result of a computation inside a text
by using strings.
First we compute the powers of 2 and comment the result with some text.
Then we do the same and give the output a nice format by computing and
adding appropriate space.

@smallexample
@c example
  // The powers of 2:
  int  n;
  for (n = 2; n <= 128; n = n * 2)
  {"n = " + string (n);}
  // The powers of 2 in a nice format
  int j;
  string space = "";
  for (n = 2; n <= 128; n = n * 2)
  {
    space = "";
    for (j = 1; j <= 5 - size (string (n)); j = j+1)
    { space = space + " "; }
    "n =" + space + string (n);
  }
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Cyclic roots, Parallelization with ssi links , Formatting output, Programming
@subsection Cyclic roots
@cindex Cyclic roots

We write a procedure returning a string that enables us to create
automatically the ideal of cyclic roots over the basering with n
variables. The procedure assumes that the variables consist of a single
letter each (hence no indexed variables are allowed; the procedure
@code{cyclic} in @code{polylib.lib} does not have this restriction). Then
we compute a standard basis of this ideal and some numerical
information.  (This ideal is used as a classical benchmark for standard
basis computations).

@smallexample
// We call the procedure 'cyclic':
proc cyclic (int n)
@{
   string vs = varstr(basering)+varstr(basering);
   int c=find(vs,",");
   while ( c!=0 )
   @{
      vs=vs[1,c-1]+vs[c+1,size(vs)];
      c=find(vs,",");
   @}
   string t,s;
   int i,j;
   for ( j=1; j<=n-1; j=j+1 )
   @{
      t="";
      for ( i=1; i <=n; i=i+1 )
      @{
         t = t + vs[i,j] + "+";
      @}
      t = t[1,size(t)-1] + ","+newline;
      s=s+t;
   @}
   s=s+vs[1,n]+"-1";
   return (s);
@}

ring r=0,(a,b,c,d,e),lp;         // basering, char 0, lex ordering
string sc=cyclic(nvars(basering));
sc;                              // the string of the ideal
@expansion{} a+b+c+d+e,
@expansion{} ab+bc+cd+de+ea,
@expansion{} abc+bcd+cde+dea+eab,
@expansion{} abcd+bcde+cdea+deab+eabc,
@expansion{} abcde-1
execute("ideal i="+sc+";");      // this defines the ideal of cyclic roots
i;
@expansion{} i[1]=a+b+c+d+e
@expansion{} i[2]=ab+bc+cd+ae+de
@expansion{} i[3]=abc+bcd+abe+ade+cde
@expansion{} i[4]=abcd+abce+abde+acde+bcde
@expansion{} i[5]=abcde-1
timer=1;
ideal j=std(i);
@expansion{} //used time: 7.5 sec
size(j);                         // number of elements in the std basis
@expansion{} 11
degree(j);
@expansion{} // codimension = 5
@expansion{} // dimension   = 0
@expansion{} // degree      = 70
@end smallexample
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node  Parallelization with ssi links, Dynamic modules, Cyclic roots, Programming
@subsection Parallelization with ssi links
@cindex Parallelization
@cindex ssi
@cindex link

In this example, we demonstrate how ssi links can be used to
parallelize computations.

To compute a standard basis for a zero-dimensional ideal in the
lexicographical ordering, one of the two powerful routines
@code{stdhilb}
@ifset singularmanual
(see @ref{stdhilb})
@end ifset
and @code{stdfglm}
@ifset singularmanual
(see @ref{stdfglm})
@end ifset
should be used. However, in general one cannot predict
which one of the two commands is faster. This very much depends on the
(input) example. Therefore, we use ssi links to let both commands
work on the problem independently and in parallel, so that the one which
finishes first delivers the result.

The example we use is the so-called "omndi example". See @i{Tim
Wichmann; Der FGLM-Algorithmus: verallgemeinert und implementiert in
Singular; Diplomarbeit Fachbereich Mathematik, Universitaet
Kaiserslautern; 1997} for more details.

@smallexample
@c example
ring r=0,(a,b,c,u,v,w,x,y,z),lp;
ideal i=a+c+v+2x-1, ab+cu+2vw+2xy+2xz-2/3,  ab2+cu2+2vw2+2xy2+2xz2-2/5,
ab3+cu3+2vw3+2xy3+2xz3-2/7, ab4+cu4+2vw4+2xy4+2xz4-2/9, vw2+2xyz-1/9,
vw4+2xy2z2-1/25, vw3+xyz2+xy2z-1/15, vw4+xyz3+xy3z-1/21;

link l_hilb,l_fglm = "ssi:fork","ssi:fork";         // 1.

open(l_fglm); open(l_hilb);

write(l_hilb, quote(stdhilb(i)));                       // 2.
write(l_fglm, quote(stdfglm(eval(i))));

list L=list(l_hilb,l_fglm);                             // 3.
int l_index=waitfirst(L);

if (l_index==1)
{
  "stdhilb won !!!!"; size(read(L[1]));
  close(L[1]); close(L[2]);
}
else                                                    // 4.
{
  "stdfglm won !!!!"; size(read(L[2]));
  close(L[1]); close(L[2]);
}
@c example
@end smallexample
Some explanatory remarks are in order:
@enumerate
@item
Instead of using links of the type @code{ssi:fork}, we alternatively
could use @code{ssi:tcp} links such that the two "competing"
@sc{Singular} processes run on different machines. This has the
advantage of "true" parallel computing since no resource sharing is
involved (as it usually is with forked processes).

@item
Notice how quoting is used in order to prevent local evaluation
(i.e., local computation of results). Since we "forked" the two
competing processes, the identifier @code{i} is defined and has
identical values in both child processes. Therefore, the innermost
@code{eval} can be omitted (as is done for the @code{l_hilb} link),
and only the identifier @code{i} needs to be communicated to the
children. However, when @code{ssi:tcp} links are used, the inner
evaluation must be applied so that actual values, and not the
identifiers are communicated (as is done for the @code{l_fglm} link
in our example).

@item
We wait until one of the two children
finished the computation. The main process sleeps
(i.e., suspends its execution) in the intermediate time.

@item
The child which has won delivers the result and is terminated with the usual
@code{close} command. The other child which is still computing needs to
be terminated by an explicit (i.e., system) kill command if running on a different computer. For ssi:fork a @code{close} is sufficient.
@end enumerate
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Dynamic modules, groebner and std, Parallelization with ssi links, Programming
@end ifset
@ifclear singularmanual
@node Dynamic modules, groebner and std, Formatting output, Programming
@end ifclear
@subsection Dynamic modules
@cindex Dynamic modules

The purpose of the following example is to illustrate the use of dynamic
modules. Giving an example on how to write a dynamic module is beyond the
scope of this manual.
A technical reference is given at @url{https://www.singular.uni-kl.de/Manual/modules.pdf}.

In this example, we use a dynamic module, residing in the file @code{kstd.so},
which allows ignoring all but the first j entries of vectors when
forming the pairs in the standard basis computation.

@smallexample
@c example
ring r=0,(x,y),dp;
module mo=[x^2-y^2,1,0,0],[xy+y^2,0,1,0],[y^2,0,0,1];
print(mo);

// load dynamic module - at the same time creating package Kstd
// procedures will be available in the packages Top and Kstd
LIB("partialgb.so");
listvar(package);

// set the number of components to be considered to 1
module mostd=partialStd(mo,1);        // calling procedure in Top
                    // obviously computation ignored pairs with leading
                    // term in the second entry
print(mostd);

// now consider 2 components
module mostd2=Partialgb::partialStd(mo,2); // calling procedure in Partialgb
                    // this time the previously unconsidered pair was
                    // treated too
print(mostd2);
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@node Computing Groebner and Standard Bases, Commutative Algebra, Programming, Examples
@section Computing Groebner and Standard Bases
@cindex Computing Groebner and Standard Bases
@cindex Groebner Bases
@cindex Standard Bases

@menu
* groebner and std::
* Groebner basis conversion::
* slim Groebner bases::
@end menu

Several operations with ideals resp. modules uses Groebner bases to compute their result.
Most allow an optional string argument to select the algorithm.
The possible arguments for the algorithm are
@itemize
@item @code{default}
@item @code{std}
see @ref{std}
@item @code{slimgb}
see @ref{slimgb}
@item  @code{sba}
see @ref{sba}; not for module operations
@item @code{singmatic}
Requires @code{singmatic.so}
@item @code{groebner}
see @ref{groebner}
@item @code{modstd}
see @ref{modStd}. Requires @ref{modstd_lib}
@item @code{ffmod}
see @ref{ffmodStd}, Requires  @ref{ffmodstd_lib}
@item  @code{nfmod}
see @ref{nfmodStd}. Requires @ref{nfmodstd_lib}
@item @code{std:sat}
Uses @code{satstd} instead of @code{std}, see @ref{satstd}. Requires @ref{customstd_lib}
@end itemize


Functions with such a choice of the algorithm:
@itemize
@item @ref{eliminate}
@item @ref{intersect}
@item @ref{modulo}
@item @ref{liftstd}
@item @ref{syz}
@end itemize

@c ----------------------------------------------------------------------------
@node groebner and std, Groebner basis conversion, Dynamic modules, Computing Groebner and Standard Bases
@subsection groebner and std
@cindex groebner
@cindex std

The basic version of Buchberger's algorithm leaves a lot
of freedom in carrying out the computational process. Considerable
improvements are obtained by implementing criteria
for reducing the number of S-polynomials to be actually considered
(e.g., by applying the product criterion or the chain criterion).
We refer to Cox, Little, and O'Shea [1997], Chapter 2 for more details
and references on these criteria and on further strategies for improving the
performance of Buchberger's algorithm (see also Greuel, Pfister [2002]).

@sc{Singular}'s implementation of Buchberger's algorithm is available via
the @code{std} command ('std' referring to @code{st}an@code{d}ard basis). The
computation of reduced Groebner and standard bases may be forced by setting
@code{option(redSB)} (see @ref{option}).

However, depending on the monomial ordering of the active basering, it may be
advisable to use the @code{groebner} command instead. This command is
provided by the @sc{Singular} library @code{standard.lib} which
is automatically loaded when starting a @sc{Singular} session. Depending
on some heuristics, @code{groebner} either refers to the @code{std}
command (e.g., for rings with ordering @code{dp}), or to one of the
algorithms described in the sections @ref{Groebner basis conversion},
@ref{slim Groebner bases}.
For information on the heuristics behind @code{groebner}, see the
library file @code{standard.lib} (see also @ref{Procedures and
libraries}).

We apply the commands @code{std} and @code{groebner} to compute a
lexicographic Groebner basis for the ideal of cyclic roots over the
basering with 6 variables (see @ref{Cyclic roots}). We set
@code{option(prot)} to make @sc{Singular} display some information on
the performed computations (see @ref{option} for an interpretation of
the displayed symbols). For long running computations, it is always
recommended to set this option.


@smallexample
@c example
LIB "polylib.lib";
ring r=32003,(a,b,c,d,e,f),lp;
ideal I=cyclic(6);
option(prot);
int t=timer;
system("--ticks-per-sec", 100);         // give time in 1/100 sec
ideal sI=std(I);
timer-t;                                // used time (in 1/100 secs)
size(sI);
t=timer;
sI=groebner(I);
timer-t;                                // used time (in 1/100 secs)
size(sI);
option(noprot);
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@node Groebner basis conversion, slim Groebner bases, groebner and std, Computing Groebner and Standard Bases
@subsection Groebner basis conversion
@cindex Groebner basis conversion
@cindex Lexicographic Groebner bases, computation of
@cindex Hilbert-driven GB algorithm

The performance of Buchberger's algorithm is
sensitive to the chosen monomial order. A Groebner basis
computation with respect to a less favorable order such as
the lexicographic ordering may easily run out of time or memory even
in cases where a Groebner basis computation with respect to a more efficient
order such as the degree reverse lexicographic ordering is very well
feasible.  Groebner basis conversion algorithms and the Hilbert-driven
Buchberger algorithm are based on this observation:


@itemize @bullet
@item
Groebner basis conversion: Given an ideal @math{I\subset K[x_1,\dots,x_n]} and a
slow monomial order, compute a Groebner basis with respect to an
appropriately chosen fast order. Then convert the result to a Groebner
basis with respect to the given slow order.

@item
Hilbert-driven Buchberger algorithm: Homogenize the given
generators for @math{I} with respect to a new variable, say, @math{x_0}.
Extend the given slow ordering on @math{K[x_1,\dots,x_n]} to a global
product ordering on @math{K[x_0,\dots,x_n]}.
Compute a Groebner basis for the ideal generated by the homogenized
polynomials with respect to a fast ordering. Read the Hilbert function, and
use this information when computing a Groebner basis with
respect to the extended (slow) ordering. Finally,
dehomogenize the elements of the resulting Groebner basis.
@end itemize

@sc{Singular} provides implementations for the FGLM conversion algorithm
(which applies to zero-dimensional ideals only, see @ref{stdfglm}) and
variants of the Groebner walk conversion algorithm (which works for
arbitrary ideals, @xref{frwalk}, @ref{grwalk_lib}).
An implementation of the Hilbert-driven Buchberger
algorithm is accessible via the @code{stdhilb} command (see also @ref{std}).

For the ideal below, @code{stdfglm} is more than 100 times
and @code{stdhilb} about 10 times faster than @code{std}.

@smallexample
@c example
  ring r =32003,(a,b,c,d,e),lp;
  ideal i=a+b+c+d, ab+bc+cd+ae+de, abc+bcd+abe+ade+cde,
          abc+abce+abde+acde+bcde, abcde-1;
  int t=timer;
  option(prot);
  ideal j1=stdfglm(i);
  timer-t;
  size(j1);   // size (no. of polys) in computed GB
  t=timer;
  ideal j2=stdhilb(i);
  timer-t;
  size(j2);   // size (no. of polys) in computed GB
  // usual Groebner basis computation for lex ordering
  t=timer;
  ideal j0 =std(i);
  option(noprot);
  timer-t;
@c example
@end smallexample


@c ----------------------------------------------------------------------------
@node slim Groebner bases, Saturation, Groebner basis conversion, Computing Groebner and Standard Bases
@subsection slim Groebner bases
@cindex Groebner bases, slim
@cindex slimgb
The command @code{slimgb} calls an implementation of an algorithm to compute
Groebner bases which is designed for keeping the polynomials slim (short
with small coefficients) during a Groebner basis computation. It provides,
in particular, a fast algorithm for computing Groebner bases over function
fields or over the rational numbers, but also in several other
cases. The algorithm which is still under development was developed in
the diploma thesis of Michael Brickenstein. It has been published as
@code{https://www.singular.uni-kl.de/reports/35/paper_35_full.ps.gz}.


In the example below (Groebner basis with respect to degree reverse
lexicographic ordering over function field) @code{slimgb} is much faster than
the @code{std} command.

@smallexample
@c example
ring r=(32003,u1, u2, u3, u4),(x1, x2, x3, x4, x5, x6, x7),dp;
timer=1;
ideal i=
-x4*u3+x5*u2,
x1*u3+2*x2*u1-2*x2*u2-2*x3*u3-u1*u4+u2*u4,
-2*x1*x5+4*x4*x6+4*x5*x7+x1*u3-2*x4*u1-2*x4*u4-2*x6*u2-2*x7*u3+u1*u2+u2*u4,
-x1*x5+x1*x7-x4*u1+x4*u2-x4*u4+x5*u3+x6*u1-x6*u2+x6*u4-x7*u3,
-x1*x4+x1*u1-x5*u1+x5*u4,
-2*x1*x3+x1*u3-2*x2*u4+u1*u4+u2*u4,
x1^2*u3+x1*u1*u2-x1*u2^2-x1*u3^2-u1*u3*u4+u3*u4^2;
i=slimgb(i);

@c example
@end smallexample

For detailed information and limitations see @ref{slimgb}.
@c ----------------------------------------------------------------------------
@node Commutative Algebra, Singularity Theory, Computing Groebner and Standard Bases, Examples
@section Commutative Algebra
@cindex Commutative Algebra

@ifset singularmanual
@menu
* Saturation::
* Finite fields::
* Elimination::
* Free resolution::
* Handling graded modules::
* Computation of Ext::
* Depth::
* Factorization::
* Primary decomposition::
* Normalization::
* Kernel of module homomorphisms::
* Algebraic dependence::
@end menu
@end ifset

@ifclear singularmanual
@menu
* Saturation::
* Elimination::
* Free resolution::
* Handling graded modules::
* Factorization::
* Primary decomposition::
* Kernel of module homomorphisms::
* Algebraic dependence::
@end menu
@end ifclear

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Saturation, Finite fields, slim Groebner bases, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Saturation, Elimination, slim Groebner bases, Commutative Algebra
@end ifclear
@subsection Saturation
@cindex Saturation

For any two ideals @math{i, j} in the basering @math{R} let
@tex
$$
\hbox{sat}(i,j)=\{x\in R\;|\; \exists\;n\hbox{ s.t. }
x\cdot(j^n)\subseteq i\}
= \bigcup_{n=1}^\infty i:j^n$$
@end tex
@ifinfo
@*sat(i,j) = @{x in @math{R} | there is an n s.t. x*(j^n) contained in i@}
@*         = union_(n=1...) of i:j^n,
@end ifinfo
@*denote the saturation of @math{i} with respect to @math{j}. This defines,
geometrically, the closure of the complement of V(@math{j}) in V(@math{i})
(where V(@math{i}) denotes the variety defined by @math{i}).



The saturation is computed by the procedure @code{sat} in
@code{elim.lib} by computing iterated ideal quotients with the maximal
ideal.  @code{sat} returns a list of two elements: the saturated ideal
and the number of iterations.

We apply saturation to show that a variety has no singular points
outside the origin (see also @ref{Critical points}).
We choose @math{m} to be the homogeneous maximal ideal
(note that @code{maxideal(n)} denotes the n-th power of the maximal
ideal).
Then @math{V(i)} has no singular point outside the origin
if and only if @math{sat(j+(f),m)} is the whole ring, that is, generated by 1.

@smallexample
@c example
  LIB "elim.lib";         // loading library elim.lib
  ring r2 = 32003,(x,y,z),dp;
  poly f = x^11+y^5+z^(3*3)+x^(3+2)*y^(3-1)+x^(3-1)*y^(3-1)*z3+
    x^(3-2)*y^3*(y^2)^2;
  ideal j=jacob(f);
  sat(j+f,maxideal(1));
  // list the variables defined so far:
  listvar();
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Finite fields, Elimination, Saturation, Commutative Algebra
@subsection Finite fields
@cindex Finite fields

We define a variety in the @math{n}-space of codimension 2 defined by
polynomials of degree @math{d} with generic coefficients over the prime
field @math{Z/p} and look for zeros on the torus. First over the prime
field and then in the finite extension field with
@tex
$p^k$
@end tex
@ifinfo
p^k
@end ifinfo
elements.
In general there will be many more solutions in the second case.
(Since the @sc{Singular} language is interpreted, the evaluation of many
@code{for}-loops is not very fast):

@smallexample
@c example
  int p=3;  int n=3;  int d=5; int k=2;
  ring rp = p,(x(1..n)),dp;
  int s = size(maxideal(d));
  s;
  // create a dense homogeneous ideal m, all generators of degree d, with
  // generic (random) coefficients:
  ideal m = maxideal(d)*random(p,s,n-2);
  m;
  // look for zeros on the torus by checking all points (with no component 0)
  // of the affine n-space over the field with p elements :
  ideal mt;
  int i(1..n);                    // initialize integers i(1),...,i(n)
  int l;
  s=0;
  for (i(1)=1;i(1)<p;i(1)=i(1)+1)
  {
    for (i(2)=1;i(2)<p;i(2)=i(2)+1)
    {
      for (i(3)=1;i(3)<p;i(3)=i(3)+1)
      {
        mt=m;
        for (l=1;l<=n;l=l+1)
        {
          mt=subst(mt,x(l),i(l));
        }
        if (size(mt)==0)
        {
          "solution:",i(1..n);
          s=s+1;
        }
      }
    }
  }
  "//",s,"solutions over GF("+string(p)+")";
  // Now go to the field with p^3 elements:
  // As long as there is no map from Z/p to the field with p^3 elements
  // implemented, use the following trick: convert the ideal to be mapped
  // to the new ring to a string and then execute this string in the
  // new ring
  string ms="ideal m="+string(m)+";";
  ms;
  // define a ring rpk with p^k elements, call the primitive element z. Hence
  // 'solution exponent: 0 1 5' means that (z^0,z^1,z^5) is a solution
  ring rpk=(p^k,z),(x(1..n)),dp;
  rpk;
  execute(ms);
  s=0;
  ideal mt;
  for (i(1)=0;i(1)<p^k-1;i(1)=i(1)+1)
  {
    for (i(2)=0;i(2)<p^k-1;i(2)=i(2)+1)
    {
      for (i(3)=0;i(3)<p^k-1;i(3)=i(3)+1)
      {
        mt=m;
        for (l=1;l<=n;l=l+1)
        {
          mt=subst(mt,x(l),z^i(l));
        }
        if (size(mt)==0)
        {
          // we show only the first 7 solutions here:
          if (s<5) {"solution exponent:",i(1..n);}
          s=s+1;
        }
      }
    }
  }
  "//",s,"solutions over GF("+string(p^k)+")";
@c example
@end smallexample
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Elimination, Free resolution, Finite fields, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Elimination, Free resolution, Saturation, Commutative Algebra
@end ifclear
@subsection Elimination
@cindex Elimination

Elimination is the algebraic counterpart of the geometric concept of
projection. If
@tex
$f=(f_1,\ldots,f_n):k^r\rightarrow k^n$
@end tex
@ifinfo
f=(f1,...,fn) : k^r --> k^n
@end ifinfo
is a polynomial map,
the Zariski-closure of the image is the zero-set of the ideal
@tex
$j=J \cap k[x_1,\ldots,x_n]$,
@end tex
@ifinfo
j = (J intersected with k[x_1,\ldots,x_n]),
@end ifinfo
where
@tex
$$J=(x_1-f_1(t_1,\ldots,t_r),\ldots,x_n-f_n(t_1,\ldots,t_r))\subseteq
k[t_1,\ldots,t_r,x_1,\ldots,x_n]$$
@end tex
@ifinfo

@smallexample
J=(x1-f1(t1,...,tr),...,xn-fn(t1,...,tr)) in k[t1,...tr,x1,...,xn]
@end smallexample

@end ifinfo
that is, of the ideal j obtained from J by eliminating the variables
@tex
$t_1,\ldots,t_r$.
@end tex
@ifinfo
t1,...,tr.
@end ifinfo
This can be done by computing a Groebner basis for J with respect to a
(global) product ordering where the block of t-variables precedes the block of
x-variables, and then selecting those polynomials which do not contain
any t. Alternatively, we may use a global monomial ordering with extra
weight vector (see @ref{Extra weight vector}), assigning to the
t-variables a positive weight and to the x-variables weight 0.

Since elimination is expensive, it may be useful to use a Hilbert-driven
approach to the elimination problem (see @ref{Groebner basis conversion}):

First compute the Hilbert function of the ideal w.r.t.
a fast ordering (e.g., @code{dp}), then make use of it to speed up the
computation by a Hilbert-driven elimination which uses the @code{intvec} provided as
third argument.

In @sc{Singular} the most convenient way is to use the
@code{eliminate} command.
In contrast to the first method, with @code{eliminate} the result needs not be a
standard basis in the given ordering.
Hence, there may be cases where the first method is the preferred one.

@strong{WARNING:} In the case of a local or a mixed ordering,
elimination needs special care. f may be considered as a map of germs
@tex
$f:(k^r,0)\rightarrow(k^n,0)$,
@end tex
@ifinfo
f : (k^r,0) --> (k^n,0),
@end ifinfo
but even
if this map germ is finite, we are in general not able to compute the image
germ because for this we would need an implementation of the Weierstrass
preparation theorem. What we can compute, and what @code{eliminate}
actually does, is the following: let V(J) be the zero-set of J in
@tex
$k^r\times(k^n,0)$,
@end tex
@ifinfo
k^r x (k^n,0),
@end ifinfo
then the
closure of the image of V(J) under the projection
@tex
$$\hbox{pr}:k^r\times(k^n,0)\rightarrow(k^n,0)$$
can be computed.
@end tex
@ifinfo
@*           pr:  k^r x (k^n,0) --> (k^n,0)
@*can be computed.
@end ifinfo
(Note that this germ contains also those components
of V(J) which meet the fiber of pr outside the origin.)
This is achieved by an ordering with the block of t-variables having a
global ordering (and preceding the x-variables) and the x-variables having
a local ordering.

In any case, if the input is weighted homogeneous (=quasihomogeneous),
the weights given to the variables should be chosen accordingly.
@sc{Singular} offers a function @code{weight} which proposes,
given an ideal or module, integer weights for the variables, such that
the ideal, resp.@: module, is as homogeneous as possible with respect to these weights.
The function finds correct weights, if the input is weighted homogeneous
(but is rather slow for many variables). In order to check, whether the
input is quasihomogeneous, use the function @code{qhweight}, which returns
an @code{intvec} of correct weights if the input is quasihomogeneous and an @code{intvec}
of zeros otherwise.

Let us give three examples:
@enumerate
@item
First we compute the equations of the simple space curve 'T[7]'
 consisting of two tangential cusps given in parametric form.
@item
We compute weights for the equations such that the
   equations are quasihomogeneous w.r.t. these weights.
@item
Then we compute the tangent developable of the rational
   normal curve in
@tex
$P^4$.
@end tex
@ifinfo
P^4.
@end ifinfo
@end enumerate

@smallexample
@c example
  // 1. Compute equations of curve given in parametric form:
  // Two transversal cusps in (k^3,0):
  ring r1 = 0,(t,x,y,z),ls;
  ideal i1 = x-t2,y-t3,z;        // parametrization of the first branch
  ideal i2 = y-t2,z-t3,x;        // parametrization of the second branch
  ideal j1 = eliminate(i1,t);
  j1;                            // equations of the first branch
  ideal j2 = eliminate(i2,t);
  j2;                            // equations of the second branch
  // Now map to a ring with only x,y,z as variables and compute the
  // intersection of j1 and j2 there:
  ring r2 = 0,(x,y,z),ds;
  ideal j1= imap(r1,j1);         // imap is a convenient ringmap for
  ideal j2= imap(r1,j2);         // inclusions and projections of rings
  ideal i = intersect(j1,j2);
  i;                             // equations of both branches
  //
  // 2. Compute the weights:
  intvec v= qhweight(i);         // compute weights
  v;
  //
  // 3. Compute the tangent developable
  // The tangent developable of a projective variety given parametrically
  // by F=(f1,...,fn) : P^r --> P^n is the union of all tangent spaces
  // of the image. The tangent space at a smooth point F(t1,...,tr)
  // is given as the image of the tangent space at (t1,...,tr) under
  // the tangent map (affine coordinates)
  //   T(t1,...,tr): (y1,...,yr) --> jacob(f)*transpose((y1,...,yr))
  // where jacob(f) denotes the jacobian matrix of f with respect to the
  // t's evaluated at the point (t1,...,tr).
  // Hence we have to create the graph of this map and then to eliminate
  // the t's and y's.
  // The rational normal curve in P^4 is given as the image of
  //        F(s,t) = (s4,s3t,s2t2,st3,t4)
  // each component being homogeneous of degree 4.
  ring P = 0,(s,t,x,y,a,b,c,d,e),dp;
  ideal M = maxideal(1);
  ideal F = M[1..2];     // take the 1st two generators of M
  F=F^4;
  // simplify(...,2); deletes 0-columns
  matrix jac = simplify(jacob(F),2);
  ideal T = x,y;
  ideal J = jac*transpose(T);
  ideal H = M[5..9];
  ideal i = matrix(H)-matrix(J);// this is tricky: difference between two
                         // ideals is not defined, but between two
                         // matrices. By type conversion
                         // the ideals are converted to matrices,
                         // subtracted and afterwards converted
                         // to an ideal. Note that '+' is defined
                         // and adds (concatenates) two ideals
  i;
  // Now we define a ring with product ordering and weights 4
  // for the variables a,...,e.
  // Then we map i from P to P1 and eliminate s,t,x,y from i.
  ring P1 = 0,(s,t,x,y,a,b,c,d,e),(dp(4),wp(4,4,4,4,4));
  ideal i = fetch(P,i);
  ideal j= eliminate(i,stxy);    // equations of tangent developable
  j;
  // We can use the product ordering to eliminate s,t,x,y from i
  // by a std-basis computation.
  // We need proc 'nselect' from elim.lib.
  LIB "elim.lib";
  j = std(i);                    // compute a std basis j
  j = nselect(j,1..4);           // select generators from j not
  j;                             // containing variable 1,...,4
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Free resolution, Handling graded modules, Elimination, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Free resolution, Handling graded modules, Elimination, Commutative Algebra
@end ifclear
@subsection  Free resolution
@cindex  Free resolution
@cindex Resolution, free

In @sc{Singular} a free resolution of a module or ideal has its own type:
@code{resolution}. It is a structure that stores all information related to
free resolutions. This allows partial computations of resolutions via
the command @code{res}. After applying @code{res}, only a pre-format of the
resolution is computed which allows to determine invariants like
Betti-numbers or homological dimension. To see the differentials
of the complex, a resolution must be converted into the type list which
yields a list of modules: the k-th module in this
list is the first syzygy-module (module of relations) of the (k-1)st module.
There are the following commands to compute a resolution:
@table @code
@item res
@ifset singularmanual
@ref{res}@*
@end ifset
computes a free resolution of an ideal or module using a heuristically
chosen method.
This is the preferred method to compute free resolutions of ideals or
modules.
@item fres
@ifset singularmanual
@ref{fres}@*
@end ifset
improved version of @ref{sres}, computes a free resolution of an ideal or
module using Schreyer's method. The input has to be a standard basis.
@item lres
@ifset singularmanual
@ref{lres}@*
@end ifset
computes a free resolution of an ideal or module with LaScala's
method. The input needs to be homogeneous.
@item mres
@ifset singularmanual
@ref{mres}@*
@end ifset
computes a minimal free resolution of an ideal or module with the syzygy
method.
@item sres
@ifset singularmanual
@ref{sres}@*
@end ifset
computes a free resolution of an ideal or module with Schreyer's
method. The input has to be a standard basis.
@item nres
@ifset singularmanual
@ref{nres}@*
@end ifset
computes a free resolution of an ideal or module with the standard basis
method.
@item minres
@ifset singularmanual
@ref{minres}@*
@end ifset
minimizes a free resolution of an ideal or module.
@item syz
@ifset singularmanual
@ref{syz}@*
@end ifset
computes the first syzygy module.
@end table
@code{res(i,r)}, @code{lres(i,r)}, @code{sres(i,r)}, @code{mres(i,r)},
@code{nres(i,r)} compute the first r modules of the resolution
of i, resp.@: the full resolution if r=0 and the basering is not a qring.
See the manual for a precise description of these commands.
@*Note: The command @code{betti} does not require a minimal
resolution for the minimal Betti numbers.

Now let us take a look at an example which uses resolutions: The Hilbert-Burch
theorem says that the ideal i of a reduced curve in
@tex
$K^3$
@end tex
@ifinfo
K^3
@end ifinfo
has a free resolution of length 2 and that i is given by the 2x2 minors
of the 2nd matrix in the resolution.
We test this for two transversal cusps in
@tex
$K^3$.
@end tex
@ifinfo
K^3.
@end ifinfo
Afterwards we compute the resolution of the ideal j of the tangent developable
of the rational normal curve in
@tex
$P^4$
@end tex
@ifinfo
P^4
@end ifinfo
from above.
Finally we demonstrate the use of the type @code{resolution} in connection with
the @code{lres} command.

@smallexample
@c example
  // Two transversal cusps in (k^3,0):
  ring r2 =0,(x,y,z),ds;
  ideal i =z2-1y3+x3y,xz,-1xy2+x4,x3z;
  resolution rs=mres(i,0);   // computes a minimal resolution
  rs;                        // the standard representation of complexes
    list resi=rs;            // conversion to a list
  print(resi[1]);            // the 1st module is i minimized
  print(resi[2]);            // the 1st syzygy module of i
  resi[3];                   // the 2nd syzygy module of i
  ideal j=minor(resi[2],2);
  reduce(j,std(i));          // check whether j is contained in i
  size(reduce(i,std(j)));    // check whether i is contained in j
  // size(<ideal>) counts the non-zero generators
  // ---------------------------------------------
  // The tangent developable of the rational normal curve in P^4:
  ring P = 0,(a,b,c,d,e),dp;
  ideal j= 3c2-4bd+ae, -2bcd+3ad2+3b2e-4ace,
           8b2d2-9acd2-9b2ce+9ac2e+2abde-1a2e2;
  resolution rs=mres(j,0);
  rs;
  list L=rs;
  print(L[2]);
  // create an intmat with graded Betti numbers
  intmat B=betti(rs);
  // this gives a nice output of Betti numbers
  print(B,"betti");
  // the user has access to all Betti numbers
  // the 2-nd column of B:
  B[1..4,2];
  ring cyc5=32003,(a,b,c,d,e,h),dp;
  ideal i=
  a+b+c+d+e,
  ab+bc+cd+de+ea,
  abc+bcd+cde+dea+eab,
  abcd+bcde+cdea+deab+eabc,
  h5-abcde;
  resolution rs=lres(i,0);   //computes the resolution according LaScala
  rs;                        //the shape of the minimal resolution
  print(betti(rs),"betti");  //shows the Betti-numbers of cyclic 5
  dim(rs);                   //the homological dimension
  size(list(rs));            //gets the full (non-reduced) resolution
  minres(rs);                //minimizes the resolution
  size(list(rs));            //gets the minimized resolution
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Handling graded modules, Computation of Ext, Free resolution, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Handling graded modules, Factorization, Free resolution, Commutative Algebra
@end ifclear
@subsection  Handling graded modules
@cindex  Handling graded modules
@cindex graded modules, handling of
@cindex Free resolution, graded

How to deal with graded modules in @sc{Singular} is best explained by
looking at an example:

@smallexample
@c example
ring R = 0, (w,x,y,z), dp;
module I = [-x,0,-z2,0,y2z], [0,-x,-yz,0,y3], [-w,0,0,yz,-z3],
           [0,-w,0,y2,-yz2], [0,-1,-w,0,xz], [0,-w,0,y2,-yz2],
           [x2,-y2,-wy2+xz2];
print(I);

// (1) Check on degrees:
// =====================
attrib(I,"isHomog"); // attribute not set => empty output
homog(I);
attrib(I,"isHomog");

print(betti(I,0),"betti");  // read degrees from Betti diagram

// (2) Shift degrees:
// ==================
def J=I;
intvec DV = 0,0,-1,-1,-2;
attrib(J,"isHomog",DV);    // assign new weight vector
attrib(J,"isHomog");
print(betti(J,0),"betti");

intmat bettiI=betti(I,0);  // degree corresponding to first non-zero row
                           // of Betti diagram is accessible via
                           // attribute "rowShift"
attrib(bettiI);
intmat bettiJ=betti(J,0);
attrib(bettiJ);

// (3) Graded free resolutions:
// ============================
resolution resJ = mres(J,0);
attrib(resJ);
print(betti(resJ),"betti");
attrib(betti(resJ));
@c example
@end smallexample

A check on degrees ((1), by using the @code{homog} command) shows that this
is a graded matrix. The @code{homog} command assigns an admissible weight
vector (here: 2,2,1,1,0) to the module @code{I} which is accessible via the
attribute @code{"isHomog"}. Thus, we may think of @code{I} as a graded
submodule of the graded free @math{R}-module
@tex
$$F=R(-2)^2\oplus R(-1)^2\oplus R\,.$$
@end tex
@ifinfo
F=R(-2)^2+ R(-1)^2+ R .
@end ifinfo
We may also read the degrees from the Betti diagram as shown above.
The degree on the left of the first nonzero row of the Betti diagram is
accessible via the attribute @code{"rowShift"} of the betti matrix
(which is of type @code{intmat}):

(2) We may shift degrees by assigning another admissible degree vector.
Note that @sc{Singular} does not check whether the assigned degree
vector really is admissible. Moreover, note that all assigned
attributes are lost under a type conversion (e.g. from @code{module} to
@code{matrix}).

(3) These considerations may be applied when computing data from free resolutions (see @ref{Computation of Ext}).

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Computation of Ext, Depth, Handling graded modules, Commutative Algebra
@subsection  Computation of Ext
@cindex  Ext, computation of

We start by showing how to calculate the n-th Ext group of an ideal. The
ingredients to do this are by the definition of Ext the following:
calculate a (minimal) resolution at least up to length n, apply the Hom
functor, and calculate the n-th homology group, that is, form the
quotient ker/im in the resolution sequence.

The Hom functor is given simply by transposing (hence dualizing) the
module or the corresponding matrix with the command @code{transpose}.
The image of the (n-1)-st map is generated by the columns of the
corresponding matrix. To calculate the kernel apply the command
@code{syz} at the (n-1)-st transposed entry of the resolution.
Finally, the quotient is obtained by the command @code{modulo}, which
gives for two modules A = ker, B = Im the module of relations of
@tex
$$A/(A \cap B)$$
@end tex
@ifinfo
A/(A intersect B)
@end ifinfo
in the usual way. As we have a chain complex, this is obviously the same
as ker/Im.

We collect these statements in the following short procedure:

@smallexample
proc ext(int n, ideal I)
@{
  resolution rs = mres(I,n+1);
  module tAn    = transpose(rs[n+1]);
  module tAn_1  = transpose(rs[n]);
  module ext_n  = modulo(syz(tAn),tAn_1);
  return(ext_n);
@}
@end smallexample

Now consider the following example:

@smallexample
ring r5 = 32003,(a,b,c,d,e),dp;
ideal I = a2b2+ab2c+b2cd, a2c2+ac2d+c2de,a2d2+ad2e+bd2e,a2e2+abe2+bce2;
print(ext(2,I));
@expansion{} 1,0,0,0,0,0,0,
@expansion{} 0,1,0,0,0,0,0,
@expansion{} 0,0,1,0,0,0,0,
@expansion{} 0,0,0,1,0,0,0,
@expansion{} 0,0,0,0,1,0,0,
@expansion{} 0,0,0,0,0,1,0,
@expansion{} 0,0,0,0,0,0,1
ext(3,I);   // too big to be displayed here
@end smallexample

The library @code{homolog.lib} contains several procedures for computing
Ext-modules and related modules, which are much more general and
sophisticated than the above one. They are used in the following
example:

If @math{M} is a module, then
@tex
$\hbox{Ext}^1(M,M)$, resp.\ $\hbox{Ext}^2(M,M)$,
@end tex
@ifinfo
Ext^1(M,M), resp.@: Ext^2(M,M),
@end ifinfo
are the modules of infinitesimal deformations, respectively of obstructions,
of @math{M} (like T1 and T2 for a singularity).  Similar to the treatment of
singularities, the semiuniversal deformation of @math{M} can be computed (if
@tex
$\hbox{Ext}^1$
@end tex
@ifinfo
Ext^1
@end ifinfo
is finite dimensional) with the help of
@tex
$\hbox{Ext}^1$, $\hbox{Ext}^2$
@end tex
@ifinfo
Ext^1, Ext^2
@end ifinfo
and the cup product.
There is an extra procedure for
@tex
$\hbox{Ext}^k(R/J,R)$
@end tex
@ifinfo
Ext^k(R/J,R)
@end ifinfo
if @math{J} is an ideal in @math{R}, since this is faster than the
general Ext.

We compute
@itemize @bullet
@item
the infinitesimal deformations
@tex
 ($=\hbox{Ext}^1(K,K)$)
@end tex
@ifinfo
 (=Ext^1(K,K))
@end ifinfo
and obstructions
@tex
 ($=\hbox{Ext}^2(K,K)$)
@end tex
@ifinfo
 (=Ext^2(K,K))
@end ifinfo
of the residue field @math{K=R/m} of an ordinary cusp,
@tex
$R=K[x,y]_m/(x^2-y^3)$, $m=(x,y)$.
@end tex
@ifinfo
R=K[x,y]_m/(x^2-y^3), m=(x,y).
@end ifinfo
To compute
@tex
$\hbox{Ext}^1(m,m)$
@end tex
@ifinfo
Ext^1(m,m),
@end ifinfo
we have to apply @code{Ext(1,syz(m),syz(m))} with
@code{syz(m)} the first syzygy module of @math{m}, which is isomorphic to
@tex
$\hbox{Ext}^2(K,K)$.
@end tex
@ifinfo
Ext^2(K,K).
@end ifinfo
@item
@tex
$\hbox{Ext}^k(R/i,R)$
@end tex
@ifinfo
Ext^k(R/i,R)
@end ifinfo
for some ideal @math{i} and with an extra option.
@end itemize

@smallexample
@c example
  LIB "homolog.lib";
  ring R=0,(x,y),ds;
  ideal i=x2-y3;
  qring q = std(i);      // defines the quotient ring k[x,y]_m/(x2-y3)
  ideal m = maxideal(1);
  module T1K = Ext(1,m,m);  // computes Ext^1(R/m,R/m)
  print(T1K);
  printlevel=2;             // gives more explanation
  module T2K=Ext(2,m,m);    // computes Ext^2(R/m,R/m)
  print(std(T2K));
  printlevel=0;
  module E = Ext(1,syz(m),syz(m));
  print(std(E));
  //The matrices which we have just computed are presentation matrices
  //of the modules T2K and E. Hence we may ignore those columns
  //containing 1 as an entry and see that T2K and E are isomorphic
  //as expected, but differently presented.
  //-------------------------------------------
  ring S=0,(x,y,z),dp;
  ideal  i = x2y,y2z,z3x;
  module E = Ext_R(2,i);
  print(E);
  // if a 3-rd argument of type int is given,
  // a list of Ext^k(R/i,R), a SB of Ext^k(R/i,R) and a vector space basis
  // is returned:
  list LE = Ext_R(3,i,0);
  LE;
  print(LE[2]);
  print(kbase(LE[2]));
@c example
@end smallexample
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Depth, Factorization, Computation of Ext, Commutative Algebra
@subsection Depth
@cindex Depth

We compute the depth of the module of Kaehler differentials
@tex
D$_k$(R)
@end tex
@ifinfo
D_k(R)
@end ifinfo
of the variety defined by the @math{(m+1)}-minors of a generic symmetric
@tex
$(n \times n)$-matrix.
@end tex
@ifinfo
(n x n)-matrix.
@end ifinfo
We do this by computing the resolution over the polynomial
ring.  Then, by the Auslander-Buchsbaum formula, the depth is equal to
the number of variables minus the length of a minimal resolution.  This
example was suggested by U.@: Vetter in order to check whether his bound
@tex
$\hbox{depth}(\hbox{D}_k(R))\geq m(m+1)/2 + m-1$
@end tex
@ifinfo
depth(D_k(R)) >= m(m+1)/2 + m-1
@end ifinfo
could be improved.

@smallexample
@c example
  LIB "matrix.lib"; LIB "sing.lib";
  int n = 4;
  int m = 3;
  int N = n*(n+1) div 2;       // will become number of variables
  ring R = 32003,x(1..N),dp;
  matrix X = symmat(n);        // proc from matrix.lib
                               // creates the symmetric generic nxn matrix
  print(X);
  ideal J = minor(X,m);
  J=std(J);
  // Kaehler differentials D_k(R)
  // of R=k[x1..xn]/J:
  module D = J*freemodule(N)+transpose(jacob(J));
  ncols(D);
  nrows(D);
  //
  // Note: D is a submodule with 110 generators of a free module
  // of rank 10 over a polynomial ring in 10 variables.
  // Compute a full resolution of D with sres.
  // This takes about 17 sec on a Mac PB 520c and 2 sec an a HP 735
  int time = timer;
  module sD = std(D);
  list Dres = sres(sD,0);                // the full resolution
  timer-time;                            // time used for std + sres
  intmat B = betti(Dres);
  print(B,"betti");
  N-ncols(B)+1;                          // the desired depth
@c example
@end smallexample
@c  killall();
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Factorization, Primary decomposition, Depth, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Factorization, Primary decomposition, Handling graded modules, Commutative Algebra
@end ifclear
@subsection Factorization
@cindex Factorization

The factorization of polynomials is implemented in the C++ libraries
Factory (written mainly by Ruediger Stobbe) and libfac (written by
Michael Messollen) which are part of the @sc{Singular} system. For
the factorization of univariate polynomials these libraries make use of
the library NTL written by Victor Shoup.

@smallexample
@c example
  ring r = 0,(x,y),dp;
  poly f = 9x16-18x13y2-9x12y3+9x10y4-18x11y2+36x8y4
         +18x7y5-18x5y6+9x6y4-18x3y6-9x2y7+9y8;
  // = 9 * (x5-1y2)^2 * (x6-2x3y2-1x2y3+y4)
  factorize(f);
  // returns factors and multiplicities,
  // first factor is a constant.
  poly g = (y4+x8)*(x2+y2);
  factorize(g);
  // The same in characteristic 2:
  ring s = 2,(x,y),dp;
  poly g = (y4+x8)*(x2+y2);
  factorize(g);
  // factorization over algebraic extension fields
  ring rext = (0,i),(x,y),dp;
  minpoly = i2+1;
  poly g = (y4+x8)*(x2+y2);
  factorize(g);
@c example
@end smallexample
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Primary decomposition, Normalization, Factorization, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Primary decomposition, Kernel of module homomorphisms, Factorization, Commutative Algebra
@end ifclear
@subsection Primary decomposition
@cindex Primary decomposition

There are two algorithms implemented in @sc{Singular} which provide
primary decomposition: @code{primdecGTZ}, based on
Gianni/Trager/Zacharias (written by Gerhard Pfister) and
@code{primdecSY}, based on Shimoyama/Yokoyama (written by Wolfram Decker
and Hans Schoenemann).

The result of @code{primdecGTZ} and @code{primdecSY} is returned as
a list of pairs of ideals,
where the second ideal is the prime ideal and the first
ideal the corresponding primary ideal.

@smallexample
@c example
  LIB "primdec.lib";
  ring r = 0,(a,b,c,d,e,f),dp;
  ideal i= f3, ef2, e2f, bcf-adf, de+cf, be+af, e3;
  primdecGTZ(i);
  // We consider now the ideal J of the base space of the
  // miniversal deformation of the cone over the rational
  // normal curve computed in section *8* and compute
  // its primary decomposition.
  ring R = 0,(A,B,C,D),dp;
  ideal J = CD, BD+D2, AD;
  primdecGTZ(J);
  // We see that there are two components which are both
  // prime, even linear subspaces, one 3-dimensional,
  // the other 1-dimensional.
  // (This is Pinkhams example and was the first known
  // surface singularity with two components of
  // different dimensions)
  //
  // Let us now produce an embedded component in the last
  // example, compute the minimal associated primes and
  // the radical. We use the Characteristic set methods
  // from primdec.lib.
  J = intersect(J,maxideal(3));
  // The following shows that the maximal ideal defines an embedded
  // (prime) component.
  primdecSY(J);
  minAssChar(J);
  radical(J);
@c example
@end smallexample
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Normalization, Kernel of module homomorphisms, Primary decomposition, Commutative Algebra
@subsection Normalization
@cindex Normalization

The normalization will be computed for a reduced ring @math{R/I}. The
result is a list of rings; ideals are always called @code{norid} in the
rings of this list. The normalization of @math{R/I} is the product of
the factor rings of the rings in the list divided out by the ideals
@code{norid}.

@smallexample
@c example
  LIB "normal.lib";
  // ----- first example: rational quadruple point -----
  ring R=32003,(x,y,z),wp(3,5,15);
  ideal I=z*(y3-x5)+x10;
  list pr=normal(I);
  def S=pr[1][1];
  setring S;
  norid;
  // ----- second example: union of straight lines -----
  ring R1=0,(x,y,z),dp;
  ideal I=(x-y)*(x-z)*(y-z);
  list qr=normal(I);
  def S1=qr[1][1]; def S2=qr[1][2];
  setring S1; norid;
  setring S2; norid;
@c example
@end smallexample
@end ifset
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Kernel of module homomorphisms, Algebraic dependence, Normalization, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Kernel of module homomorphisms, Algebraic dependence, Primary decomposition, Commutative Algebra
@end ifclear
@subsection Kernel of module homomorphisms
@cindex Kernel of module homomorphisms

Let @math{A}, @math{B} be two matrices of size
@tex
$m\times r$ and $m\times s$
@end tex
@ifinfo
m x r and m x s
@end ifinfo
over the ring @math{R} and consider the corresponding maps
@tex
$$
R^r \buildrel{A}\over{\longrightarrow}
R^m \buildrel{B}\over{\longleftarrow} R^s\;.
$$
@end tex
@ifinfo

@smallexample
   r   A     m
  R  -----> R
            ^
            |
            |
             s
            R  .
@end smallexample

@end ifinfo
We want to compute the kernel of the map
@tex
$R^r \buildrel{A}\over{\longrightarrow}
R^m\longrightarrow
R^m/\hbox{Im}(B) \;.$
@end tex
@ifinfo

@smallexample
   r   A     m         m
  R  -----> R  -----> R /Im(B) .
@end smallexample

@end ifinfo
This can be done using the @code{modulo} command:
@tex
$$
\hbox{\tt modulo}(A,B)=\hbox{ker}(R^r
\buildrel{A}\over{\longrightarrow}R^m/\hbox{Im}(B)) \; .
$$
@end tex
@ifinfo

@smallexample
                   r   A     m
  modulo(A,B)=ker(R  -----> R /Im(B)) .
@end smallexample

@end ifinfo
More precisely, the output of @code{modulo(A,B)} is a @code{module}
such that the given generating @code{vector}s span the kernel on the
right-hand side.

@smallexample
@c example
  ring r=0,(x,y,z),(c,dp);
  matrix A[2][2]=x,y,z,1;
  matrix B[2][2]=x2,y2,z2,xz;
  print(B);
  def C=modulo(A,B);
  print(C);            // matrix of generators for the kernel
  print(A*matrix(C));  // should be in Im(B)
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Algebraic dependence, Milnor and Tjurina number, Kernel of module homomorphisms, Commutative Algebra
@end ifset
@ifclear singularmanual
@node Algebraic dependence, Milnor and Tjurina number, Kernel of module homomorphisms, Commutative Algebra
@end ifclear
@subsection Algebraic dependence
@cindex Algebraic dependence

Let
@tex
$g$, $f_1$, \dots, $f_r\in K[x_1,\ldots,x_n]$.
@end tex
@ifinfo
g, f_1, @dots{}, f_r in K[x1,@dots{},xn].
@end ifinfo
We want to check whether
@enumerate
@item
@tex
$f_1$, \dots, $f_r$
@end tex
@ifinfo
f_1, @dots{}, f_r
@end ifinfo
are algebraically dependent.

Let
@tex
$I=\langle Y_1-f_1,\ldots,Y_r-f_r \rangle \subseteq
K[x_1,\ldots,x_n,Y_1,\ldots,Y_r]$.
@end tex
@ifinfo

@smallexample
I=<Y_1-f_1,@dots{},Y_r-f_r> subset K[x1,@dots{},xn,Y_1,@dots{},Y_r].
@end smallexample

@end ifinfo
Then
@tex
$I \cap K[Y_1,\ldots,Y_r]$
@end tex
@ifinfo
I intersected with K[Y_1,@dots{},Y_r]
@end ifinfo
are the algebraic relations between
@tex
$f_1$, \dots, $f_r$.
@end tex
@ifinfo
f_1, @dots{}, f_r.
@end ifinfo

@item
@tex
$g \in K [f_1,\ldots,f_r]$.
@end tex
@ifinfo
g in K[f_1,@dots{},f_r].
@end ifinfo

@tex
$g \in K[f_1,\ldots,f_r]$
@end tex
@ifinfo
g in K[f_1,@dots{},f_r]
@end ifinfo
if and only if the normal form of @math{g} with respect to @math{I} and a
block ordering with respect to
@tex
$X=(x_1,\ldots,x_n)$ and $Y=(Y_1,\ldots,Y_r)$ with $X>Y$
@end tex
@ifinfo
X=(x1,@dots{},xn) and Y=(Y_1,@dots{},Y_r) with X>Y
@end ifinfo
is in @math{K[Y]}.
@end enumerate

Both questions can be answered using the following procedure. If the
second argument is zero, it checks for algebraic dependence and returns
the ideal of relations between the generators of the given ideal.
Otherwise it checks for subring membership and returns the normal form
of the second argument with respect to the ideal I.

@smallexample
@c example
  proc algebraicDep(ideal J, poly g)
  {
    def R=basering;         // give a name to the basering
    int n=size(J);
    int k=nvars(R);
    int i;
    intvec v;

    // construction of the new ring:

    // construct a weight vector
    v[n+k]=0;         // gives a zero vector of length n+k
    for(i=1;i<=k;i++)
    {
      v[i]=1;
    }
    string orde="(a("+string(v)+"),dp);";
    string ri="ring Rhelp=("+charstr(R)+"),
                          ("+varstr(R)+",Y(1.."+string(n)+")),"+orde;
                            // ring definition as a string
    execute(ri);            // execution of the string

    // construction of the new ideal I=(J[1]-Y(1),...,J[n]-Y(n))
    ideal I=imap(R,J);
    for(i=1;i<=n;i++)
    {
      I[i]=I[i]-var(k+i);
    }
    poly g=imap(R,g);
    if(g==0)
    {
      // construction of the ideal of relations by elimination
      poly el=var(1);
      for(i=2;i<=k;i++)
      {
        el=el*var(i);
      }
      ideal KK=eliminate(I,el);
      keepring(Rhelp);
      return(KK);
    }
    // reduction of g with respect to I
    ideal KK=reduce(g,std(I));
    keepring(Rhelp);
    return(KK);
  }

  // applications of the procedure
  ring r=0,(x,y,z),dp;
  ideal i=xz,yz;
  algebraicDep(i,0);
  // Note: after call of algebraicDep(), the basering is Rhelp.
  setring r; kill Rhelp;
  ideal j=xy+z2,z2+y2,x2y2-2xy3+y4;
  algebraicDep(j,0);
  setring r; kill Rhelp;
  poly g=y2z2-xz;
  algebraicDep(i,g);
  // this shows that g is contained in i.
  setring r; kill Rhelp;
  algebraicDep(j,g);
  // this shows that g is contained in j.
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@node Singularity Theory, Invariant Theory, Commutative Algebra, Examples
@section Singularity Theory
@cindex Singularity Theory

@ifset singularmanual
@menu
* Milnor and Tjurina number::
* Critical points::
* Polar curves::
* T1 and T2::
* Deformations::
* Invariants of plane curve singularities::
* Branches of space curve singularities::
* Classification of hypersurface singularities::
* Resolution of singularities::
@end menu
@end ifset

@ifclear singularmanual
@menu
* Milnor and Tjurina number::
* Critical points::
* Deformations::
* Invariants of plane curve singularities::
* Resolution of singularities::
@end menu
@end ifclear

@c ----------------------------------------------------------------------------
@node Milnor and Tjurina number, Critical points, Algebraic dependence, Singularity Theory
@subsection Milnor and Tjurina number
@cindex Milnor number
@cindex Tjurina number

The Milnor number, resp.@: the Tjurina number, of a power
series f in
@tex
$K[[x_1,\ldots,x_n]]$
@end tex
@ifinfo
K[[x1,...,xn]]
@end ifinfo
is
@ifinfo
@*      milnor(f) = dim_K(K[[x1,...,xn]]/jacob(f))
@*resp.@:
@*      tjurina(f) = dim_K(K[[x1,...,xn]]/((f)+jacob(f)))
@*where
@end ifinfo
@tex
$$
\hbox{milnor}(f) = \hbox{dim}_K(K[[x_1,\ldots,x_n]]/\hbox{jacob}(f)),
$$
respectively
$$
\hbox{tjurina}(f) = \hbox{dim}_K(K[[x_1,\ldots,x_n]]/((f)+\hbox{jacob}(f)))
$$
where
@end tex
@code{jacob(f)} is the ideal generated by the partials
of @code{f}. @code{tjurina(f)} is finite, if and only if @code{f} has an
isolated singularity. The same holds for @code{milnor(f)} if
K has characteristic 0.
@sc{Singular} displays -1 if the dimension is infinite.

@sc{Singular} cannot compute with infinite power series. But it can
work in
@tex
$\hbox{Loc}_{(x)}K[x_1,\ldots,x_n]$,
@end tex
@ifinfo
Loc_(x)K[x1,...,xn],
@end ifinfo
the localization of
@tex
$K[x_1,\ldots,x_n]$
@end tex
@ifinfo
K[x1,...,xn]
@end ifinfo
at the maximal ideal
@tex
$(x_1,\ldots,x_n)$.
@end tex
@ifinfo
(x1,...,xn).
@end ifinfo
To do this, one has to define a ring with a local monomial ordering
such as ds, Ds, ls, ws, Ws (the second letter 's' referring to power
's'eries), or an appropriate matrix ordering.
@ifset singularmanual
See @ref{Monomial orderings} for a menu of possible orderings.
@end ifset
@ifclear singularmanual
Look at the manual to get information about the possible
monomial orderings in @sc{Singular}, or type @code{help Monomial orderings;}
to get a menu of possible orderings. For further help type, e.g.,
@code{help local orderings;}).
@end ifclear

For theoretical reasons, the vector space dimension computed over the
localization ring coincides with the Milnor (resp. Tjurina) number as
defined above (in the power series ring).

We show in the example below the following:
@itemize @bullet
@item
set option @code{prot} to have a short protocol during standard basis
computation
@item
define the ring @code{r1} of characteristic 32003 with variables @code{x,y,z},
  monomial ordering @code{ds}, series ring (i.e., K[x,y,z] localized at (x,y,z))
@item
list the information about @code{r1} by typing its name
@item
define the integers @code{a,b,c,t}
@item
define a polynomial @code{f} (depending on @code{a,b,c,t}) and display it
@item
define the jacobian ideal @code{i} of @code{f}
@item
compute a standard basis of @code{i}
@item
compute the Milnor number (=250) with @code{vdim} and create and display
  a string in order to comment the result
  (text between quotes "  "; is a 'string')
@item
compute a standard basis of @code{i+(f)}
@item
compute the Tjurina number (=195) with @code{vdim}
@item
then compute the Milnor number (=248) and the Tjurina number
(=195) for @code{t}=1
@item
reset the option to @code{noprot}
@end itemize
See also @ref{sing_lib} for the library commands for the computation of
the Milnor and Tjurina number.

@smallexample
@c example
  option(prot);
  ring r1 = 32003,(x,y,z),ds;
  r1;
  int a,b,c,t=11,5,3,0;
  poly f = x^a+y^b+z^(3*c)+x^(c+2)*y^(c-1)+x^(c-1)*y^(c-1)*z3+
           x^(c-2)*y^c*(y^2+t*x)^2;
  f;
  ideal i=jacob(f);
  i;
  ideal j=std(i);
  "The Milnor number of f(11,5,3) for t=0 is", vdim(j);
  j=i+f;    // override j
  j=std(j);
  vdim(j);  // compute the Tjurina number for t=0
  t=1;
  f=x^a+y^b+z^(3*c)+x^(c+2)*y^(c-1)+x^(c-1)*y^(c-1)*z3
    +x^(c-2)*y^c*(y^2+t*x)^2;
  ideal i1=jacob(f);
  ideal j1=std(i1);
  "The Milnor number of f(11,5,3) for t=1:",vdim(j1);
  vdim(std(j1+f));   // compute the Tjurina number for t=1
  option(noprot);
@c example
@end smallexample
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Critical points, Polar curves, Milnor and Tjurina number, Singularity Theory
@end ifset
@ifclear singularmanual
@node Critical points, Deformations, Milnor and Tjurina number, Singularity Theory
@end ifclear
@subsection Critical points
@cindex Critical points

The same computation which computes the Milnor, resp.@: the Tjurina,
number, but with ordering @code{dp} instead of @code{ds} (i.e., in
@tex
$K[x_1,\ldots,x_n]$
@end tex
@ifinfo
K[x1,...,xn]
@end ifinfo
instead of
@tex
$\hbox{Loc}_{(x)}K[x_1,\ldots,x_n])$
@end tex
@ifinfo
Loc_(x)K[x1,...,xn])
@end ifinfo
gives:
@itemize @bullet
@item
the number of critical points of @code{f} in the affine space
(counted with multiplicities)
@item
the number of singular points of @code{f} on the affine hypersurface @code{f}=0
(counted with multiplicities).
@end itemize

We start with the ring @code{r1} from section @ref{Milnor and Tjurina number} and its elements.

The following will be implemented below:
@itemize @bullet
@item
reset the protocol option and activate the timer
@item
define the ring @code{r2} of characteristic 32003 with variables @code{x,y,z} and monomial
  ordering @code{dp} (= degrevlex) (i.e., the polynomial ring = K[x,y,z]).
@item
Note that polynomials, ideals, matrices (of polys), vectors,
  modules belong to a ring, hence we have to define @code{f} and @code{jacob(f)}
  again in @code{r2}. Since these objects are local to a ring, we may use
  the same names.
  Instead of defining @code{f} again we map it from ring @code{r1} to @code{r2}
  by using the @code{imap} command
  (@code{imap} is a convenient way to map variables
  from some ring identically to variables with the same name in the
  basering, even if the ground field is different. Compare with @code{fetch}
  which works for almost identical rings,
  e.g., if the rings differ only by the ordering or by the names of the
  variables and which may be used to rename variables).
  Integers and strings, however, do not belong to any ring. Once
  defined they are globally known.
@item
The result of the computation here (together with the previous one in
 @ref{Milnor and Tjurina number}) shows that (for @code{t}=0)
@tex
$\hbox{dim}_K(\hbox{Loc}_{(x,y,z)}K[x,y,z]/\hbox{jacob}(f))$
@end tex
@ifinfo
  dim_K(Loc_(x,y,z)K[x,y,z]/jacob(f))
@end ifinfo
= 250 (previously computed) while
@tex
$\hbox{dim}_K(K[x,y,z]/\hbox{jacob}(f))$
@end tex
@ifinfo
  dim_K(K[x,y,z]/jacob(f))
@end ifinfo
= 536. Hence @code{f} has 286 critical points,
  counted with multiplicity, outside the origin.
  Moreover, since
@tex
$\hbox{dim}_K(\hbox{Loc}_{(x,y,z)}K[x,y,z]/(\hbox{jacob}(f)+(f)))$
@end tex
@ifinfo
dim_K(Loc_(x,y,z)K[x,y,z]/(jacob(f)+(f)))
@end ifinfo
= 195 =
@tex
$\hbox{dim}_K(K[x,y,z]/(\hbox{jacob}(f)+(f)))$,
@end tex
@ifinfo
dim_K(K[x,y,z]/(jacob(f)+(f))),
@end ifinfo
the affine surface @code{f}=0 is smooth outside the origin.
@end itemize

@smallexample
@c example
  ring r1 = 32003,(x,y,z),ds;
  int a,b,c,t=11,5,3,0;
  poly f = x^a+y^b+z^(3*c)+x^(c+2)*y^(c-1)+x^(c-1)*y^(c-1)*z3+
           x^(c-2)*y^c*(y^2+t*x)^2;
  option(noprot);
  timer=1;
  ring r2 = 32003,(x,y,z),dp;
  poly f=imap(r1,f);
  ideal j=jacob(f);
  vdim(std(j));
  vdim(std(j+f));
  timer=0;  // reset timer
@c example
@end smallexample
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Polar curves, T1 and T2, Critical points, Singularity Theory
@subsection Polar curves
@cindex Polar curves

The polar curve of a hypersurface given by a polynomial
@tex
$f\in k[x_1,\ldots,x_n,t]$
@end tex
@ifinfo
f in k[x1,...,xn,t]
@end ifinfo
with respect to @math{t} (we may consider @math{f=0} as a family of
hypersurfaces parametrized by @math{t}) is defined as the Zariski
closure of
@tex
$V(\partial f/\partial x_1,\ldots,\partial f/\partial x_n) \setminus V(f)$
@end tex
@ifinfo
V(diff(f,x1),...,diff(f,xn)) \ V(f)
@end ifinfo
if this happens to be a curve.  Some authors consider
@tex
$V(\partial f/\partial x_1,\ldots,\partial f/\partial x_n)$
@end tex
@ifinfo
V(diff(f,x1),...,diff(f,xn))
@end ifinfo
itself as polar curve.

We may consider projective hypersurfaces
@tex
(in $P^n$),
@end tex
@ifinfo
(in P^n),
@end ifinfo
affine hypersurfaces
@tex
(in $k^n$)
@end tex
@ifinfo
(in k^n)
@end ifinfo
or germs of hypersurfaces
@tex
(in $(k^n,0)$),
@end tex
@ifinfo
(in (k^n,0)),
@end ifinfo
getting in this way
projective, affine or local polar curves.

Now let us compute this for a family of curves.  We need the library
@code{elim.lib} for saturation and @code{sing.lib} for the singular
locus.

@smallexample
@c example
  LIB "elim.lib";
  LIB "sing.lib";
  // Affine polar curve:
  ring R = 0,(x,z,t),dp;              // global ordering dp
  poly f = z5+xz3+x2-tz6;
  dim_slocus(f);                      // dimension of singular locus
  ideal j = diff(f,x),diff(f,z);
  dim(std(j));                        // dim V(j)
  dim(std(j+ideal(f)));               // V(j,f) also 1-dimensional
  // j defines a curve, but to get the polar curve we must remove the
  // branches contained in f=0 (they exist since dim V(j,f) = 1). This
  // gives the polar curve set theoretically. But for the structure we
  // may take either j:f or j:f^k for k sufficiently large. The first is
  // just the ideal quotient, the second the iterated ideal quotient
  // or saturation. In our case both coincide.
  ideal q = quotient(j,ideal(f));     // ideal quotient
  ideal qsat = sat(j,f);              // saturation, proc from elim.lib
  ideal sq = std(q);
  dim(sq);
  // 1-dimensional, hence q defines the affine polar curve
  //
  // to check that q and qsat are the same, we show both inclusions, i.e.,
  // both reductions must give the 0-ideal
  size(reduce(qsat,sq));
  size(reduce(q,std(qsat)));
  qsat;
  // We see that the affine polar curve does not pass through the origin,
  // hence we expect the local polar "curve" to be empty
  // ------------------------------------------------
  // Local polar curve:
  ring r = 0,(x,z,t),ds;              // local ordering ds
  poly f = z5+xz3+x2-tz6;
  ideal j = diff(f,x),diff(f,z);
  dim(std(j));                        // V(j) 1-dimensional
  dim(std(j+ideal(f)));               // V(j,f) also 1-dimensional
  ideal q = quotient(j,ideal(f));     // ideal quotient
  q;
  // The local polar "curve" is empty, i.e., V(j) is contained in V(f)
  // ------------------------------------------------
  // Projective polar curve: (we need "sing.lib" and "elim.lib")
  ring P = 0,(x,z,t,y),dp;            // global ordering dp
  poly f = z5y+xz3y2+x2y4-tz6;
                                      // but consider t as parameter
  dim_slocus(f);              // projective 1-dimensional singular locus
  ideal j = diff(f,x),diff(f,z);
  dim(std(j));                        // V(j), projective 1-dimensional
  dim(std(j+ideal(f)));               // V(j,f) also projective 1-dimensional
  ideal q = quotient(j,ideal(f));
  ideal qsat = sat(j,f);              // saturation, proc from elim.lib
  dim(std(qsat));
  // projective 1-dimensional, hence q and/or qsat define the projective
  // polar curve. In this case, q and qsat are not the same, we needed
  // 2 quotients.
  // Let us check both reductions:
  size(reduce(qsat,std(q)));
  size(reduce(q,std(qsat)));
  // Hence q is contained in qsat but not conversely
  q;
  qsat;
  //
  // Now consider again the affine polar curve,
  // homogenize it with respect to y (deg t=0) and compare:
  // affine polar curve:
  ideal qa = 12zt+3z-10,5z2+12xt+3x,-144xt2-72xt-9x-50z;
  // homogenized:
  ideal qh = 12zt+3z-10y,5z2+12xyt+3xy,-144xt2-72xt-9x-50z;
  size(reduce(qh,std(qsat)));
  size(reduce(qsat,std(qh)));
  // both ideals coincide
@c example
@end smallexample
@end ifset
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node T1 and T2, Deformations, Polar curves, Singularity Theory
@subsection T1 and T2
@cindex T1
@cindex T2
@cindex Deformations, T1 and T2

@math{T^1}, resp.@: @math{T^2}, of an ideal @math{j} usually denote the modules of
infinitesimal deformations, resp.@: of obstructions.
In @sc{Singular} there are procedures @code{T_1} and @code{T_2} in
@code{sing.lib} such that
@code{T_1(j)} and @code{T_2(j)} compute a standard basis of
a presentation of these modules.
If @math{T^1,T_2} are finite dimensional K-vector spaces (e.g., for isolated
singularities), a basis can be computed by applying
@code{kbase(T_1(j));}, resp.@: @code{kbase(T_2(j));}, the dimensions by
applying @code{vdim}.
For a complete intersection j the procedure @code{Tjurina} also
computes @math{T^1}, but faster (@math{T^2=0} in this case).
For a non complete intersection, it is faster to use the procedure @code{T_12}
instead of @code{T_1} and @code{T_2}.
Type @code{help T_1;} (or @code{help T_2;} or @code{help T_12;}) to obtain
more detailed information about these procedures.

We give three examples, the first being a hypersurface, the second a complete
intersection, the third not a complete intersection:
@itemize @bullet
@item
load @code{sing.lib}
@item
check whether the ideal j is a complete intersection. It is, if
     number of variables = dimension + minimal number of generators
@item
compute the Tjurina number
@item
compute a vector space basis (kbase) of @math{T^1}
@item
compute the Hilbert function of @math{T^1}
@item
create a polynomial encoding the Hilbert series
@item
compute the dimension of @math{T^2}
@end itemize

@smallexample
@c example
  LIB "sing.lib";
  ring R=32003,(x,y,z),ds;
  // ---------------------------------------
  // hypersurface case (from series T[p,q,r]):
  int p,q,r = 3,3,4;
  poly f = x^p+y^q+z^r+xyz;
  tjurina(f);
  // Tjurina number = 8
  kbase(Tjurina(f));
  // ---------------------------------------
  // complete intersection case (from series P[k,l]):
  int k,l =3,2;
  ideal j=xy,x^k+y^l+z2;
  dim(std(j));          // Krull dimension
  size(minbase(j));     // minimal number of generators
  tjurina(j);           // Tjurina number
  module T=Tjurina(j);
  kbase(T);             // a sparse output of the k-basis of T_1
  print(kbase(T));      // columns of matrix are a k-basis of T_1
  // ---------------------------------------
  // general case (cone over rational normal curve of degree 4):
  ring r1=0,(x,y,z,u,v),ds;
  matrix m[2][4]=x,y,z,u,y,z,u,v;
  ideal i=minor(m,2);   // 2x2 minors of matrix m
  module M=T_1(i);       // a presentation matrix of T_1
  vdim(M);              // Tjurina number
  hilb(M);              // display of both Hilbert series
  bigintvec v1=hilb(M,1);  // first Hilbert series as intvec
  bigintvec v2=hilb(M,2);  // second Hilbert series as intvec
  v1;
  v2;
  v1[3];                // 3rd coefficient of the 1st Hilbert series
  module N=T_2(i);
@c example
@end smallexample
@smallexample
// In some cases it might be useful to have a polynomial in some ring
// encoding the Hilbert series. This polynomial can then be
// differentiated, evaluated etc. It can be done as follows:
ring H = 0,t,ls;
poly h1;
int ii;
for (ii=1; ii<=size(v1); ii=ii+1)
@{
   h1=h1+v1[ii]*t^(ii-1);
@}
h1;                   // 1st Hilbert series
@expansion{} 4-20t+40t2-40t3+20t4-4t5
diff(h1,t);           // differentiate  h1
@expansion{} -20+80t-120t2+80t3-20t4
subst(h1,t,1);        // substitute t by 1
@expansion{} 0

// The procedures T_1, T_2, T_12 may be called with two arguments and then
// they return a list with more information (type help T_1; etc.)
// e.g., T_12(i,<any>); returns a list with 9 nonempty objects where
// _[1] = std basis of T_1-module, _[2] = std basis of T_2-module,
// _[3]= vdim of T_1, _[4]= vdim of T_2
setring r1;           // make r1 again the basering
list L = T_12(i,1);
@expansion{} // dim T_1  =  4
@expansion{} // dim T_2  =  3
kbase(L[1]);          // kbase of T_1
@expansion{} _[1]=1*gen(2)
@expansion{} _[2]=1*gen(3)
@expansion{} _[3]=1*gen(6)
@expansion{} _[4]=1*gen(7)
kbase(L[2]);          // kbase of T_2
@expansion{} _[1]=1*gen(6)
@expansion{} _[2]=1*gen(8)
@expansion{} _[3]=1*gen(9)
L[3];                 // vdim of T_1
@expansion{} 4
L[4];                 // vdim of T_2
@expansion{} 3
@end smallexample
@c killall();            // a procedure from general.lib
@c @expansion{} // ** killing the basering for level 0
@end ifset
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Deformations, Invariants of plane curve singularities, T1 and T2, Singularity Theory
@end ifset
@ifclear singularmanual
@node Deformations, Invariants of plane curve singularities, Critical points, Singularity Theory
@end ifclear
@subsection Deformations
@cindex Deformations

@itemize @bullet
@item
The libraries @code{sing.lib}, respextively @code{deform.lib}, contain
procedures to compute total and base space of the miniversal (= semiuniversal)
deformation of an isolated complete intersection singularity, respectively
of an arbitrary isolated singularity.
@item
The procedure @code{deform} in @code{sing.lib} returns a matrix whose columns
@ifinfo
@code{h_1,..., h_r}
@end ifinfo
@tex
$h_1,\ldots,h_r$
@end tex
represent all 1st order deformations. More precisely, if
@ifinfo
I in R is the ideal generated by @code{f_1,...,f_s}, then any infinitesimal
deformation of R/I over K[e]/(e^2) is given by @code{f+eg},
where f=(f_1,...,f_s), and where g is a K-linear combination of the h_i.
@end ifinfo
@tex
$I \subset R$ is the ideal generated by $f_1,...,f_s$, then any infinitesimal
deformation of $R/I$ over $K[\varepsilon]/(\varepsilon^2)$ is given
by $f+\varepsilon g$,
where $f=(f_1,...,f_s)$, and where $g$ is a $K$-linear combination of
the $h_i$.
@end tex
@item
The procedure @code{versal} in @code{deform.lib} computes a formal
miniversal deformation up to a certain order which can be
prescribed by the user. For a complete intersection the 1st order
part is already miniversal.
@item
The procedure @code{versal} extends the basering to a new ring with
additional deformation parameters which contains the equations for the
miniversal base space and the miniversal total space.
@item
There are default names for the objects created, but the user may also
choose their own names.
@item
If the user sets @code{printlevel=2;} before running @code{versal}, some
intermediate results are shown. This is useful since @code{versal}
is already complicated and might run for some time on more
complicated examples. (type @code{help versal;})
@end itemize

@ifset singularmanual
We compute for the same examples as in the section @ref{T1 and T2}
the miniversal deformations:
@end ifset
@ifclear singularmanual
We give three examples, the first being a hypersurface, the second a
complete intersection, the third no complete intersection and compute
in each of the cases the miniversal deformation:
@end ifclear

@smallexample
@c example
  LIB "deform.lib";
  ring R=32003,(x,y,z),ds;
  //----------------------------------------------------
  // hypersurface case (from series T[p,q,r]):
  int p,q,r = 3,3,4;
  poly f = x^p+y^q+z^r+xyz;
  print(deform(f));
  // the miniversal deformation of f=0 is the projection from the
  // miniversal total space to the miniversal base space:
  // { (A,B,C,D,E,F,G,H,x,y,z) | x3+y3+xyz+z4+A+Bx+Cxz+Dy+Eyz+Fz+Gz2+Hz3 =0 }
  //  --> { (A,B,C,D,E,F,G,H) }
  //----------------------------------------------------
  // complete intersection case (from series P[k,l]):
  int k,l =3,2;
  ideal j=xy,x^k+y^l+z2;
  print(deform(j));
  def L=versal(j);            // using default names
  def Px=L[1]; setring Px;
  show(Px);                   // show is a procedure from inout.lib
  listvar(matrix);
  // ___ Equations of miniversal base space ___:
  Js;
  // ___ Equations of miniversal total space ___:
  Fs;
  // the miniversal deformation of V(j) is the projection from the
  // miniversal total space to the miniversal base space:
  // { (A,B,C,D,E,F,x,y,z) | xy+F+Ez=0, y2+z2+x3+D+Cx+Bx2+Ay=0 }
  //  --> { (A,B,C,D,E,F) }
  //----------------------------------------------------
  // general case (cone over rational normal curve of degree 4):
  kill L;
  ring r1=0,(x,y,z,u,v),ds;
  matrix m[2][4]=x,y,z,u,y,z,u,v;
  ideal i=minor(m,2);                 // 2x2 minors of matrix m
  int time=timer;
  // Call parameters of the miniversal base A(1),A(2),...:
  def L=versal(i,0,"","A(");
  "// used time:",timer-time,"sec";   // time of last command
  def Def_rPx=L[1]; setring Def_rPx;
  Fs;
  Js;
  // the miniversal deformation of V(i) is the projection from the
  // miniversal total space to the miniversal base space:
  // { (A(1..4),x,y,z,u,v) |
  //         -u^2+x*v+A(2)*u+A(4)*v=0, -z*u+y*v-A(1)*u+A(3)*u=0,
  //         -y*u+x*v+A(3)*u+A(4)*z=0,  z^2-y*u+A(1)*z+A(2)*y=0,
  //          y*z-x*u+A(2)*x-A(3)*z=0, -y^2+x*z+A(1)*x+A(3)*y=0 }
  //  --> { A(1..4) |
  //          A(2)*A(4) = -A(3)*A(4) = -A(1)*A(4)+A(4)^2 = 0 }
  //----------------------------------------------------
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Invariants of plane curve singularities, Branches of space curve singularities, Deformations, Singularity Theory
@end ifset
@ifclear singularmanual
@node Invariants of plane curve singularities, Resolution of singularities, Deformations, Singularity Theory
@end ifclear
@subsection Invariants of plane curve singularities
@cindex Puiseux pairs
@cindex Invariants of plane curve singularities

The Puiseux pairs of an irreducible and reduced plane curve singularity are
probably its most important invariants.  They can be computed from its
Hamburger-Noether expansion (which is the analogue of the Puiseux
expansion in characteristic 0 for fields of arbitrary characteristic).

The library @code{hnoether.lib} (see @ref{hnoether_lib}) uses the
algorithm of Antonio Campillo in "Algebroid curves in
positive characteristic" SLN 813, 1980.  This algorithm has the
advantage that it needs least possible field extensions and, moreover,
works in any characteristic. This fact can be used to compute the
invariants over a field of finite characteristic, say 32003, which will
most probably be the same as in characteristic 0.

We compute the Hamburger-Noether expansion of a plane curve
singularity given by a polynomial @code{f} in two variables. This
expansion is given by a matrix, and it allows us to compute a
primitive parametrization (up to a given order) for the curve
singularity defined by @code{f} and numerical invariants such as the
@itemize @bullet
@item
    characteristic exponents,
@item
    Puiseux pairs (of a complex model),
@item
    degree of the conductor,
@item
    delta invariant,
@item
    generators of the semigroup.
@end itemize
Besides commands for computing a parametrization and the invariants
mentioned above, the library @code{hnoether.lib} provides commands for the
computation of the Newton polygon of @code{f}, the square-free part of
@code{f} and a procedure to convert one set of invariants to another.

@smallexample
@c example
  LIB "hnoether.lib";
  // ======== The irreducible case ========
  ring s = 0,(x,y),ds;
  poly f = y4-2x3y2-4x5y+x6-x7;
  list hn = develop(f);
  show(hn[1]);     // Hamburger-Noether matrix
  displayHNE(hn);  // Hamburger-Noether development
  setring s;
  displayInvariants(hn);
  // invariants(hn);  returns the invariants as list
  // partial parametrization of f: param takes the first variable
  // as infinite except the ring has more than 2 variables. Then
  // the 3rd variable is chosen.
  param(hn);
  ring extring=0,(x,y,t),ds;
  poly f=x3+2xy2+y2;
  list hn=develop(f,-1);
  param(hn);       // partial parametrization of f
  list hn1=develop(f,6);
  param(hn1);     // a better parametrization
  // instead of recomputing you may extend the development:
  list hn2=extdevelop(hn,12);
  param(hn2);     // a still better parametrization
  //
  // ======== The reducible case ========
  ring r = 0,(x,y),dp;
  poly f=x11-2y2x8-y3x7-y2x6+y4x5+2y4x3+y5x2-y6;
  // = (x5-1y2) * (x6-2x3y2-1x2y3+y4)
  list L=hnexpansion(f);
  show(L[1][1]);     // Hamburger-Noether matrix of 1st branch
  displayInvariants(L);
  param(L[2]);       // parametrization of 2nd branch
@c example
@end smallexample
@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Branches of space curve singularities, Classification of hypersurface singularities, Invariants of plane curve singularities, Singularity Theory
@subsection Branches of space curve singularities
@cindex Branches of space curve singularities
@cindex Space curve singularities, branches of

In this example, the number of branches of a given quasihomogeneous isolated
space curve singularity will be computed as an example of the pitfalls
appearing in the use of primary decomposition. When dealing with singularities,
two situations are possible in which the primary decomposition algorithm
might not lead to a complete decomposition: first of all, one of the computed
components could be globally irreducible, but analytically reducible
(this is impossible for quasihomogeneous singularities) and,
as a second possibility, a component might be irreducible over the rational
numbers, but reducible over the complex numbers.
@smallexample
@c example
  ring r=0,(x,y,z),ds;
  ideal i=x^4-y*z^2,x*y-z^3,y^2-x^3*z;  // the space curve singularity
  qhweight(i);
  // The given space curve singularity is quasihomogeneous. Hence we can pass
  // to the polynomial ring.
  ring rr=0,(x,y,z),dp;
  ideal i=imap(r,i);
  resolution ires=mres(i,0);
  ires;
  // From the structure of the resolution, we see that the Cohen-Macaulay
  // type of the given singularity is 2
  //
  // Let us now look for the branches using the primdec library.
  LIB "primdec.lib";
  primdecSY(i);
  def li=_[1];
  ideal i2=li[2];       // call the first ideal i1
  // The curve seems to have 2 branches by what we computed using the
  // algorithm of Shimoyama-Yokoyama.
  // Now the same computation by the Gianni-Trager-Zacharias algorithm:
  primdecGTZ(i);
  // Having computed the primary decomposition in 2 different ways and
  // having obtained the same number of branches, we might expect that the
  // number of branches is really 2, but we can check this by formulae
  // for the invariants of space curve singularities:
  //
  // mu = tau - t + 1 (for quasihomogeneous curve singularities)
  // where mu denotes the Milnor number, tau the Tjurina number and
  // t the Cohen-Macaulay type
  //
  // mu = 2 delta - r + 1
  // where delta denotes the delta-Invariant and r the number of branches
  //
  // tau can be computed by using the corresponding procedure T1 from
  // sing.lib.
  setring r;
  LIB "sing.lib";
  T_1(i);
  setring rr;
  // Hence tau is 13 and therefore mu is 12. But then it is impossible that
  // the singularity has two branches, since mu is even and delta is an
  // integer!
  // So obviously, we did not decompose completely. Because the second branch
  // is smooth, only the first ideal can be the one which can be decomposed
  // further.
  // Let us now consider the normalization of this first ideal i1.
  LIB "normal.lib";
  normal(i2);
  def rno=_[1][1];
  setring rno;
  norid;
  // The ideal is generated by a polynomial in one variable of degree 4 which
  // factors completely into 4 polynomials of type T(2)+a.
  // From this, we know that the ring of the normalization is the direct sum of
  // 4 polynomial rings in one variable.
  // Hence our original curve has these 4 branches plus a smooth one
  // which we already determined by primary decomposition.
  // Our final result is therefore: 5 branches.
@c example
@end smallexample
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Classification of hypersurface singularities, Resolution of singularities, Branches of space curve singularities, Singularity Theory
@subsection Classification of hypersurface singularities
@cindex Classification of hypersurface singularities
@cindex Hypersurface singularities, classification of

Classification of isolated hypersurface singularities with respect to
right equivalence is provided by the command @code{classify} of the
library @code{classify.lib}. The classification is done by using the
algorithm of Arnold. Before entering this algorithm, a first guess based
on the Hilbert polynomial of the Milnor algebra is made.

@smallexample
@c example unix_only
  LIB "classify.lib";
  ring r=0,(x,y,z),ds;
  poly p=singularity("E[6k+2]",2)[1];
  p=p+z^2;
  p;
  // We received an E_14 singularity in normal form
  // from the database of normal forms. Since only the residual
  // part is saved in the database, we added z^2 to get an E_14
  // of embedding dimension 3.
  //
  // Now we apply a coordinate change in order to deal with a
  // singularity which is not in normal form:
  map phi=r,x+y,y+z,x;
  poly q=phi(p);
  // Yes, q really looks ugly, now:
  q;
  // Classification
  classify(q);
  // The library also provides routines to determine the corank of q
  // and its residual part without going through the whole
  // classification algorithm.
  corank(q);
  morsesplit(q);
@c example
@end smallexample
@end ifset
@c ----------------------------------------------------------------------------
@node Resolution of singularities
@subsection Resolution of singularities
@cindex Resolution of singularities
@cindex singularities, resolution of

Resolution of singularities and applications thereof are provided by the
libraries @code{resolve.lib} and @code{reszeta.lib}; graphical output may be
generated automatically by using external programs @code{surf} and @code{dot}
respectively to which a specialized interface is provided by the library
@code{resgraph.lib}. In this example, the basic functionality of the
resolution of singularities package is illustrated by the computation of
the intersection matrix and genera of the exceptional curves on a surface
obtained from resolving the A6 surface singularity. A separate tutorial,
which introduces the complete functionality of the package and explains
the rather complicated data structures appearing in intermediate results,
can be found at @url{https://www.singular.uni-kl.de/tutor_resol.pdf}.

@smallexample
@c example
LIB"resolve.lib";                   // load the resolution algorithm
LIB"reszeta.lib";                   // load its application algorithms

ring R=0,(x,y,z),dp;              // define the ring Q[x,y,z]
ideal I=x7+y2-z2;                 // an A6 surface singularity
list L=resolve(I);                // compute the resolution
list iD=intersectionDiv(L);       // compute intersection properties
iD;                               // show the output
// The output is a list whose first entry contains the intersection matrix
// of the exceptional divisors. The second entry is the list of genera
// of these divisors. The third and fourth entry contain the information
// how to find the corresponding divisors in the respective charts.
@c example
@end smallexample

@c ----------------------------------------------------------------------------
@node Invariant Theory, Geometric Invariant Theory, Singularity Theory, Examples
@section Invariant Theory
@cindex Invariant Theory

@ifset singularmanual
@menu
* G_a -Invariants::
* Invariants of a finite group::
@end menu
@end ifset

@c ----------------------------------------------------------------------------
@ifset singularmanual
@node G_a -Invariants, Invariants of a finite group, , Invariant Theory
@subsection G_a -Invariants
@cindex G_a -Invariants

We work in characteristic 0 and use the Lie algebra generated by one
vectorfield of the form
@tex
$\sum x_i \partial /\partial x_{i+1}$.
@end tex
@ifinfo
sum x(i)*d/dx(i+1).
@end ifinfo
@smallexample
@c example
  LIB "ainvar.lib";
  int n=5;
  int i;
  ring s=32003,(x(1..n)),wp(1,2,3,4,5);
  // definition of the vectorfield m=sum m[i,1]*d/dx(i)
  matrix m[n][1];
  for (i=1;i<=n-1;i=i+1)
  {
     m[i+1,1]=x(i);
  }
  // computation of the ring of invariants
  ideal in=invariantRing(m,x(2),x(1),0);
  in;   //invariant ring is generated by 5 invariants
  ring q=32003,(x,y,z,u,v,w),dp;
  matrix m[6][1];
  m[2,1]=x;
  m[3,1]=y;
  m[5,1]=u;
  m[6,1]=v;
  // the vectorfield is: xd/dy+yd/dz+ud/dv+vd/dw
  ideal in=invariantRing(m,y,x,0);
  in; //invariant ring is generated by 6 invariants
@c example
@end smallexample
@c kill n,i,s,q;
@end ifset


@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Invariants of a finite group, Solving systems of polynomial equations, G_a -Invariants, Invariant Theory
@subsection Invariants of a finite group
@cindex Invariants of a finite group

Two algorithms to compute the invariant ring are implemented in
@sc{Singular}, @code{invariant_ring} and @code{invariant_ring_random},
both by Agnes E. Heydtmann (@code{agnes@@math.uni-sb.de}).

Bases of homogeneous invariants are generated successively and those are
chosen as primary invariants that lower the dimension of the ideal
generated by the previously found invariants (see paper "Generating a
Noetherian Normalization of the Invariant Ring of a Finite Group" by
Decker, Heydtmann, Schreyer (J.Symb.Comput.@: 25, No.6, 727-731, 1998).
In the non-modular case secondary invariants are calculated by finding a basis
(in terms of monomials) of the basering modulo the primary invariants,
mapping to invariants with the Reynolds operator and using those or
their power products such that they are linearly independent modulo the
primary invariants (see paper "Some Algorithms in Invariant Theory of
Finite Groups" by Kemper and Steel
(In: Proceedings of the Euroconference in Essen 1997, Birkh@"auser Prog.@: Math.@: 173, 267-285, 1999)).
In the modular case they
are generated according to "Calculating Invariant Rings of Finite Groups
over Arbitrary Fields" by Kemper (J.Symb.Comput.@: 21, No.3, 351-366, 1996).

We calculate now an example from Sturmfels: "Algorithms in Invariant
Theory 2.3.7":

@smallexample
@c example
  LIB "finvar.lib";
  ring R=0,(x,y,z),dp;
  matrix A[3][3]=0,1,0,-1,0,0,0,0,-1;
  // the group G is generated by A in Gl(3,Q);
  print(A);
  print(A*A*A*A); // the fourth power of A is 1
  // Use the first method to compute the invariants of G:
  matrix B(1..3);
  B(1..3)=invariant_ring(A);
  // SINGULAR returns 2 matrices, the first containing
  // primary invariants and the second secondary
  // invariants, i.e., module generators over a Noetherian
  // normalization
  // the third result are the irreducible secondary invariants
  // if the Molien series was available
  print(B(1));
  print(B(2));
  print(B(3));
  // Use the second method,
  // with random numbers between -1 and 1:
  B(1..3)=invariant_ring_random(A,1);
  print(B(1..3));
@c example
@end smallexample
@end ifset


@c ----------------------------------------------------------------------------
@node Geometric Invariant Theory, Non-commutative Algebra,Invariant Theory, Examples
@section Geometric Invariant Theory
@cindex Geometric Invariant Theory

@ifset singularmanual
@menu
* GIT-Fans::
@end menu
@end ifset


@c ----------------------------------------------------------------------------
@ifset singularmanual
@node GIT-Fans, , , Geometric Invariant Theory
@subsection GIT-Fans
@cindex GIT-Fans


Dolgachev/Hu and Thaddeus assigned to an algebraic variety with the action of
an algebraic group the GIT-fan, a polyhedral fan enumerating the GIT-quotients in the sense
of Mumford. The case of the action of an algebraic torus H on an ane variety X has
been treated by Berchtold/Hausen. Based on their construction, an algorithm to compute
the GIT-fan in this setting has been proposed be Keicher. Note that this setting is essential for
many applications, since the torus case can be used to investigate the GIT-variation of the
action of a connected reductive group G.
In many important examples, X is symmetric under the action of a fnite group which
either is known directly from its geometry or can be computed. A prominent
instance is the Deligne-Mumford compactification M06bar of the moduli space of 6-pointed
stable curves of genus zero, which has a natural action of the symmetric group S6. The library gitfan.lib
implements an efficient algorithm for computing GIT-fans, which makes use of symmetries. We have
applied this algorithm to determine the Mori chamber decomposition of the cone of movable divisor classes of M06bar.
Each cone is represented by a single integer. The algorithm relies on Groebner basis techniques, convex geometry and actions of finite symmetry groups. It demonstrates the strength of cross-boarder methods in computer algebra, and the
efficiency of the algorithms implemented in all involved areas. The algorithm is also suitable
for parallel computations.

As an example we address in the following the computation of the GIT-Fan of M05bar.

We first compute the GIT-fan using the single line command provided by the library:

@smallexample
  LIB "gitfan.lib";
  setcores(4);
  ring R = 0,T(1..10),wp(1,1,1,1,1,1,1,1,1,1);
  ideal J =
    T(5)*T(10)-T(6)*T(9)+T(7)*T(8),
    T(1)*T(9)-T(2)*T(7)+T(4)*T(5),
    T(1)*T(8)-T(2)*T(6)+T(3)*T(5),
    T(1)*T(10)-T(3)*T(7)+T(4)*T(6),
    T(2)*T(10)-T(3)*T(9)+T(4)*T(8);
  intmat Q[5][10] =
    1, 1, 1, 1, 0, 0, 0, 0, 0, 0,
    1, 0, 0, 0, 1, 1, 1, 0, 0, 0,
    0, 1, 1, 0, 0, 0, -1, 1, 0, 0,
    0, 1, 0, 1, 0, -1, 0, 0, 1, 0,
    0, 0, 1, 1, -1, 0, 0, 0, 0, 1;
  fan GIT = GITfan(J,Q);
  fVector(GIT);
@end smallexample

The GIT-Fan can be computed using symmetries as follows:

@smallexample
  LIB "gitfan.lib";
  setcores(4);
  ring R = 0,T(1..10),wp(1,1,1,1,1,1,1,1,1,1);
  ideal J =
    T(5)*T(10)-T(6)*T(9)+T(7)*T(8),
    T(1)*T(9)-T(2)*T(7)+T(4)*T(5),
    T(1)*T(8)-T(2)*T(6)+T(3)*T(5),
    T(1)*T(10)-T(3)*T(7)+T(4)*T(6),
    T(2)*T(10)-T(3)*T(9)+T(4)*T(8);
  intmat Q[5][10] =
    1, 1, 1, 1, 0, 0, 0, 0, 0, 0,
    1, 0, 0, 0, 1, 1, 1, 0, 0, 0,
    0, 1, 1, 0, 0, 0, -1, 1, 0, 0,
    0, 1, 0, 1, 0, -1, 0, 0, 1, 0,
    0, 0, 1, 1, -1, 0, 0, 0, 0, 1;
  list simplexSymmetryGroup = G25Action();
  fan GIT2 = GITfan(J,Q,simplexSymmetryGroup);
  GIT2;
@end smallexample

Although we provide a procedure to compute the orbit decomposition of the group action on the simplex of variables
this is not fast in Singular. In the following we describe how to use GAP to obtain the orbit decomposition and then continue
with this data in Singular. This is particularly useful for more complicated examples.

The file orbits.gp in the directory doc of the Singular source tree contains GAP code to do this computation.
This result is provided in the file doc/simplexOrbitRepresentativesG25.sing.

The file doc/simplexSymmetryGroupG25.sing contains the symmetry group (which here is S5).

Moreover the file doc/elementsInTermsOfGeneratorsG25.sing contains a representation of the elements of the symmetry group
in terms of generators.

@smallexample
LIB "gitfan.lib";
setcores(4);
ring R = 0,T(1..10),wp(1,1,1,1,1,1,1,1,1,1);
ideal J =
  T(5)*T(10)-T(6)*T(9)+T(7)*T(8),
  T(1)*T(9)-T(2)*T(7)+T(4)*T(5),
  T(1)*T(8)-T(2)*T(6)+T(3)*T(5),
  T(1)*T(10)-T(3)*T(7)+T(4)*T(6),
  T(2)*T(10)-T(3)*T(9)+T(4)*T(8);
intmat Q[5][10] =
  1, 1, 1, 1, 0, 0, 0, 0, 0, 0,
  1, 0, 0, 0, 1, 1, 1, 0, 0, 0,
  0, 1, 1, 0, 0, 0, -1, 1, 0, 0,
  0, 1, 0, 1, 0, -1, 0, 0, 1, 0,
  0, 0, 1, 1, -1, 0, 0, 0, 0, 1;
intmat Qt = transpose(Q);
<"doc/simplexOrbitRepresentativesG25.sing";
list afaceOrbitRepresentatives=afaces(J,simplexOrbitRepresentatives);
<"doc/simplexSymmetryGroupG25.sing";
list fulldimAfaceOrbitRepresentatives=fullDimImages(afaceOrbitRepresentatives,Q);
list afaceOrbits=computeAfaceOrbits(fulldimAfaceOrbitRepresentatives,simplexSymmetryGroup);
apply(afaceOrbits,size);
list minAfaceOrbits = minimalAfaceOrbits(afaceOrbits);
apply(minAfaceOrbits,size);
list listOfOrbitConeOrbits = orbitConeOrbits(minAfaceOrbits,Q);
apply(listOfOrbitConeOrbits,size);
list listOfMinimalOrbitConeOrbits = minimalOrbitConeOrbits(listOfOrbitConeOrbits);
size(listOfMinimalOrbitConeOrbits);
cone mov = coneViaPoints(transpose(Q));
mov = canonicalizeCone(mov);
list OC =  listOfOrbitConeOrbits;
<"doc/elementsInTermsOfGeneratorsG25.sing";
list Asigmagens = groupActionOnQImage(generatorsG,Q);
list actionOnOrbitconeIndicesForGenerators = groupActionOnHashes(Asigmagens,OC);
list actionOnOrbitconeIndices;
for (int i =1; i<=size(elementsInTermsOfGenerators);i++)
@{
  actionOnOrbitconeIndices[i]=evaluateProduct(actionOnOrbitconeIndicesForGenerators, elementsInTermsOfGenerators[i]);
@}
list OClist = OC[1];
for (i =2;i<=size(OC);i++)
@{
OClist = OClist + OC[i];
@}
list SigmaHashes = GITfanParallelSymmetric(OClist, Q, mov, actionOnOrbitconeIndices);
SigmaHashes;
fan Sigma = hashesToFan(SigmaHashes,OClist);

@end smallexample

Note that the result is not the complete fan but only the fan generated by a minimal set of representatives of maximal cones
for the group action (by the group generated by Asigmagens).

@end ifset




@c ----------------------------------------------------------------------------
@node Non-commutative Algebra, Applications, Geometric Invariant Theory, Examples
@section Non-commutative Algebra
@cindex Non-commutative algebra

@menu
* Left and two-sided Groebner bases::
* Right Groebner bases and syzygies::
@end menu

@node Left and two-sided Groebner bases, Right Groebner bases and syzygies, Non-commutative Algebra, Non-commutative Algebra
@subsection Left and two-sided Groebner bases
@cindex Left and two-sided Groebner bases

For a set of polynomials (resp. vectors) @code{S} in a non-commutative G-algebra, @sc{Singular:Plural} provides two algorithms for computing Groebner bases.

The command @code{std} computes a left Groebner basis of a left module, generated by the set @code{S} (see @ref{std (plural)}).
The command @code{twostd (plural)} computes a  two-sided Groebner basis (which is in particular also a left Groebner basis) of a two-sided ideal, generated by the set @code{S} (see @ref{twostd (plural)}).

In the example below, we consider a particular set @code{S}
in the algebra @math{A:=U(sl_2)} with the degree reverse lexicographic ordering.
We compute a left Groebner basis @code{L} of the left ideal generated by @code{S}
and a two-sided Groebner basis @code{T} of the two-sided ideal generated by @code{S}.
@*Then, we read off the information on the vector space dimension of the
factor modules @code{A/L} and @code{A/T} using the command @code{vdim}
(see @ref{vdim (plural)}).
@*Further on, we use the command @code{reduce}
(see @ref{reduce (plural)}) to compare the left ideals generated
by @code{L} and @code{T}.

We set @code{option(redSB)} and @code{option(redTail)} to make @sc{Singular}
compute completely reduced minimal bases of ideals
(see @ref{option} and @ref{Groebner bases in G-algebras} for definitions
and further details).

For long running computations, it is always
recommended to set @code{option(prot)} to make @sc{Singular} display some information on the performed computations (see @ref{option} for an interpretation of
the displayed symbols).

@smallexample
@c example
// ----- 1.  setting up the algebra
  ring R = 0,(e,f,h),dp;
  matrix D[3][3];
  D[1,2]=-h; D[1,3]=2*e; D[2,3]=-2*f;
  def A=nc_algebra(1,D); setring A;
// ----- equivalently, you may use the following:
//  LIB "ncalg.lib";
//  def A = makeUsl2();
//  setring A;
// ----- 2.  defining the set S
  ideal S = e^3, f^3, h^3 - 4*h;
  option(redSB);
  option(redTail);
  option(prot);  // let us activate the protocol
  ideal L = std(S);
  L;
  vdim(L); // the vector space dimension of the module A/L
  option(noprot); // turn off the protocol
  ideal T = twostd(S);
  T;
  vdim(T);  // the vector space dimension of the module A/T
  print(matrix(reduce(L,T)));  // reduce L with respect to T
  // as we see, L is included in the left ideal generated by T
  print(matrix(reduce(T,L)));  // reduce T with respect to L
  // the non-zero elements belong to T only
  ideal LT = twostd(L); // the two-sided Groebner basis of L
  // LT and T coincide as left ideals:
  size(reduce(LT,T));
  size(reduce(T,LT));
@c example
@end smallexample

@c ----------------------------------------------------------------------------

@node Right Groebner bases and syzygies,  , Left and two-sided Groebner bases, Non-commutative Algebra
@subsection Right Groebner bases and syzygies
@cindex Right Groebner bases and syzygies

Most of the @sc{Singular:Plural} commands correspond to the @emph{left-sided} computations, that is left Groebner bases, left syzygies, left resolutions and so on. However, the @emph{right-sided} computations can be done, using the  @emph{left-sided} functionality and @emph{opposite} algebras.

In the example below, we consider the algebra @math{A:=U(sl_2)} and a set of generators
@tex
$I = \{ e^2, f \}$.
@end tex

We will compute a left Groebner basis @code{LI} and a left syzygy module @code{LS} of a left ideal, generated by the set @math{I}.

Then, we define the opposite algebra @code{Aop} of @code{A}, set it as a basering, and create opposite objects of already computed ones.

Further on, we compute a right Groebner basis @code{RI} and a right syzygy module @code{RS} of a right ideal, generated by the set @math{I} in @math{A}.

@smallexample
@c example
// ----- setting up the algebra:
LIB "ncalg.lib";
def A = makeUsl2();
setring A; A;
// ----- equivalently, you may use
// ring AA = 0,(e,f,h),dp;
// matrix D[3][3];
// D[1,2]=-h; D[1,3]=2*e; D[2,3]=-2*f;
// def A=nc_algebra(1,D); setring A;
option(redSB);
option(redTail);
matrix T;
// --- define a generating set
ideal   I = e2,f;
ideal  LI = std(I); // the left Groebner basis of I
LI;             // we see that I was not a Groebner basis
module LS = syz(I); // the left syzygy module of I
print(LS);
 // check: LS is a left syzygy, if T=0:
T  = transpose(LS)*transpose(I);
print(T);
// --- let us define the opposite algebra of A
def Aop = opposite(A);
setring Aop; Aop;         // see how Aop looks like
// --- we "oppose" (transfer) objects from A to Aop
ideal   Iop = oppose(A,I);
ideal  RIop = std(Iop);  // the left Groebner basis of Iop in Aop
module RSop = syz(Iop);  // the left syzygy module of Iop in Aop
module LSop = oppose(A,LS);
module RLS  = syz(transpose(LSop));
// RLS is the left syzygy of transposed LSop in Aop
// --- let us return to A and transfer (i.e. oppose)
// all the computed objects back
setring A;
ideal  RI = oppose(Aop,RIop); // the right Groebner basis of I
RI;              // it differs from the left Groebner basis LI
module RS = oppose(Aop,RSop); // the right syzygy module of I
print(RS);
 // check: RS is a right syzygy, if T=0:
T = matrix(I)*RS;
T;
module RLS;
RLS = transpose(oppose(Aop,RLS));
// RLS is the right syzygy of a left syzygy of I
// it is I itself ?
print(RLS);
@c example
@end smallexample


@c ----------------------------------------------------------------------------
@node Applications,  , Non-commutative Algebra, Examples
@section Applications
@cindex Applications

@ifset singularmanual
@menu
* Solving systems of polynomial equations::
* AG codes::
@end menu
@end ifset



@ifclear singularmanual
@menu
* Solving systems of polynomial equations::
* AG codes::
@end menu
@end ifclear


@c ----------------------------------------------------------------------------
@ifset singularmanual
@node Solving systems of polynomial equations, AG codes, Invariants of a finite group, Applications
@end ifset
@ifclear singularmanual
@node Solving systems of polynomial equations, AG codes, Resolution of singularities, Singularity Theory
@end ifclear
@subsection Solving systems of polynomial equations
@cindex Symbolic-numerical solving
@cindex Solving systems of polynomial equations

Here we turn our attention to the probably most popular aspect of
the solving problem: given a system of complex polynomial equations with only
finitely many solutions, compute floating point approximations for these
solutions. This is widely considered as a task for
numerical analysis. However, due to rounding errors, purely numerical
methods are often unstable in an unpredictable way.

Therefore, in many cases, it is worth investing more computing power to derive
additional knowledge on the geometric structure of the set of solutions (not
to mention the question of how to decide whether the set of solutions is
finite or not). The symbolic-numerical approach to the solving
problem combines numerical methods with a symbolic preprocessing.

Depending on whether we want to preserve the multiplicities of the solutions
or not, possible goals for a symbolic preprocessing are

@itemize @bullet
@item to find another system of generators (for instance, a reduced Groebner
  basis) for the ideal I generated by the polynomial equations. Alternatively,
 find a system of polynomials defining an ideal which has the same radical
as I (see @ref{Computing Groebner and Standard Bases},
resp.@: @ref{radical}).
@end itemize
In any case, the goal should be to find a system for which a
numerical solution can be found more easily and in a more stable way.
For systems with a large number of generators, the first step in a @sc{Singular}
computation could be to reduce the number of generators by applying
the @code{interred} command (see @ref{interred}). Another goal might be

@itemize @bullet
@item to decompose the system into several smaller (or, at least, more
  accessible) systems of polynomial equations. Then, the set of solutions of
  the original system is obtained by taking the union of the sets of solutions
  of the new systems.
@end itemize
Such a decomposition can be obtained in several ways: for instance, by
computing a triangular decomposition (see @ref{triang_lib}) for the ideal
I, or by applying the factorizing Buchberger algorithm (see
@ref{facstd}),  or by computing a primary decomposition of I
(see @ref{primdec_lib}).

Moreover, the equational modelling of a problem frequently causes unwanted
solutions, for instance, zero as a multiple solution. Not only for stability
reasons, one is frequently interested to get rid of those.
This can be done by computing the saturation of I with respect to
an ideal having the excess components as set of solutions (see @ref{sat}).

The @sc{Singular} libraries @code{solve.lib} and @code{triang.lib} provide
several commands for solving systems of polynomial equations (based on a
symbolic-numerical approach via Groebner bases, resp.@: resultants). In
the example below, we show some of these commands at work.

@smallexample
@c example
LIB "solve.lib";
ring r=0,x(1..5),dp;
poly f0= x(1)^3+x(2)^2+x(3)^2+x(4)^2-x(5)^2;
poly f1= x(2)^3+x(1)^2+x(3)^2+x(4)^2-x(5)^2;
poly f2=x(3)^3+x(1)^2+x(2)^2+x(4)^2-x(5)^2;
poly f3=x(4)^2+x(1)^2+x(2)^2+x(3)^2-x(5)^2;
poly f4=x(5)^2+x(1)^2+x(2)^2+x(3)^2;
ideal i=f0,f1,f2,f3,f4;
ideal si=std(i);
//
// dimension of a solution set (here: 0) can be read from a Groebner bases
// (with respect to any global monomial ordering)
dim(si);
//
// the number of complex solutions (counted with multiplicities) is:
vdim(si);
//
// The given system has a multiple solution at the origin. We use facstd
// to compute equations for the non-zero solutions:
option(redSB);
ideal maxI=maxideal(1);
ideal j=sat(si,maxI);     // output is Groebner basis
vdim(j);                  // number of non-zero solutions (with mult's)
//
// We compute a triangular decomposition for the ideal I. This requires first
// the computation of a lexicographic Groebner basis (we use the FGLM
// conversion algorithm):
ring R=0,x(1..5),lp;
ideal j=fglm(r,j);
list L=triangMH(j);
size(L);                 // number of triangular components
L[1];                    // the first component
//
// We compute floating point approximations for the solutions (with 30 digits)
def S=triang_solve(L,30);
setring S;
size(rlist);             // number of different non-zero solutions
rlist[1];                // the first solution
//
// Alternatively, we could have applied directly the solve command:
setring r;
def T=solve(i,30,1,"nodisplay");  // compute all solutions with mult's
setring T;
size(SOL);               // number of different solutions
SOL[1][1]; SOL[1][2];    // first solution and its multiplicity
SOL[size(SOL)];          // solutions of highest multiplicity
//
// Or, we could remove the multiplicities first, by computing the
// radical:
setring r;
ideal k=std(radical(i));
vdim(k);                 // number of different complex solutions
def T1=solve(k,30,"nodisplay");  // compute all solutions with mult's
setring T1;
size(SOL);               // number of different solutions
SOL[1];
@c example
@end smallexample
@c //
@c // We could also use a resultant-based method to compute the complex
@c // solutions of i (this requires that the system is quadratic, that is,
@c // number of eqns = number of vars):
@c setring r;
@c def T2=ures_solve(i,0,30);  // use sparse resultant
@c setring T1;
@c size(SOL);                  // number of different solutions
@c SOL[1];



@c ----------------------------------------------------------------------------
@ifset singularmanual
@node AG codes, , Solving systems of polynomial equations, Applications
@end ifset
@ifclear singularmanual
@node AG codes, , Solving systems of polynomial equations, Applications
@end ifclear
@subsection AG codes
@cindex AG codes
@cindex coding theory

The library @code{brnoeth.lib} provides an implementation of the Brill-Noether
algorithm for solving the Riemann-Roch problem and applications to Algebraic
Geometry codes. The procedures can be applied to plane (singular) curves
defined over a prime field of positive characteristic.

@smallexample
@c example
LIB "brnoeth.lib";
ring s=2,(x,y),lp;            // characteristic 2
poly f=x3y+y3+x;              // the Klein quartic
list KLEIN=Adj_div(f);        // compute the conductor
KLEIN=NSplaces(1..3,KLEIN);   // computes places up to degree 3
KLEIN=extcurve(3,KLEIN);      // construct Klein quartic over F_8
KLEIN[3];                     // display places (degree, number)
// We define a divisor G of degree 14=6*1+4*2:
intvec G=6,0,0,4,0,0,0,0,0,0,0;   // 6 * place #1 + 4 * place #4
// We compute an evaluation code which evaluates at all rational places
// outside the support of G (place #4 is not rational)
intvec D=2..24;
// in D, the number i refers to the i-th element of the list POINTS in
// the ring KLEIN[1][5].
def RR=KLEIN[1][5];
setring RR; POINTS[1];        // the place in the support of G (not in supp(D))
setring s;
def RR=KLEIN[1][4];
setring RR;
matrix C=AGcode_L(G,D,KLEIN); // generator matrix for the evaluation AG code
nrows(C);
ncols(C);
//
// We can also compute a generator matrix for the residual AG code
matrix CO=AGcode_Omega(G,D,KLEIN);
//
// Preparation for decoding:
// We need a divisor of degree at least 6 whose support is disjoint with the
// support of D:
intvec F=6;                   // F = 6*point #1
// in F, the i-th entry refers to the i-th element of the list POINTS in
// the ring KLEIN[1][5]
list K=prepSV(G,D,F,KLEIN);
K[size(K)][1];                // error-correcting capacity
//
//  Encoding and Decoding:
matrix word[1][11];           // a word of length 11 is encoded
word = 1,1,1,1,1,1,1,1,1,1,1;
def y=word*CO;                // the code word (length: 23)
matrix disturb[1][23];
disturb[1,1]=1;
disturb[1,10]=a;
disturb[1,12]=1+a;
y=y+disturb;                  // disturb the code word (3 errors)
def yy=decodeSV(y,K);         // error correction
yy-y;                         // display the error
@c example
@end smallexample

@c --------------------------------------------------------------------
@ifclear singularmanual
@section Further smallexamples

The example section of the @sc{Singular} manual contains further examples,
e.g.:
@itemize @bullet
@item Long coefficients
@*how they arise in innocent smallexamples
@item T1 and T2
@*compute first order deformations and obstructions
@item Finite fields
@*compute in fields with
@tex
$q=p^n$
@end tex
@ifinfo
q=p^n
@end ifinfo
elements
@item Ext
@*compute Ext groups, derived from the Hom functor
@item Polar curves
@*compute local and global polar curves
@item Depth
@*various ways to compute the depth of a module
@item Cyclic roots
@*create and compute with this standard benchmark smallexample
@item Invariants of finite group
@*compute invariant rings for finite group
@item Normalization
@*compute the normalization of a ring
@item Classification of hypersurface singularities
@*determine type and normal form of a hypersurface singularity
after Arnold
@item Parallelization with ssi links
@*use ssi for distributed and parallel computation
@end itemize

In this list the names of the items are the names of the
examples in the online help system. So by the command
@code{help T1 and T2} the example about the computation of
first order deformations and obstructions is displayed.
@end ifclear