Index: Singular/LIB/powers.lib
===================================================================
--- Singular/LIB/powers.lib	(revision 1c026142c2c815bbafb6ac5b7a32ca944f9e6be3)
+++ Singular/LIB/powers.lib	(revision 1c026142c2c815bbafb6ac5b7a32ca944f9e6be3)
@@ -0,0 +1,390 @@
+// version string automatically expanded by CVS
+version="$Id: powers.lib,v 1.1 2008-05-29 16:57:53 motsak Exp $";
+category="Linear algebra";
+
+info="
+LIBRARY:   powers.lib  Symmetric and Exterior k-powers
+AUTHOR:    Oleksandr Motsak, email: {U@D}, where U={motsak}, D={mathematik.uni-kl.de}.
+
+KEYWORDS: symmetric power, exterior power
+
+PROCEDURES:
+    SymmetricPower(module,int)     return k-th symmetric power of a matrix
+    SymmetricBasis(int, int)       return a basis of a k-th symmetric power
+
+    ExteriorPower(module,int)      return k-th exterior power of a matrix
+    ExteriorPower(int, int)        return a basis of a k-th exteriorPower power
+       
+";
+    
+
+LIB "general.lib"; // for sort
+LIB "nctools.lib"; // for SuperCommutative
+
+
+// thanks to Oleksandr Iena for his persistence ;-)
+
+proc SymmetricBasis(int n, int k)
+"USAGE:    SymmetricBasis(n, k); n int, k int
+RETURN:   ring: the basis of the k-th symmetric power of a n-dim vector space
+NOTE:     The output consists of a polynomial ring in n variables, together with the ideal called I.
+          The ideal I is the basis of the k-th symmetric power of a n-dim vector space (ordered lexicographically).
+          The char. of the returned ring is the same as of current basis ring or zero.
+SEE ALSO: ExteriorBasis
+KEYWORDS: symmetric basis
+EXAMPLE:  example SymmetricBasis; shows an example"
+{
+  int p = 0;
+  if( nameof(basering) != "basering" )
+  {
+    p = char(basering);
+  }
+
+  ring T = (p), (e(1..n)), dp;
+  ideal I = simplify( NF(maxideal(k), std(0)), 1 + 2 + 8 );
+  int N = ncols(I);
+  I = sort(I)[1]; // lex
+  I = I[N..1];
+
+  export I;
+  return(T);
+}
+example
+{ "EXAMPLE:"; echo = 2;
+
+def r = SymmetricBasis(2, 3); setring r; r;
+I; // basis of 3rd sym. power of a 2-dim v.s.
+kill r;
+
+ring r = (0),(a, b, c, d), dp; r;
+def g = SymmetricBasis(3, 2); setring g; g; I;
+
+kill g, r;
+
+ring r = (32003),(a, b, c, d), dp; r;
+def g = SymmetricBasis(4, 2); setring g; g; I;
+}
+
+proc ExteriorBasis(int n, int k)
+"USAGE:    ExteriorBasis(n, k); n int, k int
+RETURN:   ring: the basis of the k-th exterior power of a n-dim vector space
+NOTE:     The output consists of a polynomial ring in n variables, together with the ideal called I.
+          The ideal I is the basis of the k-th exterior power of a n-dim vector space (ordered lexicographically).
+          The char. of the returned ring is the same as of current basis ring or zero.
+
+SEE ALSO: SymmetricBasis
+KEYWORDS: exterior basis
+EXAMPLE:  example ExteriorBasis; shows an example"
+{
+  int p = 0;
+  if( nameof(basering) != "basering" )
+  {
+    p = char(basering);
+  }
+    
+  ring S = (p), (e(1..n)), dp;
+  def T = SuperCommutative(); setring T;
+  ideal I = simplify( NF(maxideal(k), std(0)), 1 + 2 + 8 );
+  int N = ncols(I);
+  I = sort(I)[1];
+  I = I[N..1 ];
+
+  export I;
+  return(T);
+}
+example
+{ "EXAMPLE:"; echo = 2;
+
+def r = ExteriorBasis(2, 3); setring r; r;
+"Basis: ", I; // basis of 3rd sym. power of a 2-dim v.s.
+kill r;
+
+ring r = (0),(a, b, c, d), dp; r;
+def g = ExteriorBasis(3, 2); setring g; g; "Basis: ", I;
+
+kill g, r;
+
+ring r = (32003),(a, b, c, d), dp; r;
+def g = ExteriorBasis(4, 2); setring g; g; "Basis: ", I;
+}
+
+
+
+proc SymmetricPower(module A, int k)
+"USAGE:    SymmetricPower(A, k); A module, k int
+RETURN:   module: the k-th symmetric power of A
+NOTE:     the chosen bases (ordered lexicographically) and
+          temporary data  may be shown if printlevel is big enough
+SEE ALSO: ExteriorPower
+KEYWORDS: symmetric power
+EXAMPLE:  example SymmetricPower; shows an example"
+{
+  int p = printlevel - voice + 2;
+  def save = basering;
+
+  int M = nrows(A);
+  int N = ncols(A);
+
+  int i, j;
+
+  /////////////////////////////////////////////////////////////////
+  def Tn = SymmetricBasis(N, k); setring Tn;
+  ideal Ink = I;
+  dbprint(p-3, "Temporary Source Ring", basering);
+  dbprint(p, "S^k(Source Basis)", Ink);
+
+  /////////////////////////////////////////////////////////////////
+  def Tm = SymmetricBasis(M, k); setring Tm;
+  ideal Imk = I;
+  dbprint(p-3, "Temporary Image Ring", basering);
+  dbprint(p, "S^k(Image Basis)", Imk);
+
+  /////////////////////////////////////////////////////////////////
+  def Rm = save + Tm;
+  setring Rm;
+  dbprint(p-2, "Temporary Working Ring", Rm);
+
+  module A = imap(save, A);
+
+  ideal B; poly t;
+
+  for( i = N; i > 0; i-- )
+  {
+    t = 0;
+    for( j = M; j > 0; j-- )
+    {
+      t = t + A[i][j] * var( nvars(save) + j); // tensor product!!!
+    }
+
+    B[i] = t;
+  }
+
+  dbprint(p-1, "Matrix of simgle images", B);
+
+  map TMap = Tn, B; ideal C = TMap(Ink); // apply S^k A to Source basis vectors... (Ink)
+  C = NF(C, std(0));
+
+  dbprint(p-1, "Image Matrix: ", C);
+
+  // and write it in Image basis (Imk)
+
+  ideal Imk = imap(Tm, Imk);
+
+  module D; poly lm; vector tt;
+
+  for( i = ncols(C); i > 0; i-- )
+  {
+    t = C[i];
+    tt = 0;
+
+    while( t != 0 )
+    {
+      lm = leadmonom(t);
+      //    lm;
+      for( j = ncols(Imk); j > 0; j-- )
+      {
+        if( lm / Imk[j] != 0 )
+        {
+          tt = tt + (lead(t) / Imk[j]) * gen(j);
+          break;
+        }
+      }
+      t = t - lead(t);
+    }
+
+    D[i] = tt;
+//    tt;
+  }
+
+  // D; print(D);
+
+
+  setring save;
+
+  // basering;
+
+  module AA = imap(Rm, D);
+
+//  nrows(AA) = Nk; // ????
+//  ncols(AA) = Mk;
+
+//  Nk, Mk;
+
+//  AA[Mk] = AA[Mk] * 1;
+
+  return(AA);
+}
+example
+{ "EXAMPLE:"; echo = 2;
+int save = printlevel; printlevel = 1;
+
+ring r = (0),(a, b, c, d), dp; r;
+
+module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1);
+
+print(SymmetricPower(A, 2));
+print(SymmetricPower(A, 3));
+
+module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2);
+
+print(SymmetricPower(B, 2));
+print(SymmetricPower(B, 3));
+
+kill r;
+ring r = (0),(a, b, c, d), dp;def g = SuperCommutative();setring g;
+
+module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1);
+
+print(SymmetricPower(A, 2));
+print(SymmetricPower(A, 3));
+
+module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2);
+
+print(SymmetricPower(B, 2));
+print(SymmetricPower(B, 3));
+
+kill r;
+
+printlevel = save;
+}
+
+
+
+proc ExteriorPower(module A, int k)
+"USAGE:    ExteriorPower(A, k); A module, k int
+RETURN:   module: the k-th exterior power of A
+NOTE:     the chosen bases and temporary data
+          may be shown if printlevel is big enough.
+          Last rows may be invisible if zero.
+SEE ALSO: SymmetricPower
+KEYWORDS: exterior power
+EXAMPLE:  example ExteriorPower; shows an example"
+{
+  int p = printlevel - voice + 2;
+  def save = basering;
+
+  int M = nrows(A);
+  int N = ncols(A);
+
+  int i, j;
+
+  /////////////////////////////////////////////////////////////////
+  def Tn = ExteriorBasis(N, k); setring Tn;
+  ideal Ink = I;
+  dbprint(p-3, "Temporary Source Ring", basering);
+  dbprint(p, "E^k(Source Basis)", Ink);
+
+  /////////////////////////////////////////////////////////////////
+  def Tm = ExteriorBasis(M, k); setring Tm;
+  ideal Imk = I;
+  dbprint(p-3, "Temporary Image Ring", basering);
+  dbprint(p, "E^k(Image Basis)", Imk);
+
+
+  def Rm = save + Tm;
+  setring Rm;
+  dbprint(p-2, "Temporary Working Ring", Rm);
+
+  module A = imap(save, A);
+
+  ideal B; poly t;
+
+  for( i = N; i > 0; i-- )
+  {
+    t = 0;
+    for( j = M; j > 0; j-- )
+    {
+      t = t + A[i][j] * var( nvars(save) + j); // tensor product!!!
+    }
+
+    B[i] = t;
+  }
+
+  dbprint(p-1, "Matrix of simgle images", B);
+
+  map TMap = Tn, B; ideal C = TMap(Ink); // apply S^k A to Source basis vectors... (Ink)
+  C = NF(C, std(0));
+
+  dbprint(p-1, "Image Matrix: ", C);
+
+  // and write it in Image basis (Imk)
+
+  ideal Imk = imap(Tm, Imk);
+
+  module D; poly lm; vector tt;
+
+  for( i = ncols(C); i > 0; i-- )
+  {
+    t = C[i];
+    tt = 0;
+
+    while( t != 0 )
+    {
+      lm = leadmonom(t);
+      //    lm;
+      for( j = ncols(Imk); j > 0; j-- )
+      {
+        if( lm / Imk[j] != 0 )
+        {
+          tt = tt + (lead(t) / Imk[j]) * gen(j);
+          break;
+        }
+      }
+      t = t - lead(t);
+    }
+
+    D[i] = tt;
+//    tt;
+  }
+
+  // D; print(D);
+
+
+  setring save;
+
+  // basering;
+
+  module AA = imap(Rm, D);
+
+//  nrows(AA) = Nk; // ????
+//  ncols(AA) = Mk;
+
+//  Nk, Mk;
+
+//  AA[Mk] = AA[Mk] * 1;
+
+  return(AA);
+}
+example
+{ "EXAMPLE:"; echo = 2;
+  int save = printlevel; printlevel = 1;
+
+  ring r = (0),(a, b, c, d), dp; r;
+
+  module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1);
+
+  print(ExteriorPower(A, 2));
+  print(ExteriorPower(A, 3));
+
+  module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2);
+
+  print(ExteriorPower(B, 2));
+  print(ExteriorPower(B, 3));
+
+  kill r;
+  ring r = (0),(a, b, c, d), dp;def g = SuperCommutative();setring g;
+
+  module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1);
+
+  print(ExteriorPower(A, 2));
+  print(ExteriorPower(A, 3));
+
+  module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2);
+
+  print(ExteriorPower(B, 2));
+  print(ExteriorPower(B, 3));
+
+  kill r;
+
+  printlevel = save;
+}
