--OSalgebra (page 4 of lecture)
E1= matrix{{1,0,0},{0,1,0},{1,-1,0},{0,0,1}}
A3= matrix{{1,0,0},{1,1,0},{0,1,0},{1,-1,1},{0,1,-1},{1,0,1}}
ds = (M)->(P:=M**(ZZ/101[z]);
           c := rank source P;
           r := rank target P;
           out={};
           scan(c, i->(c1 = subsets(r,i+2);
                       scan(c1,j->(c2 = submatrix(P,j,);
                       if rank c2 < (i+2) then out = append(out,j)))));
           out)

--take a matrix, and returns a list of the dependencies. Rows = H_i
--example for A3

osbd = (L) -> (b1:=#L;
               b2:=reverse subsets(L,b1-1);
               b3:=0;
               b4:=vars R;
               scan(b1, i->(
                   b3=b3+((-1)^i)*
(det(diagonalMatrix transpose submatrix(b4,,b2#i),Strategy=>Cofactor))));
               b3)
--take a list, e.g. {1,2,3} and return the boundary map of the
--corresponding monomial e.g. x_2*x_3-x_1*x_3+x_1*x_2.

osbds = (M)->(I := ideal (matrix{{0}}**R);
              scan(M, i->(I = I + ideal osbd(i)));
              ideal mingens I)
--Return the relations of the OS algebra. 

osideal = (M)->(t0=min mingle M;
                t1=max mingle M;
                R=ZZ/31991[x_(t0)..x_(t1),SkewCommutative=>true];
                J = transpose mingens osbds M)
--return the OS ideal of a hyperplane arrangement, input = dep sets. 
--------------------------------------------------------------------------------------------------
--Example 3 
U=ZZ/31991[a_1..a_4]
d1 = transpose vars U
d2=matrix {{a_2+a_3 , -a_1     , -a_1 , 0},
           {a_4     , 0        , 0 , -a_1   },
           {-a_2    ,  a_1+a_3 , -a_2 , 0},
           {0       , a_4      , 0 , -a_2      },
           {0       , 0        , a_4 , -a_3}}
A = chainComplex(d2, d1)
A.dd^2 == 0
HH_1(A)     
I3 = minors(3,d2)   -- rank <=2 => ker >= 2.
radical I3          --the R^1 locus.

--Example 4
I = ideal osideal ds E1
rI =res (I, DegreeLimit=>5, LengthLimit=>5)
betti rI
rI.dd
-------------------------------------------------------------------------------------------------
osalg = (M)->(t0=min mingle M;
              t1=max mingle M;
              S1= ZZ/31991[a_(t0)..a_(t1)];
              S = frac S1;
              R=S[x_(t0)..x_(t1),SkewCommutative=>true];
              J = (osbds M);
              --print transpose gens J;
              R/J)

oshom = (M)->(A=osalg(M); 
              out={};
              i=0;
              f = map(A^1,A^{-1},(vars R * transpose vars S1)**A);
              while (hilbertFunction(i,A) =!= 0) do 
              (out=append(out, entries inducedMap(image super basis(i+1,target f), image super basis(i+1,source f), f));
              i=i+1);
              out)
--this version of oshom returns a list of lists, with the ith element
--the list version of the matrix from OS_i->OS_{i+1}
--modified to work for any arrangement, tested, works, HKS 1/13/05


EPY = (M)->(osmaps = oshom M;
             l = #osmaps-1; --number of maps in EPY complex.
             s = map(S,A,vars S);
             sprime = map(S1,S,vars S1);
             f0=sprime(s(matrix osmaps#(l-1))); --F(A) is coker of last map
             FA = coker f0**S1^{l})             --get grading right
 
Aomoto=(M)->(res EPY M)
----------------------------------------------------------------------------
--Example 8: A3
B=Ext^2(EPY(ds A3),S1)                                
Delta = transpose relations prune B                   --presentation for B
scan(primaryDecomposition annihilator B, i->print i)  --V(ann(B))=R^1(A)
hilbertSeries B                                       --hilbSeries
scan(10, i->print hilbertFunction(i,B))               --see it!
--MAGIC TRICK1!
I = ideal osideal ds A3
betti res (I, DegreeLimit=>5, LengthLimit=>5)
--MAGIC TRICK2!
AA3 = Aomoto ds A3
--make screen wide:
AA3.dd
FA = coker AA3.dd_1
betti res FA
--The Aomoto Complex is a free S-resolution for the cokernel
--of the last map!!!
hilbertSeries(Ext^2(FA,S1))  
--the module B is an Ext module!
---------------------------------------------------------------------------
--Example 13 
V=ZZ/31991[e_0,e_1, SkewCommutative=>true]
d1 = matrix{{e_0},{e_1}}
d2 = map(V^1,V^2, matrix{{0,e_1}}) --trick for degrees (else degs = 0,1).
d3=matrix{{e_1}}
B = chainComplex(d3,d2,d1)
B.dd^2 == 0
--but, why do this by hand? code!!!
---------------------------------------------------------------------------
bgg= (i,M,E)->(
     S:=ring(M);
     numvarsE := rank source vars E;
     ev := map(E,S,vars E);
     f0:= basis(i,M);
     f1:= basis(i+1,M);
     g:=((vars S)**f0)//f1;
     b:=(ev g)*((transpose vars E)**(ev source f0));
     map(E^{(rank target b):i+1},E^{(rank source b):i},b))

RM = chainComplex apply(4, i->bgg(3-i,M,V))
--this builds the chain complex we built by hand above!
scan(4, i->print prune HH_i(RM)) -- of course first is nonzero, we 
                                 --need to keep going. Interesting
                                 --one is HH_3(RM).
print hilbertFunction(-2,HH_3(RM));
print hilbertFunction(-1,HH_3(RM));  -- =1
print hilbertFunction(0,HH_3(RM));   -- =2
print hilbertFunction(1,HH_3(RM));

-- H^1_(1+1)=(E(3)/e0e1)_2 = 2 = Tor_1(M,C)_2 = 2 
-- H^1_(2+1)=(E(3)/e0e1)_3 = 1 = Tor_2(M,C)_3 = 1
-----------------------------------------------------------------------------             
--Final example: Tate resoln of twisted cubic
LastE=ZZ/31991[e_0..e_3, SkewCommutative=>true]
LastS=ZZ/31991[x_0..x_3]
I=minors(2,matrix{{x_0,x_1,x_2},{x_1,x_2,x_3}})
M=coker gens I
diffsRM = apply(3, i->bgg(3-i,M,LastE))
RM = chainComplex diffsRM
K = ker RM.dd_2
splicePart = res(K, DegreeLimit=>6, LengthLimit=>3)
diffsSpliceP = apply(3, i->splicePart.dd_(i+1))
Tdiffs = join(diffsRM, diffsSpliceP)
TateR = chainComplex(Tdiffs)
betti TateR
--i31 : betti TateR
--
--              0  1 2 3 4 5 6
--o31 = total: 13 10 7 4 3 5 8
--         -4: 13 10 7 4 1 . .
--         -3:  .  . . . 2 5 8