Groebner basis not working over symbolic ring

[rhinton]

I'm not sure if this is a problem in the multivariate polynomials (which seem to raise the actual error) or somewhere in the symbolics.

sage: R2.<a,b> = SR[]
sage: I2 = [a*b+a, a*a] * R2
sage: G2 = I2.groebner_basis()
verbose 0 (2247: multi_polynomial_ideal.py, groebner_basis) Warning: falling back to very slow toy implementation.
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
...
AttributeError: 'MPolynomialRing_polydict' object has no attribute 'monomial_divides'
sage: 

Apply:

trac_6581_enable_more_ideals_for_toy_buchberger.patch

[john_perry]

The failure is because the polynomials' parent is of type MPolynomialRing_polydict, but MPolynomialRing_polydict_domain is the type with the necessary methods monomial_divides, monomial_pairwise_prime, etc.

• Should a polynomial ring in this case always initialize as ..._domain?
• If not, is there an easy option to create the ideal of such a ring? (in this case, at most a documentation change is necessary)
• If no such method currently exists, what would be the best way to let the user do this?
  1. Have the user create R2 as above, then call a (new) function that makes it re-initialize as ..._domain? or
  2. Have a new creation mechanism for such rings, e.g., a symbolic ring SD that is a domain, and thus can be initialized as ..._domain?

[malb]

I don't think GB computations over SR make sense. How do we compare with zero, e.g., this: cos(x)^2 *X + sin(x)^2 * X - X. How do we decide divisibility, e.g. is "cos(x)" zero?

[john_perry]

I've had situations where you can still compute a Groebner basis with some indeterminate coefficients. It's annoying to have to do all the work by hand, though. I can see implementing this with a warning that it may not make sense.

[malb]

If you only want symbolic variables, wouldn't this work:

sage: P.<a,b,c> = QQ[]
sage: K = Frac(P)
sage: P.<x,y> = PolynomialRing(K)
sage: P
Multivariate Polynomial Ring in x, y over Fraction Field of Multivariate Polynomial Ring in a, b, c over Rational Field

This gets mapped to a Singular ring with parameters, i.e., to this one:

> ring r = (0,a,b,c),(x,y,z),dp;
> r;
//   characteristic : 0
//   3 parameter    : a b c 
//   minpoly        : 0
//   number of vars : 3
//        block   1 : ordering dp
//                  : names    x y z
//        block   2 : ordering C

So, this works:

sage: singular(P)
//   characteristic : 0
//   3 parameter    : a b c 
//   minpoly        : 0
//   number of vars : 2
//        block   1 : ordering dp
//                  : names    x y
//        block   2 : ordering C
sage: I = Ideal([P.random_element(), P.random_element()])
sage: %time gb = I.groebner_basis()

but it's slow.

[malb]

Here's a smaller example:

sage: P.<a> = QQ[]               
sage: K = Frac(P)
sage: P.<x,y,z> = PolynomialRing(K)
sage: I = Ideal([P.random_element(), P.random_element(), P.random_element()])
sage: I
Ideal (((a^2 - 3/4*a - 1)/(a^2 + a))*y*z + ((2/3*a^2 - 5/3*a - 1)/(-9*a^2 + 5/6*a - 9/2))*z^2 + ((3/2*a^2 - a - 5/2)/(-1/20*a - 5/57))*x + ((2/3*a^2 - a - 1)/(-1/2*a^2 + 5*a + 7/4))*y + (7/3*a^2 - 5*a - 1)/(1/3*a^2 - 1/3*a + 3), ((1/2*a^2 - a - 1/2)/(-a^2 + 7))*x*y + ((-5/6*a^2 + 3*a - 1)/(-a^2 + 2/3*a - 1/4))*x*z + ((a^2 + 3/5*a + 1)/(-6*a^2 - a + 95))*y*z + ((-3*a^2 + 21/2)/(-a^2 - 1/5*a - 1))*z^2 + ((-1/3*a - 1/2)/(-1/3*a^2 - 1/8*a - 3))*z, ((a^2 + a)/(-3/2*a^2 - 1/3*a - 1))*x*y + 1/2*a*y^2 + ((1/31*a^2 - 1/8)/(-4*a^2 + 1/5*a - 2/3))*x*z + ((a^2 - a)/(a^2 - 1/20*a - 1/2))*z^2 + (-1/15*a^2 + 1/315*a + 7/15)*y) of Multivariate Polynomial Ring in x, y, z over Fraction Field of Univariate Polynomial Ring in a over Rational Field
sage: %time gb = I.groebner_basis()
CPU times: user 0.08 s, sys: 0.02 s, total: 0.09 s
Wall time: 0.33 s
sage: gb[-1]
y*z + ((-5472*a^4 + 8208*a^3 + 21888*a^2 + 8208*a)/(73872*a^4 - 62244*a^3 - 31806*a^2 - 20862*a - 36936))*z^2 + ((-2216160*a^4 - 738720*a^3 + 5171040*a^2 + 3693600*a)/(73872*a^3 + 74196*a^2 - 171072*a - 129600))*x + ((-98496*a^4 + 49248*a^3 + 295488*a^2 + 147744*a)/(73872*a^4 - 794124*a^3 + 221616*a^2 + 932634*a + 258552))*y + (517104*a^4 - 590976*a^3 - 1329696*a^2 - 221616*a)/(73872*a^4 - 129276*a^3 + 646380*a^2 - 424764*a - 664848)

[john_perry]

Would this work when we place some assumptions on the variables; e.g., a^2+b^2==1? That was my particular situations, and I was under the impression that doing it in the way you suggest won't work for that.

[john_perry]

(Sorry, I meant assumptions on the coefficients.)

[malb]

This might be naive, but can't you add a^2+b^2 - 1 to your ideal?

[john_perry]

Replying to malb:

This might be naive, but can't you add a^2+b^2 - 1 to your ideal?

Doesn't that change the ideal? These are elements of the base ring.

[john_perry]

Hey, I got it to work.

sage: R.<a,b> = QQ[]
sage: I = R.ideal(a^2+b^2-1)
sage: R2 = QuotientRing(R,I)
sage: K = Frac(R2)
sage: P.<x,y,z> = K[]
sage: f = (a^2+b^2)*x + 1
sage: g = x + 1
sage: f - g
0

Sorry -- I'm a bit slow sometimes. :-)

Okay, so we should mark this as "invalid/wontfix" or something like that? Or do we wait to see if the original reporter has a scenario in mind that really does require SR?

[malb]

Unfortunately, your construction bombs out when we try to construct an actual ideal:

sage: R.<a,b> = QQ[]
sage: I = R.ideal(a^2+b^2-1)
sage: R2 = QuotientRing(R,I)
sage: K = Frac(R2)
sage: P.<x,y,z> = K[]
sage: f = (a^2+b^2)*x + 1
sage: g = x + 1
sage: f - g
0
sage: Ideal([f,g])
NotImplementedError

[john_perry]

It works (initially) if you do this:

sage: J = P.ideal([f,g])

...but then you can't do anything with J. The characteristic isn't set in R2. Curiously, I can still work with ideals in R2.

Likewise,

sage: Z7 = QuotientRing(ZZ,7*ZZ)
sage: Z7.characteristic()
...
NotImplementedError:

Wow. We should probably look at that one day. ;-)

Anyway, do you know offhand why one needs the characteristic, but not the other?

[malb]

Nope, unfortunately, I don't.

[john_perry]

Anyway, do you know offhand why one needs the characteristic, but not the other?

Bingo. In lines 302--307 of multi_polynomial_sequence.py we encounter:

    if k.characteristic() != 2:
        return PolynomialSequence_generic(parts, ring, immutable=immutable, cr=cr, cr_str=cr_str)
    elif k.degree() == 1:
        return PolynomialSequence_gf2(parts, ring, immutable=immutable, cr=cr, cr_str=cr_str)
    elif k.degree() > 1:
        return PolynomialSequence_gf2e(parts, ring, immutable=immutable, cr=cr, cr_str=cr_str)

I can make this work via judicious use of a try/catch. I found some other instances where it tried to compute a Singular representation of itself (the reduce function in multi_polynomial_element, for instance). That has allowed me to compute several Gröbner bases successfully, including one similar to the ideal you couldn't get to work:

sage: J = R2.ideal([(a^2+b^2)*x + y, x+y])
sage: J
Ideal (x + y, x + y) of Multivariate Polynomial Ring in x, y over Fraction Field of
Quotient of Multivariate Polynomial Ring in a, b over Rational Field by the ideal
(a^2 + b^2 - 1)
sage: J.groebner_basis()
verbose 0 (2854: multi_polynomial_ideal.py, groebner_basis) Warning: falling back to
very slow toy implementation.
[x + y]

I think the patch is worth adding. However, I'm still not sure this is what the original requester wanted. Should I open a new ticket?

[john_perry]

I decided to upload the patch here. I didn't think a doctest was appropriate, but I can add one if desired.

[john_perry]

Whoops -- the patch introduces doctest failures. I'll look into them.

[john_perry]

Okay, that was pretty dumb. Fixed now. Added a doctest while I was at it.

[malb]

Looks good.

[jdemeyer]

You should replace the line number in

verbose 0 (2854: multi_polynomial_ideal.py, groebner_basis) Warning: falling back to very slow toy implementation

by "..." because the line number will change all the time.

[john_perry]

Sorry, I should have thought of that myself.

I remembered this time to check the box to replace the file of the same name, so I've changed the directions (& link) on which patch to apply.

[jdemeyer]

Please write a proper commit message for your patch, use hg qrefresh -e for this.

[john_perry]

I should have thought of that, too. Sorry again...