Complexes for version 1-24-11, part 1.

M2/Macaulay2/m2/exports.m2

@@ -829,6 +829,7 @@ export {
 	"isModule",
 	"isMonomialIdeal",
 	"isMutable",
+        "isNormal",
 	"isOpen",
 	"isOutputFile",
 	"isPolynomialRing",
@@ -954,6 +955,7 @@ export {
 	"nextPrime",
 	"nextkey",
 	"norm",
+        "normalCone",
 	"not",
 	"notImplemented",
 	"null",

M2/Macaulay2/m2/shared.m2

@@ -24,6 +24,10 @@ isEmpty = method(TypicalValue => Boolean)
 isSmooth = method(TypicalValue => Boolean, Options => true)
 isVeryAmple = method(TypicalValue => Boolean, Options => true)
 
+isNormal = method()
+
+normalCone = method(Options => true)
+
 -- symbols
 
 protect Base

M2/Macaulay2/packages/AInfinity.m2

@@ -126,7 +126,7 @@ d1d1 := hashTable for i from min B to max B2 -2 list
 
 D := map(labeledTensorComplex{A0,B},B2,d1d1, Degree => -2);
 assert (isComplexMap D);
-m0 := nullHomotopy D;
+m0 := nullHomotopy(D, FreeToExact => true);
 for i from 4 to limit do(
     (C,K) := componentsAndIndices B2_i;
     for k in K do (
@@ -212,7 +212,7 @@ for i from lo+1 to hi list
 	 (A0**G).dd_(i-2)*D_i - D_(i-1)*BG.dd_i
 	 );
 
-m0 := nullHomotopy D;
+m0 := nullHomotopy(D, FreeToExact => true);
 for i from 2 to min(limit, 1+(concentration G)_1) do( --was just 1+(concentration G)_1
     (C,K) := componentsAndIndices BG_i;
     for k in K do (
@@ -386,7 +386,7 @@ p = new MutableHashTable from
 for i from 2 to length G list i=>(presentation R ** G_(i-2))*(BG_i)^[{2,i-2}]
 p#2
 p#3
-nullHomotopy map(G[-2],BG,p)
+nullHomotopy(map(G[-2],BG,p), FreeToExact => true)
 G[-2]
 BG
 
@@ -646,7 +646,7 @@ if n >= 2 then (
     d1d1 := hashTable for i from min B to max B list 
        i+2 => (d1**id_(B_i))*(B2_(i+2))^[{2,i}] - (id_(B_i)**d1)*(B2_(i+2))^[{i,2}];
     D := map(A0**B,B2,d1d1, Degree => -2);
-    m0 := nullHomotopy D;
+    m0 := nullHomotopy(D, FreeToExact => true);
 
     tlist := flatten for i from 4 to 1+(concentration B)_1 list(
         (C,K) := componentsAndIndices B2_i;

M2/Macaulay2/packages/Complexes.m2

@@ -32,6 +32,7 @@ export {
     "connectingExtMap",
     "connectingTorMap",
     "cylinder",
+    "epicResolutionMap",
     "freeResolution",
     "homotopyMap",
     "horseshoeResolution",
@@ -45,10 +46,14 @@ export {
     "isNullHomotopyOf",
     "isShortExactSequence",
     "liftMapAlongQuasiIsomorphism",
+--    "minimalBetti",
     "minimizingMap",
     "nullHomotopy",
+    --"nullhomotopy" => "nullHomotopy",
     "naiveTruncation",
     "randomComplexMap",
+--    "res" => "resolution",
+--    "resolution",
     "resolutionMap",
     "tensorCommutativity",
     "torSymmetry",
@@ -58,6 +63,7 @@ export {
     "yonedaMap'",
     "yonedaProduct",
     -- Option names
+    "FreeToExact", -- used in nullHomotopy
     "OverField",
     "OverZZ",
     "Homogenization",
@@ -112,6 +118,7 @@ load "Complexes/ChainComplex.m2"
 load "Complexes/FreeResolutions.m2"
 load "Complexes/ChainComplexMap.m2"
 load "Complexes/Tor.m2"
+load "Complexes/Ext.m2"
 
 --------------------------------------------------------------------
 -- interface code to legacy types ----------------------------------
@@ -125,7 +132,7 @@ chainComplex Complex := ChainComplex => (cacheValue symbol ChainComplex) (C -> (
     D
     ))
 
-complex ChainComplex := Complex => opts -> (cacheValue symbol Complex)(D -> (
+complex ChainComplex := Complex => {} >> opts -> (cacheValue symbol Complex)(D -> (
     (lo,hi) := (min D, max D);
     while lo < hi and (D_lo).numgens == 0 do lo = lo+1;
     while lo < hi and (D_hi).numgens == 0 do hi = hi-1;
@@ -146,7 +153,7 @@ chainComplex ComplexMap := ChainComplexMap => f -> (
     g
     )
 
-complex ChainComplexMap := ComplexMap => opts -> g -> (
+complex ChainComplexMap := ComplexMap => {} >> opts -> g -> (
     map(complex target g, complex source g, i -> g_i, Degree => degree g)
     )
 --------------------------------------------------------------------
@@ -395,6 +402,9 @@ doc ///
 
 
 S = ZZ/101[a..d]
-K = res coker vars S
+K = freeResolution coker vars S
 L = K ** K
+-- would be nice if these were fast(er):
 elapsedTime L**L;
+elapsedTime (oo ** K)
+elapsedTime (K ** ooo)

M2/Macaulay2/packages/Complexes/ChainComplex.m2

@@ -53,8 +53,10 @@ concentration ComplexMap := Sequence => f -> (
 max Complex := ZZ => C -> max concentration C
 min Complex := ZZ => C -> min concentration C
 
-complex = method(Options => {Base=>0})
-complex HashTable := Complex => opts -> maps -> (
+complexOptions = {Base => 0}
+--complex = method(Options => {Base=>0})
+complex = method(Options => true)
+complex HashTable := Complex => complexOptions >> opts -> maps -> (
     spots := sort keys maps;
     if #spots === 0 then
       error "expected at least one matrix";
@@ -80,7 +82,7 @@ complex HashTable := Complex => opts -> maps -> (
     C.dd = map(C,C,maps,Degree=>-1);
     C
     )
-complex List := Complex => opts -> L -> (
+complex List := Complex => complexOptions >> opts -> L -> (
     -- L is a list of matrices or a list of modules
     if not instance(opts.Base, ZZ) then
       error "expected Base to be an integer";
@@ -104,7 +106,10 @@ complex List := Complex => opts -> L -> (
         );
     error "expected a list of matrices or a list of modules";
     )
-complex Module := Complex => opts -> (M) -> (
+complex Matrix := Complex => complexOptions >> opts -> M -> (
+    complex({M}, opts)
+    )
+complex Module := Complex => complexOptions >> opts -> (M) -> (
     if not instance(opts.Base, ZZ) then
       error "complex: expected base to be an integer";
     if M.cache.?Complex and opts.Base === 0 then return M.cache.Complex;
@@ -118,9 +123,9 @@ complex Module := Complex => opts -> (M) -> (
     C.dd = map(C,C,0,Degree=>-1);
     C
     )
-complex Ring := Complex => opts -> R -> complex(R^1, opts)
-complex Ideal := Complex => opts -> I -> complex(module I, opts)
-complex Complex := Complex => opts -> C -> (
+complex Ring := Complex => complexOptions >> opts -> R -> complex(R^1, opts)
+complex Ideal := Complex => complexOptions >> opts -> I -> complex(module I, opts)
+complex Complex := Complex => complexOptions >> opts -> C -> (
     -- all this does is change the homological degrees 
     -- so the concentration begins at opts.Base
     (lo,hi) := concentration C;
@@ -133,7 +138,7 @@ complex Complex := Complex => opts -> C -> (
         complex(L, Base=>opts.Base)
         )
     )
-complex ComplexMap := Complex => opts -> f -> (
+complex ComplexMap := Complex => complexOptions >> opts -> f -> (
     if degree f === -1 then (
         if source f =!= target f then error "expected a differential";
         (lo,hi) := concentration source f;
@@ -483,6 +488,8 @@ defaultLengthLimit = (R, baselen, len) -> (
       len
     )
 
+-- MES: note, this list of options is all of the ones from resolution, in the Core,
+-- except FastNonminimal is not present (use instead: Strategy => Nonminimal).
 freeResolution = method(Options => {
 	StopBeforeComputation	=> false,
 	LengthLimit		=> infinity,	-- (infinity means numgens R)
@@ -533,6 +540,29 @@ freeResolution Matrix := ComplexMap => opts -> f -> extend(
     matrix f
     )
 
+-- TODO: reinstate these once we remove all uses of ChainComplex...
+-- resolution Module := Complex => opts  -> M -> (
+--     o := pairs opts;
+--     o2 := new OptionTable from select(pairs opts, x -> x#0 =!= FastNonminimal);
+--     if opts.FastNonminimal then (
+--         o2 = o2 ++ {Strategy => Nonminimal};
+--         << "warning: `FastNonminimal => true` is deprecated.  Use: res(..., Strategy => Nonminimal) instead" << endl;
+--         );
+--     freeResolution(M, o2)
+--     )
+-- resolution Ideal := Complex => opts -> I -> resolution(comodule I, opts)
+-- resolution MonomialIdeal := Complex => opts -> I -> resolution(comodule ideal I, opts)
+-- resolution Matrix := ComplexMap => opts -> f -> extend(
+--     resolution(target f, opts), 
+--     resolution(source f, opts),
+--     matrix f
+--     )
+
+complete Complex := C -> C
+complete ComplexMap := F -> F
+nullhomotopy ComplexMap := F -> nullHomotopy F
+status Complex := C -> << "resolution status of a Complex needs to be implemented" << endl;
+
 isHomogeneous Complex := (C) -> isHomogeneous dd^C
 
 -- These next two local functions are lifted from previous code in m2/chaincomplexes.m2
@@ -592,6 +622,11 @@ poincareN Complex := C -> (
     f
     )
 
+rank Complex := ZZ => C -> (
+    (lo, hi) := concentration C;
+    sum for i from lo to hi list (-1)^i * rank C_i
+    )
+
 minimalPresentation Complex := 
 prune Complex := Complex => opts -> (cacheValue symbol minimalPresentation)(C -> (
     -- opts is ignored here
@@ -783,15 +818,15 @@ homomorphism(ZZ, Matrix, Complex) := ComplexMap => (i, f, E) -> (
     if not E.cache.?homomorphism then error "expected target of map to be of the form 'Hom(C,D)'";
     if not isFreeModule source f
     or not rank source f == 1 then error "expected source of map to be free of rank 1";
-    if E_i =!= target f then (
-        -- if f arises from a kernel computation, then the target is not E_i
-        -- it is instead a submodule of E_i.  The next line provides the 'f'
-        -- that maps directly to E_i.
-        -- BUT: if you just use 'ambient f', which seems like it should
-        -- work, the problem is that the target of the map 'ambient f'
-        -- doesn't retain the information about the components of E_i
-        f = map(E_i, source f, super f);
-        );
+    -- we redefine f for two reasons:
+    -- (1) there might be a bug in 'super Matrix' where it gives a module which is === to E_i
+    --     but is missing the direct sum component information.
+    -- (2) If f is a map to the kernel of the differential then we use 'super'
+    --     to make the target exactly (this) E_i (with component info).
+    f = if E_i == target f then
+            map(E_i, source f, f)
+        else
+            map(E_i, source f, super f);
     (C,D) := E.cache.homomorphism;
     (lo,hi) := concentration C;
     H := hashTable for j from lo to hi list j => 
@@ -899,26 +934,65 @@ Complex ** RingMap := Complex => (C, phi) -> tensor(phi, C)
 --------------------------------------------------------------------
 -- resolutions -----------------------------------------------------
 --------------------------------------------------------------------
+-- private function
 nextLambda = method()
-nextLambda ComplexMap := ComplexMap => (lambda) -> (
+-- nextLambda ComplexMap := ComplexMap => (lambda) -> (
+--     C := target lambda;
+--     L0 := source lambda;
+--     (lo,hi) := concentration L0;
+--     D := cone naiveTruncation(lambda, (hi,hi+2), (hi-1, hi));
+--     HC1 := HH_(hi+1) D;
+--     pHC1 := prune HC1;
+--     if pHC1 == 0 then return null;
+--     a1 := inducedMap(pHC1, cover pHC1);
+--     a2 := pHC1.cache.pruningMap;
+--     g1 := map(D_(hi+1), source gens HC1, (gens HC1) // (gens D_(hi+1)));
+--     g2 := map(HC1, source gens HC1, 1);
+--     h := g1 * ((a2 * a1)//g2);
+--     L1 := complex(append(for i from lo+1 to hi list dd^L0_i, h^[0]), Base=>lo);
+--     map(C,L1,i -> if i === hi+1 then -h^[1] else lambda_i)
+--     )
+-- -- private function
+-- nextLambdaEpi = method()
+-- nextLambdaEpi(ComplexMap) := ComplexMap => (lambda) -> (
+--     -- This version is for creating an epimorphism
+--     C := target lambda;
+--     L0 := source lambda;
+--     (lo,hi) := concentration L0;
+--     D := cone naiveTruncation(lambda, (hi,hi+2), (hi-1, hi));
+--     ZC1 := ker dd^D_(hi+1);
+--     --HC1 := HH_(hi+1) D;
+--     pZC1 := prune ZC1;
+--     if pZC1 == 0 then return null;
+--     a1 := inducedMap(pZC1, cover pZC1);
+--     a2 := pZC1.cache.pruningMap;
+--     g1 := map(D_(hi+1), source gens ZC1, (gens ZC1) // (gens D_(hi+1)));
+--     g2 := map(ZC1, source gens ZC1, 1);
+--     h := g1 * ((a2 * a1)//g2);
+--     L1 := complex(append(for i from lo+1 to hi list dd^L0_i, h^[0]), Base=>lo);
+--     map(C,L1,i -> if i === hi+1 then -h^[1] else lambda_i)
+--     )
+-- private function
+nextLambda(ComplexMap, Boolean) := ComplexMap => (lambda, isEpi) -> (
+    -- This version is for creating an epimorphism
     C := target lambda;
     L0 := source lambda;
     (lo,hi) := concentration L0;
     D := cone naiveTruncation(lambda, (hi,hi+2), (hi-1, hi));
-    HC1 := HH_(hi+1) D;
-    pHC1 := prune HC1;
-    if pHC1 == 0 then return null;
-    a1 := inducedMap(pHC1, cover pHC1);
-    a2 := pHC1.cache.pruningMap;
-    g1 := map(D_(hi+1), source gens HC1, (gens HC1) // (gens D_(hi+1)));
-    g2 := map(HC1, source gens HC1, 1);
+    keyModule := if isEpi then ker dd^D_(hi+1) else HH_(hi+1) D;
+    pkeyModule := prune keyModule;
+    if pkeyModule == 0 then return null;
+    a1 := inducedMap(pkeyModule, cover pkeyModule);
+    a2 := pkeyModule.cache.pruningMap;
+    g1 := map(D_(hi+1), source gens keyModule, (gens keyModule) // (gens D_(hi+1)));
+    g2 := map(keyModule, source gens keyModule, 1);
     h := g1 * ((a2 * a1)//g2);
     L1 := complex(append(for i from lo+1 to hi list dd^L0_i, h^[0]), Base=>lo);
     map(C,L1,i -> if i === hi+1 then -h^[1] else lambda_i)
     )
 
-resolutionMap = method(Options => options freeResolution)
-resolutionMap Complex := ComplexMap => opts -> C -> (
+resolutionMapPrivate = method(Options => options freeResolution)
+resolutionMapPrivate(Complex, Boolean) := ComplexMap => opts -> (C, isEpi) -> (
     if opts.LengthLimit < 0 then error "expected a non-negative value for LengthLimit";
     if not C.cache.?resolutionMap
       or C.cache.resolutionMap.cache.LengthLimit < opts.LengthLimit then (
@@ -937,7 +1011,7 @@ resolutionMap Complex := ComplexMap => opts -> C -> (
             local g;
             -- how to implement length limit here.  What does length limit mean?
             while (
-                g = nextLambda f;
+                g = nextLambda(f, isEpi);
                 (len <= hi - lo or g =!= null) and len <= lengthlimit
                 ) do (
                 if g === null then (
@@ -964,6 +1038,12 @@ resolutionMap Complex := ComplexMap => opts -> C -> (
     else fC
     )
 
+resolutionMap = method(Options => options freeResolution)
+resolutionMap Complex := ComplexMap => opts -> C -> resolutionMapPrivate(C, false, opts)
+
+epicResolutionMap = method(Options => options freeResolution)
+epicResolutionMap Complex := ComplexMap => opts -> C -> resolutionMapPrivate(C, true, opts)
+
 resolution Complex := opts -> C -> (
     -- TODO: remove this hack once resolution doesn't have FastNonminimal anymore and is defined in Complexes).
     opts1 := new OptionTable from for k in keys opts list if k === FastNonminimal then continue else k => opts#k;