Bug in reduction of exterior algebra elements modulo ideal

[trevorkarn]

Steps To Reproduce

In Sage 10.2 and 10.2.rc4, reduction of elements of exterior algebra modulo an ideal does not give a reduced form.

sage: N = 4
sage: E = ExteriorAlgebra(QQ,binomial(N,2))
sage: e = E.gens()
sage: K = matroids.CompleteGraphic(N);
sage: C = K.circuits()
sage: ideal_gens = [sum([(-1)^j*E.prod(e[i] for i in list(c)[:j] + list(c)[j+1:]) for j in range(len(c))]) for c in C]
sage: OS = E.quo(ideal_gens)
sage: [OS(b).lift() for b in E.homogeneous_component_basis(4)]
[0,
 0,
 0,
 0,
 0,
 -e0*e1*e3*e4 + e0*e1*e3*e5,
 0,
 0,
 0,
 0,
 0,
 0,
 e0*e1*e3*e4 - e0*e1*e3*e5,
 e0*e1*e3*e4 - e0*e1*e3*e5,
 0]
sage: [OS(OS(b).lift()) for b in E.homogeneous_component_basis(4)]
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]

Expected Behavior

There should be no nonzero degree-4 elements in the previous example.

Actual Behavior

One needs to do the reduction step manually twice in order to reduce elements modulo the ideal.

Additional Information

A guess at the bug is that it could be caused by a loop terminating early during the Groebner basis reduction step.

Environment

- **OS**: MacOS Ventura 13.4.1
- **Sage Version**: 10.2, 10.2.rc4

Checklist

• I have searched the existing issues for a bug report that matches the one I want to file, without success.
• I have read the documentation and troubleshoot guide

[trevorkarn]

FYSA @tscrim

[tscrim]

Strictly speaking, your example doesn’t necessarily show a bug: it is only required to lift to some element (which goes to 0 under the natural quotient map). That being said, yes, it should reduce everything to 0, which uniquely lifts to 0.

Some additional useful information for here would be the GB for the ideal. It might also be an issue in the quotient with not using a GB.

[trevorkarn]

Some additional useful information for here would be the GB for the ideal.

It seems like it is using a computed GB. The computed GB is not reduced, but after manually reducing the GB, it still has this same discrepancy.

It might also be an issue in the quotient with not using a GB.

Could you say more about this? Do you mean something like the Groebner strategy is trying to do the reduction with an empty GB?

[tscrim]

Some additional useful information for here would be the GB for the ideal.

It seems like it is using a computed GB. The computed GB is not reduced, but after manually reducing the GB, it still has this same discrepancy.

I am not surprised. A reduced GB is not necessary to do the reduction and make it unique.

It might also be an issue in the quotient with not using a GB.

Could you say more about this? Do you mean something like the Groebner strategy is trying to do the reduction with an empty GB?

I was thinking it was possible the quotient wasn't constructing the GB before doing the reduction, but that is not the case.

Your comment about needed to do the reduction twice I think is on the money for the bug and the fix. We need to reduce as many times until we no longer change the element. This should be easy to fix:

    def reduce():
        if not f:
            return f
        # Make a copy to mutate
        f = type(f)(f._parent, copy(f._monomial_coefficients))
        was_reduced = True
        while was_reduced:
            was_reduced = False
            for g in self.groebner_basis:
                was_reduced = self.reduce_single(f, g)
                if was_reduced:
                    break
        return f

We could potentially speed this part up by keeping better track of the degree of the GB leading terms, but that is a more invasive change and optimization that can be done later.

We are likely making the same mistake in the reduced GB computation code.

[trevorkarn]

I am trying to go through and verify the GroebnerStrategy.reduce_single method, but for some reason I can't access the reduce_single method:

sage: N = 4
....: E = ExteriorAlgebra(QQ,binomial(N,2))
....: e = E.gens()
....: E.inject_variables()
....: K = matroids.CompleteGraphic(N);
....: C = K.circuits()
....: ideal_gens = [sum([(-1)^j*E.prod(e[i] for i in list(c)[:j] + list(c)[j+1:]) for j in range(len(c))]) for c in C]
....: OS = E.quo(ideal_gens)
Defining e0, e1, e2, e3, e4, e5
sage: I = OS.defining_ideal()
sage: I.groebner_basis()
(e0*e1 - e0*e3 + e1*e3,
 e0*e2 - e0*e4 + e2*e4,
 e1*e2 - e1*e5 + e2*e5,
 -e0*e1*e2 + e0*e1*e5 - e0*e2*e3 - e0*e3*e5 + e1*e2*e3 + e1*e3*e5,
 e3*e4 - e3*e5 + e4*e5)
sage: GS = I._groebner_strategy
sage: GS.reduce_single
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In [16], line 1
----> 1 GS.reduce_single

AttributeError: 'sage.algebras.exterior_algebra_groebner.GroebnerSt' object has no attribute 'reduce_single'

Is this some sort of Cython issue?

Here's the source code for reduce_single:

    cdef bint reduce_single(self, CliffordAlgebraElement f, CliffordAlgebraElement g) except -1:
        r"""
        Reduce ``f`` by ``g``.

        .. WARNING::

            This modifies the element ``f``.
        """
        cdef FrozenBitset lm = self.leading_support(g), s, t
        cdef bint did_reduction = True, was_reduced=False
        cdef CliffordAlgebraElement gp

        one = self.E._base.one()
        while did_reduction:
            did_reduction = False
            for s in f._monomial_coefficients:
                if lm.issubset(s):
                    t = s
                    did_reduction = True
                    was_reduced = True
                    break
            if did_reduction:
                if self.side == 0:
                    gp = g._mul_term_self(t - lm, one)
                else:
                    gp = g._mul_self_term(t - lm, one)
                coeff = f[t] / gp._monomial_coefficients[t]
                iaxpy(-coeff, gp._monomial_coefficients, f._monomial_coefficients)
        return was_reduced

Can you help me to understand the line for s in f._monomial_coefficients? If I understand what is going on, this is the part of the division of f by g where we should be looking at the leading term. Is there a reason we are iterating over all ._monomial_coefficients instead of just the leading support of f?

[trevorkarn]

I am trying to go through and verify the GroebnerStrategy.reduce_single method, but for some reason I can't access the reduce_single method

I also have the same question with g._mul_term_self on line 511.

[tscrim]

Those are cdef methods, so they are not accessible from Python. For testing purposes, you can make them a cpdef (but you have to also modify the pxd file accordingly).

When performing a reduction of an arbitrary element $f$ by a GB element $g$, we want to see if there is any term of $f$ that is divisible by the leading term of $g$. So we need to look at every monomial, which is what that line is doing. It is only $g$ that provides the leading term. For example, you want to reduce $x + y$ by the GB $(y)$ (with $x < y$), which in the quotient is equivalent to $x$.

[trevorkarn]

Those are cdef methods, so they are not accessible from Python. For testing purposes, you can make them a cpdef (but you have to also modify the pxd file accordingly).

Ah that makes sense, thanks!

When performing a reduction of an arbitrary element f by a GB element g, we want to see if there is any term of f that is divisible by the leading term of g. So we need to look at every monomial, which is what that line is doing. It is only g that provides the leading term. For example, you want to reduce x+y by the GB (y) (with x<y), which in the quotient is equivalent to x.

Ah duh - thanks for breaking it down for me.

[tscrim]

No problem. Thanks for pushing the fix.