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

Commit

Permalink
Bringing more to the GB strategy classes.
Browse files Browse the repository at this point in the history
  • Loading branch information
tscrim authored and trevorkarn committed Jul 26, 2022
1 parent f2a25e3 commit dba5089
Show file tree
Hide file tree
Showing 4 changed files with 76 additions and 33 deletions.
30 changes: 23 additions & 7 deletions src/sage/algebras/clifford_algebra.py
Original file line number Diff line number Diff line change
Expand Up @@ -2824,7 +2824,6 @@ def __init__(self, ring, gens, coerce=True, side="twosided"):
"""
Initialize ``self``.
"""
self._groebner_basis = None
self._groebner_strategy = None
self._homogeneous = all(x.is_super_homogeneous() for x in gens)
if self._homogeneous:
Expand All @@ -2835,7 +2834,24 @@ def reduce(self, f):
"""
Reduce ``f`` modulo ``self``.
"""
return f.reduce(self)
if self._groebner_strategy is None:
self.groebner_basis()
R = self.ring()
return self._groebner_strategy.reduce(R(f))

def _contains_(self, f):
r"""
Return ``True`` if ``f`` is in this ideal,
``False`` otherwise.
EXAMPLES::
.. NOTE::
Requires computation of a Groebner basis, which can be a very
expensive operation.
"""
return self.reduce(f).is_zero()

def groebner_basis(self, term_order="neglex"):
r"""
Expand Down Expand Up @@ -2903,9 +2919,9 @@ def groebner_basis(self, term_order="neglex"):
from sage.algebras.exterior_algebra_groebner import GroebnerStrategyDegLex as strategy
else:
raise ValueError("invalid term order")
if strategy == self._groebner_strategy:
return self._groebner_basis
self._groebner_strategy = strategy
self._groebner_basis = self._groebner_strategy(self).compute_groebner()
return self._groebner_basis
if strategy == type(self._groebner_strategy):
return self._groebner_strategy.groebner_basis
self._groebner_strategy = strategy(self)
self._groebner_strategy.compute_groebner()
return self._groebner_strategy.groebner_basis

11 changes: 6 additions & 5 deletions src/sage/algebras/clifford_algebra_element.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -390,18 +390,19 @@ cdef class ExteriorAlgebraElement(CliffordAlgebraElement):
INPUT:
- ``I`` -- a list of exterior algebra elements or an ideal
- ``side`` -- the side, ignored if ``I`` is an ideal
- ``left`` -- boolean; if reduce as a left ideal (``True``)
or right ideal (``False``), ignored if ``I`` is an ideal
"""
from sage.algebras.clifford_algebra import ExteriorAlgebraIdeal
if isinstance(I, ExteriorAlgebraIdeal):
left = (I.side() == "left")
I = I.groebner_basis()
return I.reduce(self)

f = self
E = self._parent
from sage.algebras.exterior_algebra_cython import leading_support

cdef FrozenBitset lm, s
for g in I:
lm = leading_support(g)
lm = g.leading_support()
reduction = True
while reduction:
supp = f.support()
Expand Down
4 changes: 4 additions & 0 deletions src/sage/algebras/exterior_algebra_groebner.pxd
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@ cdef class GroebnerStrategy:
cdef int side
cdef MonoidElement ideal
cdef bint homogeneous
cdef public tuple groebner_basis

cdef inline bint build_S_poly(self, CliffordAlgebraElement f, CliffordAlgebraElement g)

Expand All @@ -26,6 +27,9 @@ cdef class GroebnerStrategy:
cdef set preprocessing(self, list P, list G)
cdef list reduction(self, list P, list G)

cpdef CliffordAlgebraElement reduce(self, CliffordAlgebraElement f)
cdef CliffordAlgebraElement reduce_single(self, CliffordAlgebraElement f, CliffordAlgebraElement g)

# These are the methods that determine the ordering of the monomials.
# These must be implemented in subclasses. Declare them as "inline" there.
cdef Integer bitset_to_int(self, FrozenBitset X)
Expand Down
64 changes: 43 additions & 21 deletions src/sage/algebras/exterior_algebra_groebner.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -51,6 +51,7 @@ cdef class GroebnerStrategy:
Initialize ``self``.
"""
self.ideal = I
self.groebner_basis = (None,)
self.E = <Parent> I.ring()
self.homogeneous = all(x.is_super_homogeneous() for x in I.gens())
if self.homogeneous or I.side() == "left":
Expand Down Expand Up @@ -221,39 +222,60 @@ cdef class GroebnerStrategy:

# Now that we have a Gröbner basis, we make this into a reduced Gröbner basis
cdef set pairs = set((i, j) for i in range(n) for j in range(n) if i != j)
cdef list supp
cdef tuple supp
cdef bint did_reduction
cdef FrozenBitset lm, s
while pairs:
i,j = pairs.pop()
# We perform the classical reduction algorithm here on each pair
# TODO: Make this faster by using the previous technique
f = G[i]
g = G[j]
lm = self.leading_supp(g)
did_reduction = True
while did_reduction:
supp = f.support()
did_reduction = False
for s in supp:
if lm <= s:
did_reduction = True
mon = self.E.monomial(s - lm)
if self.side == 0:
gp = mon * g
f = f - f[s] / gp[s] * gp
else:
gp = g * mon
f = f - f[s] / gp[s] * gp
break
# TODO: Make this faster by using the previous technique?
f = self.reduce_single(G[i], G[j])
if G[i] != f:
G[i] = f
if not f:
pairs.difference_update((k, i) for k in range(n))
else:
pairs.update((k, i) for k in range(n) if k != i)

return tuple([~f[self.leading_supp(f)] * f for f in G if f])
self.groebner_basis = tuple([~f[self.leading_supp(f)] * f for f in G if f])

cpdef CliffordAlgebraElement reduce(self, CliffordAlgebraElement f):
"""
Reduce ``f`` modulo the ideal with Gröbner basis ``G``.
"""
for g in self.groebner_basis:
f = self.reduce_single(f, <CliffordAlgebraElement> g)
return f

cdef CliffordAlgebraElement reduce_single(self, CliffordAlgebraElement f, CliffordAlgebraElement g):
"""
Reduce ``f`` by ``g``.
.. TODO::
Optimize this by doing it in-place and changing the underlying dict of ``f``.
"""
cdef FrozenBitset lm, s
cdef tuple supp

lm = self.leading_supp(g)
did_reduction = True
while did_reduction:
supp = tuple(f._monomial_coefficients)
did_reduction = False
for s in supp:
if lm <= s:
did_reduction = True
mon = self.E.monomial(s - lm)
if self.side == 0:
gp = mon * g
f -= f[s] / gp[s] * gp
else:
gp = g * mon
f -= f[s] / gp[s] * gp
break
return f


cdef Integer bitset_to_int(self, FrozenBitset X):
raise NotImplementedError
Expand Down

0 comments on commit dba5089

Please sign in to comment.