newPackage(
    	  "PieriMaps",
   	  Version => "0.4.0",
	  Date => "May 5, 2008",
	  Authors => {{
		    Name => "Steven V Sam",
		    Email => "ssam@berkeley.edu",
		    HomePage => "http://www.ocf.berkeley.edu/~ssam/"
		    }},
	  Headline => "computations of Schur modules and maps between them",
	  DebuggingMode => false
	  )

export({straighten, standardTableaux, dualPart, pieri})

--------------------------------
-- subroutines (not exported) --
--------------------------------

-- All matrices computed are in terms of the ordered bases given by 
-- standardTableaux

-- Input: 
-- List T: a tableau, e.g., {{0,1,2},{2,3}}
-- Output:
-- If T not standard (weakly increasing rows, increasing columns), return 
-- first violating entry (starting from bottom to top, left to right);
-- otherwise return null
isStandard = T -> (
     i := #T-2;
     while i >= 0 do (
	  a := T#i;
	  b := T#(i+1);
	  n := #b;
	  for j from 0 to n-1 do if a#j >= b#j then return (i,j);
	  i = i-1;
	  );
     null
     )

-- Input:
-- List T: a tableau
-- Integer col: specify a column
-- Integer row{1,2}, specify two rows
-- Output:
-- Take all entries in row1 from column col to the end and all entries in
-- row2 from beginning to col, and return all possible permutations of these
-- entries (ignoring that some entries might be equal). The output is given in 
-- the form a hash table where the keys are the resulting tableau and the 
-- values are -1
shuffle = (T, col, row1, row2) -> (
     len1 := #(T#row1);
     len2 := #(T#row2);
     truncatedrow1 := (T#row1)_{0..col-2}; -- grab row1 entries
     truncatedrow2 := (T#row2)_{col..len2-1}; -- grab row2 entries
     L := join((T#row1)_{col-1..len1-1}, (T#row2)_{0..col-1});
     P := permutations L;
     output := {};
     P = apply(P, i-> (for j from 0 to #T-1 list (
		    if j == row1 then sort join(truncatedrow1, i_{0..len1-col})
	       	    else if j == row2 then sort join(i_{len1-col+1..#i-1}, truncatedrow2)
	       	    else T#j)));
     coeff := 0;
     for i in P do if i == T then coeff = coeff + 1;
     for i in P do if i != T then output = append(output, (i, -1 / coeff));
     return hashTable(plus, output);
     )

-- Input:
-- List T: a tableau
-- Output:
-- Writes T as a linear combination of other tableaux T' s.t. T'<T if T is not 
-- standard, otherwise writes T = T. The equalities are expressed in a hash 
-- table which contains tableaux T_i as keys and their values c_i which 
-- represent coefficients: T = c_1T_1 + ... + c_nT_n
towardStandard = T -> (
     x := isStandard T;
     if x === null then return new HashTable from {T=>1};
     H := new MutableHashTable from shuffle(T, x#1+1, x#0, x#0+1);
     if H #? T then (
	  coeff := -(H#T) + 1;
	  remove(H,T);
	  prehash := {};
	  for i in keys H do 
	  prehash = append(prehash, (i, H#i / coeff));
	  return hashTable(prehash)
     	  ) 
     else return new HashTable from H
     )

-- The following functions are based on the Olver map defined in 
-- Andrzej Daszkiewicz, ``On the Invariant Ideals of the Symmetric
-- Algebra S. (V \oplus \wedge^2 V)''
-- The notation used is taken from there

-- Input:
-- List mu: a partition
-- ZZ k: an index
-- Output:
-- Subtract one from the kth (in human count, not computer count) row of mu
subtractOne = (mu, k) -> (
     result := {};
     for i from 0 to #mu-1 when true do
     if i == k-1 then result = append(result, mu_i - 1)
     else result = append(result, mu_i);
     return result
     )

-------------------
-- Main routines --
-------------------

-- Input:
-- List mu: A partition
-- Output:
-- Returns the dual (conjugate) partition of mu
dualPart = method()
dualPart(List) := mu -> (
     result := {};
     for i from 1 to mu#0 when true do (
	  counter := 0;
     	  for j from 0 to #mu-1 when true do if mu#j >= i then counter = counter + 1;
	  result = append(result, counter);
	  );
     return result;
     )

-- Input:
-- Integer dim: number of labels to be used
-- List mu: a partition
-- Output: 
-- All standard tableaux of shape mu and with labels from 0..dim-1
standardTableaux = method()
standardTableaux(ZZ, List) := (dim, mu) -> (
     if #mu == 0 then return {{}};
     output := {};
     otherrows := standardTableaux(dim, remove(mu, 0));
     firstrow := rsort compositions(dim, mu#0);
     for i in firstrow do (
     	  temp := {};
     	  for j from 0 to #i-1 do
     	  for k from 1 to i#j when true do
     	  temp = append(temp,j);
     	  for j in otherrows do (
	       temp2 := prepend(temp,j);
     	       if isStandard(temp2) === null then
     	       output = append(output, temp2);
	       );
     	  );
     return output
     )

-- Input:
-- List t: a tableau
-- MutableHashTable h: a hash table which memoizes recursive calls and 
-- stores the coefficients of the straightening of t into standard tableaux
straighten = method()
straighten(List, MutableHashTable) := (t, h) -> (
     t = apply(t, i -> sort i);
     if h #? t then return null;
     if isStandard(t) === null then 
     h#t = new HashTable from {t => 1};
     
     firstIter := towardStandard(t);
     H := hashTable({}); -- straightening of t
     for i in keys firstIter do (
	  coeff := firstIter#i;
     	  straighten(i, h);
	  temp := {};
	  for j in keys h#i do temp = append(temp, (j, coeff * (h#i)#j));
	  H = merge(H, hashTable(temp), plus);
     	  temp = {};
	  for j in keys H do if H#j != 0 then temp = append(temp,(j,H#j));
	  H = hashTable(temp);
	  );
     h#t = H;
     return null;
     )

-- Input:
-- List mu: A partition
-- ZZ k: a number specifying which row to remove a box from mu
-- Ring P: a polynomial ring k[x_0, ..., x_(n-1)], k a field
-- Output:
-- An GL(n,k)-equivariant map P \otimes S_mu k^n -> P \otimes S_lambda k^n 
-- where lambda is mu with a box removed from the kth row, and S_mu denotes 
-- the Schur functor associated to mu
pieriHelper = method()
pieriHelper(List, ZZ, Ring) := (mu, k, P) -> (
     d := numgens P;
     output := {};
     Sbasis := standardTableaux(d, mu);
     Tbasis := standardTableaux(d, subtractOne(mu, k));
     mu = prepend(0,mu);
     A := {};
     for p from 0 to k when true do A = join(A, apply(subsets(1..(k-1), p), s -> prepend(0, append(s, k))));
     for s in apply(Sbasis, i->prepend({},i)) do (
	  row := {};
     	  H := new HashTable from {};
     	  for J in A do (
	       cJ := 1;
	       for q from 1 to #J-2 when true do cJ = cJ * (mu_(J_q) - mu_k + k - J_q);
     	       h := hashTable({(s, (-1)^#J / cJ)});
	       for i from 0 to #J-2 when true do (
	       	    temp := {};
	       	    for T in keys h do (
		    	 for b from 0 to #(T_(J_(i+1)))-1 when true do (
		    	      U := new MutableList from T;
     	       	    	      U#(J_(i+1)) = remove(U#(J_(i+1)), b);
			      U#(J_i) = append(U#(J_i), (T_(J_(i+1)))_b);
			      temp = append(temp, (new List from U, h#T));
			      );
		    	 );
	       	    h = hashTable(plus, temp);
	       	    );
	       H = merge(H, h, plus);
	       );
     	  H = new MutableHashTable from H;
     	  memo := new MutableHashTable from {};
     	  for T in keys H do (
	       U := apply(remove(T,0), i->sort(i));
	       coeff := H#T;
	       remove(H, T);
	       straighten(U, memo);
	       for i in keys memo#U do 
	       if H #? i then H#i = H#i + coeff * (memo#U)#i * x_((T_0)_0)
	       else H#i = coeff * (memo#U)#i * x_((T_0)_0);
	       );
     	  for t in Tbasis do if H #? t then row = append(row, H#t) else row = append(row, 0);
	  output = append(output, row);
	  );
     return map(P^(#Tbasis), P^{#Sbasis:-1}, transpose output);
     )

-- Input:
-- List mu: A partition
-- List boxes: A list of numbers which specifies the rows that boxes are to be 
-- removed from mu
-- Module V: An n-dimensional vector space over a field k of characteristic 0
-- Output:
-- A GL(V)-equivariant map P \otimes S_mu V -> P \otimes S_lambda V where
-- P = k[x_0, ..., x_(n-1)], S_mu denotes the Schur functor associated to mu, 
-- and lambda is the partition obtained from mu by deleting boxes_0, boxes_1,
-- ..., in order.
pieri = method()
pieri(List, List, Module) := (mu, boxes, V) -> (
     R := ring V;
     d := rank V;
     P := R[x_0 .. x_(d-1)];
     result := pieriHelper(mu, boxes#0, P);
     for b from 1 to #boxes-1 do (
	  mu = subtractOne(mu, boxes#(b-1));
     	  result = pieriHelper(mu, boxes#b, P) * result;
	  );
     return result;
     )

-------------------
-- Documentation --
-------------------

beginDocumentation()
document {
     Key => PieriMaps,
     Headline => "Methods for computing Schur functors and Pieri inclusions",
     PARA {
	  "Some references:"
	  },
     UL {
	  TEX "Andrzej Daszkiewicz, On the Invariant Ideals of the Symmetric Algebra $S.(V \\oplus \\wedge^2 V)$, J. Algebra 125, 1989, 444-473.",
	  "William Fulton, Young Tableaux: With Applications to Representation Theory and Geometry, London Math. Society Student Texts 35, 1997.",
	  "Jerzy Weyman, Cohomology of Vector Bundles and Syzygies, Cambridge University Press, 2002."
	  }
     }

document {
     Key => {standardTableaux, (standardTableaux, ZZ, List)},
     Headline => "list all standard tableaux of a certain shape with bounded labels",
     Usage => "standardTableaux(dim, mu)",
     Inputs => { 
	  "dim" => {ofClass ZZ, ", number of labels to be used"},
	  "mu" => {ofClass List, ", a partition (mu_1, ..., mu_r) which gives the shape, where mu_i is the number of boxes in the ith row"}
	  },
     Outputs => {
	  List => {"list of all standard tableaux"}
	  },
     EXAMPLE lines ///
     	  standardTableaux(3, {2,2}) -- lists all standard tableaux on the square with entries 0,1,2
	  ///
     }

document {
     Key => {dualPart, (dualPart, List)},
     Headline => "give the dual (conjugate) of a partition",
     Usage => "dualPart(mu)",
     Inputs => {
	  "mu" => {ofClass List, ", a partition (mu_1, ..., m_r) where mu_i is the number of boxes in the ith row"}
	  },
     Outputs => {
	  List => {"the dual partition to mu"}
	  },
     EXAMPLE lines ///
     	  dualPart({4,3,1}) -- returns {3,2,2,1}
	  dualPart({5,5,2,2}) -- returns {4,4,2,2,2}
	  dualPart(dualPart({7,5,1})) -- returns {7,5,1}
	  ///
     }

document {
     Key => {straighten, (straighten, List, MutableHashTable)},
     Headline => "computers straightening of a tableau",
     Usage => "straighten(List, MutableHashTable)",
     Inputs => {
	  "t" => {ofClass List, ", a tableau to straighten; a tableau looks like {{3,4}, {1,2}} for example, where we list the entries from left to right, top to bottom"},
	  "h" => {ofClass MutableHashTable, ", where the answers are to be stored (any data already in h may be used)"}
	  },
     "Outputs nothing, just modifies h. When looking up values, remember that the keys are stored with rows weakly increasing.",
     EXAMPLE lines ///
     	  h = new MutableHashTable from {}
	  straighten({{3,4}, {1,2}}, h)
	  h#{{3,4}, {1,2}} -- get the coefficients
	  ///
     }

document {
     Key => {pieri, (pieri, List, List, Module)},
     Headline => "computes a matrix representation for a Pieri inclusion of representations of a general linear group",
     Usage => "pieri(List, List, Module)",
     Inputs => {
	  "mu" => {ofClass List, "a partition (mu_1, ..., m_r) where mu_i is the number of boxes in the ith row"},
	  "boxes" => {ofClass List, "" },
	  "V" => {ofClass Module, "an n-dimensional vector space over a field k of characteristic 0" }
	  },
     Outputs => {
	  Matrix => {"A GL(V)-equivariant map P otimes S_mu V -> P otimes S_lambda V where P = k[x_0, ..., x_(n-1)], S_mu denotes the Schur functor associated to mu, and lambda is the partition obtained from mu by deleting boxes_0, boxes_1, ..., in order."}
	  },
     EXAMPLE lines ///
     	  pieri({3,1}, {1}, QQ^3) -- removes the first box from the partition {3,1}
	  res coker oo -- resolve this map
	  betti oo -- check that the resolution is pure
	  ///
     }

TEST ///
assert ( #(standardTableaux(3,{3,2})) == 15 )
///

TEST ///
h = new MutableHashTable from {}
straighten({{1,2,3},{1,1}},h)
h = new MutableHashTable from {}
straighten({{3,4},{1,2}},h)
h = new MutableHashTable from {}
straighten({{1,4},{2,3}},h)
///

TEST ///
t = new BettiTally from {(0,{0},0)=>8, (1,{1},1) =>15, (2,{3},3)=>10, (3,{5},5)=>3};
assert ( t == (betti res coker pieri({3,1},{1},QQ^3)) )
t = new BettiTally from {(0,{0},0)=>15, (1,{1},1) =>24, (2,{4},4)=>15, (3,{6},6)=>6};
assert ( t == (betti res coker pieri({4,1},{1},QQ^3)) )
t = new BettiTally from {(0,{0},0)=>20, (1,{2},2) =>60, (2,{3},3)=>64, (3,{4},4)=>20};
assert ( t == (betti res coker pieri({3,2},{2,2},QQ^4)) )
///

end
loadPackage "PieriMaps"
installPackage PieriMaps