Skip to content
This repository has been archived by the owner on Jan 30, 2023. It is now read-only.

Commit

Permalink
#34519 Gröbner bases using msolve
Browse files Browse the repository at this point in the history
  • Loading branch information
mezzarobba committed Sep 16, 2022
1 parent 7db5678 commit e43bc6e
Show file tree
Hide file tree
Showing 2 changed files with 76 additions and 3 deletions.
48 changes: 46 additions & 2 deletions src/sage/rings/polynomial/msolve.py
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,7 @@
from sage.rings.real_double import RealDoubleField_class
from sage.rings.real_mpfr import RealField_class
from sage.rings.real_mpfi import RealIntervalField_class, RealIntervalField
from sage.structure.sequence import Sequence

def _run_msolve(ideal, options):
r"""
Expand Down Expand Up @@ -67,8 +68,51 @@ def _run_msolve(ideal, options):

return msolve_out.stdout

def groebner_basis_degrevlex(ideal, *, proof=True):
pass
def groebner_basis_degrevlex(ideal):
r"""
Compute a degrevlex Gröbner basis using msolve
EXAMPLES::
sage: from sage.rings.polynomial.msolve import groebner_basis_degrevlex
sage: R.<a,b,c> = PolynomialRing(GF(101), 3, order='lex')
sage: I = sage.rings.ideal.Katsura(R,3)
sage: gb = groebner_basis_degrevlex(I); gb # optional - msolve
[a + 2*b + 2*c - 1, b*c - 19*c^2 + 10*b + 40*c,
b^2 - 41*c^2 + 20*b - 20*c, c^3 + 28*c^2 - 37*b + 13*c]
sage: gb.universe() is R # optional - msolve
False
sage: gb.universe().term_order() # optional - msolve
Degree reverse lexicographic term order
sage: ideal(gb).transformed_basis(other_ring=R) # optional - msolve
[c^4 + 38*c^3 - 6*c^2 - 6*c, 30*c^3 + 32*c^2 + b - 14*c,
a + 2*b + 2*c - 1]
TESTS::
sage: R.<foo, bar> = PolynomialRing(GF(536870909), 2)
sage: I = Ideal([ foo^2 - 1, bar^2 - 1 ])
sage: I.groebner_basis(algorithm='msolve') # optional - msolve
[bar^2 - 1, foo^2 - 1]
sage: R.<x, y> = PolynomialRing(QQ, 2)
sage: I = Ideal([ x*y - 1, (x-2)^2 + (y-1)^2 - 1])
sage: I.groebner_basis(algorithm='msolve') # optional - msolve
Traceback (most recent call last):
...
NotImplementedError: unsupported base field: Rational Field
"""

base = ideal.base_ring()
if not (isinstance(base, FiniteField) and base.is_prime_field() and
base.characteristic() < 2**31):
raise NotImplementedError(f"unsupported base field: {base}")

drlpolring = ideal.ring().change_ring(order='degrevlex')
msolve_out = _run_msolve(ideal, ["-g", "2"])
gbasis = sage_eval(msolve_out[:-2], locals=drlpolring.gens_dict())
return Sequence(gbasis)

def variety(ideal, ring, *, proof=True):
r"""
Expand Down
31 changes: 30 additions & 1 deletion src/sage/rings/polynomial/multi_polynomial_ideal.py
Original file line number Diff line number Diff line change
Expand Up @@ -4012,6 +4012,10 @@ def groebner_basis(self, algorithm='', deg_bound=None, mult_bound=None, prot=Fal
'macaulay2:mgb'
Macaulay2's ``GroebnerBasis`` command with the strategy "MGB" (if available)
'msolve'
`msolve <https://msolve.lip6.fr/>`_ (if available, degrevlex order,
prime fields)
'magma:GroebnerBasis'
Magma's ``Groebnerbasis`` command (if available)
Expand Down Expand Up @@ -4135,6 +4139,15 @@ def groebner_basis(self, algorithm='', deg_bound=None, mult_bound=None, prot=Fal
sage: I.groebner_basis('macaulay2:mgb') # optional - macaulay2
[c^3 + 28*c^2 - 37*b + 13*c, b^2 - 41*c^2 + 20*b - 20*c, b*c - 19*c^2 + 10*b + 40*c, a + 2*b + 2*c - 1]
Over prime fields of small characteristic, we can also use
`msolve <https://msolve.lip6.fr/>`_ (if available in the system program
search path)::
sage: R.<a,b,c> = PolynomialRing(GF(101), 3)
sage: I = sage.rings.ideal.Katsura(R,3) # regenerate to prevent caching
sage: I.groebner_basis('msolve') # optional - msolve
[a + 2*b + 2*c - 1, b*c - 19*c^2 + 10*b + 40*c, b^2 - 41*c^2 + 20*b - 20*c, c^3 + 28*c^2 - 37*b + 13*c]
::
sage: I = sage.rings.ideal.Katsura(P,3) # regenerate to prevent caching
Expand Down Expand Up @@ -4353,6 +4366,15 @@ def groebner_basis(self, algorithm='', deg_bound=None, mult_bound=None, prot=Fal
sage: I = sage.rings.ideal.Katsura(P,3) # regenerate to prevent caching
sage: I.groebner_basis('magma:GroebnerBasis') # optional - magma
[a + (-60)*c^3 + 158/7*c^2 + 8/7*c - 1, b + 30*c^3 + (-79/7)*c^2 + 3/7*c, c^4 + (-10/21)*c^3 + 1/84*c^2 + 1/84*c]
msolve currently supports the degrevlex order only::
sage: R.<a,b,c> = PolynomialRing(GF(101), 3, order='lex')
sage: I = sage.rings.ideal.Katsura(R,3) # regenerate to prevent caching
sage: I.groebner_basis('msolve') # optional - msolve
Traceback (most recent call last):
...
NotImplementedError: msolve only supports the degrevlex order (use transformed_basis())
"""
from sage.rings.polynomial.multi_polynomial_sequence import PolynomialSequence
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
Expand Down Expand Up @@ -4438,7 +4460,14 @@ def groebner_basis(self, algorithm='', deg_bound=None, mult_bound=None, prot=Fal
elif algorithm == 'giac:gbasis':
from sage.libs.giac import groebner_basis as groebner_basis_libgiac
gb = groebner_basis_libgiac(self, prot=prot, *args, **kwds)

elif algorithm == 'msolve':
if self.ring().term_order() != 'degrevlex':
raise NotImplementedError("msolve only supports the degrevlex order "
"(use transformed_basis())")
if not (deg_bound is mult_bound is None) or prot:
raise NotImplementedError("unsupported options for msolve")
from . import msolve
return msolve.groebner_basis_degrevlex(self, *args, **kwds)
else:
raise NameError("Algorithm '%s' unknown."%algorithm)

Expand Down

0 comments on commit e43bc6e

Please sign in to comment.