Discriminant of multivariate polynomial over RR or RDF gives an error

[rburing]

Inspired by Ask SageMath question #50317:

sage: R.<x,z> = RR[]
sage: f = x**3 + x*z + 1
sage: f.discriminant(x)

results in an error. Traceback:

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
<ipython-input-1-e8d3020693f4> in <module>()
      1 R = RR['x, z']; (x, z,) = R._first_ngens(2)
      2 f = x**Integer(3) + x*z + Integer(1)
----> 3 f.discriminant(x)

/home/sc_serv/sage/local/lib/python3.7/site-packages/sage/rings/polynomial/multi_polynomial.pyx in sage.rings.polynomial.multi_polynomial.MPolynomial.discriminant (build/cythonized/sage/rings/polynomial/multi_polynomial.c:16730)()
   1464             u = 1
   1465         an = self.coefficient(variable**n)**(n - k - 2)
-> 1466         return self.parent()(u * self.resultant(d, variable) * an)
   1467 
   1468     def macaulay_resultant(self, *args):

/home/sc_serv/sage/local/lib/python3.7/site-packages/sage/rings/qqbar_decorators.py in wrapper(*args, **kwds)
     68         if not any([isinstance(a, (Polynomial, MPolynomial, Ideal_generic))
     69                     and isinstance(a.base_ring(), AlgebraicField_common) for a in args]):
---> 70             return func(*args, **kwds)
     71 
     72         polynomials = []

/home/sc_serv/sage/local/lib/python3.7/site-packages/sage/rings/polynomial/multi_polynomial_element.py in resultant(self, other, variable)
   2020             r = rt.sage_poly(R)
   2021         else:
-> 2022             r = self.sylvester_matrix(other, variable).det()
   2023         if R.ngens() <= 1 and r.degree() <= 0:
   2024             return R.base_ring()(r[0])

/home/sc_serv/sage/local/lib/python3.7/site-packages/sage/matrix/matrix2.pyx in sage.matrix.matrix2.Matrix.det (build/cythonized/sage/matrix/matrix2.c:14214)()
   1519             6
   1520         """
-> 1521         return self.determinant(*args, **kwds)
   1522 
   1523     def determinant(self, algorithm=None):

/home/sc_serv/sage/local/lib/python3.7/site-packages/sage/matrix/matrix_mpolynomial_dense.pyx in sage.matrix.matrix_mpolynomial_dense.Matrix_mpolynomial_dense.determinant (build/cythonized/sage/matrix/matrix_mpolynomial_dense.cpp:7144)()
    550 
    551             elif can_convert_to_singular(self.base_ring()):
--> 552                 d = R(self._singular_().det())
    553 
    554             else:

/home/sc_serv/sage/local/lib/python3.7/site-packages/sage/rings/polynomial/multi_polynomial_ring.py in __call__(self, x, check)
    499             self._singular_().set_ring()
    500             try:
--> 501                 return x.sage_poly(self)
    502             except TypeError:
    503                 raise TypeError("unable to coerce singular object")

/home/sc_serv/sage/local/lib/python3.7/site-packages/sage/interfaces/singular.py in sage_poly(self, R, kcache)
   1775                         else:
   1776                             power=1
-> 1777                         exp[var_dict[var]]=power
   1778 
   1779                 if kcache is None:

KeyError: '(1.000e+00)'

By contrast, the discriminant of the same polynomial over RIF or RBF or CC or CDF works fine.

Component: number theory

Keywords: discriminant, multivariate, RR, RDF

Author: Dave Morris

Branch/Commit: fcd0354

Reviewer: Michael Orlitzky

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

[DaveWitteMorris]

Branch: public/29396

[DaveWitteMorris]

comment:2

This one-line patch to SingularElement.sage_poly() seems to fix the problem, but I do not fully understand the code for this function yet, so the patch may need further work.

sage: R.<x,z> = RR[]
sage: f = x**3 + x*z + 1
sage: f.discriminant(x)
-4.00000000000000*z^3 - 27.0000000000000

---

New commits:

5adab40

fix bug in SingularElement.sage_poly

[DaveWitteMorris]

Author: Dave Morris

[DaveWitteMorris]

Commit: 5adab40

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

687208d

same correction in another place

[sagetrac-git]

Changed commit from 5adab40 to 687208d

[DaveWitteMorris]

comment:4

I now understand the code well enough to be confident that the patch will not have any bad side-effects. But I discovered that exactly the same bug is in another part of the code, so I patched that too (and added another doctest).

[DaveWitteMorris]

comment:5

Here is an example to show the effect of the latest patch. Before the patch, the bug can cause incorrect results:

sage: Sing = singular.ring('real', '(x,y)', 'lp') # initialize singular
sage: S.<x> = RR[]
sage: g = singular.new('x + 5') # define g to be the polynomial x + 5
sage: g.sage_poly(S) # should be the same polynomial (in sage), but it is not
5.00000000000000*x
sage: g.sage_poly(S) == x + 5
False

We get the correct result after applying the patch:

sage: Sing = singular.ring('real', '(x,y)', 'lp')
sage: S.<x> = RR[]
sage: g = singular.new('x + 5')
sage: g.sage_poly(S)
x + 5.00000000000000
sage: g.sage_poly(S) == x + 5
True

[orlitzky]

Reviewer: Michael Orlitzky

[orlitzky]

comment:6

This pattern already appears in the code once,

# Directly treat constants
if singular_poly_list[0] in ['1', '(1.000e+00)']:
    return R(singular_poly_list[1])

so it can't be too far off the mark. You could use "not in" instead of "not( .. in .. )", but that's just a matter of taste. Otherwise, it looks reasonable and fixes the problem.

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1 and set ticket back to needs_review. New commits:

fcd0354

better to use "not in"

[sagetrac-git]

Changed commit from 687208d to fcd0354

[DaveWitteMorris]

comment:8

Thanks for the review and the suggestion. (I am not fluent in python.) I believe "not in" is much better, so I made the change.

[orlitzky]

comment:9

No problem, thanks for the fix!

[vbraun]

Changed branch from public/29396 to fcd0354