Center of skew rings and kernel of RingMap to skew rings

M2/Macaulay2/m2/center.m2

@@ -0,0 +1,52 @@
+-*
+computes the center of a ring R and returns the map from center(R) -> R
+
+maybe this should just be called "center" but the term is ubiquitous across
+mathematics and appears in several different. trying to merge them all into a
+single method runs into the awkwardness of method options leaking out to all
+the different overrides.
+*-
+
+center = method(Options => {Strategy => null})
+center(Ring) := RingMap => opts -> (R) -> R.cache.center ??= (
+    z := runHooks((center, Ring), (opts, R), Strategy => opts.Strategy);
+    if z =!= null then z else error "center: no strategy for computing center of the given ring"
+)
+
+addHook((center, Ring), Strategy => "isCommutative", (opts, R) -> if isCommutative R then map(R, R) else null)
+addHook((center, Ring), Strategy => "isSkewAffineRing", (opts, R) -> (
+    -- David Eisenbud observed this strategy for computing the center of a skew polynomial ring
+    if not isPolynomialRing R then return null;
+    if not isSkewAffineRing R then return null;
+
+    (R', phi) := flattenRing R;
+    -- get the natural prime coefficient ring
+    kk := coefficientRing R';
+
+    x := symbol x; -- for pairs of skew variables
+    y := symbol y; -- for commuting variables
+    skews := R'.SkewCommutative;
+    n := #skews;
+    xvars := flatten for i from 0 to n - 1 list
+        for j from i + 1 to n - 1 list
+            x_{i, j};
+    ximages := flatten for i from 0 to n - 1 list
+        for j from i + 1 to n - 1 list
+            R'_(skews_i) * R'_(skews_j);
+
+    commuters := toList (0..numgens R' - 1) - set skews;
+    yvars := toList(y_0..y_(#commuters - 1));
+    yimages := apply(commuters, i -> R'_i);
+    S := kk[xvars, yvars, Degrees => apply(ximages| yimages, degree)];
+
+    -- hack to lazy-load PushForward package only when this hook needs it
+    PushForward := needsPackage "PushForward";
+    pushFwd := value PushForward.Dictionary#"pushFwd";
+
+    -- the ideal defining the center is the annihilator of R^1 as a module over the polynomial ring S
+    K := annihilator pushFwd(map(R', S, ximages | yimages), R'^1);
+    Z := S/K;
+    Z.formation = FunctionApplication {center, R};
+
+    map(R, Z, phi^-1 matrix {ximages | yimages})
+))

M2/Macaulay2/m2/loadsequence

@@ -80,6 +80,7 @@ monomcurve.m2
 fano.m2
 schubert.m2
 local.m2
+center.m2
 
 packages.m2
 document.m2

M2/Macaulay2/m2/ringmap.m2

@@ -238,126 +238,144 @@ kernel RingMap := Ideal => opts -> (cacheValue (symbol kernel => opts)) (f -> (
 -- This is a map from method keys to strategy hash tables
 algorithms := new MutableHashTable from {}
 algorithms#(kernel, RingMap) = new MutableHashTable from {
+    "Center" => (opts, f) -> (
+        -- exploit the fact that any map from a commutative ring to a
+        -- non-commutative ring factors through the center.
+        (R, S) := (target f, source f);
+        if not isCommutative S then return null;
+        if isCommutative R then return null;
+        inc := try(center R) else return null;
+
+        C := source inc;
+        G := graphRing inc;
+
+        -- include R into G
+        includeR := map(G, R, submatrix(vars G, toList(0..(numgens R) - 1)));
+        -- project from G to C
+        projectToCenter := map(C, G, map(C^1, C^(numgens R), 0) | vars C);
+        kernel map(C, S, projectToCenter * (includeR * f))
+    ),
+
     FractionField => (opts, f) -> (
-		(F, R) := (target f, source f);
-		C := last F.baseRings;
-		if not instance(F, FractionField)
-			or not coefficientRing R === (k := coefficientRing C)
-			or not(isPolynomialRing C or isQuotientOf(PolynomialRing, C))
-			or not(isPolynomialRing R or isQuotientOf(PolynomialRing, R))
-		then return null;
-
-		prs := presentation C;
-		B := ring prs;
-		images := apply(generators R, x -> (
-			w := f x;
-			new Divide from {numerator w, denominator w}
-		));
-
-		-- now make a common denominator for all images
-		images = new MutableList from images;
-		i := 1;
-		while i < #images do (
-			z := syz(
-				matrix{{denominator images#0,denominator images#i}},
-				SyzygyLimit => 1
-			);
-			a := -z_(0,0);
-			b := z_(1,0);
-
-			j := 0;
-			while j < i do (
-				images#j = apply(images#j, s -> s*a);
-				j += 1;
-			);
-			images#i = apply(images#i, s -> s*b);
-			i += 1;
-		);
-		images = toList images;
-		commonDenominator := images#0#1;
-
-		d := symbol d;
-		h := symbol h;
-		x := symbol x;
-		y := symbol y;
-		(n1, n2) := (numgens F, numgens R);
-		S := k[x_1 .. x_n1, d, y_1 .. y_n2, h,
-			MonomialOrder => Eliminate (n1 + 1),
-			Degrees => join(degrees C, {{1}}, degrees R, {{1}})
-		];
-
-		in1 := map(S,C,matrix {take (generators S, n1)});
-		in2 := map(S,B,matrix {take (generators S, n1)});
-		in3 := map(S,R,matrix {take (generators S, {n1 + 1, n1 + n2})});
-		back := map(R,S,map(R^1,R^(n1 + 1),0) | vars R | 1 );
-		ideal back selectInSubring(
-			1,
-			generators gb(
-				in2 prs | homogenize (
-					in3 vars source in3 - d * in1 matrix {apply(images, first)}
-					| d * in1 commonDenominator - 1,
-					h
-				),
-				Strategy => LongPolynomial,
-				opts
-			)
-		)
-	),
+        (F, R) := (target f, source f);
+        C := last F.baseRings;
+        if not instance(F, FractionField)
+            or not coefficientRing R === (k := coefficientRing C)
+            or not(isPolynomialRing C or isQuotientOf(PolynomialRing, C))
+            or not(isPolynomialRing R or isQuotientOf(PolynomialRing, R))
+        then return null;
+
+        prs := presentation C;
+        B := ring prs;
+        images := apply(generators R, x -> (
+            w := f x;
+            new Divide from {numerator w, denominator w}
+        ));
+
+        -- now make a common denominator for all images
+        images = new MutableList from images;
+        i := 1;
+        while i < #images do (
+            z := syz(
+                    matrix{{denominator images#0,denominator images#i}},
+                    SyzygyLimit => 1
+            );
+            a := -z_(0,0);
+            b := z_(1,0);
+
+            j := 0;
+            while j < i do (
+                    images#j = apply(images#j, s -> s*a);
+                    j += 1;
+            );
+            images#i = apply(images#i, s -> s*b);
+            i += 1;
+        );
+        images = toList images;
+        commonDenominator := images#0#1;
+
+        d := symbol d;
+        h := symbol h;
+        x := symbol x;
+        y := symbol y;
+        (n1, n2) := (numgens F, numgens R);
+        S := k[x_1 .. x_n1, d, y_1 .. y_n2, h,
+            MonomialOrder => Eliminate (n1 + 1),
+            Degrees => join(degrees C, {{1}}, degrees R, {{1}})
+        ];
+
+        in1 := map(S,C,matrix {take (generators S, n1)});
+        in2 := map(S,B,matrix {take (generators S, n1)});
+        in3 := map(S,R,matrix {take (generators S, {n1 + 1, n1 + n2})});
+        back := map(R,S,map(R^1,R^(n1 + 1),0) | vars R | 1 );
+        ideal back selectInSubring(
+            1,
+            generators gb(
+                in2 prs | homogenize (
+                    in3 vars source in3 - d * in1 matrix {apply(images, first)}
+                    | d * in1 commonDenominator - 1,
+                    h
+                ),
+                Strategy => LongPolynomial,
+                opts
+            )
+        )
+    ),
 
     "AffineRing" => (opts, f) -> (
-		(F, R) := (target f, source f);
-		if not isAffineRing R
-			or not isAffineRing F
-			or not instance(ambient R, PolynomialRing)
-			or not instance(ambient F, PolynomialRing)
-			or not isField coefficientRing R
-			or not coefficientRing R === coefficientRing F
-		then return null;
-
-		graph := generators graphIdeal f;
-		assert( not isHomogeneous f or isHomogeneous graph );
-		SS := ring graph;
-		chh := canUseHilbertHint graph;
-		if chh then (
-			-- compare with pushNonLinear
-			hf := poincare module target f;
-			T := degreesRing SS;
-			hf = hf * product(degrees source graph, d -> 1 - T_d);
-			-- cache poincare
-			poincare cokernel graph = hf;
-		);
-		n1 := numgens F;
-		mapback := map(R, ring graph, map(R^1, R^n1, 0) | vars R);
-		G := gb(graph,opts);
-		assert (not chh or G#?"rawGBSetHilbertFunction log"); -- ensure the Hilbert function hint was actually used in gb.m2
-		ideal mapback selectInSubring(1,generators G)
-	),
+        (F, R) := (target f, source f);
+        if not isAffineRing R
+            or not isAffineRing F
+            or not instance(ambient R, PolynomialRing)
+            or not instance(ambient F, PolynomialRing)
+            or not isField coefficientRing R
+            or not coefficientRing R === coefficientRing F
+        then return null;
+
+        graph := generators graphIdeal f;
+        assert( not isHomogeneous f or isHomogeneous graph );
+        SS := ring graph;
+        chh := canUseHilbertHint graph;
+        if chh then (
+            -- compare with pushNonLinear
+            hf := poincare module target f;
+            T := degreesRing SS;
+            hf = hf * product(degrees source graph, d -> 1 - T_d);
+            -- cache poincare
+            poincare cokernel graph = hf;
+        );
+        n1 := numgens F;
+        mapback := map(R, ring graph, map(R^1, R^n1, 0) | vars R);
+        G := gb(graph,opts);
+        assert (not chh or G#?"rawGBSetHilbertFunction log"); -- ensure the Hilbert function hint was actually used in gb.m2
+        ideal mapback selectInSubring(1,generators G)
+    ),
 
     ZZ => (opts, f) -> if source f === ZZ then ideal char target f,
 
     Default => (opts, f) -> (
-		(F, R) := (target f, source f);
-		numsame := 0;
-		while (
-		    R.baseRings#?numsame and
-		    F.baseRings#?numsame and
-		    R.baseRings#numsame === F.baseRings#numsame
-		) do numsame = numsame + 1;
-		while not (
-		    isField F.baseRings#(numsame-1)
-		    or
-		    F.baseRings#(numsame-1).?isBasic
-		) do numsame = numsame - 1;
-		k := F.baseRings#(numsame-1);
-		(R',p) := flattenRing(R, CoefficientRing => k);
-		(F',r) := flattenRing(F, CoefficientRing => k);
-		if R' === R and F' === F then return null;
-		p^-1 kernel (r * f * p^-1)
-	),
+        (F, R) := (target f, source f);
+        numsame := 0;
+        while (
+            R.baseRings#?numsame and
+            F.baseRings#?numsame and
+            R.baseRings#numsame === F.baseRings#numsame
+        ) do numsame = numsame + 1;
+        while not (
+            isField F.baseRings#(numsame-1)
+            or
+            F.baseRings#(numsame-1).?isBasic
+        ) do numsame = numsame - 1;
+        k := F.baseRings#(numsame-1);
+        (R',p) := flattenRing(R, CoefficientRing => k);
+        (F',r) := flattenRing(F, CoefficientRing => k);
+        if R' === R and F' === F then return null;
+        p^-1 kernel (r * f * p^-1)
+    )
 }
 
 -- Installing hooks for kernel RingMap
-scan({Default, ZZ, "AffineRing", FractionField}, strategy ->
+scan({Default, ZZ, "AffineRing", FractionField, "Center"}, strategy ->
     addHook(key := (kernel, RingMap), algorithms#key#strategy, Strategy => strategy))
 
 -----------------------------------------------------------------------------

M2/Macaulay2/m2/shared.m2

@@ -4,7 +4,6 @@ needs "methods.m2"
 
 -- methods
 
-center = method(Options => true)
 chi = method()
 euler  = method()
 eulers = method()

M2/Macaulay2/packages/GradedLieAlgebras.m2

@@ -563,7 +563,7 @@ annihilator(FGLieSubAlgebra):=opts->(S)->(
     quotient(zeroIdeal L,S)
     );
 
-center(LieAlgebra):= {} >> o -> L->(
+center(LieAlgebra):= opts -> L->(
     c:=annihilator(fullLieSubAlgebra L);
     new LieIdeal from c
     );

M2/Macaulay2/packages/Graphs.m2

@@ -726,7 +726,7 @@ breadthFirstSearch (Digraph, Thing) := List => (G, v) -> (
 
 BFS = breadthFirstSearch
 
-center Graph := List => {} >> o -> G -> select(vertexSet G, i -> eccentricity(G, i) == radius G)
+center Graph := List => opts -> G -> select(vertexSet G, i -> eccentricity(G, i) == radius G)
 
 children = method()
 children (Digraph, Thing) := Set => (G, v) -> (

M2/Macaulay2/packages/Macaulay2Doc/functions.m2

@@ -23,6 +23,7 @@ load "./functions/between-doc.m2"
 load "./functions/binomial-doc.m2"
 load "./functions/canUseHilbertHint-doc.m2"
 load "./functions/ceiling-doc.m2"
+load "./functions/center-doc.m2"
 load "./functions/changeBase-doc.m2"
 load "./functions/char-doc.m2"
 load "./functions/characters-doc.m2"