sagemath / sage

Main repository of SageMath
https://www.sagemath.org
Other
1.31k stars 449 forks source link

Groebner bases for exterior algebras native to Sage #34138

Closed trevorkarn closed 2 years ago

trevorkarn commented 2 years ago

Add a native Sage implementation to compute a Groebner basis for an ideal of a Sage ExteriorAlgebra

Depends on #32369

CC: @tscrim @simon-king-jena

Component: algebra

Keywords: groebner basis, gsoc2022, f4

Author: Trevor K. Karn, Travis Scrimshaw

Branch/Commit: 36b6b20

Reviewer: Travis Scrimshaw, Trevor K. Karn

Issue created by migration from https://trac.sagemath.org/ticket/34138

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 2 years ago

Changed commit from e7e84b7 to 36b6b20

tscrim commented 2 years ago
comment:28

Here is an update for the GBs that does the reduction in-place, which required me to remove duplicate entries during the GB. For "small" examples, it is slower:

sage: E.<x0,x1,x2,x3,x4,x5> = ExteriorAlgebra(QQ)
sage: I = E.ideal(5*x0^2*x3+37*x1*x3*x4+32*x1*x3*x5+21*x3*x5+55*x4*x5,
....:                     39*x0*x1*x5+23*x1^2*x4+57*x1*x2*x4+56*x1*x4^2+10*x2^2+52*x3*x4*x5,
....:                     33*x0^2*x3+51*x0^2+42*x0*x3*x5+51*x1^2*x4+32*x1*x3^2+x5^3,
....:                     44*x0*x3^2+42*x1*x3+47*x1*x4^2+12*x2*x3+2*x2*x4*x5+43*x3*x4^2,
....:                     49*x0^2*x2+11*x0*x1*x2+39*x0*x3*x4+44*x0*x3*x4+54*x0*x3+45*x1^2*x4,
....:                     48*x0*x2*x3+2*x2^2*x3+59*x2^2*x5+17*x2+36*x3^3+45*x4)
sage: %time GB = I.groebner_basis()
CPU times: user 116 ms, sys: 0 ns, total: 116 ms
Wall time: 115 ms

versus before

CPU times: user 76.7 ms, sys: 299 µs, total: 77 ms
Wall time: 76.7 ms

However, I think that is just due to the extra amount of information we have to store and manipulate. I strongly believe that will get dominated (and become insignificant) once the examples get larger. The data structure for that is not very efficient and could be improved. Yet, I still haven't gotten it to be anywhere close to what M2 can do; in particular this example that took you ~7s in M2:

sage: E.<p0,p1,p2,p3,p4,p5,p6,p7> = ExteriorAlgebra(QQ)
sage: I = E.ideal(3821*p2*p3*p4*p5*p6*p7-7730*p2*p3*p4*p5*p6-164*p2*p3*p4*p5*p7-2536*p2*p3*p4*p6*p7-4321*p2*p3*p5*p6*p7-2161*p2*p4*p5*p6*p7-2188*p3*p4*p5*p6*p7-486*p2*p3*p4*p5+3491*p2*p3*p4*p6+4247*p2*p3*p5*p6+3528*p2*p4*p5*p6+2616*p3*p4*p5*p6-101*p2*p3*p4*p7+1765*p2*p3*p5*p7+258*p2*p4*p5*p7-378*p3*p4*p5*p7+1246*p2*p3*p6*p7+2320*p2*p4*p6*p7+1776*p3*p4*p6*p7+1715*p2*p5*p6*p7+728*p3*p5*p6*p7+842*p4*p5*p6*p7+69*p2*p3*p4-1660*p2*p3*p5+1863*p2*p4*p5+1520*p3*p4*p5-245*p2*p3*p6-804*p2*p4*p6-2552*p3*p4*p6-3152*p2*p5*p6+40*p3*p5*p6-1213*p4*p5*p6+270*p2*p3*p7-851*p2*p4*p7+327*p3*p4*p7-1151*p2*p5*p7+1035*p3*p5*p7-161*p4*p5*p7-230*p2*p6*p7-294*p3*p6*p7-973*p4*p6*p7-264*p5*p6*p7+874*p2*p3-2212*p2*p4+168*p3*p4+511*p2*p5-918*p3*p5-2017*p4*p5-76*p2*p6+465*p3*p6+1629*p4*p6+856*p5*p6-54*p2*p7-1355*p3*p7+227*p4*p7+77*p5*p7-220*p6*p7-696*p2+458*p3+486*p4+661*p5-650*p6+671*p7-439,
....:       -6157*p1*p3*p4*p5*p6*p7+13318*p1*p3*p4*p5*p6+5928*p1*p3*p4*p5*p7+1904*p1*p3*p4*p6*p7+2109*p1*p3*p5*p6*p7+8475*p1*p4*p5*p6*p7+2878*p3*p4*p5*p6*p7-8339*p1*p3*p4*p5-2800*p1*p3*p4*p6-9649*p1*p3*p5*p6-10964*p1*p4*p5*p6-4481*p3*p4*p5*p6+251*p1*p3*p4*p7-4245*p1*p3*p5*p7-7707*p1*p4*p5*p7-2448*p3*p4*p5*p7+1057*p1*p3*p6*p7-3605*p1*p4*p6*p7+546*p3*p4*p6*p7-3633*p1*p5*p6*p7-699*p3*p5*p6*p7-4126*p4*p5*p6*p7-730*p1*p3*p4+5519*p1*p3*p5+8168*p1*p4*p5+4366*p3*p4*p5+2847*p1*p3*p6+2058*p1*p4*p6-1416*p3*p4*p6+8004*p1*p5*p6+4740*p3*p5*p6+5361*p4*p5*p6-677*p1*p3*p7+1755*p1*p4*p7-760*p3*p4*p7+3384*p1*p5*p7+2038*p3*p5*p7+4119*p4*p5*p7+812*p1*p6*p7+11*p3*p6*p7+2022*p4*p6*p7+2642*p5*p6*p7+1276*p1*p3-1723*p1*p4+121*p3*p4-6456*p1*p5-3710*p3*p5-4525*p4*p5-2187*p1*p6-1559*p3*p6-848*p4*p6-4041*p5*p6-83*p1*p7-12*p3*p7-1180*p4*p7-2747*p5*p7-1970*p6*p7+2575*p1-161*p3+2149*p4+4294*p5+1687*p6+958*p7-1950,
....:       182*p1*p2*p4*p5*p6*p7-2824*p1*p2*p4*p5*p6-3513*p1*p2*p4*p5*p7-3386*p1*p2*p4*p6*p7-2330*p1*p2*p5*p6*p7-2838*p1*p4*p5*p6*p7+1294*p2*p4*p5*p6*p7+4764*p1*p2*p4*p5+1647*p1*p2*p4*p6+4221*p1*p2*p5*p6+814*p1*p4*p5*p6+2738*p2*p4*p5*p6+4057*p1*p2*p4*p7+2403*p1*p2*p5*p7+2552*p1*p4*p5*p7+471*p2*p4*p5*p7+448*p1*p2*p6*p7+2336*p1*p4*p6*p7+1617*p2*p4*p6*p7+2220*p1*p5*p6*p7-1543*p2*p5*p6*p7+402*p4*p5*p6*p7-5184*p1*p2*p4-3983*p1*p2*p5+44*p1*p4*p5-1327*p2*p4*p5-581*p1*p2*p6-389*p1*p4*p6-2722*p2*p4*p6+443*p1*p5*p6-2893*p2*p5*p6-154*p4*p5*p6-1277*p1*p2*p7-2018*p1*p4*p7-509*p2*p4*p7-1254*p1*p5*p7+602*p2*p5*p7-464*p4*p5*p7-647*p1*p6*p7+922*p2*p6*p7-1463*p4*p6*p7+729*p5*p6*p7+2665*p1*p2+591*p1*p4+981*p2*p4-444*p1*p5+1818*p2*p5-1985*p4*p5-1818*p1*p6+197*p2*p6+1038*p4*p6+340*p5*p6+399*p1*p7-835*p2*p7+787*p4*p7-753*p5*p7-221*p6*p7+481*p1+260*p2+1713*p4+1219*p5+794*p6+762*p7-1231,
....:       2923*p1*p2*p3*p5*p6*p7-4328*p1*p2*p3*p5*p6-3674*p1*p2*p3*p5*p7-3291*p1*p2*p3*p6*p7-4955*p1*p2*p5*p6*p7-8*p1*p3*p5*p6*p7-135*p2*p3*p5*p6*p7+7784*p1*p2*p3*p5+3471*p1*p2*p3*p6+1544*p1*p2*p5*p6+1607*p1*p3*p5*p6+1710*p2*p3*p5*p6+2434*p1*p2*p3*p7+1408*p1*p2*p5*p7-215*p1*p3*p5*p7+507*p2*p3*p5*p7+2208*p1*p2*p6*p7+1920*p1*p3*p6*p7-389*p2*p3*p6*p7+1304*p1*p5*p6*p7+2480*p2*p5*p6*p7+102*p3*p5*p6*p7-2683*p1*p2*p3-3508*p1*p2*p5-3505*p1*p3*p5-2400*p2*p3*p5-2236*p1*p2*p6-1727*p1*p3*p6-1354*p2*p3*p6-1022*p1*p5*p6-2099*p2*p5*p6-918*p3*p5*p6-495*p1*p2*p7-109*p1*p3*p7+474*p2*p3*p7+268*p1*p5*p7+1084*p2*p5*p7-190*p3*p5*p7-666*p1*p6*p7-497*p2*p6*p7+615*p3*p6*p7-912*p5*p6*p7+473*p1*p2+742*p1*p3+186*p2*p3+1021*p1*p5+2556*p2*p5+2312*p3*p5+1075*p1*p6+920*p2*p6+164*p3*p6+80*p5*p6-199*p1*p7-1270*p2*p7-1050*p3*p7-724*p5*p7+136*p6*p7+740*p1-474*p2+37*p3-1056*p5+303*p6+833*p7+736,
....:       4990*p1*p2*p3*p4*p6*p7-3067*p1*p2*p3*p4*p6-1661*p1*p2*p3*p4*p7-4064*p1*p2*p3*p6*p7-223*p1*p2*p4*p6*p7-5229*p1*p3*p4*p6*p7-4636*p2*p3*p4*p6*p7+5720*p1*p2*p3*p4+4872*p1*p2*p3*p6+1643*p1*p2*p4*p6+4536*p1*p3*p4*p6+2451*p2*p3*p4*p6+1264*p1*p2*p3*p7+70*p1*p2*p4*p7+2213*p1*p3*p4*p7+1734*p2*p3*p4*p7+1698*p1*p2*p6*p7+3799*p1*p3*p6*p7+1622*p2*p3*p6*p7+901*p1*p4*p6*p7-496*p2*p4*p6*p7+3782*p3*p4*p6*p7-5591*p1*p2*p3-1303*p1*p2*p4-6383*p1*p3*p4-2332*p2*p3*p4-3179*p1*p2*p6-6257*p1*p3*p6-3654*p2*p3*p6-1830*p1*p4*p6-1473*p2*p4*p6-3278*p3*p4*p6-1462*p1*p2*p7-1495*p1*p3*p7-468*p2*p3*p7-400*p1*p4*p7+431*p2*p4*p7-1907*p3*p4*p7-1547*p1*p6*p7-214*p2*p6*p7-1423*p3*p6*p7-1625*p4*p6*p7+5708*p1*p2+3809*p1*p3+2053*p2*p3+2824*p1*p4+1122*p2*p4+3653*p3*p4+3658*p1*p6+3001*p2*p6+3890*p3*p6+2371*p4*p6+602*p1*p7+185*p2*p7+899*p3*p7+963*p4*p7+560*p6*p7-4557*p1-3536*p2-1635*p3-2552*p4-2595*p6-207*p7+2740,
....:       -1407*p1*p2*p3*p4*p5*p7+4444*p1*p2*p3*p4*p5+2350*p1*p2*p3*p4*p7+5424*p1*p2*p3*p5*p7-2524*p1*p2*p4*p5*p7-985*p1*p3*p4*p5*p7-431*p2*p3*p4*p5*p7-2662*p1*p2*p3*p4-5342*p1*p2*p3*p5-39*p1*p2*p4*p5-2525*p1*p3*p4*p5-2650*p2*p3*p4*p5-3553*p1*p2*p3*p7-71*p1*p2*p4*p7-3268*p1*p3*p4*p7-1140*p2*p3*p4*p7-702*p1*p2*p5*p7-924*p1*p3*p5*p7-2198*p2*p3*p5*p7+4087*p1*p4*p5*p7+2709*p2*p4*p5*p7+587*p3*p4*p5*p7+968*p1*p2*p3-150*p1*p2*p4+909*p1*p3*p4+4587*p2*p3*p4+929*p1*p2*p5+1804*p1*p3*p5+2226*p2*p3*p5-916*p1*p4*p5+906*p2*p4*p5+2735*p3*p4*p5+1894*p1*p2*p7+2998*p1*p3*p7+1611*p2*p3*p7+304*p1*p4*p7-1601*p2*p4*p7+2066*p3*p4*p7-1971*p1*p5*p7-480*p2*p5*p7-500*p3*p5*p7-2617*p4*p5*p7-532*p1*p2+2016*p1*p3-2574*p2*p3+529*p1*p4-1251*p2*p4-2082*p3*p4+280*p1*p5-852*p2*p5-476*p3*p5-340*p4*p5-924*p1*p7+253*p2*p7-1090*p3*p7+170*p4*p7+1204*p5*p7-869*p1+1394*p2-264*p3+719*p4+219*p5-128*p7+506,
....:       -901*p1*p2*p3*p4*p5*p6+1805*p1*p2*p3*p4*p5-1103*p1*p2*p3*p4*p6-1746*p1*p2*p3*p5*p6-1968*p1*p2*p4*p5*p6+3957*p1*p3*p4*p5*p6+1293*p2*p3*p4*p5*p6-523*p1*p2*p3*p4-2498*p1*p2*p3*p5+693*p1*p2*p4*p5-2805*p1*p3*p4*p5-722*p2*p3*p4*p5-770*p1*p2*p3*p6+1088*p1*p2*p4*p6-232*p1*p3*p4*p6+2657*p2*p3*p4*p6+3281*p1*p2*p5*p6-1066*p1*p3*p5*p6+240*p2*p3*p5*p6-1174*p1*p4*p5*p6+1304*p2*p4*p5*p6-2070*p3*p4*p5*p6+2571*p1*p2*p3+115*p1*p2*p4+3899*p1*p3*p4-4641*p2*p3*p4-752*p1*p2*p5+1531*p1*p3*p5+1178*p2*p3*p5+11*p1*p4*p5-1144*p2*p4*p5-1701*p3*p4*p5+592*p1*p2*p6+1140*p1*p3*p6+130*p2*p3*p6+304*p1*p4*p6-2273*p2*p4*p6-1224*p3*p4*p6-2*p1*p5*p6-1090*p2*p5*p6+585*p3*p5*p6+670*p4*p5*p6-1867*p1*p2-4780*p1*p3+1079*p2*p3-2435*p1*p4+2901*p2*p4+2073*p3*p4+499*p1*p5+908*p2*p5+323*p3*p5+1631*p4*p5-966*p1*p6-315*p2*p6-481*p3*p6+759*p4*p6-595*p5*p6+3233*p1-1978*p2+729*p3-1184*p4-40*p5+446*p6+282, side='left')
sage: %time GB = I.groebner_basis(reduced=False)
tscrim commented 2 years ago
comment:29

My guess is M2 is doing something very special and particular with the linear algebra (and maybe some stuff in parallel). Anyways, this is as best as I can get for now.

tscrim commented 2 years ago
comment:30

Yay, green patchbot. :).

tscrim commented 2 years ago
comment:32

Thank you.

Next steps for followup tickets:

cc-ing Simon to let him know about this as he might be interested in it.

simon-king-jena commented 2 years ago
comment:33

Replying to @tscrim:

Thank you.

Next steps for followup tickets:

  • Generalizing this implementation to work with any PBW-like basis.
  • The improvements noted at the end of comment:26.
  • Finding other speed improvements.
  • Implement conversion to/from M2 and/or Singular (via plural) to have that as an alternative algorithm.

cc-ing Simon to let him know about this as he might be interested in it.

Indeed I'm interested. General question (before heaving looked at the code): What was the rationale for a "native" implementation (aka "reinventing the wheel" instead of using what (e.g.) Singular has to offer?

tscrim commented 2 years ago
comment:34

Replying to @simon-king-jena:

Replying to @tscrim:

cc-ing Simon to let him know about this as he might be interested in it.

Indeed I'm interested. General question (before heaving looked at the code): What was the rationale for a "native" implementation (aka "reinventing the wheel" instead of using what (e.g.) Singular has to offer?

I couldn't find a dedicated implementation of GB and the exterior algebra through some (not super extensive) Google searching. With Singular using the generic g-algebra and Plural code, I didn't think it would be particularly fast (although I have no timings to back this up) and it seems quite cumbersome to go back-and-forth with. Not everyone has Macaulay2 installed. So I figured this would be the easiest way to get this functionality and, hopefully (but unfortunately does not (yet) seem to be done), the fastest way.

simon-king-jena commented 2 years ago
comment:35

I may be wrong, but I thought that the main feature of "exterior algebra" is being graded-commutative. And that's provided by what Singular calls "super-commutative". And that belongs to the parts of Plural that are wrapped in Sage already.

simon-king-jena commented 2 years ago
comment:36

Looking at SCA (from sage.rings.polynomial.plural import SCA), it seems that ideals over super-commutative algebras are mainly dealt with generically. On the other hand, I see uses of libsingular's groeber_basis functions in the SCA-code. Strange. Anyway, it is available via libsingular, and I do use it in my cohomology code. Perhaps all what's missing is an interface, so that it can be used without directly invoking libsingular.

tscrim commented 2 years ago
comment:37

To build quotients using GradedCommutativeAlgebra (which uses Plural), the ideals are currently mandated to be homogeneous. Although if I pass check=False, it seems like it is doing the GB computation:

sage: E.<p0,p1,p2,p3,p4,p5,p6,p7> = GradedCommutativeAlgebra(QQ)
sage: I = E.ideal(5*x0^2*x3+37*x1*x3*x4+32*x1*x3*x5+21*x3*x5+55*x4*x5,
....:                     39*x0*x1*x5+23*x1^2*x4+57*x1*x2*x4+56*x1*x4^2+10*x2^2+52*x3*x4*x5,
....:                     33*x0^2*x3+51*x0^2+42*x0*x3*x5+51*x1^2*x4+32*x1*x3^2+x5^3,
....:                     44*x0*x3^2+42*x1*x3+47*x1*x4^2+12*x2*x3+2*x2*x4*x5+43*x3*x4^2,
....:                     49*x0^2*x2+11*x0*x1*x2+39*x0*x3*x4+44*x0*x3*x4+54*x0*x3+45*x1^2*x4,
....:                     48*x0*x2*x3+2*x2^2*x3+59*x2^2*x5+17*x2+36*x3^3+45*x4)
sage: Q = E.quotient(I, check=False)
sage: Q.gens()
(x0, x1, -45/17*x4, x3, x4, x5)

Indeed, the interface for GBs does seem to be missing from these ideals. So I cannot run timings.

The specialized multiplication for the exterior algebra is apparently faster:

sage: E.<x0,x1,x2,x3,x4,x5> = ExteriorAlgebra(QQ)
sage: r = sum(E.basis())
sage: %timeit r * r
527 µs ± 6.81 µs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)

sage: E.<x0,x1,x2,x3,x4,x5> = GradedCommutativeAlgebra(QQ)
sage: r = sum(sum(E.basis(n)) for n in range(7))
sage: %timeit r * r
712 µs ± 5.39 µs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)

So it should be possible to make the version of GBs here faster. It should just require being more clever about how we do the linear algebra.

Another reason why it is good to have a Sage version is this can handle any field that Sage has, which I believe is a lot more than what Singular can do (and take advantage of the corresponding specialized linear algebra code). With some hopefully minor future tweaks, this could work over more general (skew?) commutative rings.

tscrim commented 2 years ago
comment:38

Also, do you know if Plural handles left/right and two-sided ideals?

vbraun commented 2 years ago

Changed branch from public/algebras/exterior_groebner-34138 to 36b6b20