sympy.polys.polyerrors.HeuristicGCDFailed: no luck

[oscarbenjamin]

from sympy import *

x, y = symbols('x, y')

expr = ((2*y**3*sin(x/y)**2 + x)**2*(y*(-6*y**2*sin(x/y)**2 +
    4*y*x*sin(x/y)*cos(x/y))/(2*y**3*sin(x/y)**2 + x)**2 +
    1/(2*y**3*sin(x/y)**2 + x))/(4*y*(2*y**2*(3*y*sin(x/y) -
        2*x*cos(x/y))**2*sin(x/y)**2/(2*y**3*sin(x/y)**2 + x) - 3*y*sin(x/y)**2
        + 4*x*sin(x/y)*cos(x/y) - (3*y*sin(x/y) - 2*x*cos(x/y))*sin(x/y) +
        x**2*sin(x/y)**2/y - x**2*cos(x/y)**2/y)))

cancel(expr)

This gives:

In [2]: cancel(expr)
---------------------------------------------------------------------------
HeuristicGCDFailed                        Traceback (most recent call last)
<ipython-input-2-f5a0e740cc3b> in <module>
----> 1 cancel(expr)

~/current/sympy/sympy/sympy/polys/polytools.py in cancel(f, *gens, **args)
   6729             return f.xreplace(dict(reps))
   6730 
-> 6731     c, (P, Q) = 1, F.cancel(G)
   6732     if opt.get('polys', False) and not 'gens' in opt:
   6733         opt['gens'] = R.symbols

~/current/sympy/sympy/sympy/polys/rings.py in cancel(self, g)
   2222 
   2223         if not (domain.is_Field and domain.has_assoc_Ring):
-> 2224             _, p, q = f.cofactors(g)
   2225         else:
   2226             new_ring = ring.clone(domain=domain.get_ring())

~/current/sympy/sympy/sympy/polys/rings.py in cofactors(f, g)
   2138 
   2139         J, (f, g) = f.deflate(g)
-> 2140         h, cff, cfg = f._gcd(g)
   2141 
   2142         return (h.inflate(J), cff.inflate(J), cfg.inflate(J))

~/current/sympy/sympy/sympy/polys/rings.py in _gcd(f, g)
   2171             return f._gcd_QQ(g)
   2172         elif ring.domain.is_ZZ:
-> 2173             return f._gcd_ZZ(g)
   2174         else: # TODO: don't use dense representation (port PRS algorithms)
   2175             return ring.dmp_inner_gcd(f, g)

~/current/sympy/sympy/sympy/polys/rings.py in _gcd_ZZ(f, g)
   2176 
   2177     def _gcd_ZZ(f, g):
-> 2178         return heugcd(f, g)
   2179 
   2180     def _gcd_QQ(self, g):

~/current/sympy/sympy/sympy/polys/heuristicgcd.py in heugcd(f, g)
     78                 h, cff, cfg = domain.cofactors(ff, gg)
     79             else:
---> 80                 h, cff, cfg = heugcd(ff, gg)
     81 
     82             h = _gcd_interpolate(h, x, ring)

~/current/sympy/sympy/sympy/polys/heuristicgcd.py in heugcd(f, g)
     78                 h, cff, cfg = domain.cofactors(ff, gg)
     79             else:
---> 80                 h, cff, cfg = heugcd(ff, gg)
     81 
     82             h = _gcd_interpolate(h, x, ring)

~/current/sympy/sympy/sympy/polys/heuristicgcd.py in heugcd(f, g)
     78                 h, cff, cfg = domain.cofactors(ff, gg)
     79             else:
---> 80                 h, cff, cfg = heugcd(ff, gg)
     81 
     82             h = _gcd_interpolate(h, x, ring)

~/current/sympy/sympy/sympy/polys/heuristicgcd.py in heugcd(f, g)
    116         x = 73794*x * domain.sqrt(domain.sqrt(x)) // 27011
    117 
--> 118     raise HeuristicGCDFailed('no luck')
    119 
    120 def _gcd_interpolate(h, x, ring):

HeuristicGCDFailed: no luck

[oscarbenjamin]

This is a regression since 1.8. Bisected to 657df3d from #21016

[oscarbenjamin]

A more direct reproducer:

from sympy import *

_, z = ring([Symbol('z')], ZZ)

f = (-213375767148027074625732447889785229321393181328*z**6 +
    29940160751184796497659306841761998431157572554609743759334191398912*z**5
    - 4116497187077889158951734821171228*z**4 +
    770150482473545951341668749547430109371574445215318016*z**3 +
    4952643462600347854340917482003231670272*z + 226981)

g = (4800247059247192345130702637232498966798272*z**6 -
    1347108535294359924305833697408922113390365686018750703441281024*z**5 +
    94510840277438599115334503936503817716223778927494169931732233144176687623854617688*z**4
    + 810367161391734969018769459197659802927544712552024943757613312084716*z**2
    + 222836087494110272180194707970719744*z -
    5211266779035792276688740773745655154995394429159735296)

f._gcd(g)

[oscarbenjamin]

The bug shows up only if you have gmpy (and not gmpy2 installed).

[oscarbenjamin]

The bug is caused by mpmath using gmpy when sympy isn't.

[oscarbenjamin]

With gmpy installed:

In [1]: ZZ.sqrt(2)
Out[1]: mpz(1)

However without

In [1]: ZZ.sqrt(2)
Out[1]: 1

That changes the type of x here:

sympy/sympy/polys/heuristicgcd.py

Line 68 in cc44619

x = max(min(B, 99*domain.sqrt(B)),

In turn that means that here convert has an mpz instead of an int:

sympy/sympy/polys/rings.py

Line 2307 in cc44619

a = ring.domain.convert(a)

Then because sympy isn't using gmpy HAS_GMPY is false so this line skips checking for mpz:

sympy/sympy/polys/domains/domain.py

Lines 429 to 432 in cc44619

	if HAS_GMPY:
	integers = ZZ
	if isinstance(element, integers.tp):
	return self.convert_from(element, integers)

But then we get here:

sympy/sympy/polys/domains/domain.py

Lines 458 to 465 in cc44619

	else: # TODO: remove this branch
	if not is_sequence(element):
	try:
	element = sympify(element, strict=True)
	if isinstance(element, Basic):
	return self.from_sympy(element)
	except (TypeError, ValueError):
	pass

And at this point sympify converts the mpz to a Float which ZZ.from_sympy then rounds back to an int here:

sympy/sympy/polys/domains/integerring.py

Lines 57 to 62 in cc44619

	def from_sympy(self, a):
	"""Convert SymPy's Integer to ``dtype``. """
	if a.is_Integer:
	return MPZ(a.p)
	elif a.is_Float and int(a) == a:
	return MPZ(int(a))

That means that mpz is not exactly converted to an int:

In [2]: from gmpy import mpz

In [3]: ZZ.convert(mpz(12345678901234567890))
Out[3]: 12345678901234567168

The inexact conversion then messes up heugcd.

The root cause is the return type of ZZ.sqrt. There are a number of bug-magnets here though:

1. The sympify branch in Domain.convert should be removed.
2. ZZ.from_sympy should not convert a Float.

These kinds of problems would show up much quicker if a TypeError was raised immediately.

[oscarbenjamin]

A possible fix:

diff --git a/sympy/external/gmpy.py b/sympy/external/gmpy.py
index dd893f6ceb..e387a48d9b 100644
--- a/sympy/external/gmpy.py
+++ b/sympy/external/gmpy.py
@@ -98,4 +98,4 @@
     MPQ = PythonMPQ
 
     factorial = mlib.ifac
-    sqrt = mlib.isqrt
+    sqrt = lambda x: int(mlib.isqrt(x))

[oscarbenjamin]

There are lots of test failures if gmpy is installed (without the diff above):

$ pytest sympy/core/ sympy/polys -n4 -d -m 'not slow'
...
                                                                       DO *NOT* COMMIT!                                                                       
================================================================== short test summary info ===================================================================
FAILED sympy/core/tests/test_arit.py::test_issue_8247_8354 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/core/tests/test_arit.py::test_Add_is_negative_positive - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/core/tests/test_expr.py::test_equals - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_euclidtools.py::test_dup_gcd - assert 3654318780237...73459383222272 == 3654318780237...30321821976049
FAILED sympy/polys/tests/test_factortools.py::test_dup_qq_i_factor - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_factortools.py::test_dmp_qq_i_factor - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_factortools.py::test_dup_zz_i_factor - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_factortools.py::test_dmp_zz_i_factor - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_factortools.py::test_dup_ext_factor - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_factortools.py::test_dmp_ext_factor - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_factortools.py::test_dup_factor_list - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_factortools.py::test_dmp_factor_list - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_issue_19760 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_partfrac.py::test_apart_extension - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_partfrac.py::test_apart_list - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_partfrac.py::test_assemble_partfrac_list - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polyroots.py::test_issue_14522 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polytools.py::test_issue_20427 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_minimal_polynomial - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polytools.py::test_factor_large - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_minimal_polynomial_issue_19732 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_minimal_polynomial_hi_prec - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polytools.py::test_sturm - sympy.polys.polyerrors.HeuristicGCDFailed: no luck
FAILED sympy/polys/tests/test_polyroots.py::test_issue_21287 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_minpoly_compose - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polytools.py::test_cancel - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polytools.py::test_factor - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polyroots.py::test_roots_cyclotomic - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_primitive_element - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polyroots.py::test_roots0 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_rings.py::test_PolyElement_sturm - sympy.polys.polyerrors.HeuristicGCDFailed: no luck
FAILED sympy/polys/tests/test_numberfields.py::test_field_isomorphism - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polyroots.py::test_roots_mixed - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_minpoly_fraction_field - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_polytools.py::test_Poly__unify - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_numberfields.py::test_issue_14831 - TypeError: PythonMPQ() requires numeric or string argument
FAILED sympy/polys/tests/test_solvers.py::test_solve_lin_sys_6x6_2 - sympy.polys.polyerrors.HeuristicGCDFailed: no luck
======================================= 37 failed, 2922 passed, 69 skipped, 43 xfailed, 1 xpassed in 165.32s (0:02:45) =======================================

[oscarbenjamin]

This fixes all test failures in core and polys:

diff --git a/sympy/external/gmpy.py b/sympy/external/gmpy.py
index dd893f6ceb..bdf4b31d33 100644
--- a/sympy/external/gmpy.py
+++ b/sympy/external/gmpy.py
@@ -97,5 +97,5 @@
     MPZ = int
     MPQ = PythonMPQ
 
-    factorial = mlib.ifac
-    sqrt = mlib.isqrt
+    factorial = lambda x: int(mlib.ifac(x))
+    sqrt = lambda x: int(mlib.isqrt(x))