"""OO layer for several polynomial representations. """ from __future__ import annotations from sympy.external.gmpy import GROUND_TYPES from sympy.utilities.exceptions import sympy_deprecation_warning from sympy.core.numbers import oo from sympy.core.sympify import CantSympify from sympy.polys.polyutils import PicklableWithSlots, _sort_factors from sympy.polys.domains import Domain, ZZ, QQ from sympy.polys.polyerrors import ( CoercionFailed, ExactQuotientFailed, DomainError, NotInvertible, ) from sympy.polys.densebasic import ( ninf, dmp_validate, dup_normal, dmp_normal, dup_convert, dmp_convert, dmp_from_sympy, dup_strip, dmp_degree_in, dmp_degree_list, dmp_negative_p, dmp_ground_LC, dmp_ground_TC, dmp_ground_nth, dmp_one, dmp_ground, dmp_zero, dmp_zero_p, dmp_one_p, dmp_ground_p, dup_from_dict, dmp_from_dict, dmp_to_dict, dmp_deflate, dmp_inject, dmp_eject, dmp_terms_gcd, dmp_list_terms, dmp_exclude, dup_slice, dmp_slice_in, dmp_permute, dmp_to_tuple,) from sympy.polys.densearith import ( dmp_add_ground, dmp_sub_ground, dmp_mul_ground, dmp_quo_ground, dmp_exquo_ground, dmp_abs, dmp_neg, dmp_add, dmp_sub, dmp_mul, dmp_sqr, dmp_pow, dmp_pdiv, dmp_prem, dmp_pquo, dmp_pexquo, dmp_div, dmp_rem, dmp_quo, dmp_exquo, dmp_add_mul, dmp_sub_mul, dmp_max_norm, dmp_l1_norm, dmp_l2_norm_squared) from sympy.polys.densetools import ( dmp_clear_denoms, dmp_integrate_in, dmp_diff_in, dmp_eval_in, dup_revert, dmp_ground_trunc, dmp_ground_content, dmp_ground_primitive, dmp_ground_monic, dmp_compose, dup_decompose, dup_shift, dmp_shift, dup_transform, dmp_lift) from sympy.polys.euclidtools import ( dup_half_gcdex, dup_gcdex, dup_invert, dmp_subresultants, dmp_resultant, dmp_discriminant, dmp_inner_gcd, dmp_gcd, dmp_lcm, dmp_cancel) from sympy.polys.sqfreetools import ( dup_gff_list, dmp_norm, dmp_sqf_p, dmp_sqf_norm, dmp_sqf_part, dmp_sqf_list, dmp_sqf_list_include) from sympy.polys.factortools import ( dup_cyclotomic_p, dmp_irreducible_p, dmp_factor_list, dmp_factor_list_include) from sympy.polys.rootisolation import ( dup_isolate_real_roots_sqf, dup_isolate_real_roots, dup_isolate_all_roots_sqf, dup_isolate_all_roots, dup_refine_real_root, dup_count_real_roots, dup_count_complex_roots, dup_sturm, dup_cauchy_upper_bound, dup_cauchy_lower_bound, dup_mignotte_sep_bound_squared) from sympy.polys.polyerrors import ( UnificationFailed, PolynomialError) _flint_domains: tuple[Domain, ...] if GROUND_TYPES == 'flint': import flint _flint_domains = (ZZ, QQ) else: flint = None _flint_domains = () class DMP(CantSympify): """Dense Multivariate Polynomials over `K`. """ __slots__ = () def __new__(cls, rep, dom, lev=None): if lev is None: rep, lev = dmp_validate(rep) elif not isinstance(rep, list): raise CoercionFailed("expected list, got %s" % type(rep)) return cls.new(rep, dom, lev) @classmethod def new(cls, rep, dom, lev): # It would be too slow to call _validate_args always at runtime. # Ideally this checking would be handled by a static type checker. # #cls._validate_args(rep, dom, lev) if flint is not None: if lev == 0 and dom in _flint_domains: return DUP_Flint._new(rep, dom, lev) return DMP_Python._new(rep, dom, lev) @property def rep(f): """Get the representation of ``f``. """ sympy_deprecation_warning(""" Accessing the ``DMP.rep`` attribute is deprecated. The internal representation of ``DMP`` instances can now be ``DUP_Flint`` when the ground types are ``flint``. In this case the ``DMP`` instance does not have a ``rep`` attribute. Use ``DMP.to_list()`` instead. Using ``DMP.to_list()`` also works in previous versions of SymPy. """, deprecated_since_version="1.13", active_deprecations_target="dmp-rep", ) return f.to_list() def to_best(f): """Convert to DUP_Flint if possible. This method should be used when the domain or level is changed and it potentially becomes possible to convert from DMP_Python to DUP_Flint. """ if flint is not None: if isinstance(f, DMP_Python) and f.lev == 0 and f.dom in _flint_domains: return DUP_Flint.new(f._rep, f.dom, f.lev) return f @classmethod def _validate_args(cls, rep, dom, lev): assert isinstance(dom, Domain) assert isinstance(lev, int) and lev >= 0 def validate_rep(rep, lev): assert isinstance(rep, list) if lev == 0: assert all(dom.of_type(c) for c in rep) else: for r in rep: validate_rep(r, lev - 1) validate_rep(rep, lev) @classmethod def from_dict(cls, rep, lev, dom): rep = dmp_from_dict(rep, lev, dom) return cls.new(rep, dom, lev) @classmethod def from_list(cls, rep, lev, dom): """Create an instance of ``cls`` given a list of native coefficients. """ return cls.new(dmp_convert(rep, lev, None, dom), dom, lev) @classmethod def from_sympy_list(cls, rep, lev, dom): """Create an instance of ``cls`` given a list of SymPy coefficients. """ return cls.new(dmp_from_sympy(rep, lev, dom), dom, lev) @classmethod def from_monoms_coeffs(cls, monoms, coeffs, lev, dom): return cls(dict(list(zip(monoms, coeffs))), dom, lev) def convert(f, dom): """Convert ``f`` to a ``DMP`` over the new domain. """ if f.dom == dom: return f elif f.lev or flint is None: return f._convert(dom) elif isinstance(f, DUP_Flint): if dom in _flint_domains: return f._convert(dom) else: return f.to_DMP_Python()._convert(dom) elif isinstance(f, DMP_Python): if dom in _flint_domains: return f._convert(dom).to_DUP_Flint() else: return f._convert(dom) else: raise RuntimeError("unreachable code") def _convert(f, dom): raise NotImplementedError @classmethod def zero(cls, lev, dom): return DMP(dmp_zero(lev), dom, lev) @classmethod def one(cls, lev, dom): return DMP(dmp_one(lev, dom), dom, lev) def _one(f): raise NotImplementedError def __repr__(f): return "%s(%s, %s)" % (f.__class__.__name__, f.to_list(), f.dom) def __hash__(f): return hash((f.__class__.__name__, f.to_tuple(), f.lev, f.dom)) def __getnewargs__(self): return self.to_list(), self.dom, self.lev def ground_new(f, coeff): """Construct a new ground instance of ``f``. """ raise NotImplementedError def unify_DMP(f, g): """Unify and return ``DMP`` instances of ``f`` and ``g``. """ if not isinstance(g, DMP) or f.lev != g.lev: raise UnificationFailed("Cannot unify %s with %s" % (f, g)) if f.dom != g.dom: dom = f.dom.unify(g.dom) f = f.convert(dom) g = g.convert(dom) return f, g def to_dict(f, zero=False): """Convert ``f`` to a dict representation with native coefficients. """ return dmp_to_dict(f.to_list(), f.lev, f.dom, zero=zero) def to_sympy_dict(f, zero=False): """Convert ``f`` to a dict representation with SymPy coefficients. """ rep = f.to_dict(zero=zero) for k, v in rep.items(): rep[k] = f.dom.to_sympy(v) return rep def to_sympy_list(f): """Convert ``f`` to a list representation with SymPy coefficients. """ def sympify_nested_list(rep): out = [] for val in rep: if isinstance(val, list): out.append(sympify_nested_list(val)) else: out.append(f.dom.to_sympy(val)) return out return sympify_nested_list(f.to_list()) def to_list(f): """Convert ``f`` to a list representation with native coefficients. """ raise NotImplementedError def to_tuple(f): """ Convert ``f`` to a tuple representation with native coefficients. This is needed for hashing. """ raise NotImplementedError def to_ring(f): """Make the ground domain a ring. """ return f.convert(f.dom.get_ring()) def to_field(f): """Make the ground domain a field. """ return f.convert(f.dom.get_field()) def to_exact(f): """Make the ground domain exact. """ return f.convert(f.dom.get_exact()) def slice(f, m, n, j=0): """Take a continuous subsequence of terms of ``f``. """ if not f.lev and not j: return f._slice(m, n) else: return f._slice_lev(m, n, j) def _slice(f, m, n): raise NotImplementedError def _slice_lev(f, m, n, j): raise NotImplementedError def coeffs(f, order=None): """Returns all non-zero coefficients from ``f`` in lex order. """ return [ c for _, c in f.terms(order=order) ] def monoms(f, order=None): """Returns all non-zero monomials from ``f`` in lex order. """ return [ m for m, _ in f.terms(order=order) ] def terms(f, order=None): """Returns all non-zero terms from ``f`` in lex order. """ if f.is_zero: zero_monom = (0,)*(f.lev + 1) return [(zero_monom, f.dom.zero)] else: return f._terms(order=order) def _terms(f, order=None): raise NotImplementedError def all_coeffs(f): """Returns all coefficients from ``f``. """ if f.lev: raise PolynomialError('multivariate polynomials not supported') if not f: return [f.dom.zero] else: return list(f.to_list()) def all_monoms(f): """Returns all monomials from ``f``. """ if f.lev: raise PolynomialError('multivariate polynomials not supported') n = f.degree() if n < 0: return [(0,)] else: return [ (n - i,) for i, c in enumerate(f.to_list()) ] def all_terms(f): """Returns all terms from a ``f``. """ if f.lev: raise PolynomialError('multivariate polynomials not supported') n = f.degree() if n < 0: return [((0,), f.dom.zero)] else: return [ ((n - i,), c) for i, c in enumerate(f.to_list()) ] def lift(f): """Convert algebraic coefficients to rationals. """ return f._lift().to_best() def _lift(f): raise NotImplementedError def deflate(f): """Reduce degree of `f` by mapping `x_i^m` to `y_i`. """ raise NotImplementedError def inject(f, front=False): """Inject ground domain generators into ``f``. """ raise NotImplementedError def eject(f, dom, front=False): """Eject selected generators into the ground domain. """ raise NotImplementedError def exclude(f): r""" Remove useless generators from ``f``. Returns the removed generators and the new excluded ``f``. Examples ======== >>> from sympy.polys.polyclasses import DMP >>> from sympy.polys.domains import ZZ >>> DMP([[[ZZ(1)]], [[ZZ(1)], [ZZ(2)]]], ZZ).exclude() ([2], DMP_Python([[1], [1, 2]], ZZ)) """ J, F = f._exclude() return J, F.to_best() def _exclude(f): raise NotImplementedError def permute(f, P): r""" Returns a polynomial in `K[x_{P(1)}, ..., x_{P(n)}]`. Examples ======== >>> from sympy.polys.polyclasses import DMP >>> from sympy.polys.domains import ZZ >>> DMP([[[ZZ(2)], [ZZ(1), ZZ(0)]], [[]]], ZZ).permute([1, 0, 2]) DMP_Python([[[2], []], [[1, 0], []]], ZZ) >>> DMP([[[ZZ(2)], [ZZ(1), ZZ(0)]], [[]]], ZZ).permute([1, 2, 0]) DMP_Python([[[1], []], [[2, 0], []]], ZZ) """ return f._permute(P) def _permute(f, P): raise NotImplementedError def terms_gcd(f): """Remove GCD of terms from the polynomial ``f``. """ raise NotImplementedError def abs(f): """Make all coefficients in ``f`` positive. """ raise NotImplementedError def neg(f): """Negate all coefficients in ``f``. """ raise NotImplementedError def add_ground(f, c): """Add an element of the ground domain to ``f``. """ return f._add_ground(f.dom.convert(c)) def sub_ground(f, c): """Subtract an element of the ground domain from ``f``. """ return f._sub_ground(f.dom.convert(c)) def mul_ground(f, c): """Multiply ``f`` by a an element of the ground domain. """ return f._mul_ground(f.dom.convert(c)) def quo_ground(f, c): """Quotient of ``f`` by a an element of the ground domain. """ return f._quo_ground(f.dom.convert(c)) def exquo_ground(f, c): """Exact quotient of ``f`` by a an element of the ground domain. """ return f._exquo_ground(f.dom.convert(c)) def add(f, g): """Add two multivariate polynomials ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._add(G) def sub(f, g): """Subtract two multivariate polynomials ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._sub(G) def mul(f, g): """Multiply two multivariate polynomials ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._mul(G) def sqr(f): """Square a multivariate polynomial ``f``. """ return f._sqr() def pow(f, n): """Raise ``f`` to a non-negative power ``n``. """ if not isinstance(n, int): raise TypeError("``int`` expected, got %s" % type(n)) return f._pow(n) def pdiv(f, g): """Polynomial pseudo-division of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._pdiv(G) def prem(f, g): """Polynomial pseudo-remainder of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._prem(G) def pquo(f, g): """Polynomial pseudo-quotient of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._pquo(G) def pexquo(f, g): """Polynomial exact pseudo-quotient of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._pexquo(G) def div(f, g): """Polynomial division with remainder of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._div(G) def rem(f, g): """Computes polynomial remainder of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._rem(G) def quo(f, g): """Computes polynomial quotient of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._quo(G) def exquo(f, g): """Computes polynomial exact quotient of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._exquo(G) def _add_ground(f, c): raise NotImplementedError def _sub_ground(f, c): raise NotImplementedError def _mul_ground(f, c): raise NotImplementedError def _quo_ground(f, c): raise NotImplementedError def _exquo_ground(f, c): raise NotImplementedError def _add(f, g): raise NotImplementedError def _sub(f, g): raise NotImplementedError def _mul(f, g): raise NotImplementedError def _sqr(f): raise NotImplementedError def _pow(f, n): raise NotImplementedError def _pdiv(f, g): raise NotImplementedError def _prem(f, g): raise NotImplementedError def _pquo(f, g): raise NotImplementedError def _pexquo(f, g): raise NotImplementedError def _div(f, g): raise NotImplementedError def _rem(f, g): raise NotImplementedError def _quo(f, g): raise NotImplementedError def _exquo(f, g): raise NotImplementedError def degree(f, j=0): """Returns the leading degree of ``f`` in ``x_j``. """ if not isinstance(j, int): raise TypeError("``int`` expected, got %s" % type(j)) return f._degree(j) def _degree(f, j): raise NotImplementedError def degree_list(f): """Returns a list of degrees of ``f``. """ raise NotImplementedError def total_degree(f): """Returns the total degree of ``f``. """ raise NotImplementedError def homogenize(f, s): """Return homogeneous polynomial of ``f``""" td = f.total_degree() result = {} new_symbol = (s == len(f.terms()[0][0])) for term in f.terms(): d = sum(term[0]) if d < td: i = td - d else: i = 0 if new_symbol: result[term[0] + (i,)] = term[1] else: l = list(term[0]) l[s] += i result[tuple(l)] = term[1] return DMP.from_dict(result, f.lev + int(new_symbol), f.dom) def homogeneous_order(f): """Returns the homogeneous order of ``f``. """ if f.is_zero: return -oo monoms = f.monoms() tdeg = sum(monoms[0]) for monom in monoms: _tdeg = sum(monom) if _tdeg != tdeg: return None return tdeg def LC(f): """Returns the leading coefficient of ``f``. """ raise NotImplementedError def TC(f): """Returns the trailing coefficient of ``f``. """ raise NotImplementedError def nth(f, *N): """Returns the ``n``-th coefficient of ``f``. """ if all(isinstance(n, int) for n in N): return f._nth(N) else: raise TypeError("a sequence of integers expected") def _nth(f, N): raise NotImplementedError def max_norm(f): """Returns maximum norm of ``f``. """ raise NotImplementedError def l1_norm(f): """Returns l1 norm of ``f``. """ raise NotImplementedError def l2_norm_squared(f): """Return squared l2 norm of ``f``. """ raise NotImplementedError def clear_denoms(f): """Clear denominators, but keep the ground domain. """ raise NotImplementedError def integrate(f, m=1, j=0): """Computes the ``m``-th order indefinite integral of ``f`` in ``x_j``. """ if not isinstance(m, int): raise TypeError("``int`` expected, got %s" % type(m)) if not isinstance(j, int): raise TypeError("``int`` expected, got %s" % type(j)) return f._integrate(m, j) def _integrate(f, m, j): raise NotImplementedError def diff(f, m=1, j=0): """Computes the ``m``-th order derivative of ``f`` in ``x_j``. """ if not isinstance(m, int): raise TypeError("``int`` expected, got %s" % type(m)) if not isinstance(j, int): raise TypeError("``int`` expected, got %s" % type(j)) return f._diff(m, j) def _diff(f, m, j): raise NotImplementedError def eval(f, a, j=0): """Evaluates ``f`` at the given point ``a`` in ``x_j``. """ if not isinstance(j, int): raise TypeError("``int`` expected, got %s" % type(j)) elif not (0 <= j <= f.lev): raise ValueError("invalid variable index %s" % j) if f.lev: return f._eval_lev(a, j) else: return f._eval(a) def _eval(f, a): raise NotImplementedError def _eval_lev(f, a, j): raise NotImplementedError def half_gcdex(f, g): """Half extended Euclidean algorithm, if univariate. """ F, G = f.unify_DMP(g) if F.lev: raise ValueError('univariate polynomial expected') return F._half_gcdex(G) def _half_gcdex(f, g): raise NotImplementedError def gcdex(f, g): """Extended Euclidean algorithm, if univariate. """ F, G = f.unify_DMP(g) if F.lev: raise ValueError('univariate polynomial expected') if not F.dom.is_Field: raise DomainError('ground domain must be a field') return F._gcdex(G) def _gcdex(f, g): raise NotImplementedError def invert(f, g): """Invert ``f`` modulo ``g``, if possible. """ F, G = f.unify_DMP(g) if F.lev: raise ValueError('univariate polynomial expected') return F._invert(G) def _invert(f, g): raise NotImplementedError def revert(f, n): """Compute ``f**(-1)`` mod ``x**n``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._revert(n) def _revert(f, n): raise NotImplementedError def subresultants(f, g): """Computes subresultant PRS sequence of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._subresultants(G) def _subresultants(f, g): raise NotImplementedError def resultant(f, g, includePRS=False): """Computes resultant of ``f`` and ``g`` via PRS. """ F, G = f.unify_DMP(g) if includePRS: return F._resultant_includePRS(G) else: return F._resultant(G) def _resultant(f, g, includePRS=False): raise NotImplementedError def discriminant(f): """Computes discriminant of ``f``. """ raise NotImplementedError def cofactors(f, g): """Returns GCD of ``f`` and ``g`` and their cofactors. """ F, G = f.unify_DMP(g) return F._cofactors(G) def _cofactors(f, g): raise NotImplementedError def gcd(f, g): """Returns polynomial GCD of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._gcd(G) def _gcd(f, g): raise NotImplementedError def lcm(f, g): """Returns polynomial LCM of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._lcm(G) def _lcm(f, g): raise NotImplementedError def cancel(f, g, include=True): """Cancel common factors in a rational function ``f/g``. """ F, G = f.unify_DMP(g) if include: return F._cancel_include(G) else: return F._cancel(G) def _cancel(f, g): raise NotImplementedError def _cancel_include(f, g): raise NotImplementedError def trunc(f, p): """Reduce ``f`` modulo a constant ``p``. """ return f._trunc(f.dom.convert(p)) def _trunc(f, p): raise NotImplementedError def monic(f): """Divides all coefficients by ``LC(f)``. """ raise NotImplementedError def content(f): """Returns GCD of polynomial coefficients. """ raise NotImplementedError def primitive(f): """Returns content and a primitive form of ``f``. """ raise NotImplementedError def compose(f, g): """Computes functional composition of ``f`` and ``g``. """ F, G = f.unify_DMP(g) return F._compose(G) def _compose(f, g): raise NotImplementedError def decompose(f): """Computes functional decomposition of ``f``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._decompose() def _decompose(f): raise NotImplementedError def shift(f, a): """Efficiently compute Taylor shift ``f(x + a)``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._shift(f.dom.convert(a)) def shift_list(f, a): """Efficiently compute Taylor shift ``f(X + A)``. """ a = [f.dom.convert(ai) for ai in a] return f._shift_list(a) def _shift(f, a): raise NotImplementedError def transform(f, p, q): """Evaluate functional transformation ``q**n * f(p/q)``.""" if f.lev: raise ValueError('univariate polynomial expected') P, Q = p.unify_DMP(q) F, P = f.unify_DMP(P) F, Q = F.unify_DMP(Q) return F._transform(P, Q) def _transform(f, p, q): raise NotImplementedError def sturm(f): """Computes the Sturm sequence of ``f``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._sturm() def _sturm(f): raise NotImplementedError def cauchy_upper_bound(f): """Computes the Cauchy upper bound on the roots of ``f``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._cauchy_upper_bound() def _cauchy_upper_bound(f): raise NotImplementedError def cauchy_lower_bound(f): """Computes the Cauchy lower bound on the nonzero roots of ``f``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._cauchy_lower_bound() def _cauchy_lower_bound(f): raise NotImplementedError def mignotte_sep_bound_squared(f): """Computes the squared Mignotte bound on root separations of ``f``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._mignotte_sep_bound_squared() def _mignotte_sep_bound_squared(f): raise NotImplementedError def gff_list(f): """Computes greatest factorial factorization of ``f``. """ if f.lev: raise ValueError('univariate polynomial expected') return f._gff_list() def _gff_list(f): raise NotImplementedError def norm(f): """Computes ``Norm(f)``.""" raise NotImplementedError def sqf_norm(f): """Computes square-free norm of ``f``. """ raise NotImplementedError def sqf_part(f): """Computes square-free part of ``f``. """ raise NotImplementedError def sqf_list(f, all=False): """Returns a list of square-free factors of ``f``. """ raise NotImplementedError def sqf_list_include(f, all=False): """Returns a list of square-free factors of ``f``. """ raise NotImplementedError def factor_list(f): """Returns a list of irreducible factors of ``f``. """ raise NotImplementedError def factor_list_include(f): """Returns a list of irreducible factors of ``f``. """ raise NotImplementedError def intervals(f, all=False, eps=None, inf=None, sup=None, fast=False, sqf=False): """Compute isolating intervals for roots of ``f``. """ if f.lev: raise PolynomialError("Cannot isolate roots of a multivariate polynomial") if all and sqf: return f._isolate_all_roots_sqf(eps=eps, inf=inf, sup=sup, fast=fast) elif all and not sqf: return f._isolate_all_roots(eps=eps, inf=inf, sup=sup, fast=fast) elif not all and sqf: return f._isolate_real_roots_sqf(eps=eps, inf=inf, sup=sup, fast=fast) else: return f._isolate_real_roots(eps=eps, inf=inf, sup=sup, fast=fast) def _isolate_all_roots(f, eps, inf, sup, fast): raise NotImplementedError def _isolate_all_roots_sqf(f, eps, inf, sup, fast): raise NotImplementedError def _isolate_real_roots(f, eps, inf, sup, fast): raise NotImplementedError def _isolate_real_roots_sqf(f, eps, inf, sup, fast): raise NotImplementedError def refine_root(f, s, t, eps=None, steps=None, fast=False): """ Refine an isolating interval to the given precision. ``eps`` should be a rational number. """ if f.lev: raise PolynomialError( "Cannot refine a root of a multivariate polynomial") return f._refine_real_root(s, t, eps=eps, steps=steps, fast=fast) def _refine_real_root(f, s, t, eps, steps, fast): raise NotImplementedError def count_real_roots(f, inf=None, sup=None): """Return the number of real roots of ``f`` in ``[inf, sup]``. """ raise NotImplementedError def count_complex_roots(f, inf=None, sup=None): """Return the number of complex roots of ``f`` in ``[inf, sup]``. """ raise NotImplementedError @property def is_zero(f): """Returns ``True`` if ``f`` is a zero polynomial. """ raise NotImplementedError @property def is_one(f): """Returns ``True`` if ``f`` is a unit polynomial. """ raise NotImplementedError @property def is_ground(f): """Returns ``True`` if ``f`` is an element of the ground domain. """ raise NotImplementedError @property def is_sqf(f): """Returns ``True`` if ``f`` is a square-free polynomial. """ raise NotImplementedError @property def is_monic(f): """Returns ``True`` if the leading coefficient of ``f`` is one. """ raise NotImplementedError @property def is_primitive(f): """Returns ``True`` if the GCD of the coefficients of ``f`` is one. """ raise NotImplementedError @property def is_linear(f): """Returns ``True`` if ``f`` is linear in all its variables. """ raise NotImplementedError @property def is_quadratic(f): """Returns ``True`` if ``f`` is quadratic in all its variables. """ raise NotImplementedError @property def is_monomial(f): """Returns ``True`` if ``f`` is zero or has only one term. """ raise NotImplementedError @property def is_homogeneous(f): """Returns ``True`` if ``f`` is a homogeneous polynomial. """ raise NotImplementedError @property def is_irreducible(f): """Returns ``True`` if ``f`` has no factors over its domain. """ raise NotImplementedError @property def is_cyclotomic(f): """Returns ``True`` if ``f`` is a cyclotomic polynomial. """ raise NotImplementedError def __abs__(f): return f.abs() def __neg__(f): return f.neg() def __add__(f, g): if isinstance(g, DMP): return f.add(g) else: try: return f.add_ground(g) except CoercionFailed: return NotImplemented def __radd__(f, g): return f.__add__(g) def __sub__(f, g): if isinstance(g, DMP): return f.sub(g) else: try: return f.sub_ground(g) except CoercionFailed: return NotImplemented def __rsub__(f, g): return (-f).__add__(g) def __mul__(f, g): if isinstance(g, DMP): return f.mul(g) else: try: return f.mul_ground(g) except CoercionFailed: return NotImplemented def __rmul__(f, g): return f.__mul__(g) def __truediv__(f, g): if isinstance(g, DMP): return f.exquo(g) else: try: return f.mul_ground(g) except CoercionFailed: return NotImplemented def __rtruediv__(f, g): if isinstance(g, DMP): return g.exquo(f) else: try: return f._one().mul_ground(g).exquo(f) except CoercionFailed: return NotImplemented def __pow__(f, n): return f.pow(n) def __divmod__(f, g): return f.div(g) def __mod__(f, g): return f.rem(g) def __floordiv__(f, g): if isinstance(g, DMP): return f.quo(g) else: try: return f.quo_ground(g) except TypeError: return NotImplemented def __eq__(f, g): if f is g: return True if not isinstance(g, DMP): return NotImplemented try: F, G = f.unify_DMP(g) except UnificationFailed: return False else: return F._strict_eq(G) def _strict_eq(f, g): raise NotImplementedError def eq(f, g, strict=False): if not strict: return f == g else: return f._strict_eq(g) def ne(f, g, strict=False): return not f.eq(g, strict=strict) def __lt__(f, g): F, G = f.unify_DMP(g) return F.to_list() < G.to_list() def __le__(f, g): F, G = f.unify_DMP(g) return F.to_list() <= G.to_list() def __gt__(f, g): F, G = f.unify_DMP(g) return F.to_list() > G.to_list() def __ge__(f, g): F, G = f.unify_DMP(g) return F.to_list() >= G.to_list() def __bool__(f): return not f.is_zero class DMP_Python(DMP): """Dense Multivariate Polynomials over `K`. """ __slots__ = ('_rep', 'dom', 'lev') @classmethod def _new(cls, rep, dom, lev): obj = object.__new__(cls) obj._rep = rep obj.lev = lev obj.dom = dom return obj def _strict_eq(f, g): if type(f) != type(g): return False return f.lev == g.lev and f.dom == g.dom and f._rep == g._rep def per(f, rep): """Create a DMP out of the given representation. """ return f._new(rep, f.dom, f.lev) def ground_new(f, coeff): """Construct a new ground instance of ``f``. """ return f._new(dmp_ground(coeff, f.lev), f.dom, f.lev) def _one(f): return f.one(f.lev, f.dom) def unify(f, g): """Unify representations of two multivariate polynomials. """ # XXX: This function is not really used any more since there is # unify_DMP now. if not isinstance(g, DMP) or f.lev != g.lev: raise UnificationFailed("Cannot unify %s with %s" % (f, g)) if f.dom == g.dom: return f.lev, f.dom, f.per, f._rep, g._rep else: lev, dom = f.lev, f.dom.unify(g.dom) F = dmp_convert(f._rep, lev, f.dom, dom) G = dmp_convert(g._rep, lev, g.dom, dom) def per(rep): return f._new(rep, dom, lev) return lev, dom, per, F, G def to_DUP_Flint(f): """Convert ``f`` to a Flint representation. """ return DUP_Flint._new(f._rep, f.dom, f.lev) def to_list(f): """Convert ``f`` to a list representation with native coefficients. """ return list(f._rep) def to_tuple(f): """Convert ``f`` to a tuple representation with native coefficients. """ return dmp_to_tuple(f._rep, f.lev) def _convert(f, dom): """Convert the ground domain of ``f``. """ return f._new(dmp_convert(f._rep, f.lev, f.dom, dom), dom, f.lev) def _slice(f, m, n): """Take a continuous subsequence of terms of ``f``. """ rep = dup_slice(f._rep, m, n, f.dom) return f._new(rep, f.dom, f.lev) def _slice_lev(f, m, n, j): """Take a continuous subsequence of terms of ``f``. """ rep = dmp_slice_in(f._rep, m, n, j, f.lev, f.dom) return f._new(rep, f.dom, f.lev) def _terms(f, order=None): """Returns all non-zero terms from ``f`` in lex order. """ return dmp_list_terms(f._rep, f.lev, f.dom, order=order) def _lift(f): """Convert algebraic coefficients to rationals. """ r = dmp_lift(f._rep, f.lev, f.dom) return f._new(r, f.dom.dom, f.lev) def deflate(f): """Reduce degree of `f` by mapping `x_i^m` to `y_i`. """ J, F = dmp_deflate(f._rep, f.lev, f.dom) return J, f.per(F) def inject(f, front=False): """Inject ground domain generators into ``f``. """ F, lev = dmp_inject(f._rep, f.lev, f.dom, front=front) # XXX: domain and level changed here return f._new(F, f.dom.dom, lev) def eject(f, dom, front=False): """Eject selected generators into the ground domain. """ F = dmp_eject(f._rep, f.lev, dom, front=front) # XXX: domain and level changed here return f._new(F, dom, f.lev - len(dom.symbols)) def _exclude(f): """Remove useless generators from ``f``. """ J, F, u = dmp_exclude(f._rep, f.lev, f.dom) # XXX: level changed here return J, f._new(F, f.dom, u) def _permute(f, P): """Returns a polynomial in `K[x_{P(1)}, ..., x_{P(n)}]`. """ return f.per(dmp_permute(f._rep, P, f.lev, f.dom)) def terms_gcd(f): """Remove GCD of terms from the polynomial ``f``. """ J, F = dmp_terms_gcd(f._rep, f.lev, f.dom) return J, f.per(F) def _add_ground(f, c): """Add an element of the ground domain to ``f``. """ return f.per(dmp_add_ground(f._rep, c, f.lev, f.dom)) def _sub_ground(f, c): """Subtract an element of the ground domain from ``f``. """ return f.per(dmp_sub_ground(f._rep, c, f.lev, f.dom)) def _mul_ground(f, c): """Multiply ``f`` by a an element of the ground domain. """ return f.per(dmp_mul_ground(f._rep, c, f.lev, f.dom)) def _quo_ground(f, c): """Quotient of ``f`` by a an element of the ground domain. """ return f.per(dmp_quo_ground(f._rep, c, f.lev, f.dom)) def _exquo_ground(f, c): """Exact quotient of ``f`` by a an element of the ground domain. """ return f.per(dmp_exquo_ground(f._rep, c, f.lev, f.dom)) def abs(f): """Make all coefficients in ``f`` positive. """ return f.per(dmp_abs(f._rep, f.lev, f.dom)) def neg(f): """Negate all coefficients in ``f``. """ return f.per(dmp_neg(f._rep, f.lev, f.dom)) def _add(f, g): """Add two multivariate polynomials ``f`` and ``g``. """ return f.per(dmp_add(f._rep, g._rep, f.lev, f.dom)) def _sub(f, g): """Subtract two multivariate polynomials ``f`` and ``g``. """ return f.per(dmp_sub(f._rep, g._rep, f.lev, f.dom)) def _mul(f, g): """Multiply two multivariate polynomials ``f`` and ``g``. """ return f.per(dmp_mul(f._rep, g._rep, f.lev, f.dom)) def sqr(f): """Square a multivariate polynomial ``f``. """ return f.per(dmp_sqr(f._rep, f.lev, f.dom)) def _pow(f, n): """Raise ``f`` to a non-negative power ``n``. """ return f.per(dmp_pow(f._rep, n, f.lev, f.dom)) def _pdiv(f, g): """Polynomial pseudo-division of ``f`` and ``g``. """ q, r = dmp_pdiv(f._rep, g._rep, f.lev, f.dom) return f.per(q), f.per(r) def _prem(f, g): """Polynomial pseudo-remainder of ``f`` and ``g``. """ return f.per(dmp_prem(f._rep, g._rep, f.lev, f.dom)) def _pquo(f, g): """Polynomial pseudo-quotient of ``f`` and ``g``. """ return f.per(dmp_pquo(f._rep, g._rep, f.lev, f.dom)) def _pexquo(f, g): """Polynomial exact pseudo-quotient of ``f`` and ``g``. """ return f.per(dmp_pexquo(f._rep, g._rep, f.lev, f.dom)) def _div(f, g): """Polynomial division with remainder of ``f`` and ``g``. """ q, r = dmp_div(f._rep, g._rep, f.lev, f.dom) return f.per(q), f.per(r) def _rem(f, g): """Computes polynomial remainder of ``f`` and ``g``. """ return f.per(dmp_rem(f._rep, g._rep, f.lev, f.dom)) def _quo(f, g): """Computes polynomial quotient of ``f`` and ``g``. """ return f.per(dmp_quo(f._rep, g._rep, f.lev, f.dom)) def _exquo(f, g): """Computes polynomial exact quotient of ``f`` and ``g``. """ return f.per(dmp_exquo(f._rep, g._rep, f.lev, f.dom)) def _degree(f, j=0): """Returns the leading degree of ``f`` in ``x_j``. """ return dmp_degree_in(f._rep, j, f.lev) def degree_list(f): """Returns a list of degrees of ``f``. """ return dmp_degree_list(f._rep, f.lev) def total_degree(f): """Returns the total degree of ``f``. """ return max(sum(m) for m in f.monoms()) def LC(f): """Returns the leading coefficient of ``f``. """ return dmp_ground_LC(f._rep, f.lev, f.dom) def TC(f): """Returns the trailing coefficient of ``f``. """ return dmp_ground_TC(f._rep, f.lev, f.dom) def _nth(f, N): """Returns the ``n``-th coefficient of ``f``. """ return dmp_ground_nth(f._rep, N, f.lev, f.dom) def max_norm(f): """Returns maximum norm of ``f``. """ return dmp_max_norm(f._rep, f.lev, f.dom) def l1_norm(f): """Returns l1 norm of ``f``. """ return dmp_l1_norm(f._rep, f.lev, f.dom) def l2_norm_squared(f): """Return squared l2 norm of ``f``. """ return dmp_l2_norm_squared(f._rep, f.lev, f.dom) def clear_denoms(f): """Clear denominators, but keep the ground domain. """ coeff, F = dmp_clear_denoms(f._rep, f.lev, f.dom) return coeff, f.per(F) def _integrate(f, m=1, j=0): """Computes the ``m``-th order indefinite integral of ``f`` in ``x_j``. """ return f.per(dmp_integrate_in(f._rep, m, j, f.lev, f.dom)) def _diff(f, m=1, j=0): """Computes the ``m``-th order derivative of ``f`` in ``x_j``. """ return f.per(dmp_diff_in(f._rep, m, j, f.lev, f.dom)) def _eval(f, a): return dmp_eval_in(f._rep, f.dom.convert(a), 0, f.lev, f.dom) def _eval_lev(f, a, j): rep = dmp_eval_in(f._rep, f.dom.convert(a), j, f.lev, f.dom) return f.new(rep, f.dom, f.lev - 1) def _half_gcdex(f, g): """Half extended Euclidean algorithm, if univariate. """ s, h = dup_half_gcdex(f._rep, g._rep, f.dom) return f.per(s), f.per(h) def _gcdex(f, g): """Extended Euclidean algorithm, if univariate. """ s, t, h = dup_gcdex(f._rep, g._rep, f.dom) return f.per(s), f.per(t), f.per(h) def _invert(f, g): """Invert ``f`` modulo ``g``, if possible. """ s = dup_invert(f._rep, g._rep, f.dom) return f.per(s) def _revert(f, n): """Compute ``f**(-1)`` mod ``x**n``. """ return f.per(dup_revert(f._rep, n, f.dom)) def _subresultants(f, g): """Computes subresultant PRS sequence of ``f`` and ``g``. """ R = dmp_subresultants(f._rep, g._rep, f.lev, f.dom) return list(map(f.per, R)) def _resultant_includePRS(f, g): """Computes resultant of ``f`` and ``g`` via PRS. """ res, R = dmp_resultant(f._rep, g._rep, f.lev, f.dom, includePRS=True) if f.lev: res = f.new(res, f.dom, f.lev - 1) return res, list(map(f.per, R)) def _resultant(f, g): res = dmp_resultant(f._rep, g._rep, f.lev, f.dom) if f.lev: res = f.new(res, f.dom, f.lev - 1) return res def discriminant(f): """Computes discriminant of ``f``. """ res = dmp_discriminant(f._rep, f.lev, f.dom) if f.lev: res = f.new(res, f.dom, f.lev - 1) return res def _cofactors(f, g): """Returns GCD of ``f`` and ``g`` and their cofactors. """ h, cff, cfg = dmp_inner_gcd(f._rep, g._rep, f.lev, f.dom) return f.per(h), f.per(cff), f.per(cfg) def _gcd(f, g): """Returns polynomial GCD of ``f`` and ``g``. """ return f.per(dmp_gcd(f._rep, g._rep, f.lev, f.dom)) def _lcm(f, g): """Returns polynomial LCM of ``f`` and ``g``. """ return f.per(dmp_lcm(f._rep, g._rep, f.lev, f.dom)) def _cancel(f, g): """Cancel common factors in a rational function ``f/g``. """ cF, cG, F, G = dmp_cancel(f._rep, g._rep, f.lev, f.dom, include=False) return cF, cG, f.per(F), f.per(G) def _cancel_include(f, g): """Cancel common factors in a rational function ``f/g``. """ F, G = dmp_cancel(f._rep, g._rep, f.lev, f.dom, include=True) return f.per(F), f.per(G) def _trunc(f, p): """Reduce ``f`` modulo a constant ``p``. """ return f.per(dmp_ground_trunc(f._rep, p, f.lev, f.dom)) def monic(f): """Divides all coefficients by ``LC(f)``. """ return f.per(dmp_ground_monic(f._rep, f.lev, f.dom)) def content(f): """Returns GCD of polynomial coefficients. """ return dmp_ground_content(f._rep, f.lev, f.dom) def primitive(f): """Returns content and a primitive form of ``f``. """ cont, F = dmp_ground_primitive(f._rep, f.lev, f.dom) return cont, f.per(F) def _compose(f, g): """Computes functional composition of ``f`` and ``g``. """ return f.per(dmp_compose(f._rep, g._rep, f.lev, f.dom)) def _decompose(f): """Computes functional decomposition of ``f``. """ return list(map(f.per, dup_decompose(f._rep, f.dom))) def _shift(f, a): """Efficiently compute Taylor shift ``f(x + a)``. """ return f.per(dup_shift(f._rep, a, f.dom)) def _shift_list(f, a): """Efficiently compute Taylor shift ``f(X + A)``. """ return f.per(dmp_shift(f._rep, a, f.lev, f.dom)) def _transform(f, p, q): """Evaluate functional transformation ``q**n * f(p/q)``.""" return f.per(dup_transform(f._rep, p._rep, q._rep, f.dom)) def _sturm(f): """Computes the Sturm sequence of ``f``. """ return list(map(f.per, dup_sturm(f._rep, f.dom))) def _cauchy_upper_bound(f): """Computes the Cauchy upper bound on the roots of ``f``. """ return dup_cauchy_upper_bound(f._rep, f.dom) def _cauchy_lower_bound(f): """Computes the Cauchy lower bound on the nonzero roots of ``f``. """ return dup_cauchy_lower_bound(f._rep, f.dom) def _mignotte_sep_bound_squared(f): """Computes the squared Mignotte bound on root separations of ``f``. """ return dup_mignotte_sep_bound_squared(f._rep, f.dom) def _gff_list(f): """Computes greatest factorial factorization of ``f``. """ return [ (f.per(g), k) for g, k in dup_gff_list(f._rep, f.dom) ] def norm(f): """Computes ``Norm(f)``.""" r = dmp_norm(f._rep, f.lev, f.dom) return f.new(r, f.dom.dom, f.lev) def sqf_norm(f): """Computes square-free norm of ``f``. """ s, g, r = dmp_sqf_norm(f._rep, f.lev, f.dom) return s, f.per(g), f.new(r, f.dom.dom, f.lev) def sqf_part(f): """Computes square-free part of ``f``. """ return f.per(dmp_sqf_part(f._rep, f.lev, f.dom)) def sqf_list(f, all=False): """Returns a list of square-free factors of ``f``. """ coeff, factors = dmp_sqf_list(f._rep, f.lev, f.dom, all) return coeff, [ (f.per(g), k) for g, k in factors ] def sqf_list_include(f, all=False): """Returns a list of square-free factors of ``f``. """ factors = dmp_sqf_list_include(f._rep, f.lev, f.dom, all) return [ (f.per(g), k) for g, k in factors ] def factor_list(f): """Returns a list of irreducible factors of ``f``. """ coeff, factors = dmp_factor_list(f._rep, f.lev, f.dom) return coeff, [ (f.per(g), k) for g, k in factors ] def factor_list_include(f): """Returns a list of irreducible factors of ``f``. """ factors = dmp_factor_list_include(f._rep, f.lev, f.dom) return [ (f.per(g), k) for g, k in factors ] def _isolate_real_roots(f, eps, inf, sup, fast): return dup_isolate_real_roots(f._rep, f.dom, eps=eps, inf=inf, sup=sup, fast=fast) def _isolate_real_roots_sqf(f, eps, inf, sup, fast): return dup_isolate_real_roots_sqf(f._rep, f.dom, eps=eps, inf=inf, sup=sup, fast=fast) def _isolate_all_roots(f, eps, inf, sup, fast): return dup_isolate_all_roots(f._rep, f.dom, eps=eps, inf=inf, sup=sup, fast=fast) def _isolate_all_roots_sqf(f, eps, inf, sup, fast): return dup_isolate_all_roots_sqf(f._rep, f.dom, eps=eps, inf=inf, sup=sup, fast=fast) def _refine_real_root(f, s, t, eps, steps, fast): return dup_refine_real_root(f._rep, s, t, f.dom, eps=eps, steps=steps, fast=fast) def count_real_roots(f, inf=None, sup=None): """Return the number of real roots of ``f`` in ``[inf, sup]``. """ return dup_count_real_roots(f._rep, f.dom, inf=inf, sup=sup) def count_complex_roots(f, inf=None, sup=None): """Return the number of complex roots of ``f`` in ``[inf, sup]``. """ return dup_count_complex_roots(f._rep, f.dom, inf=inf, sup=sup) @property def is_zero(f): """Returns ``True`` if ``f`` is a zero polynomial. """ return dmp_zero_p(f._rep, f.lev) @property def is_one(f): """Returns ``True`` if ``f`` is a unit polynomial. """ return dmp_one_p(f._rep, f.lev, f.dom) @property def is_ground(f): """Returns ``True`` if ``f`` is an element of the ground domain. """ return dmp_ground_p(f._rep, None, f.lev) @property def is_sqf(f): """Returns ``True`` if ``f`` is a square-free polynomial. """ return dmp_sqf_p(f._rep, f.lev, f.dom) @property def is_monic(f): """Returns ``True`` if the leading coefficient of ``f`` is one. """ return f.dom.is_one(dmp_ground_LC(f._rep, f.lev, f.dom)) @property def is_primitive(f): """Returns ``True`` if the GCD of the coefficients of ``f`` is one. """ return f.dom.is_one(dmp_ground_content(f._rep, f.lev, f.dom)) @property def is_linear(f): """Returns ``True`` if ``f`` is linear in all its variables. """ return all(sum(monom) <= 1 for monom in dmp_to_dict(f._rep, f.lev, f.dom).keys()) @property def is_quadratic(f): """Returns ``True`` if ``f`` is quadratic in all its variables. """ return all(sum(monom) <= 2 for monom in dmp_to_dict(f._rep, f.lev, f.dom).keys()) @property def is_monomial(f): """Returns ``True`` if ``f`` is zero or has only one term. """ return len(f.to_dict()) <= 1 @property def is_homogeneous(f): """Returns ``True`` if ``f`` is a homogeneous polynomial. """ return f.homogeneous_order() is not None @property def is_irreducible(f): """Returns ``True`` if ``f`` has no factors over its domain. """ return dmp_irreducible_p(f._rep, f.lev, f.dom) @property def is_cyclotomic(f): """Returns ``True`` if ``f`` is a cyclotomic polynomial. """ if not f.lev: return dup_cyclotomic_p(f._rep, f.dom) else: return False class DUP_Flint(DMP): """Dense Multivariate Polynomials over `K`. """ lev = 0 __slots__ = ('_rep', 'dom', '_cls') def __reduce__(self): return self.__class__, (self.to_list(), self.dom, self.lev) @classmethod def _new(cls, rep, dom, lev): rep = cls._flint_poly(rep[::-1], dom, lev) return cls.from_rep(rep, dom) def to_list(f): """Convert ``f`` to a list representation with native coefficients. """ return f._rep.coeffs()[::-1] @classmethod def _flint_poly(cls, rep, dom, lev): assert dom in _flint_domains assert lev == 0 flint_cls = cls._get_flint_poly_cls(dom) return flint_cls(rep) @classmethod def _get_flint_poly_cls(cls, dom): if dom.is_ZZ: return flint.fmpz_poly elif dom.is_QQ: return flint.fmpq_poly else: raise RuntimeError("Domain %s is not supported with flint" % dom) @classmethod def from_rep(cls, rep, dom): """Create a DMP from the given representation. """ if dom.is_ZZ: assert isinstance(rep, flint.fmpz_poly) _cls = flint.fmpz_poly elif dom.is_QQ: assert isinstance(rep, flint.fmpq_poly) _cls = flint.fmpq_poly else: raise RuntimeError("Domain %s is not supported with flint" % dom) obj = object.__new__(cls) obj.dom = dom obj._rep = rep obj._cls = _cls return obj def _strict_eq(f, g): if type(f) != type(g): return False return f.dom == g.dom and f._rep == g._rep def ground_new(f, coeff): """Construct a new ground instance of ``f``. """ return f.from_rep(f._cls([coeff]), f.dom) def _one(f): return f.ground_new(f.dom.one) def unify(f, g): """Unify representations of two polynomials. """ raise RuntimeError def to_DMP_Python(f): """Convert ``f`` to a Python native representation. """ return DMP_Python._new(f.to_list(), f.dom, f.lev) def to_tuple(f): """Convert ``f`` to a tuple representation with native coefficients. """ return tuple(f.to_list()) def _convert(f, dom): """Convert the ground domain of ``f``. """ if dom == QQ and f.dom == ZZ: return f.from_rep(flint.fmpq_poly(f._rep), dom) elif dom == ZZ and f.dom == QQ: # XXX: python-flint should provide a faster way to do this. return f.to_DMP_Python()._convert(dom).to_DUP_Flint() else: raise RuntimeError(f"DUP_Flint: Cannot convert {f.dom} to {dom}") def _slice(f, m, n): """Take a continuous subsequence of terms of ``f``. """ coeffs = f._rep.coeffs()[m:n] return f.from_rep(f._cls(coeffs), f.dom) def _slice_lev(f, m, n, j): """Take a continuous subsequence of terms of ``f``. """ # Only makes sense for multivariate polynomials raise NotImplementedError def _terms(f, order=None): """Returns all non-zero terms from ``f`` in lex order. """ if order is None or order.alias == 'lex': terms = [ ((n,), c) for n, c in enumerate(f._rep.coeffs()) if c ] return terms[::-1] else: # XXX: InverseOrder (ilex) comes here. We could handle that case # efficiently by reversing the coefficients but it is not clear # how to test if the order is InverseOrder. # # Otherwise why would the order ever be different for univariate # polynomials? return f.to_DMP_Python()._terms(order=order) def _lift(f): """Convert algebraic coefficients to rationals. """ # This is for algebraic number fields which DUP_Flint does not support raise NotImplementedError def deflate(f): """Reduce degree of `f` by mapping `x_i^m` to `y_i`. """ # XXX: Check because otherwise this segfaults with python-flint: # # >>> flint.fmpz_poly([]).deflation() # Exception (fmpz_poly_deflate). Division by zero. # Aborted (core dumped # if f.is_zero: return (1,), f g, n = f._rep.deflation() return (n,), f.from_rep(g, f.dom) def inject(f, front=False): """Inject ground domain generators into ``f``. """ # Ground domain would need to be a poly ring raise NotImplementedError def eject(f, dom, front=False): """Eject selected generators into the ground domain. """ # Only makes sense for multivariate polynomials raise NotImplementedError def _exclude(f): """Remove useless generators from ``f``. """ # Only makes sense for multivariate polynomials raise NotImplementedError def _permute(f, P): """Returns a polynomial in `K[x_{P(1)}, ..., x_{P(n)}]`. """ # Only makes sense for multivariate polynomials raise NotImplementedError def terms_gcd(f): """Remove GCD of terms from the polynomial ``f``. """ # XXX: python-flint should have primitive, content, etc methods. J, F = f.to_DMP_Python().terms_gcd() return J, F.to_DUP_Flint() def _add_ground(f, c): """Add an element of the ground domain to ``f``. """ return f.from_rep(f._rep + c, f.dom) def _sub_ground(f, c): """Subtract an element of the ground domain from ``f``. """ return f.from_rep(f._rep - c, f.dom) def _mul_ground(f, c): """Multiply ``f`` by a an element of the ground domain. """ return f.from_rep(f._rep * c, f.dom) def _quo_ground(f, c): """Quotient of ``f`` by a an element of the ground domain. """ return f.from_rep(f._rep // c, f.dom) def _exquo_ground(f, c): """Exact quotient of ``f`` by a an element of the ground domain. """ q, r = divmod(f._rep, c) if r: raise ExactQuotientFailed(f, c) return f.from_rep(q, f.dom) def abs(f): """Make all coefficients in ``f`` positive. """ return f.to_DMP_Python().abs().to_DUP_Flint() def neg(f): """Negate all coefficients in ``f``. """ return f.from_rep(-f._rep, f.dom) def _add(f, g): """Add two multivariate polynomials ``f`` and ``g``. """ return f.from_rep(f._rep + g._rep, f.dom) def _sub(f, g): """Subtract two multivariate polynomials ``f`` and ``g``. """ return f.from_rep(f._rep - g._rep, f.dom) def _mul(f, g): """Multiply two multivariate polynomials ``f`` and ``g``. """ return f.from_rep(f._rep * g._rep, f.dom) def sqr(f): """Square a multivariate polynomial ``f``. """ return f.from_rep(f._rep ** 2, f.dom) def _pow(f, n): """Raise ``f`` to a non-negative power ``n``. """ return f.from_rep(f._rep ** n, f.dom) def _pdiv(f, g): """Polynomial pseudo-division of ``f`` and ``g``. """ d = f.degree() - g.degree() + 1 q, r = divmod(g.LC()**d * f._rep, g._rep) return f.from_rep(q, f.dom), f.from_rep(r, f.dom) def _prem(f, g): """Polynomial pseudo-remainder of ``f`` and ``g``. """ d = f.degree() - g.degree() + 1 q = (g.LC()**d * f._rep) % g._rep return f.from_rep(q, f.dom) def _pquo(f, g): """Polynomial pseudo-quotient of ``f`` and ``g``. """ d = f.degree() - g.degree() + 1 r = (g.LC()**d * f._rep) // g._rep return f.from_rep(r, f.dom) def _pexquo(f, g): """Polynomial exact pseudo-quotient of ``f`` and ``g``. """ d = f.degree() - g.degree() + 1 q, r = divmod(g.LC()**d * f._rep, g._rep) if r: raise ExactQuotientFailed(f, g) return f.from_rep(q, f.dom) def _div(f, g): """Polynomial division with remainder of ``f`` and ``g``. """ if f.dom.is_Field: q, r = divmod(f._rep, g._rep) return f.from_rep(q, f.dom), f.from_rep(r, f.dom) else: # XXX: python-flint defines division in ZZ[x] differently q, r = f.to_DMP_Python()._div(g.to_DMP_Python()) return q.to_DUP_Flint(), r.to_DUP_Flint() def _rem(f, g): """Computes polynomial remainder of ``f`` and ``g``. """ return f.from_rep(f._rep % g._rep, f.dom) def _quo(f, g): """Computes polynomial quotient of ``f`` and ``g``. """ return f.from_rep(f._rep // g._rep, f.dom) def _exquo(f, g): """Computes polynomial exact quotient of ``f`` and ``g``. """ q, r = f._div(g) if r: raise ExactQuotientFailed(f, g) return q def _degree(f, j=0): """Returns the leading degree of ``f`` in ``x_j``. """ d = f._rep.degree() if d == -1: d = ninf return d def degree_list(f): """Returns a list of degrees of ``f``. """ return ( f._degree() ,) def total_degree(f): """Returns the total degree of ``f``. """ return f._degree() def LC(f): """Returns the leading coefficient of ``f``. """ return f._rep[f._rep.degree()] def TC(f): """Returns the trailing coefficient of ``f``. """ return f._rep[0] def _nth(f, N): """Returns the ``n``-th coefficient of ``f``. """ [n] = N return f._rep[n] def max_norm(f): """Returns maximum norm of ``f``. """ return f.to_DMP_Python().max_norm() def l1_norm(f): """Returns l1 norm of ``f``. """ return f.to_DMP_Python().l1_norm() def l2_norm_squared(f): """Return squared l2 norm of ``f``. """ return f.to_DMP_Python().l2_norm_squared() def clear_denoms(f): """Clear denominators, but keep the ground domain. """ denom = f._rep.denom() numer = f.from_rep(f._cls(f._rep.numer()), f.dom) return denom, numer def _integrate(f, m=1, j=0): """Computes the ``m``-th order indefinite integral of ``f`` in ``x_j``. """ assert j == 0 if f.dom.is_QQ: rep = f._rep for i in range(m): rep = rep.integral() return f.from_rep(rep, f.dom) else: return f.to_DMP_Python()._integrate(m=m, j=j).to_DUP_Flint() def _diff(f, m=1, j=0): """Computes the ``m``-th order derivative of ``f``. """ assert j == 0 rep = f._rep for i in range(m): rep = rep.derivative() return f.from_rep(rep, f.dom) def _eval(f, a): return f.to_DMP_Python()._eval(a) def _eval_lev(f, a, j): # Only makes sense for multivariate polynomials raise NotImplementedError def _half_gcdex(f, g): """Half extended Euclidean algorithm. """ s, h = f.to_DMP_Python()._half_gcdex(g.to_DMP_Python()) return s.to_DUP_Flint(), h.to_DUP_Flint() def _gcdex(f, g): """Extended Euclidean algorithm. """ h, s, t = f._rep.xgcd(g._rep) return f.from_rep(s, f.dom), f.from_rep(t, f.dom), f.from_rep(h, f.dom) def _invert(f, g): """Invert ``f`` modulo ``g``, if possible. """ if f.dom.is_QQ: gcd, F_inv, _ = f._rep.xgcd(g._rep) if gcd != 1: raise NotInvertible("zero divisor") return f.from_rep(F_inv, f.dom) else: return f.to_DMP_Python()._invert(g.to_DMP_Python()).to_DUP_Flint() def _revert(f, n): """Compute ``f**(-1)`` mod ``x**n``. """ return f.to_DMP_Python()._revert(n).to_DUP_Flint() def _subresultants(f, g): """Computes subresultant PRS sequence of ``f`` and ``g``. """ R = f.to_DMP_Python()._subresultants(g.to_DMP_Python()) return [ g.to_DUP_Flint() for g in R ] def _resultant_includePRS(f, g): """Computes resultant of ``f`` and ``g`` via PRS. """ res, R = f.to_DMP_Python()._resultant_includePRS(g.to_DMP_Python()) return res, [ g.to_DUP_Flint() for g in R ] def _resultant(f, g): """Computes resultant of ``f`` and ``g``. """ return f.to_DMP_Python()._resultant(g.to_DMP_Python()) def discriminant(f): """Computes discriminant of ``f``. """ return f.to_DMP_Python().discriminant() def _cofactors(f, g): """Returns GCD of ``f`` and ``g`` and their cofactors. """ h = f.gcd(g) return h, f.exquo(h), g.exquo(h) def _gcd(f, g): """Returns polynomial GCD of ``f`` and ``g``. """ return f.from_rep(f._rep.gcd(g._rep), f.dom) def _lcm(f, g): """Returns polynomial LCM of ``f`` and ``g``. """ # XXX: python-flint should have a lcm method if not (f and g): return f.ground_new(f.dom.zero) l = f._mul(g)._exquo(f._gcd(g)) if l.dom.is_Field: l = l.monic() elif l.LC() < 0: l = l.neg() return l def _cancel(f, g): """Cancel common factors in a rational function ``f/g``. """ # Think carefully about how to handle denominators and coefficient # canonicalisation if more domains are permitted... assert f.dom == g.dom in (ZZ, QQ) if f.dom.is_QQ: cG, F = f.clear_denoms() cF, G = g.clear_denoms() else: cG, F = f.dom.one, f cF, G = g.dom.one, g cH = cF.gcd(cG) cF, cG = cF // cH, cG // cH H = F._gcd(G) F, G = F.exquo(H), G.exquo(H) f_neg = F.LC() < 0 g_neg = G.LC() < 0 if f_neg and g_neg: F, G = F.neg(), G.neg() elif f_neg: cF, F = -cF, F.neg() elif g_neg: cF, G = -cF, G.neg() return cF, cG, F, G def _cancel_include(f, g): """Cancel common factors in a rational function ``f/g``. """ cF, cG, F, G = f._cancel(g) return F._mul_ground(cF), G._mul_ground(cG) def _trunc(f, p): """Reduce ``f`` modulo a constant ``p``. """ return f.to_DMP_Python()._trunc(p).to_DUP_Flint() def monic(f): """Divides all coefficients by ``LC(f)``. """ return f._exquo_ground(f.LC()) def content(f): """Returns GCD of polynomial coefficients. """ # XXX: python-flint should have a content method return f.to_DMP_Python().content() def primitive(f): """Returns content and a primitive form of ``f``. """ cont = f.content() prim = f._exquo_ground(cont) return cont, prim def _compose(f, g): """Computes functional composition of ``f`` and ``g``. """ return f.from_rep(f._rep(g._rep), f.dom) def _decompose(f): """Computes functional decomposition of ``f``. """ return [ g.to_DUP_Flint() for g in f.to_DMP_Python()._decompose() ] def _shift(f, a): """Efficiently compute Taylor shift ``f(x + a)``. """ x_plus_a = f._cls([a, f.dom.one]) return f.from_rep(f._rep(x_plus_a), f.dom) def _transform(f, p, q): """Evaluate functional transformation ``q**n * f(p/q)``.""" F, P, Q = f.to_DMP_Python(), p.to_DMP_Python(), q.to_DMP_Python() return F.transform(P, Q).to_DUP_Flint() def _sturm(f): """Computes the Sturm sequence of ``f``. """ return [ g.to_DUP_Flint() for g in f.to_DMP_Python()._sturm() ] def _cauchy_upper_bound(f): """Computes the Cauchy upper bound on the roots of ``f``. """ return f.to_DMP_Python()._cauchy_upper_bound() def _cauchy_lower_bound(f): """Computes the Cauchy lower bound on the nonzero roots of ``f``. """ return f.to_DMP_Python()._cauchy_lower_bound() def _mignotte_sep_bound_squared(f): """Computes the squared Mignotte bound on root separations of ``f``. """ return f.to_DMP_Python()._mignotte_sep_bound_squared() def _gff_list(f): """Computes greatest factorial factorization of ``f``. """ F = f.to_DMP_Python() return [ (g.to_DUP_Flint(), k) for g, k in F.gff_list() ] def norm(f): """Computes ``Norm(f)``.""" # This is for algebraic number fields which DUP_Flint does not support raise NotImplementedError def sqf_norm(f): """Computes square-free norm of ``f``. """ # This is for algebraic number fields which DUP_Flint does not support raise NotImplementedError def sqf_part(f): """Computes square-free part of ``f``. """ return f._exquo(f._gcd(f._diff())) def sqf_list(f, all=False): """Returns a list of square-free factors of ``f``. """ coeff, factors = f.to_DMP_Python().sqf_list(all=all) return coeff, [ (g.to_DUP_Flint(), k) for g, k in factors ] def sqf_list_include(f, all=False): """Returns a list of square-free factors of ``f``. """ factors = f.to_DMP_Python().sqf_list_include(all=all) return [ (g.to_DUP_Flint(), k) for g, k in factors ] def factor_list(f): """Returns a list of irreducible factors of ``f``. """ if f.dom.is_ZZ: # python-flint matches polys here coeff, factors = f._rep.factor() factors = [ (f.from_rep(g, f.dom), k) for g, k in factors ] elif f.dom.is_QQ: # python-flint returns monic factors over QQ whereas polys returns # denominator free factors. coeff, factors = f._rep.factor() factors_monic = [ (f.from_rep(g, f.dom), k) for g, k in factors ] # Absorb the denominators into coeff factors = [] for g, k in factors_monic: d, g = g.clear_denoms() coeff /= d**k factors.append((g, k)) else: # Check carefully when adding more domains here... raise RuntimeError("Domain %s is not supported with flint" % f.dom) # We need to match the way that polys orders the factors factors = f._sort_factors(factors) return coeff, factors def factor_list_include(f): """Returns a list of irreducible factors of ``f``. """ # XXX: factor_list_include seems to be broken in general: # # >>> Poly(2*(x - 1)**3, x).factor_list_include() # [(Poly(2*x - 2, x, domain='ZZ'), 3)] # # Let's not try to implement it here. factors = f.to_DMP_Python().factor_list_include() return [ (g.to_DUP_Flint(), k) for g, k in factors ] def _sort_factors(f, factors): """Sort a list of factors to canonical order. """ # Convert the factors to lists and use _sort_factors from polys factors = [ (g.to_list(), k) for g, k in factors ] factors = _sort_factors(factors, multiple=True) to_dup_flint = lambda g: f.from_rep(f._cls(g[::-1]), f.dom) return [ (to_dup_flint(g), k) for g, k in factors ] def _isolate_real_roots(f, eps, inf, sup, fast): return f.to_DMP_Python()._isolate_real_roots(eps, inf, sup, fast) def _isolate_real_roots_sqf(f, eps, inf, sup, fast): return f.to_DMP_Python()._isolate_real_roots_sqf(eps, inf, sup, fast) def _isolate_all_roots(f, eps, inf, sup, fast): return f.to_DMP_Python()._isolate_all_roots(eps, inf, sup, fast) def _isolate_all_roots_sqf(f, eps, inf, sup, fast): return f.to_DMP_Python()._isolate_all_roots_sqf(eps, inf, sup, fast) def _refine_real_root(f, s, t, eps, steps, fast): return f.to_DMP_Python()._refine_real_root(s, t, eps, steps, fast) def count_real_roots(f, inf=None, sup=None): """Return the number of real roots of ``f`` in ``[inf, sup]``. """ return f.to_DMP_Python().count_real_roots(inf=inf, sup=sup) def count_complex_roots(f, inf=None, sup=None): """Return the number of complex roots of ``f`` in ``[inf, sup]``. """ return f.to_DMP_Python().count_complex_roots(inf=inf, sup=sup) @property def is_zero(f): """Returns ``True`` if ``f`` is a zero polynomial. """ return not f._rep @property def is_one(f): """Returns ``True`` if ``f`` is a unit polynomial. """ return f._rep == f.dom.one @property def is_ground(f): """Returns ``True`` if ``f`` is an element of the ground domain. """ return f._rep.degree() <= 0 @property def is_linear(f): """Returns ``True`` if ``f`` is linear in all its variables. """ return f._rep.degree() <= 1 @property def is_quadratic(f): """Returns ``True`` if ``f`` is quadratic in all its variables. """ return f._rep.degree() <= 2 @property def is_monomial(f): """Returns ``True`` if ``f`` is zero or has only one term. """ return f.to_DMP_Python().is_monomial @property def is_monic(f): """Returns ``True`` if the leading coefficient of ``f`` is one. """ return f.LC() == f.dom.one @property def is_primitive(f): """Returns ``True`` if the GCD of the coefficients of ``f`` is one. """ return f.to_DMP_Python().is_primitive @property def is_homogeneous(f): """Returns ``True`` if ``f`` is a homogeneous polynomial. """ return f.to_DMP_Python().is_homogeneous @property def is_sqf(f): """Returns ``True`` if ``f`` is a square-free polynomial. """ return f.to_DMP_Python().is_sqf @property def is_irreducible(f): """Returns ``True`` if ``f`` has no factors over its domain. """ return f.to_DMP_Python().is_irreducible @property def is_cyclotomic(f): """Returns ``True`` if ``f`` is a cyclotomic polynomial. """ if f.dom.is_ZZ: return bool(f._rep.is_cyclotomic()) else: return f.to_DMP_Python().is_cyclotomic def init_normal_DMF(num, den, lev, dom): return DMF(dmp_normal(num, lev, dom), dmp_normal(den, lev, dom), dom, lev) class DMF(PicklableWithSlots, CantSympify): """Dense Multivariate Fractions over `K`. """ __slots__ = ('num', 'den', 'lev', 'dom') def __init__(self, rep, dom, lev=None): num, den, lev = self._parse(rep, dom, lev) num, den = dmp_cancel(num, den, lev, dom) self.num = num self.den = den self.lev = lev self.dom = dom @classmethod def new(cls, rep, dom, lev=None): num, den, lev = cls._parse(rep, dom, lev) obj = object.__new__(cls) obj.num = num obj.den = den obj.lev = lev obj.dom = dom return obj def ground_new(self, rep): return self.new(rep, self.dom, self.lev) @classmethod def _parse(cls, rep, dom, lev=None): if isinstance(rep, tuple): num, den = rep if lev is not None: if isinstance(num, dict): num = dmp_from_dict(num, lev, dom) if isinstance(den, dict): den = dmp_from_dict(den, lev, dom) else: num, num_lev = dmp_validate(num) den, den_lev = dmp_validate(den) if num_lev == den_lev: lev = num_lev else: raise ValueError('inconsistent number of levels') if dmp_zero_p(den, lev): raise ZeroDivisionError('fraction denominator') if dmp_zero_p(num, lev): den = dmp_one(lev, dom) else: if dmp_negative_p(den, lev, dom): num = dmp_neg(num, lev, dom) den = dmp_neg(den, lev, dom) else: num = rep if lev is not None: if isinstance(num, dict): num = dmp_from_dict(num, lev, dom) elif not isinstance(num, list): num = dmp_ground(dom.convert(num), lev) else: num, lev = dmp_validate(num) den = dmp_one(lev, dom) return num, den, lev def __repr__(f): return "%s((%s, %s), %s)" % (f.__class__.__name__, f.num, f.den, f.dom) def __hash__(f): return hash((f.__class__.__name__, dmp_to_tuple(f.num, f.lev), dmp_to_tuple(f.den, f.lev), f.lev, f.dom)) def poly_unify(f, g): """Unify a multivariate fraction and a polynomial. """ if not isinstance(g, DMP) or f.lev != g.lev: raise UnificationFailed("Cannot unify %s with %s" % (f, g)) if f.dom == g.dom: return (f.lev, f.dom, f.per, (f.num, f.den), g._rep) else: lev, dom = f.lev, f.dom.unify(g.dom) F = (dmp_convert(f.num, lev, f.dom, dom), dmp_convert(f.den, lev, f.dom, dom)) G = dmp_convert(g._rep, lev, g.dom, dom) def per(num, den, cancel=True, kill=False, lev=lev): if kill: if not lev: return num/den else: lev = lev - 1 if cancel: num, den = dmp_cancel(num, den, lev, dom) return f.__class__.new((num, den), dom, lev) return lev, dom, per, F, G def frac_unify(f, g): """Unify representations of two multivariate fractions. """ if not isinstance(g, DMF) or f.lev != g.lev: raise UnificationFailed("Cannot unify %s with %s" % (f, g)) if f.dom == g.dom: return (f.lev, f.dom, f.per, (f.num, f.den), (g.num, g.den)) else: lev, dom = f.lev, f.dom.unify(g.dom) F = (dmp_convert(f.num, lev, f.dom, dom), dmp_convert(f.den, lev, f.dom, dom)) G = (dmp_convert(g.num, lev, g.dom, dom), dmp_convert(g.den, lev, g.dom, dom)) def per(num, den, cancel=True, kill=False, lev=lev): if kill: if not lev: return num/den else: lev = lev - 1 if cancel: num, den = dmp_cancel(num, den, lev, dom) return f.__class__.new((num, den), dom, lev) return lev, dom, per, F, G def per(f, num, den, cancel=True, kill=False): """Create a DMF out of the given representation. """ lev, dom = f.lev, f.dom if kill: if not lev: return num/den else: lev -= 1 if cancel: num, den = dmp_cancel(num, den, lev, dom) return f.__class__.new((num, den), dom, lev) def half_per(f, rep, kill=False): """Create a DMP out of the given representation. """ lev = f.lev if kill: if not lev: return rep else: lev -= 1 return DMP(rep, f.dom, lev) @classmethod def zero(cls, lev, dom): return cls.new(0, dom, lev) @classmethod def one(cls, lev, dom): return cls.new(1, dom, lev) def numer(f): """Returns the numerator of ``f``. """ return f.half_per(f.num) def denom(f): """Returns the denominator of ``f``. """ return f.half_per(f.den) def cancel(f): """Remove common factors from ``f.num`` and ``f.den``. """ return f.per(f.num, f.den) def neg(f): """Negate all coefficients in ``f``. """ return f.per(dmp_neg(f.num, f.lev, f.dom), f.den, cancel=False) def add_ground(f, c): """Add an element of the ground domain to ``f``. """ return f + f.ground_new(c) def add(f, g): """Add two multivariate fractions ``f`` and ``g``. """ if isinstance(g, DMP): lev, dom, per, (F_num, F_den), G = f.poly_unify(g) num, den = dmp_add_mul(F_num, F_den, G, lev, dom), F_den else: lev, dom, per, F, G = f.frac_unify(g) (F_num, F_den), (G_num, G_den) = F, G num = dmp_add(dmp_mul(F_num, G_den, lev, dom), dmp_mul(F_den, G_num, lev, dom), lev, dom) den = dmp_mul(F_den, G_den, lev, dom) return per(num, den) def sub(f, g): """Subtract two multivariate fractions ``f`` and ``g``. """ if isinstance(g, DMP): lev, dom, per, (F_num, F_den), G = f.poly_unify(g) num, den = dmp_sub_mul(F_num, F_den, G, lev, dom), F_den else: lev, dom, per, F, G = f.frac_unify(g) (F_num, F_den), (G_num, G_den) = F, G num = dmp_sub(dmp_mul(F_num, G_den, lev, dom), dmp_mul(F_den, G_num, lev, dom), lev, dom) den = dmp_mul(F_den, G_den, lev, dom) return per(num, den) def mul(f, g): """Multiply two multivariate fractions ``f`` and ``g``. """ if isinstance(g, DMP): lev, dom, per, (F_num, F_den), G = f.poly_unify(g) num, den = dmp_mul(F_num, G, lev, dom), F_den else: lev, dom, per, F, G = f.frac_unify(g) (F_num, F_den), (G_num, G_den) = F, G num = dmp_mul(F_num, G_num, lev, dom) den = dmp_mul(F_den, G_den, lev, dom) return per(num, den) def pow(f, n): """Raise ``f`` to a non-negative power ``n``. """ if isinstance(n, int): num, den = f.num, f.den if n < 0: num, den, n = den, num, -n return f.per(dmp_pow(num, n, f.lev, f.dom), dmp_pow(den, n, f.lev, f.dom), cancel=False) else: raise TypeError("``int`` expected, got %s" % type(n)) def quo(f, g): """Computes quotient of fractions ``f`` and ``g``. """ if isinstance(g, DMP): lev, dom, per, (F_num, F_den), G = f.poly_unify(g) num, den = F_num, dmp_mul(F_den, G, lev, dom) else: lev, dom, per, F, G = f.frac_unify(g) (F_num, F_den), (G_num, G_den) = F, G num = dmp_mul(F_num, G_den, lev, dom) den = dmp_mul(F_den, G_num, lev, dom) return per(num, den) exquo = quo def invert(f, check=True): """Computes inverse of a fraction ``f``. """ return f.per(f.den, f.num, cancel=False) @property def is_zero(f): """Returns ``True`` if ``f`` is a zero fraction. """ return dmp_zero_p(f.num, f.lev) @property def is_one(f): """Returns ``True`` if ``f`` is a unit fraction. """ return dmp_one_p(f.num, f.lev, f.dom) and \ dmp_one_p(f.den, f.lev, f.dom) def __neg__(f): return f.neg() def __add__(f, g): if isinstance(g, (DMP, DMF)): return f.add(g) elif g in f.dom: return f.add_ground(f.dom.convert(g)) try: return f.add(f.half_per(g)) except (TypeError, CoercionFailed, NotImplementedError): return NotImplemented def __radd__(f, g): return f.__add__(g) def __sub__(f, g): if isinstance(g, (DMP, DMF)): return f.sub(g) try: return f.sub(f.half_per(g)) except (TypeError, CoercionFailed, NotImplementedError): return NotImplemented def __rsub__(f, g): return (-f).__add__(g) def __mul__(f, g): if isinstance(g, (DMP, DMF)): return f.mul(g) try: return f.mul(f.half_per(g)) except (TypeError, CoercionFailed, NotImplementedError): return NotImplemented def __rmul__(f, g): return f.__mul__(g) def __pow__(f, n): return f.pow(n) def __truediv__(f, g): if isinstance(g, (DMP, DMF)): return f.quo(g) try: return f.quo(f.half_per(g)) except (TypeError, CoercionFailed, NotImplementedError): return NotImplemented def __rtruediv__(self, g): return self.invert(check=False)*g def __eq__(f, g): try: if isinstance(g, DMP): _, _, _, (F_num, F_den), G = f.poly_unify(g) if f.lev == g.lev: return dmp_one_p(F_den, f.lev, f.dom) and F_num == G else: _, _, _, F, G = f.frac_unify(g) if f.lev == g.lev: return F == G except UnificationFailed: pass return False def __ne__(f, g): try: if isinstance(g, DMP): _, _, _, (F_num, F_den), G = f.poly_unify(g) if f.lev == g.lev: return not (dmp_one_p(F_den, f.lev, f.dom) and F_num == G) else: _, _, _, F, G = f.frac_unify(g) if f.lev == g.lev: return F != G except UnificationFailed: pass return True def __lt__(f, g): _, _, _, F, G = f.frac_unify(g) return F < G def __le__(f, g): _, _, _, F, G = f.frac_unify(g) return F <= G def __gt__(f, g): _, _, _, F, G = f.frac_unify(g) return F > G def __ge__(f, g): _, _, _, F, G = f.frac_unify(g) return F >= G def __bool__(f): return not dmp_zero_p(f.num, f.lev) def init_normal_ANP(rep, mod, dom): return ANP(dup_normal(rep, dom), dup_normal(mod, dom), dom) class ANP(CantSympify): """Dense Algebraic Number Polynomials over a field. """ __slots__ = ('_rep', '_mod', 'dom') def __new__(cls, rep, mod, dom): if isinstance(rep, DMP): pass elif type(rep) is dict: # don't use isinstance rep = DMP(dup_from_dict(rep, dom), dom, 0) else: if isinstance(rep, list): rep = [dom.convert(a) for a in rep] else: rep = [dom.convert(rep)] rep = DMP(dup_strip(rep), dom, 0) if isinstance(mod, DMP): pass elif isinstance(mod, dict): mod = DMP(dup_from_dict(mod, dom), dom, 0) else: mod = DMP(dup_strip(mod), dom, 0) return cls.new(rep, mod, dom) @classmethod def new(cls, rep, mod, dom): if not (rep.dom == mod.dom == dom): raise RuntimeError("Inconsistent domain") obj = super().__new__(cls) obj._rep = rep obj._mod = mod obj.dom = dom return obj # XXX: It should be possible to use __getnewargs__ rather than __reduce__ # but it doesn't work for some reason. Probably this would be easier if # python-flint supported pickling for polynomial types. def __reduce__(self): return ANP, (self.rep, self.mod, self.dom) @property def rep(self): return self._rep.to_list() @property def mod(self): return self.mod_to_list() def to_DMP(self): return self._rep def mod_to_DMP(self): return self._mod def per(f, rep): return f.new(rep, f._mod, f.dom) def __repr__(f): return "%s(%s, %s, %s)" % (f.__class__.__name__, f._rep.to_list(), f._mod.to_list(), f.dom) def __hash__(f): return hash((f.__class__.__name__, f.to_tuple(), f._mod.to_tuple(), f.dom)) def convert(f, dom): """Convert ``f`` to a ``ANP`` over a new domain. """ if f.dom == dom: return f else: return f.new(f._rep.convert(dom), f._mod.convert(dom), dom) def unify(f, g): """Unify representations of two algebraic numbers. """ # XXX: This unify method is not used any more because unify_ANP is used # instead. if not isinstance(g, ANP) or f.mod != g.mod: raise UnificationFailed("Cannot unify %s with %s" % (f, g)) if f.dom == g.dom: return f.dom, f.per, f.rep, g.rep, f.mod else: dom = f.dom.unify(g.dom) F = dup_convert(f.rep, f.dom, dom) G = dup_convert(g.rep, g.dom, dom) if dom != f.dom and dom != g.dom: mod = dup_convert(f.mod, f.dom, dom) else: if dom == f.dom: mod = f.mod else: mod = g.mod per = lambda rep: ANP(rep, mod, dom) return dom, per, F, G, mod def unify_ANP(f, g): """Unify and return ``DMP`` instances of ``f`` and ``g``. """ if not isinstance(g, ANP) or f._mod != g._mod: raise UnificationFailed("Cannot unify %s with %s" % (f, g)) # The domain is almost always QQ but there are some tests involving ZZ if f.dom != g.dom: dom = f.dom.unify(g.dom) f = f.convert(dom) g = g.convert(dom) return f._rep, g._rep, f._mod, f.dom @classmethod def zero(cls, mod, dom): return ANP(0, mod, dom) @classmethod def one(cls, mod, dom): return ANP(1, mod, dom) def to_dict(f): """Convert ``f`` to a dict representation with native coefficients. """ return f._rep.to_dict() def to_sympy_dict(f): """Convert ``f`` to a dict representation with SymPy coefficients. """ rep = dmp_to_dict(f.rep, 0, f.dom) for k, v in rep.items(): rep[k] = f.dom.to_sympy(v) return rep def to_list(f): """Convert ``f`` to a list representation with native coefficients. """ return f._rep.to_list() def mod_to_list(f): """Return ``f.mod`` as a list with native coefficients. """ return f._mod.to_list() def to_sympy_list(f): """Convert ``f`` to a list representation with SymPy coefficients. """ return [ f.dom.to_sympy(c) for c in f.to_list() ] def to_tuple(f): """ Convert ``f`` to a tuple representation with native coefficients. This is needed for hashing. """ return f._rep.to_tuple() @classmethod def from_list(cls, rep, mod, dom): return ANP(dup_strip(list(map(dom.convert, rep))), mod, dom) def add_ground(f, c): """Add an element of the ground domain to ``f``. """ return f.per(f._rep.add_ground(c)) def sub_ground(f, c): """Subtract an element of the ground domain from ``f``. """ return f.per(f._rep.sub_ground(c)) def mul_ground(f, c): """Multiply ``f`` by an element of the ground domain. """ return f.per(f._rep.mul_ground(c)) def quo_ground(f, c): """Quotient of ``f`` by an element of the ground domain. """ return f.per(f._rep.quo_ground(c)) def neg(f): return f.per(f._rep.neg()) def add(f, g): F, G, mod, dom = f.unify_ANP(g) return f.new(F.add(G), mod, dom) def sub(f, g): F, G, mod, dom = f.unify_ANP(g) return f.new(F.sub(G), mod, dom) def mul(f, g): F, G, mod, dom = f.unify_ANP(g) return f.new(F.mul(G).rem(mod), mod, dom) def pow(f, n): """Raise ``f`` to a non-negative power ``n``. """ if not isinstance(n, int): raise TypeError("``int`` expected, got %s" % type(n)) mod = f._mod F = f._rep if n < 0: F, n = F.invert(mod), -n # XXX: Need a pow_mod method for DMP return f.new(F.pow(n).rem(f._mod), mod, f.dom) def exquo(f, g): F, G, mod, dom = f.unify_ANP(g) return f.new(F.mul(G.invert(mod)).rem(mod), mod, dom) def div(f, g): return f.exquo(g), f.zero(f._mod, f.dom) def quo(f, g): return f.exquo(g) def rem(f, g): F, G, mod, dom = f.unify_ANP(g) s, h = F.half_gcdex(G) if h.is_one: return f.zero(mod, dom) else: raise NotInvertible("zero divisor") def LC(f): """Returns the leading coefficient of ``f``. """ return f._rep.LC() def TC(f): """Returns the trailing coefficient of ``f``. """ return f._rep.TC() @property def is_zero(f): """Returns ``True`` if ``f`` is a zero algebraic number. """ return f._rep.is_zero @property def is_one(f): """Returns ``True`` if ``f`` is a unit algebraic number. """ return f._rep.is_one @property def is_ground(f): """Returns ``True`` if ``f`` is an element of the ground domain. """ return f._rep.is_ground def __pos__(f): return f def __neg__(f): return f.neg() def __add__(f, g): if isinstance(g, ANP): return f.add(g) try: g = f.dom.convert(g) except CoercionFailed: return NotImplemented else: return f.add_ground(g) def __radd__(f, g): return f.__add__(g) def __sub__(f, g): if isinstance(g, ANP): return f.sub(g) try: g = f.dom.convert(g) except CoercionFailed: return NotImplemented else: return f.sub_ground(g) def __rsub__(f, g): return (-f).__add__(g) def __mul__(f, g): if isinstance(g, ANP): return f.mul(g) try: g = f.dom.convert(g) except CoercionFailed: return NotImplemented else: return f.mul_ground(g) def __rmul__(f, g): return f.__mul__(g) def __pow__(f, n): return f.pow(n) def __divmod__(f, g): return f.div(g) def __mod__(f, g): return f.rem(g) def __truediv__(f, g): if isinstance(g, ANP): return f.quo(g) try: g = f.dom.convert(g) except CoercionFailed: return NotImplemented else: return f.quo_ground(g) def __eq__(f, g): try: F, G, _, _ = f.unify_ANP(g) except UnificationFailed: return NotImplemented return F == G def __ne__(f, g): try: F, G, _, _ = f.unify_ANP(g) except UnificationFailed: return NotImplemented return F != G def __lt__(f, g): F, G, _, _ = f.unify_ANP(g) return F < G def __le__(f, g): F, G, _, _ = f.unify_ANP(g) return F <= G def __gt__(f, g): F, G, _, _ = f.unify_ANP(g) return F > G def __ge__(f, g): F, G, _, _ = f.unify_ANP(g) return F >= G def __bool__(f): return bool(f._rep)