diff --git a/src/doc/en/reference/coding/index.rst b/src/doc/en/reference/coding/index.rst index ae6b46587d5..eacbbf46b0f 100644 --- a/src/doc/en/reference/coding/index.rst +++ b/src/doc/en/reference/coding/index.rst @@ -30,7 +30,6 @@ Linear codes and related constructions sage/coding/hamming_code sage/coding/guruswami_sudan/gs_decoder sage/coding/guruswami_sudan/interpolation - sage/coding/guruswami_sudan/rootfinding sage/coding/guruswami_sudan/utils sage/coding/subfield_subcode sage/coding/code_constructions diff --git a/src/sage/coding/guruswami_sudan/gs_decoder.py b/src/sage/coding/guruswami_sudan/gs_decoder.py index 2f1a039071d..6185bf0f38d 100644 --- a/src/sage/coding/guruswami_sudan/gs_decoder.py +++ b/src/sage/coding/guruswami_sudan/gs_decoder.py @@ -31,7 +31,6 @@ from sage.rings.integer_ring import ZZ from sage.coding.decoder import Decoder from sage.coding.guruswami_sudan.interpolation import gs_interpolation_linalg, gs_interpolation_lee_osullivan -from sage.coding.guruswami_sudan.rootfinding import rootfind_roth_ruckenstein from sage.coding.guruswami_sudan.utils import (johnson_radius, gilt, solve_degree2_to_integer_range) @@ -85,6 +84,24 @@ def n_k_params(C, n_k): elif n_k is not None: return n_k +def roth_ruckenstein_root_finder(p, maxd=None, precision=None): + """ + Wrapper for Roth-Ruckenstein algorithm to compute the roots of a polynomial + with coefficients in ``F[x]``. + + TESTS:: + + sage: from sage.coding.guruswami_sudan.gs_decoder import roth_ruckenstein_root_finder + sage: R. = GF(13)[] + sage: S. = R[] + sage: p = (y - x^2 - x - 1) * (y + x + 1) + sage: roth_ruckenstein_root_finder(p, maxd = 2) + [12*x + 12, x^2 + x + 1] + """ + gens = p.parent().gens() + if len(gens) == 2: + p = p.polynomial(gens[1]) + return p.roots(multiplicities=False, degree_bound=maxd, algorithm="Roth-Ruckenstein") class GRSGuruswamiSudanDecoder(Decoder): r""" @@ -155,8 +172,7 @@ class GRSGuruswamiSudanDecoder(Decoder): ``my_rootfinder(Q, maxd=default_value, precision=default_value)``. `Q` will be given as an element of `F[x][y]`. The function must return the roots as a list of polynomials over a univariate polynomial ring. See - :meth:`sage.coding.guruswami_sudan.rootfinding.rootfind_roth_ruckenstein` - for an example. + :meth:`roth_ruckenstein_root_finder` for an example. If one provides a function as ``interpolation_alg``, its signature has to be: ``my_inter(interpolation_points, tau, s_and_l, wy)``. See @@ -178,8 +194,8 @@ class GRSGuruswamiSudanDecoder(Decoder): One can pass a method as ``root_finder`` (works also for ``interpolation_alg``):: - sage: from sage.coding.guruswami_sudan.rootfinding import rootfind_roth_ruckenstein - sage: rf = rootfind_roth_ruckenstein + sage: from sage.coding.guruswami_sudan.gs_decoder import roth_ruckenstein_root_finder + sage: rf = roth_ruckenstein_root_finder sage: D = codes.decoders.GRSGuruswamiSudanDecoder(C, parameters = (1,2), root_finder = rf) sage: D Guruswami-Sudan decoder for [250, 70, 181] Generalized Reed-Solomon Code over Finite Field of size 251 decoding 97 errors with parameters (1, 2) @@ -189,8 +205,6 @@ class GRSGuruswamiSudanDecoder(Decoder): This works for ``interpolation_alg`` as well:: - sage: from sage.coding.guruswami_sudan.rootfinding import rootfind_roth_ruckenstein - sage: rf = rootfind_roth_ruckenstein sage: D = codes.decoders.GRSGuruswamiSudanDecoder(C, parameters = (1,2), root_finder="RothRuckenstein") sage: D Guruswami-Sudan decoder for [250, 70, 181] Generalized Reed-Solomon Code over Finite Field of size 251 decoding 97 errors with parameters (1, 2) @@ -576,7 +590,7 @@ def __init__(self, code, tau = None, parameters = None, interpolation_alg = None if hasattr(root_finder, '__call__'): self._root_finder = root_finder elif root_finder == None or root_finder == "RothRuckenstein": - self._root_finder = rootfind_roth_ruckenstein + self._root_finder = roth_ruckenstein_root_finder else: raise ValueError("Please provide a method or one of the allowed strings for root_finder") super(GRSGuruswamiSudanDecoder, self).__init__(code, code.ambient_space(), "EvaluationPolynomial") @@ -640,8 +654,8 @@ def interpolation_algorithm(self): sage: C = codes.GeneralizedReedSolomonCode(GF(251).list()[:250], 70) sage: D = C.decoder("GuruswamiSudan", tau = 97) - sage: D.interpolation_algorithm() #random - + sage: D.interpolation_algorithm() + """ return self._interpolation_alg @@ -651,15 +665,16 @@ def rootfinding_algorithm(self): Remember that its signature has to be: ``my_rootfinder(Q, maxd=default_value, precision=default_value)``. - See :meth:`sage.coding.guruswami_sudan.rootfinding.rootfind_roth_ruckenstein` + See :meth:`roth_ruckenstein_root_finder` for an example. EXAMPLES:: + sage: from sage.coding.guruswami_sudan.gs_decoder import roth_ruckenstein_root_finder sage: C = codes.GeneralizedReedSolomonCode(GF(251).list()[:250], 70) sage: D = C.decoder("GuruswamiSudan", tau = 97) - sage: D.rootfinding_algorithm() #random - + sage: D.rootfinding_algorithm() + """ return self._root_finder diff --git a/src/sage/coding/guruswami_sudan/rootfinding.py b/src/sage/coding/guruswami_sudan/rootfinding.py deleted file mode 100644 index bbd6bf5ec2d..00000000000 --- a/src/sage/coding/guruswami_sudan/rootfinding.py +++ /dev/null @@ -1,373 +0,0 @@ -r""" -Finding `F[x]`-roots for polynomials over `F[x][y]`, with`F` is a (finite) field, as used in the Guruswami-Sudan decoding. - -This module contains functions for finding two types of `F[x]` roots in a -polynomial over `F[x][y]`, where `F` is a field. Note that if `F` is an infinite -field, then these functions should work, but no measures are taken to limit -coefficient growth. The functions also assume the existence of a root finding -procedure in `F[x]`. - -Given a `Q(x,y) \in F[x,y]`, the first type of root are actual `F[x]` roots, -i.e. polynomials `f(x) \in F[x]` such that `Q(x, f(x)) = 0`. - -The second type of root are modular roots: given `Q(x,y) \in F[x,y]` and a -precision `d \in \ZZ_+`, then we find all `f(x) \in F[x]` such that `Q(x, f(x)) -\equiv 0 \mod x^d`. Since this set is infinite, we return a succinct description -of all such polynomials: this is given as pairs `(f(x), h)` such that `f(x) + -g(x)x^h` is a modular root of `Q` for any `g \in F[x]`. - - -AUTHORS: - -- Johan S. R. Nielsen, original implementation (see [Nielsen]_ for details) -- David Lucas, ported the original implementation in Sage -""" - -#***************************************************************************** -# Copyright (C) 2015 David Lucas -# 2015 Johan S. R. Nielsen -# -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 2 of the License, or -# (at your option) any later version. -# http://www.gnu.org/licenses/ -#***************************************************************************** - - -from sage.rings.infinity import infinity -from sage.functions.other import binomial, floor - -def _convert_Q_representation(Q): - r""" - Converts the bivariate polynomial ``Q(x,y)`` from any accepted type to the - internally used `F[x]` list, for a field `F`. - - INPUT: - - - ``Q`` -- a bivariate polynomial, represented either over `F[x,y]`, `F[x][y]` or `F[x]` list. - - EXAMPLES:: - - sage: from sage.coding.guruswami_sudan.rootfinding import _convert_Q_representation - sage: F = GF(17) - sage: Px. = F[] - sage: Pxy. = Px[] - sage: Q1 = (y - (x**2 + x + 1)) * (y**2 - x + 1) * (y - (x**3 + 4*x + 16)) - sage: _convert_Q_representation(Q1) - ([16*x^6 + 13*x^4 + 2*x^3 + 4*x + 16, - x^4 + 4*x^2 + 12*x, - x^5 + x^4 + 5*x^3 + 3*x^2 + 2*x, - 16*x^3 + 16*x^2 + 12*x, - 1], Univariate Polynomial Ring in x over Finite Field of size 17) - """ - if isinstance(Q, list): - if Q == []: - return ([], None) - Rx = Q[0].parent() - if not hasattr(Rx,'gen'): - raise ValueError("Q must be given as F[x][y], F[x,y] or as F[x] list.") - return (Q, Rx) - else: - # Find out if Q is in F[x,y] or F[x][y] - Qorig = Q - Ryx = Q.parent() - #TODO: Check Ryx is a polynomial ring over a field - if len(Ryx.gens())==1: - # Ok, Q is in F[x][y] - Rx = Q.base_ring() - elif len(Ryx.gens())==2: - F = Ryx.base_ring() - (xs,ys) = Ryx.variable_names() - Rx = F[xs] - Ryx = Rx[ys] - x, y = Rx.gen(), Ryx.gen() - Q = Ryx(Q) - else: - raise ValueError("Q must be given as F[x][y], F[x,y] or as F[x] list.") - # Then make sure Q is a list of F[x] elements - return (Q.list(), Rx) - -def _sanitise_rootfinding_input(Q, maxd, precision): - r""" - Verifies, converts and sanitises legal input to the root-finding procedures, - as well as returning relevant helper variables. - - INPUT: - - - ``Q`` -- a bivariate polynomial, represented either over `F[x,y]`, `F[x][y]` or `F[x]` list. - - - ``maxd``, an integer, the maximal degree of a root of ``Q`` that we're - interested in, possibly ``None``. - - - ``precision``, an integer, the precision asked for modular roots of ``Q``, possibly ``None``. - - OUTPUT: - - - ``Q``, a modified version of ``Q``, where all monomials have been - truncated to ``precision``. Represented as an `F[x]` list. - - - ``Qinp``, the original ``Q`` passed in input, represented as an `F[x]` list. - - - ``F``, the base ring of the coefficients in ``Q``'s first variable. - - - ``Rx``, the polynomial ring `F[x]`. - - - ``x``, the generator of ``Rx``. - - - ``maxd``, the maximal degree of a root of ``Q`` that we're interested in, - possibly inferred according ``precision``. - - EXAMPLES:: - - sage: from sage.coding.guruswami_sudan.rootfinding import _sanitise_rootfinding_input - sage: F = GF(17) - sage: Px. = F[] - sage: Pxy. = Px[] - sage: Q = (y - (x**2 + x + 1)) * (y**2 - x + 1) * (y - (x**3 + 4*x + 16)) - sage: _sanitise_rootfinding_input(Q, None, None) - ([16*x^6 + 13*x^4 + 2*x^3 + 4*x + 16, - x^4 + 4*x^2 + 12*x, - x^5 + x^4 + 5*x^3 + 3*x^2 + 2*x, - 16*x^3 + 16*x^2 + 12*x, - 1], - [16*x^6 + 13*x^4 + 2*x^3 + 4*x + 16, - x^4 + 4*x^2 + 12*x, - x^5 + x^4 + 5*x^3 + 3*x^2 + 2*x, - 16*x^3 + 16*x^2 + 12*x, - 1], - Finite Field of size 17, - Univariate Polynomial Ring in x over Finite Field of size 17, - x, - 3) - """ - (Q, Rx) = _convert_Q_representation(Q) - if Q == []: - return ([],[],None,Rx,None,0) # Q == 0 so just bail - Qinp = Q - F = Rx.base_ring() - x = Rx.gen() - - if not maxd: - if precision: - maxd = precision-1 - else: - #The maximal degree of a root is at most - #(di-dl)/(l-i) d being the degree of a monomial - maxd = 0 - l = len(Q) -1 - dl = Q[l].degree() - for i in range(l): - qi = Q[i] - if not qi.is_zero(): - tmp = floor((qi.degree() - dl) / (l - i)) - if tmp > maxd: - maxd = tmp - if precision: - for t in range(len(Q)): - if Q[t].degree >= precision: - Q[t] = Q[t].truncate(precision) - return (Q, Qinp, F, Rx, x, maxd) - -def _strip_x_pows(Q): - r""" - Returns ``(Q', s)`` where ``Q'`` is ``Q`` whose elements have all been - divided by the largest power ``s`` of ``x`` possible such that all these - elements remain polynomials. - - INPUT: - - - ``Q`` -- a bivariate polynomial in `F[x][y]` as a list of `F[x]` polynomials. - - OUTPUT: - - - ``(Q', s)`` a list of two elements: - - - ``Q'``, the reduced bivariate polynomial, as a list of univariate ones. - - ``s``, an integer, the power of `x` that was stripped from `Q`. - - EXAMPLES:: - - sage: from sage.coding.guruswami_sudan.rootfinding import _strip_x_pows - sage: F = GF(17) - sage: Px. = F[] - sage: Q = [ Px(3*x^2 + 2*x), Px(5*x^7 + x^6) ] - sage: _strip_x_pows(Q) - ([3*x + 2, 5*x^6 + x^5], 1) - """ - def lead_zeroes(p): - if p.is_zero(): - return infinity - i = 0 - while p[i].is_zero(): - i+=1 - return i - strip = min([lead_zeroes(p) for p in Q]) - if strip == 0: - return (Q, 0) - if strip == infinity: - return ([ Q[0].parent().zero() ], infinity) - return ([ p.shift(-strip) for p in Q ] , strip) - -def _roth_ruckenstein_i(Q, F, Rx, x, maxd, precision): - r""" - Returns all polynomials which are a solution to the, possibly modular, - root-finding problem. - - This is the core of Roth-Ruckenstein's algorithm where all conversions, - checks and parent-extraction have been done. Most of the inputs corresponds - to the output of ``_sanitise_rootfinding_input``. - - INPUT: - - - ``Q``, a modified version of ``Q``, where all monomials have been - truncated to ``precision``. Represented as an `F[x]` list. - - - ``Qinp``, the original ``Q`` passed in input, represented as an `F[x]` list. - - - ``F``, the base ring of the coefficients in ``Q``'s first variable. - - - ``Rx``, the polynomial ring `F[x]`. - - - ``x``, the generator of ``Rx``. - - - ``maxd``, the maximal degree of a root of ``Q`` that we're interested in, - possibly inferred according ``precision``. - - - ``precision``, a non-negative integer or `None`. If given, it is the - sought precision for modular roots of `Q`. Otherwise, we will find - unconditional roots. - - OUTPUT: - - - a list, containing all `F[x]` roots of `Q(x,y)`, possibly modular. If - ``precision`` is given, we return a list of pairs `(f, h)`, where `f \in - F[x]` and `h` is a non-negative integer, and where `f + h*g \equiv 0 \mod - x^{d}` for any `g \in F[x]`, and where `d` is ``precision``. - - EXAMPLES:: - - sage: from sage.coding.guruswami_sudan.rootfinding import _sanitise_rootfinding_input - sage: from sage.coding.guruswami_sudan.rootfinding import _roth_ruckenstein_i - sage: F = GF(17) - sage: Px. = F[] - sage: Pxy. = Px[] - sage: Q = (y - (x**2 + x + 1)) * (y**2 - x + 1) * (y - (x**3 + 4*x + 16)) - sage: (Q, Qinp, F, Rx, x, maxd) = _sanitise_rootfinding_input(Q, None, None) - sage: set(_roth_ruckenstein_i(Q, F, Rx, x, maxd, None)) - {x^2 + x + 1, x^3 + 4*x + 16} - sage: set(_roth_ruckenstein_i(Q, F, Rx, x, maxd, precision=2)) - {(x + 1, 2), (2*x + 13, 2), (4*x + 16, 2), (15*x + 4, 2)} - """ - solutions = [] - g = [F.zero()] * (maxd+1) - - def roth_rec(Q, lam, k): - r""" - Recursion of the root finding: - Q is the remaining poly, lam is the power of x whose coefficient we are - to determine now, and k is the remaining precision to handle (if ``precision`` is given) - """ - if precision and k <= 0: - solutions.append((Rx(g[:lam]), lam)) - return - (T, strip) = _strip_x_pows(Q) - if precision: - k = k - strip - Ty = Rx([ p[0] for p in T ]) - if Ty.is_zero() or (precision and k <= 0): - if precision: - solutions.append((Rx(g[:lam]), lam)) - else: - solutions.append(Rx(g[:lam])) - return - roots = Ty.roots(multiplicities=False) - for gamma in roots: - g[lam] = gamma - if lam = F[] - sage: Py. = Px[] - sage: Q = (y - (x**2 + x + 1)) * (y**2 - x + 1) * (y - (x**3 + 4*x + 16)) - sage: roots = rootfind_roth_ruckenstein(Q); set(roots) - {x^2 + x + 1, x^3 + 4*x + 16} - sage: Q(roots[0]) - 0 - sage: set(rootfind_roth_ruckenstein(Q, maxd = 2)) - {x^2 + x + 1} - sage: modroots = rootfind_roth_ruckenstein(Q, precision = 3); set(modroots) - {(4*x + 16, 3), (x^2 + x + 1, 3), (8*x^2 + 15*x + 4, 3), (9*x^2 + 2*x + 13, 3)} - sage: (f,h) = modroots[0] - sage: Q(f + x^h * Px.random_element()) % x^3 - 0 - sage: modroots2 = rootfind_roth_ruckenstein(Q, maxd=1, precision = 3); set(modroots2) - {(x + 1, 2), (2*x + 13, 2), (4*x + 16, 2), (15*x + 4, 2)} - - TESTS: - - Test that if `Q = 0`, then the appropriate response is given - - sage: F = GF(17) - sage: R. = F[] - sage: rootfind_roth_ruckenstein(R.zero()) - ValueError('The zero polynomial has infinitely many roots.',) - sage: rootfind_roth_ruckenstein(R.zero(), precision=1) - [(0, 0)] - - """ - (Q, Qinp, F, Rx, x, maxd) = _sanitise_rootfinding_input(Q, maxd, precision) - if all(p.is_zero() for p in Q): - if precision: - return [(Rx.zero() if hasattr(Rx,'zero') else 0, 0)] - else: - return ValueError("The zero polynomial has infinitely many roots.") - return _roth_ruckenstein_i(Q, F, Rx, x, maxd, precision) diff --git a/src/sage/rings/polynomial/polynomial_element.pyx b/src/sage/rings/polynomial/polynomial_element.pyx index e20766979a9..c02702960c3 100644 --- a/src/sage/rings/polynomial/polynomial_element.pyx +++ b/src/sage/rings/polynomial/polynomial_element.pyx @@ -3918,6 +3918,19 @@ cdef class Polynomial(CommutativeAlgebraElement): if self.degree() == 0: return Factorization([], unit=self[0]) + # Use multivariate implementations for polynomials over polynomial rings + variables = self._parent.variable_names_recursive() + if len(variables) > 1: + base = self._parent._mpoly_base_ring() + ring = PolynomialRing(base, variables) + if ring._has_singular: + try: + d = self._mpoly_dict_recursive() + F = ring(d).factor(**kwargs) + return Factorization([(self._parent(f),m) for (f,m) in F], unit = F.unit()) + except NotImplementedError: + pass + R = self.parent().base_ring() if hasattr(R, '_factor_univariate_polynomial'): return R._factor_univariate_polynomial(self, **kwargs) @@ -6291,7 +6304,7 @@ cdef class Polynomial(CommutativeAlgebraElement): return self.parent()(v) - def roots(self, ring=None, multiplicities=True, algorithm=None): + def roots(self, ring=None, multiplicities=True, algorithm=None, **kwds): """ Return the roots of this polynomial (by default, in the base ring of this polynomial). @@ -6873,7 +6886,10 @@ cdef class Polynomial(CommutativeAlgebraElement): """ K = self.parent().base_ring() if hasattr(K, '_roots_univariate_polynomial'): - return K._roots_univariate_polynomial(self, ring=ring, multiplicities=multiplicities, algorithm=algorithm) + return K._roots_univariate_polynomial(self, ring=ring, multiplicities=multiplicities, algorithm=algorithm, **kwds) + + if kwds: + raise TypeError("roots() got unexpected keyword argument(s): {}".format(kwds.keys())) L = K if ring is None else ring @@ -8621,7 +8637,6 @@ cdef class Polynomial(CommutativeAlgebraElement): raise ValueError("not a %s power"%m.ordinal_str()) raise ValueError("not a %s power"%m.ordinal_str()) - # ----------------- inner functions ------------- # Cython can't handle function definitions inside other function diff --git a/src/sage/rings/polynomial/polynomial_ring.py b/src/sage/rings/polynomial/polynomial_ring.py index 43cc6512660..126dac037ac 100644 --- a/src/sage/rings/polynomial/polynomial_ring.py +++ b/src/sage/rings/polynomial/polynomial_ring.py @@ -1509,6 +1509,42 @@ def weyl_algebra(self): from sage.algebras.weyl_algebra import DifferentialWeylAlgebra return DifferentialWeylAlgebra(self) + def _roots_univariate_polynomial(self, p, ring=None, multiplicities=True, algorithm=None, degree_bound=None): + """ + Return the list of roots of ``p``. + + INPUT: + + - ``p`` -- the polynomial whose roots are computed + - ``ring`` -- the ring to find roots (default is the base ring of ``p``) + - ``multiplicities`` -- bool (default: True): if ``True``, return a list of pairs ``(root, multiplicity)``; if ``False`` return a list of roots + - ``algorithm`` -- ignored (TODO: remove) + - ``degree_bound``-- if not ``None``, return only roots of degree at most ``degree_bound`` + + EXAMPLES:: + + sage: R. = QQ[] + sage: S. = R[] + sage: p = y^3 + (-x^2 - 3)*y^2 + (2*x^3 - x^2 + 3)*y - x^4 + 2*x^2 - 1 + sage: p.roots() + [(x^2 - 2*x + 1, 1), (x + 1, 2)] + sage: p.roots(multiplicities=False) + [x^2 - 2*x + 1, x + 1] + sage: p.roots(degree_bound=1) + [(x + 1, 2)] + """ + if ring is not None and ring is not self: + p = p.change_ring(ring) + return p.roots(multiplicities, algorithm, degree_bound) + roots = p._roots_from_factorization(p.factor(), multiplicities) + if degree_bound is not None: + if multiplicities: + roots = [(r,m) for (r,m) in roots if r.degree() <= degree_bound] + else: + roots = [r for r in roots if r.degree() <= degree_bound] + return roots + + class PolynomialRing_integral_domain(PolynomialRing_commutative, ring.IntegralDomain): def __init__(self, base_ring, name="x", sparse=False, implementation=None, element_class=None, category=None): @@ -2054,6 +2090,142 @@ def irreducible_element(self, n, algorithm=None): else: raise ValueError("no such algorithm for finding an irreducible polynomial: %s" % algorithm) + def _roth_ruckenstein(self, p, degree_bound, precision): + r""" + Returns all polynomials which are a solution to the, possibly modular, + root-finding problem. + + This is the core of Roth-Ruckenstein's algorithm where all conversions, + checks and parent-extraction have been done. + + INPUT: + + - ``p`` -- a nonzero polynomial over ``F[x][y]``. The polynomial ``p`` + should be first truncated to ``precision`` + + - ``degree_bound`` -- a bound on the degree of the roots of ``p`` that + the algorithm computes + + - ``precision`` -- a non-negative integer or `None`. If given, it is the + sought precision for modular roots of `p`. Otherwise, the algorithm + computes unconditional roots. + + OUTPUT: + + - a list containing all `F[x]` roots of `p(y)`, possibly modular. If + ``precision`` is given, the algorithm returns a list of pairs `(f, h)`, + where `f` is a polynomial and `h` is a non-negative integer such that + `p(f + x^{h}*g) \equiv 0 \mod x^{d}` for any `g \in F[x]`, where `d` is + ``precision``. + + EXAMPLES:: + + sage: F = GF(17) + sage: Px. = F[] + sage: Pxy. = Px[] + sage: p = (y - (x**2 + x + 1)) * (y**2 - x + 1) * (y - (x**3 + 4*x + 16)) + sage: Px._roth_ruckenstein(p, 3, None) + [x^3 + 4*x + 16, x^2 + x + 1] + sage: Px._roth_ruckenstein(p, 2, None) + [x^2 + x + 1] + sage: Px._roth_ruckenstein(p, 1, 2) + [(4*x + 16, 2), (2*x + 13, 2), (15*x + 4, 2), (x + 1, 2)] + """ + def roth_rec(p, lam, k, g): + r""" + Recursive core method for Roth-Ruckenstein algorithm. + + INPUT: + + - ``p`` -- the current value of the polynomial + - ``lam`` -- is the power of x whose coefficient is being computed + - ``k`` -- the remaining precision to handle (if ``precision`` is given) + - ``g`` -- the root being computed + """ + if precision and k <= 0: + solutions.append((g, lam)) + return + val = min(c.valuation() for c in p) + if precision: + k = k - val + T = p.map_coefficients(lambda c:c.shift(-val)) + Ty = T.map_coefficients(lambda c:c[0]).change_ring(F) + if Ty.is_zero() or (precision and k <= 0): + if precision: + solutions.append((g, lam)) + else: + solutions.append(g) + return + roots = Ty.roots(multiplicities=False) + for gamma in roots: + g_new = g + gamma*x**lam + if lam < degree_bound: + Tg = T(x*y + gamma) + roth_rec(Tg , lam+1, k, g_new) + else: + if precision: + solutions.append((g_new, lam+1)) + elif p(gamma).is_zero(): + solutions.append(g_new) + return + + x = self.gen() + y = p.parent().gen() + F = self.base_ring() + solutions = [] + g = self.zero() + + roth_rec(p, 0, precision, g) + return solutions + + def _roots_univariate_polynomial(self, p, ring=None, multiplicities=False, algorithm=None, degree_bound=None): + """ + Return the list of roots of ``p``. + + INPUT: + + - ``p`` -- the polynomial whose roots are computed + - ``ring`` -- the ring to find roots (default is the base ring of ``p``) + - ``multiplicities`` -- bool (default: True): currently, only roots + without multiplicities are computed. + - ``algorithm`` -- the algorithm to use; currently, the only supported + algorithm is ``"Roth-Ruckenstein"`` + - ``degree_bound``-- if not ``None``, return only roots of degree at + most ``degree_bound`` + + EXAMPLES:: + + sage: R. = GF(13)[] + sage: S. = R[] + sage: p = y^2 + (12*x^2 + x + 11)*y + x^3 + 12*x^2 + 12*x + 1 + sage: p.roots(multiplicities=False) + [x^2 + 11*x + 1, x + 1] + sage: p.roots(multiplicities=False, degree_bound=1) + [x + 1] + """ + if multiplicities: + raise NotImplementedError("Use multiplicities=False") + +# Q = p.list() + + if degree_bound is None: + l = p.degree() + dl = p[l].degree() + degree_bound = max((p[i].degree() - dl)//(l - i) for i in range(l) if p[i]) + + if algorithm is None: + algorithm = "Roth-Ruckenstein" + + if algorithm == "Roth-Ruckenstein": + return self._roth_ruckenstein(p, degree_bound, None) + + elif algorithm == "Alekhnovich": + raise NotImplementedError("Alekhnovich algorithm is not implemented yet") + + else: + raise ValueError("unknown algorithm '{}'".format(algorithm)) + + class PolynomialRing_cdvr(PolynomialRing_integral_domain): r""" @@ -2484,7 +2656,7 @@ def irreducible_element(self, n, algorithm=None): x^33 + x^10 + 1 In degree 1:: - + sage: GF(97)['x'].irreducible_element(1) x + 96 sage: GF(97)['x'].irreducible_element(1, algorithm="conway")