Enable algebraic extensions by default (stage 2)

.github/workflows/test.yml

@@ -49,12 +49,12 @@ jobs:
     - name: Tests
       if: matrix.coverage == false
       env:
-        PYTEST_ADDOPTS: --doctest-modules --color yes
+        PYTEST_ADDOPTS: -n0 --doctest-modules --color yes
       run: pytest
     - name: Coverage tests
       if: matrix.coverage
       env:
-        PYTEST_ADDOPTS: --color yes --cov diofant --cov-append -m "not slow and not xfail and not regression" --ignore docs
+        PYTEST_ADDOPTS: -n0 --color yes --cov diofant --cov-append -m "not slow and not xfail and not regression" --ignore docs
       run: |
         pytest
         pip uninstall -y gmpy2 numpy IPython

diofant/domains/expressiondomain.py

@@ -108,7 +108,7 @@ def __bool__(self):
 
         def gcd(self, other):
             from ..polys import gcd
-            return self.__class__(gcd(self.ex, self.__class__(other).ex))
+            return self.__class__(gcd(self.ex, self.__class__(other).ex, extension=False))
 
         def lcm(self, other):
             from ..polys import lcm

diofant/domains/rationalfield.py

@@ -31,7 +31,7 @@ def from_expr(self, expr):
         if expr.is_Float:
             from . import RR
             return self.dtype(*RR.to_rational(expr))
-        if expr.is_algebraic:
+        if expr.is_number and expr.is_algebraic:
             from ..polys import primitive_element
             try:
                 _, _, [[rep]] = primitive_element([expr], domain=self)

diofant/integrals/meijerint.py

@@ -1335,7 +1335,7 @@ def _rewrite_single(f, x, recursive=True):
 
     if isinstance(f, meijerg):
         from ..polys import factor
-        coeff, m = factor(f.argument, x).as_coeff_mul(x)
+        coeff, m = factor(f.argument, x, extension=False).as_coeff_mul(x)
         if len(m) > 1:
             return
         m = m[0]

diofant/polys/constructor.py

@@ -171,7 +171,7 @@ def _construct_composite(coeffs, opt):
             coeffs.update(list(numer.values()))
             coeffs.update(list(denom.values()))
 
-    rationals, reals = False, False
+    algebraics, rationals, reals = False, False, False
 
     for coeff in coeffs:
         if coeff.is_Rational:
@@ -180,11 +180,15 @@ def _construct_composite(coeffs, opt):
         elif coeff.is_Float:
             reals = True
             break
+        elif coeff.is_number and coeff.is_algebraic:
+            algebraics = True
         else:
             raise NotImplementedError
 
     if reals:
         ground = RR
+    elif algebraics:
+        ground, _ = _construct_algebraic(coeffs, opt)
     elif rationals:
         ground = QQ
     else:

diofant/polys/factorization_alg_field.py

@@ -661,7 +661,7 @@ def _diophantine(F, c, A, d, minpoly, p):
                 C = C.quo_ground(ring.domain.factorial(k + 1))
                 T = _diophantine(G, C, A, d, minpoly, p)
                 if T is None:
-                    raise NotImplementedError
+                    return
 
                 for i, t in enumerate(T):
                     T[i] = t.set_ring(ring) * M
@@ -673,8 +673,6 @@ def _diophantine(F, c, A, d, minpoly, p):
                     c = c - t * b
 
                 c = _trunc(c, minpoly, p)
-            else:
-                raise NotImplementedError
 
         S = [_trunc(s, minpoly, p) for s in S]
 

diofant/polys/polyroots.py

@@ -533,7 +533,7 @@ def roots_quintic(f):
     a, b = symbols('a b', cls=Dummy)
     _sol = _solve(sol**5 - a - I*b, sol)
     for i in range(5):
-        _sol[i] = factor(_sol[i][sol])
+        _sol[i] = factor(_sol[i][sol], extension=False)
     R1 = R1.as_real_imag()
     R2 = R2.as_real_imag()
     R3 = R3.as_real_imag()
@@ -594,7 +594,7 @@ def roots_quintic(f):
 
 def _quintic_simplify(expr):
     expr = powsimp(expr)
-    expr = cancel(expr)
+    expr = cancel(expr, extension=False)
     return together(expr)
 
 
@@ -859,7 +859,7 @@ def _try_heuristics(f):
         n = f.degree()
 
         if n == 2:
-            result += list(map(cancel, roots_quadratic(f)))
+            result += roots_quadratic(f)
         elif f.is_cyclotomic:
             result += roots_cyclotomic(f)
         elif n == 3 and cubics:
@@ -899,7 +899,7 @@ def _try_heuristics(f):
             for r in roots_fun(f):
                 _update_dict(result, r, 1)
         else:
-            _, factors = f.as_expr().as_poly().factor_list()
+            _, factors = f.as_expr().as_poly(extension=False).factor_list()
             if len(factors) == 1 and f.degree() == 2:
                 for r in roots_quadratic(f):
                     _update_dict(result, r, 1)

diofant/polys/polyutils.py

@@ -131,9 +131,10 @@ def _is_coeff(expr, opt):
     """Check if expr is a coefficient."""
     if opt.domain is not None:
         return expr in opt.domain
-    if opt.extension:
-        return expr.is_number and expr.is_algebraic
     if opt.greedy is not False:
+        if opt.extension is not False:
+            return ((expr.is_Number and expr.is_finite is not False) or
+                    (expr.is_number and expr.is_algebraic))
         return expr.is_Number and expr.is_finite is not False
     return expr.is_number and expr.is_finite is not False
 

diofant/simplify/hyperexpand.py

@@ -2022,7 +2022,7 @@ def carryout_plan(f, ops):
     if res is not None:
         p = apply_operators(res, ops, lambda f: z0*f.diff(z0))
         p = apply_operators(p*premult, ops0, lambda f: z0*f.diff(z0))
-        return unpolarify(simplify(p).subs({z0: z}))
+        return unpolarify(simplify(unpolarify(p)).subs({z0: z}))
 
     # Try to recognise a shifted sum.
     p = Integer(0)

diofant/simplify/ratsimp.py

@@ -17,7 +17,7 @@ def ratsimp(expr):
     (x + y)/(x*y)
 
     """
-    f, g = cancel(expr).as_numer_denom()
+    f, g = cancel(expr, extension=False).as_numer_denom()
     try:
         Q, r = reduced(f, [g], field=True, expand=False)
     except ComputationFailedError:

diofant/simplify/simplify.py

@@ -618,8 +618,8 @@ def shorter(*choices):
     expr = bottom_up(lambda e: e if isinstance(e, (Integral, Product, Sum)) else e.doit(deep=False))(expr)
     expr = bottom_up(lambda w: w.normal() if hasattr(w, 'normal') else w)(expr)
     expr = Mul(*powsimp(expr).as_content_primitive())
-    _e = cancel(expr)
-    expr1 = shorter(_e, _mexpand(_e).cancel())  # issue sympy/sympy#6829
+    _e = cancel(expr, extension=False)
+    expr1 = shorter(_e, _mexpand(_e).cancel(extension=False))  # issue sympy/sympy#6829
     expr2 = shorter(together(expr, deep=True), together(expr1, deep=True))
 
     if ratio is oo:

diofant/simplify/sqrtdenest.py

@@ -679,7 +679,7 @@ def _canonical(eq, cov):
         # make the sign canonical
         free = eq.free_symbols
         if len(free) == 1:
-            if eq.coeff(free.pop()**degree(eq)).could_extract_minus_sign():
+            if eq.coeff(free.pop()**degree(eq, extension=False)).could_extract_minus_sign():
                 eq = -eq
         elif eq.could_extract_minus_sign():
             eq = -eq
@@ -721,7 +721,7 @@ def _take(d, take_int_pow):
         return
 
     syms = set(syms) or eq.free_symbols
-    poly = eq.as_poly()
+    poly = eq.as_poly(extension=False)
     gens = [g for g in poly.gens if _take(g, True)]
     if not gens:
         return

diofant/simplify/trigsimp.py

@@ -1123,16 +1123,16 @@ def trigs(x):
                 TR3,  # canonical angles
                 TR1,  # sec-csc -> cos-sin
                 TR12,  # expand tan of sum
-                lambda x: _eapply(factor, x, trigs),
+                lambda x: _eapply(lambda e: factor(e, extension=False), x, trigs),
                 TR2,  # tan-cot -> sin-cos
                 [identity, lambda x: _eapply(_mexpand, x, trigs)],
                 TR2i,  # sin-cos ratio -> tan
-                lambda x: _eapply(lambda i: factor(i.normal()), x, trigs),
+                lambda x: _eapply(lambda i: factor(i.normal(), extension=False), x, trigs),
                 TR14,  # factored identities
                 TR5,  # sin-pow -> cos_pow
                 TR10,  # sin-cos of sums -> sin-cos prod
                 TR11, TR6,  # reduce double angles and rewrite cos pows
-                lambda x: _eapply(factor, x, trigs),
+                lambda x: _eapply(lambda e: factor(e, extension=False), x, trigs),
                 TR14,  # factored powers of identities
                 [identity, lambda x: _eapply(_mexpand, x, trigs)],
                 TRmorrie,

diofant/solvers/ode.py

@@ -2564,7 +2564,7 @@ def constantsimp(expr, constants):
     expr = __remove_linear_redundancies(expr, Cs)
 
     def _conditional_term_factoring(expr):
-        new_expr = terms_gcd(expr, clear=False, deep=True, expand=False)
+        new_expr = terms_gcd(expr, clear=False, deep=True, expand=False, extension=False)
 
         # we do not want to factor exponentials, so handle this separately
         if new_expr.is_Mul:

diofant/solvers/solvers.py

@@ -514,7 +514,7 @@ def _solve(f, symbol, **flags):
         # as a polynomial, followed (perhaps) by a change of variables if the
         # generator is not a symbol
 
-        poly = Poly(f_num.expand(power_base=False, log=False))
+        poly = Poly(f_num.expand(power_base=False, log=False), extension=False)
         gens = [g for g in poly.gens if g.has(symbol)]
 
         def _as_base_q(x):

diofant/tests/calculus/test_demidovich.py

@@ -1,8 +1,8 @@
 # Numbers listed with the tests refer to problem numbers in books
 # :cite:`demidovich1970problems` and :cite:`anti-demidovich2001` (with '*').
 
-from diofant import (Rational, asin, cos, exp, limit, log, oo, pi, root, sign,
-                     sin, sqrt, tan)
+from diofant import (I, Rational, asin, cos, exp, limit, log, oo, pi, root,
+                     sign, sin, sqrt, tan)
 from diofant.abc import a, h, k, m, n, x
 
 
@@ -68,7 +68,7 @@ def test_Limits_simple_4b():
     # issue sympy/sympy#3511
     assert limit(x - root(x**3 - 1, 3), x, oo) == 0
     assert limit(x + root(1 - x**3, 3),
-                 x, oo) == oo*sign(1 + root(-1, 3))  # 215
+                 x, oo) == oo*sign(3 + I*sqrt(3))  # 215
 
 
 def test_Limits_simple_4c():

diofant/tests/functions/test_spec_polynomials.py

@@ -124,12 +124,10 @@ def test_legendre():
     assert legendre(10, 0) != 0
     assert legendre(11, 0) == 0
 
-    assert roots(legendre(4, x), x) == {
-        sqrt(Rational(3, 7) - Rational(2, 35)*sqrt(30)): 1,
-        -sqrt(Rational(3, 7) - Rational(2, 35)*sqrt(30)): 1,
-        sqrt(Rational(3, 7) + Rational(2, 35)*sqrt(30)): 1,
-        -sqrt(Rational(3, 7) + Rational(2, 35)*sqrt(30)): 1,
-    }
+    assert roots(legendre(4, x), x) == {-sqrt(35)*sqrt(-2*sqrt(30) + 15)/35: 1,
+                                        +sqrt(35)*sqrt(-2*sqrt(30) + 15)/35: 1,
+                                        -sqrt(35)*sqrt(+2*sqrt(30) + 15)/35: 1,
+                                        +sqrt(35)*sqrt(+2*sqrt(30) + 15)/35: 1}
 
     X = legendre(n, x)
     assert isinstance(X, legendre)

diofant/tests/integrals/test_rationaltools.py

@@ -105,7 +105,7 @@ def test_ratint():
            sqrt(2)*log(x**2 + x*(-2 + sqrt(2)) - sqrt(2) + 2)/8 -
            sqrt(2)*atan(-sqrt(2)*x + 1 + sqrt(2))/4 +
            sqrt(2)*atan(sqrt(2)*x - sqrt(2) + 1)/4)
-    assert ratint(1/((x - 1)**4 + 1), x) == ans
+    assert (ratint(1/((x - 1)**4 + 1), x) - ans).equals(0) is True
 
     ans = RootSum(776887*t**7 + 27216*t**5 - 15120*t**4 + 3780*t**3 -
                   504*t**2 + 35*t - 1,

diofant/tests/polys/test_constructor.py

@@ -95,7 +95,9 @@ def test_construct_domain():
 
     assert construct_domain(x**2 + 2*x + E) == (dom, dom(x**2 + 2*x + E))
 
-    assert construct_domain(x + y + GoldenRatio) == (EX, EX(x + y + GoldenRatio))
+    dom = QQ.algebraic_field(GoldenRatio).inject(x, y)
+
+    assert construct_domain(x + y + GoldenRatio) == (dom, dom(x + y + GoldenRatio))
 
 
 def test_composite_option():

diofant/tests/polys/test_polytools.py

@@ -2619,12 +2619,12 @@ def test_cancel():
 
     f = (x**2 - 2)/(x + sqrt(2))
 
-    assert cancel(f) == f
+    assert cancel(f, extension=False) == f
     assert cancel(f, greedy=False) == x - sqrt(2)
 
     f = (x**2 - 2)/(x - sqrt(2))
 
-    assert cancel(f) == f
+    assert cancel(f, extension=False) == f
     assert cancel(f, greedy=False) == x + sqrt(2)
 
     assert cancel((x**2/4 - 1, x/2 - 1)) == (Rational(1, 2), x + 2, 1)
@@ -2651,7 +2651,7 @@ def test_cancel():
     f = x**3 + (sqrt(2) - 2)*x**2 - (2*sqrt(2) + 3)*x - 3*sqrt(2)
     g = x**2 - 2
 
-    assert cancel((f, g), extension=True) == (1, x**2 - 2*x - 3, x - sqrt(2))
+    assert cancel((f, g)) == (1, x**2 - 2*x - 3, x - sqrt(2))
 
     f = (-2*x + 3).as_poly()
     g = (-x**9 + x**8 + x**6 - x**5 + 2*x**2 - 3*x + 1).as_poly()
@@ -3286,6 +3286,10 @@ def test_sympyissue_22673():
                                   for i in range(7)])
 
 
+def test_sympyissue_18391():
+    assert gcd(x**2 + 1, x + I) == x + I
+
+
 @pytest.mark.timeout(10)
 def test_sympyissue_23766():
     assert factor(exp((x + 1)**25 + 1),

diofant/tests/polys/test_polyutils.py

@@ -243,11 +243,38 @@ def test__dict_from_expr_no_gens():
     assert _parallel_dict_from_expr([x], opt) == ([{(1,): 1}], build_options([x], {}))
     assert _parallel_dict_from_expr([y], opt) == ([{(1,): 1}], build_options([y], {}))
 
+    opt = build_options([], {'extension': False})
+    assert _parallel_dict_from_expr([sqrt(2)],
+                                    opt) == ([{(1,): 1}],
+                                             build_options([sqrt(2)],
+                                                           {'extension': False}))
+    opt = build_options([], {'greedy': False})
+    pytest.raises(GeneratorsNeededError,
+                  lambda: _parallel_dict_from_expr([sqrt(2)], opt))
+
+    opt = build_options([], {'domain': ZZ.inject(x)})
+    assert _parallel_dict_from_expr([x*y],
+                                    opt) == ([{(1,): x}], build_options([y], {'domain': ZZ.inject(x)}))
+    opt = build_options([], {'domain': ZZ.inject(y)})
+    assert _parallel_dict_from_expr([x*y],
+                                    opt) == ([{(1,): y}],
+                                             build_options([x], {'domain': ZZ.inject(y)}))
+
+    opt = build_options([], {'extension': False})
+    assert _parallel_dict_from_expr([3*sqrt(2)*pi*x*y],
+                                    opt) == ([{(1, 1, 1, 1): 3}],
+                                             build_options([x, y, pi, sqrt(2)], {'extension': False}))
+    opt = build_options([], {})
+    assert _parallel_dict_from_expr([3*sqrt(2)*pi*x*y], opt) == ([{(1, 1, 1):
+                                                                   3*sqrt(2)}],
+                                                                 build_options([x, y, pi], {}))
+
     assert _parallel_dict_from_expr([x*y], opt) == ([{(1, 1): 1}], build_options([x, y], {}))
     assert _parallel_dict_from_expr([x + y], opt) == ([{(1, 0): 1, (0, 1): 1}],
                                                       build_options([x, y], {}))
 
-    assert _parallel_dict_from_expr([sqrt(2)], opt) == ([{(1,): 1}], build_options([sqrt(2)], {}))
+    pytest.raises(GeneratorsNeededError,
+                  lambda: _parallel_dict_from_expr([sqrt(2)], opt))
 
     f = cos(x)*sin(x) + cos(x)*sin(y) + cos(y)*sin(x) + cos(y)*sin(y)
 
@@ -265,14 +292,10 @@ def test__dict_from_expr_no_gens():
     opt = build_options([], {'domain': ZZ.inject(y)})
     assert _parallel_dict_from_expr([x*y], opt) == ([{(1,): y}], build_options([x], {'domain': ZZ.inject(y)}))
 
-    opt = build_options([], {'extension': None})
-    assert _parallel_dict_from_expr([3*sqrt(2)*pi*x*y],
-                                    opt) == ([{(1, 1, 1, 1): 3}],
-                                             build_options([x, y, pi, sqrt(2)], {'extension': None}))
-    opt = build_options([], {'extension': True})
+    opt = build_options([], {})
     assert _parallel_dict_from_expr([3*sqrt(2)*pi*x*y],
                                     opt) == ([{(1, 1, 1): 3*sqrt(2)}],
-                                             build_options([x, y, pi], {'extension': True}))
+                                             build_options([x, y, pi], {}))
 
 
 def test__parallel_dict_from_expr_if_gens():

diofant/tests/test_wester.py

@@ -407,7 +407,7 @@ def test_H20():
     f = x**3 + (sqrt(2) - 2)*x**2 - (2*sqrt(2) + 3)*x - 3*sqrt(2)
     g = x**2 - 2
     r = (x**2 - 2*x - 3)/(x - sqrt(2))
-    assert cancel(f/g, extension=True) == cancel(f/g, extension=sqrt(2)) == r
+    assert cancel(f/g) == cancel(f/g, extension=sqrt(2)) == r
 
 
 def test_H22():
@@ -793,8 +793,8 @@ def test_M22():
 
 
 def test_M23():
-    assert solve(x - 1/sqrt(1 + x**2)) == [{x: -I*sqrt(Rational(1, 2) + sqrt(5)/2)},
-                                           {x: sqrt(Rational(-1, 2) + sqrt(5)/2)}]
+    assert solve(x - 1/sqrt(1 + x**2)) == [{x: sqrt(2)*sqrt(-1 + sqrt(5))/2},
+                                           {x: -sqrt(2)*I*sqrt(1 + sqrt(5))/2}]
 
 
 def test_M24():