Ticket #13215: trac_13215_skew_polynomials.patch

File trac_13215_skew_polynomials.patch, 254.4 KB (added by caruso, 7 years ago)
  • module_list.py

    # HG changeset patch
    # User Xavier Caruso <xavier.caruso@normalesup.org>
    # Date 1350912684 -7200
    # Node ID 055730e426275f1a961973a8e22fe55776088242
    # Parent  4977d9a6c35a4c78910a846057e861f3e17e248d
    imported patch skew_polynomial.patch
    
    diff --git a/module_list.py b/module_list.py
    a b  
    17951795    Extension('sage.rings.polynomial.symmetric_reduction',
    17961796              sources = ['sage/rings/polynomial/symmetric_reduction.pyx']),
    17971797
     1798    Extension('sage.rings.polynomial.skew_polynomial_element',
     1799              sources = ['sage/rings/polynomial/skew_polynomial_element.pyx'],
     1800              depends = numpy_depends),
     1801
     1802    Extension('sage.rings.polynomial.skew_polynomial_finite_field',
     1803              sources = ['sage/rings/polynomial/skew_polynomial_finite_field.pyx'],
     1804              depends = numpy_depends),
     1805
    17981806    ################################
    17991807    ##
    18001808    ## sage.schemes
  • sage/rings/polynomial/all.py

    diff --git a/sage/rings/polynomial/all.py b/sage/rings/polynomial/all.py
    a b  
    4545# Infinite Polynomial Rings
    4646from sage.rings.polynomial.infinite_polynomial_ring import InfinitePolynomialRing
    4747
     48# Skew Polynomial Rings
     49from sage.rings.polynomial.skew_polynomial_ring_constructor import SkewPolynomialRing
     50from sage.rings.polynomial.skew_polynomial_ring import is_SkewPolynomialRing
     51from sage.rings.polynomial.skew_polynomial_element import is_SkewPolynomial
    4852# Evaluation of cyclotomic polynomials
    4953from sage.rings.polynomial.cyclotomic import cyclotomic_value
  • new file sage/rings/polynomial/skew_polynomial_element.pxd

    diff --git a/sage/rings/polynomial/skew_polynomial_element.pxd b/sage/rings/polynomial/skew_polynomial_element.pxd
    new file mode 100644
    - +  
     1include "../../ext/interrupt.pxi"
     2include "../../ext/cdefs.pxi"
     3include '../../ext/stdsage.pxi'
     4
     5
     6from sage.rings.integer cimport Integer
     7
     8from sage.structure.element import Element, AlgebraElement
     9from sage.structure.element cimport Element, AlgebraElement, ModuleElement
     10from sage.structure.parent cimport Parent
     11from polynomial_compiled import CompiledPolynomialFunction
     12from polynomial_compiled cimport CompiledPolynomialFunction
     13
     14from sage.rings.morphism cimport RingHomomorphism
     15from sage.structure.element cimport RingElement
     16
     17from sage.rings.polynomial.polynomial_element cimport Polynomial_generic_dense
     18
     19
     20cdef class CenterSkewPolynomial_generic_dense(Polynomial_generic_dense):
     21    pass
     22
     23
     24cdef class SkewPolynomial(AlgebraElement):
     25    cdef _is_gen
     26
     27    cdef long _hash_c(self)
     28    cdef list _list_c(self)
     29    cdef SkewPolynomial _new_c(self,list coeffs,Parent P,char check=*)
     30    cpdef SkewPolynomial _new_constant_poly(self,RingElement a,Parent P,char check=*)
     31    cpdef ModuleElement _neg_(self)
     32
     33
     34
     35cdef class SkewPolynomial_generic_dense(SkewPolynomial):
     36    cdef list __coeffs
     37    cdef void __normalize(self)
     38    cpdef _pow_(self,right,modulus=*,side=*)
     39
     40    # Inplace functions
     41    cdef void _inplace_add(self, SkewPolynomial_generic_dense right)
     42    cdef void _inplace_sub(self, SkewPolynomial_generic_dense right)
     43    cdef void _inplace_rmul(self, SkewPolynomial_generic_dense right)
     44    cdef void _inplace_lmul(self, SkewPolynomial_generic_dense right)
     45    cdef void _inplace_pow(self, Py_ssize_t n)
     46
     47
     48cdef class SkewPolynomialBaseringInjection(RingHomomorphism):
     49    cdef RingElement _an_element
     50    cdef object _new_constant_poly_
  • new file sage/rings/polynomial/skew_polynomial_element.pyx

    diff --git a/sage/rings/polynomial/skew_polynomial_element.pyx b/sage/rings/polynomial/skew_polynomial_element.pyx
    new file mode 100644
    - +  
     1r"""
     2This module implements elements in skew polynomial rings.
     3
     4DEFINITION::
     5
     6Let `R` be a commutative ring equipped with an endomorphism `\sigma`.
     7
     8The skew polynomial ring over `(R, \sigma)` is the ring `S = `R[X,\sigma]`
     9is the usual ring of polynomials over `R` equipped with the skew
     10multiplication defined by the rule `X*a = \sigma(a)*X` for all `a`
     11in `R`.
     12
     13EXAMPLES::
     14
     15We illustrate some functionalities implemented in this class.
     16
     17We create the skew polynomial ring::
     18
     19    sage: R.<t> = ZZ[]
     20    sage: sigma = R.hom([t+1])
     21    sage: S.<x> = R['x',sigma]; S
     22    Skew Polynomial Ring in x over Univariate Polynomial Ring in t over Integer Ring twisted by t |--> t + 1
     23
     24and some elements in it::
     25
     26    sage: a = t + x + 1; a
     27    x + t + 1
     28    sage: b = S([t^2,t+1,1]); b
     29    x^2 + (t + 1)*x + t^2
     30    sage: c = S.random_element(degree=3,monic=True); c
     31    x^3 + (-95*t^2 + t + 2)*x^2 + (-t^2 + t)*x + 2*t - 8
     32
     33Ring operations are supported::
     34
     35    sage: a + b
     36    x^2 + (t + 2)*x + t^2 + t + 1
     37    sage: a - b
     38    -x^2 - t*x - t^2 + t + 1
     39
     40    sage: a * b
     41    x^3 + (2*t + 3)*x^2 + (2*t^2 + 4*t + 2)*x + t^3 + t^2
     42    sage: b * a
     43    x^3 + (2*t + 4)*x^2 + (2*t^2 + 3*t + 2)*x + t^3 + t^2
     44    sage: a * b == b * a
     45    False
     46
     47    sage: b^2
     48    x^4 + (2*t + 4)*x^3 + (3*t^2 + 7*t + 6)*x^2 + (2*t^3 + 4*t^2 + 3*t + 1)*x + t^4
     49    sage: b^2 == b*b
     50    True
     51
     52
     53Sage implements also some arithmetics over skew polynomial rings.
     54You will find below a short panorama.
     55
     56DEFINITION:
     57
     58Let `a` and `b` be two skew polynomials in the same ring `S`.
     59The *left (resp. right) euclidean division* of `a` by `b` is
     60a couple `(q,r)` of elements in `S` such that
     61
     62-  `a = q*b + r` (resp. `a = b*q + r`)
     63
     64-  the degree of `r` is less than the degree of `b`
     65
     66`q` (resp. `r`) is called the *quotient* (resp. the remainder)
     67of this euclidean division.
     68
     69PROPERTIES:
     70
     71Keeping the previous notations, if the leading coefficient of `b`
     72is a unit (e.g. if `b` is monic) then the quotient and the remainder
     73in the *right* euclidean division exist and are unique.
     74
     75The same result holds for the *left* euclidean division if in addition
     76the twist map defining the skew polynomial ring is invertible.
     77
     78EXAMPLES::
     79
     80    sage: q,r = c.quo_rem(b)   # default side is right
     81    sage: q
     82    x - 95*t^2
     83    sage: r
     84    (95*t^3 + 93*t^2 - t - 1)*x + 95*t^4 + 2*t - 8
     85    sage: c == q*b + r
     86    True
     87
     88The operators ``//`` and ``%`` give respectively the quotient
     89and the remainder of the right euclidean division::
     90
     91    sage: q == c // b
     92    True
     93    sage: r == c % b
     94    True
     95
     96If we want left euclidean division, we need to specify
     97``side=Left``. Nonetheless, it won't work over our current
     98`S` because Sage can't invert the twist map::
     99
     100    sage: q,r = c.quo_rem(b,side=Left)
     101    Traceback (most recent call last):
     102    ...
     103    NotImplementedError
     104
     105Here is a working example over a finite field::
     106
     107    sage: k.<t> = GF(5^3)
     108    sage: Frob = k.frobenius_endomorphism()
     109    sage: S.<x> = k['x',Frob]
     110    sage: a = x^4 + (4*t + 1)*x^3 + (t^2 + 3*t + 3)*x^2 + (3*t^2 + 2*t + 2)*x + 3*t^2 + 3*t + 1
     111    sage: b = (2*t^2 + 3)*x^2 + (3*t^2 + 1)*x + 4*t + 2
     112    sage: q,r = a.quo_rem(b,side=Left)
     113    sage: q
     114    (4*t^2 + t + 1)*x^2 + (2*t^2 + 2*t + 2)*x + 2*t^2 + 4*t + 3
     115    sage: r
     116    (t + 2)*x + 3*t^2 + 2*t + 4
     117    sage: a == b*q + r
     118    True
     119
     120
     121Once we have euclidean divisions, we have for free gcd and lcm
     122(at least if the base ring is a field).
     123This class provides an implementation of gcd and lcm::
     124
     125    sage: a = (x + t) * (x + t^2)^2
     126    sage: b = (x + t) * (t*x + t + 1) * (x + t^2)
     127    sage: a.gcd(b)  # default side is right
     128    x + t^2
     129    sage: a.gcd(b,side=Left)
     130    x + t
     131
     132For lcm, the default side is left but be very careful: by
     133convention, a left (resp. right) lcm is common multiple on
     134the right (resp. left)::
     135
     136    sage: c = a.lcm(b); c  # default side is left
     137    x^5 + (4*t^2 + t + 3)*x^4 + (3*t^2 + 4*t)*x^3 + 2*t^2*x^2 + (2*t^2 + t)*x + 4*t^2 + 4
     138    sage: c.is_divisible_by(a)
     139    True
     140    sage: c.is_divisible_by(b)
     141    True
     142
     143    sage: d = a.lcm(b,side=Right); d
     144    x^5 + (t^2 + 1)*x^4 + (3*t^2 + 3*t + 3)*x^3 + (3*t^2 + t + 2)*x^2 + (4*t^2 + 3*t)*x + 4*t + 4
     145    sage: d.is_divisible_by(a,side=Left)
     146    True
     147    sage: d.is_divisible_by(b,side=Left)
     148    True
     149
     150AUTHOR:
     151
     152-  Xavier Caruso (2012-06-29)
     153"""
     154
     155#############################################################################
     156#    Copyright (C) 2012 Xavier Caruso <xavier.caruso@normalesup.org>
     157#
     158#  Distributed under the terms of the GNU General Public License (GPL)
     159#
     160#                  http://www.gnu.org/licenses/
     161#****************************************************************************
     162
     163
     164include "../../ext/stdsage.pxi"
     165
     166import operator, copy, re
     167
     168import skew_polynomial_ring
     169import sage.rings.infinity as infinity
     170from sage.misc.latex import latex
     171from sage.structure.factorization import Factorization
     172
     173from sage.categories.homset import Hom
     174
     175from sage.structure.element import RingElement
     176from sage.structure.element cimport Element, RingElement, ModuleElement
     177
     178from sage.rings.ring import Field
     179
     180from sage.structure.parent_gens cimport ParentWithGens
     181
     182from sage.rings.integer cimport Integer
     183from sage.categories.map cimport Map
     184from sage.rings.morphism cimport RingHomomorphism
     185from sage.rings.polynomial.polynomial_element cimport Polynomial_generic_dense
     186
     187from sage.structure.side import Left, Right
     188
     189def is_SkewPolynomial(a):
     190    """
     191    Return True if `a` is a skew polynomial (over some base).
     192
     193    INPUT:
     194
     195    -  ``a`` -- an object
     196
     197    EXAMPLES::
     198
     199        sage: from sage.rings.polynomial.skew_polynomial_element import is_SkewPolynomial
     200        sage: R.<t> = ZZ[]
     201        sage: sigma = R.hom([t+1])
     202        sage: S.<x> = R['x',sigma]
     203        sage: a = x^2
     204        sage: is_SkewPolynomial(a)
     205        True
     206    """
     207    return PY_TYPE_CHECK(a, SkewPolynomial)
     208
     209
     210cdef class CenterSkewPolynomial_generic_dense(Polynomial_generic_dense):
     211    """
     212    A class for elements in the center of a skew polynomial ring.
     213    """
     214    pass
     215
     216
     217cdef class SkewPolynomial(AlgebraElement):
     218    """
     219    A skew polynomial.
     220    """
     221    def __init__(self,parent,is_gen=False,construct=False):
     222        """
     223        The following examples illustrate the creation of elements of
     224        skew polynomial rings.
     225
     226            sage: R.<t> = ZZ[]
     227            sage: sigma = R.hom([t+1])
     228            sage: S.<x> = R['x',sigma]
     229            sage: P = x+t; P
     230            x + t
     231            sage: Q = S([1,t,t+2]); Q
     232            (t + 2)*x^2 + t*x + 1
     233        """
     234        AlgebraElement.__init__(self, parent)
     235        self._is_gen = is_gen
     236
     237    # you may have to replicate this boilerplate code in derived classes if you override
     238    # __richcmp__.  The python documentation at  http://docs.python.org/api/type-structs.html
     239    # explains how __richcmp__, __hash__, and __cmp__ are tied together.
     240    def __hash__(self):
     241        return self._hash_c()
     242
     243    cdef long _hash_c(self):
     244        """
     245        This hash incorporates the name of the variable.
     246        """
     247        cdef long result = 0 # store it in a c-int and just let the overflowing additions wrap
     248        cdef long result_mon
     249        cdef long c_hash
     250        cdef long var_name_hash
     251        cdef int i
     252        for i from 0<= i <= self.degree():
     253            if i == 1:
     254                # we delay the hashing until now to not waste it one a constant poly
     255                var_name_hash = hash((<ParentWithGens>self._parent)._names[0])
     256            #  I'm assuming (incorrectly) that hashes of zero indicate that the element is 0.
     257            # This assumption is not true, but I think it is true enough for the purposes and it
     258            # it allows us to write fast code that omits terms with 0 coefficients.  This is
     259            # important if we want to maintain the '==' relationship with sparse polys.
     260            c_hash = hash(self[i])
     261            if c_hash != 0:
     262                if i == 0:
     263                    result += c_hash
     264                else:
     265                    # Hash (self[i], generator, i) as a tuple according to the algorithm.
     266                    result_mon = c_hash
     267                    result_mon = (1000003 * result_mon) ^ var_name_hash
     268                    result_mon = (1000003 * result_mon) ^ i
     269                    result += result_mon
     270        if result == -1:
     271            return -2
     272        return result
     273
     274
     275    # Comparison
     276    # ----------
     277
     278    def __richcmp__(left, right, int op):
     279        return (<Element>left)._richcmp(right,op)
     280
     281
     282    cdef int _cmp_c_impl(self, Element other) except -2:
     283        """
     284        Compare the two skew polynomials self and other.
     285
     286        We order polynomials first by degree, then in dictionary order
     287        starting with the coefficient of largest degree.
     288        """
     289        cdef list x = (<SkewPolynomial>self)._list_c()
     290        cdef list y = (<SkewPolynomial>other)._list_c()
     291        cdef Py_ssize_t dx = len(x), dy = len(y)
     292        cdef Py_ssize_t i
     293        cdef int c
     294        c = cmp(dx,dy)
     295        if c: return c
     296        for i from dx > i >= 0:
     297            c = cmp(x[i],y[i])
     298            if c: return c
     299        return 0
     300
     301
     302    # Some c functions
     303    # ----------------
     304
     305    cdef SkewPolynomial _new_c(self,list coeffs,Parent P,char check=0):
     306        """
     307        Fast creation of a new skew polynomial
     308
     309        .. NOTE::
     310
     311            Override this function in classes which inherit
     312            from SkewPolynomial.
     313        """
     314        return P(list)
     315
     316
     317    cpdef SkewPolynomial _new_constant_poly(self,RingElement a,Parent P,char check=0):
     318        """
     319        Fast creation of a new constant skew polynomial
     320        """
     321        if a:
     322            return self._new_c([a],P,check)
     323        else:
     324            return self._new_c([],P)
     325
     326
     327    cdef list _list_c(self):
     328        """
     329        Return the list of the underlying elements of this
     330        skew polynomial.
     331
     332        .. WARNING::
     333
     334            It is *not* a copy; do not modify this list!
     335        """
     336        raise NotImplementedError
     337
     338
     339    # A skew polynomial as a list of coefficients
     340    # -------------------------------------------
     341
     342    def list(self):
     343        """
     344        Return a new copy of the list of the underlying elements of self.
     345
     346        EXAMPLES::
     347
     348            sage: R.<t> = QQ[]
     349            sage: sigma = R.hom([t+1])
     350            sage: S.<x> = R['x',sigma]
     351            sage: a = 1 + x^4 + (t+1)*x^2 + t^2
     352            sage: l = a.list(); l
     353            [t^2 + 1, 0, t + 1, 0, 1]
     354
     355        Note that v is a list, it is mutable, and each call to the list
     356        method returns a new list::
     357
     358            sage: type(l)
     359            <type 'list'>
     360            sage: l[0] = 5
     361            sage: a.list()
     362            [t^2 + 1, 0, t + 1, 0, 1]
     363        """
     364        return list((<SkewPolynomial>self)._list_c())
     365
     366
     367    def __getitem__(self,n):
     368        """
     369        INPUT::
     370
     371        -  ``n`` -- an integer
     372
     373        OUTPUT::
     374
     375        -  the `n`-th coefficient of self
     376
     377        .. NOTE::
     378
     379            Coefficients are on the left (ie before the variable)
     380
     381        EXAMPLES::
     382
     383            sage: R.<t> = QQ[]
     384            sage: sigma = R.hom([t+1])
     385            sage: S.<x> = R['x',sigma]
     386            sage: a = t*x^2 + (t + 3/7)*x + t^2
     387            sage: a[1]
     388            t + 3/7
     389            sage: a[3]
     390            0
     391        """
     392        try:
     393            return (<SkewPolynomial>self)._list_c()[n]
     394        except IndexError:
     395            return self.base_ring()(0)
     396
     397
     398    def __getslice__(self, Py_ssize_t i, Py_ssize_t j):
     399        """
     400        EXAMPLES::
     401            sage: R.<t> = QQ[]
     402            sage: sigma = R.hom([t+1])
     403            sage: S.<x> = R['x',sigma]
     404            sage: a = t*x^2 + (t + 3/7)*x + t^2
     405            sage: a[1:]
     406            t*x^2 + (t + 3/7)*x
     407            sage: a[:1]
     408            t^2
     409            sage: a[3:]
     410            0
     411        """
     412        if i <= 0:
     413            i = 0
     414            zeros = []
     415        elif i > 0:
     416            zeros = [self._parent.base_ring()(0)] * i
     417        return self._new_c(zeros + self._list_c()[i:j], self._parent, 1)
     418
     419
     420    def __setitem__(self, n, value):
     421        """
     422        Set the `n`-th coefficient of this skew polynomial. This always
     423        raises an IndexError, since in Sage polynomials are immutable.
     424
     425        INPUT:
     426
     427        -  ``n`` - an integer
     428
     429        -  ``value`` - value to set the `n`-th coefficient to
     430
     431        OUTPUT: an IndexError is always raised.
     432
     433        EXAMPLES::
     434
     435            sage: k.<t> = GF(5^3)
     436            sage: Frob = k.frobenius_endomorphism()
     437            sage: S.<x> = k['x',Frob]
     438            sage: a = x + t
     439            sage: a[1] = t + 1
     440            Traceback (most recent call last):
     441            ...
     442            IndexError: skew polynomials are immutable
     443        """
     444        raise IndexError, "skew polynomials are immutable"
     445
     446
     447    # Basic operations
     448    # ----------------
     449
     450    def degree(self):
     451        """
     452        Return the degree of this skew polynomial. The zero
     453        skew polynomial has degree -1.
     454
     455        EXAMPLES::
     456
     457            sage: R.<t> = ZZ[]
     458            sage: sigma = R.hom([t+1])
     459            sage: S.<x> = R['x',sigma]
     460            sage: a = x^2 + t*x^3 + t^2*x + 1
     461            sage: a.degree()
     462            3
     463
     464        By convention, the degree of 0 is -1::
     465
     466            sage: S(0).degree()
     467            -1
     468        """
     469        return len((<SkewPolynomial>self)._list_c())-1
     470
     471
     472    def valuation(self):
     473        """
     474        Return the valuation of this skew polynomial.
     475        The zero skew polynomial has valuation +Infinity.
     476
     477        EXAMPLES::
     478
     479            sage: R.<t> = ZZ[]
     480            sage: sigma = R.hom([t+1])
     481            sage: S.<x> = R['x',sigma]
     482            sage: a = x^2 + t*x^3 + t^2*x
     483            sage: a.valuation()
     484            1
     485
     486        By convention, the valuation of 0 is +Infinity::
     487
     488            sage: S(0).valuation()
     489            +Infinity
     490        """
     491        cdef list x = (<SkewPolynomial>self)._list_c()
     492        if len(x) == 0:
     493            return infinity.infinity
     494        cdef Py_ssize_t v = 0
     495        while x[v].is_zero():
     496            v += 1
     497        return v
     498
     499
     500    cpdef ModuleElement _add_(self, ModuleElement right):
     501        """
     502        Compute self + right
     503
     504        INPUT:
     505
     506        -  right -- a skew polynomial over the same base
     507
     508        TESTS::
     509
     510            sage: R.<t> = QQ[]
     511            sage: sigma = R.hom([t+1])
     512            sage: S.<x> = R['x',sigma]
     513            sage: a = S.random_element(monic=True); a
     514            x^2 + (-12*t^2 + 1/2*t - 1/95)*x - 1/2*t^2 - 4
     515            sage: b = -S.random_element(monic=True); b
     516            -x^2 + (5/2*t - 2/3)*x + 1/4*t^2 - 1/2*t + 1
     517            sage: c = a+b; c
     518            (-12*t^2 + 3*t - 193/285)*x - 1/4*t^2 - 1/2*t - 3
     519            sage: c.degree()
     520            1
     521        """
     522        cdef Py_ssize_t i, min
     523        cdef list x = (<SkewPolynomial>self)._list_c()
     524        cdef list y = (<SkewPolynomial>right)._list_c()
     525        cdef Py_ssize_t dx = len(x), dy = len(y)
     526
     527        if dx > dy:
     528            return self._new_c([x[i] + y[i] for i from 0 <= i < dy] + x[dy:], self._parent, 0)
     529        elif dx < dy:
     530            return self._new_c([x[i] + y[i] for i from 0 <= i < dx] + y[dx:], self._parent, 0)
     531        else:
     532            return self._new_c([x[i] + y[i] for i from 0 <= i < dx], self._parent, 1)
     533
     534
     535    cpdef ModuleElement _sub_(self, ModuleElement right):
     536        """
     537        Compute self - right
     538
     539        INPUT:
     540
     541        -  right -- a skew polynomial over the same base
     542
     543        TESTS::
     544
     545            sage: R.<t> = QQ[]
     546            sage: sigma = R.hom([t+1])
     547            sage: S.<x> = R['x',sigma]
     548            sage: a = S.random_element(monic=True); a
     549            x^2 + (-12*t^2 + 1/2*t - 1/95)*x - 1/2*t^2 - 4
     550            sage: b = S.random_element(monic=True); b
     551            x^2 + (-5/2*t + 2/3)*x - 1/4*t^2 + 1/2*t - 1
     552            sage: c = a-b; c
     553            (-12*t^2 + 3*t - 193/285)*x - 1/4*t^2 - 1/2*t - 3
     554            sage: c.degree()
     555            1
     556        """
     557        cdef Py_ssize_t i, min
     558        cdef list x = (<SkewPolynomial>self)._list_c()
     559        cdef list y = (<SkewPolynomial>right)._list_c()
     560        cdef Py_ssize_t dx = len(x), dy = len(y)
     561        cdef RingElement c
     562
     563        if dx > dy:
     564            return self._new_c([x[i] - y[i] for i from 0 <= i < dy] + x[dy:], self._parent, 0)
     565        elif dx < dy:
     566            return self._new_c([x[i] - y[i] for i from 0 <= i < dx] + [ -c for c in y[dx:] ], self._parent, 0)
     567        else:
     568            return self._new_c([x[i] - y[i] for i from 0 <= i < dx], self._parent, 1)
     569
     570
     571    cpdef ModuleElement _neg_(self):
     572        """
     573        Return the opposite of self
     574
     575        TESTS::
     576
     577            sage: R.<t> = QQ[]
     578            sage: sigma = R.hom([t+1])
     579            sage: S.<x> = R['x',sigma]
     580            sage: a = t*x^2 + x - 3
     581            sage: -a
     582            -t*x^2 - x + 3
     583        """
     584        return self._new_c([-x for x in (<SkewPolynomial>self)._list_c()], self._parent, 0)
     585
     586
     587    cpdef ModuleElement _lmul_(self, RingElement right):
     588        """
     589        Compute self * right (in this order)
     590
     591        INPUT:
     592
     593        -  right -- an element of the base ring
     594
     595        TESTS::
     596
     597            sage: R.<t> = QQ[]
     598            sage: sigma = R.hom([t+1])
     599            sage: S.<x> = R['x',sigma]
     600            sage: a = x + t
     601            sage: b = t
     602            sage: a * b
     603            (t + 1)*x + t^2
     604            sage: a * b == b * a
     605            False
     606        """
     607        if right == 0:
     608            return self._parent(0)
     609        cdef list x = (<SkewPolynomial>self)._list_c()
     610        cdef Py_ssize_t i
     611        map = self._parent._map
     612        return self._new_c([ (map**i)(right)*x[i] for i from 0 <= i < len(x) ], self._parent, 0)
     613
     614
     615    cpdef ModuleElement _rmul_(self, RingElement left):
     616        """
     617        Compute left * self (in this order)
     618
     619        INPUT:
     620
     621        -  left -- an element of the base ring
     622
     623        TESTS::
     624
     625            sage: R.<t> = QQ[]
     626            sage: sigma = R.hom([t+1])
     627            sage: S.<x> = R['x',sigma]
     628            sage: a = t
     629            sage: b = x + t
     630            sage: a * b
     631            t*x + t^2
     632            sage: a * b == b * a
     633            False
     634        """
     635        if left == 0:
     636            return self.parent().zero_element()
     637        cdef list x = (<SkewPolynomial>self)._list_c()
     638        cdef Py_ssize_t i
     639        map = self._parent._map
     640        return self._new_c([ left*x[i] for i from 0 <= i < len(x) ], self._parent, 0)
     641
     642
     643    cpdef RingElement _mul_(self, RingElement right):
     644        """
     645        Compute self * right (in this order)
     646
     647        INPUT:
     648
     649        -  right -- a skew polynomial in the same ring
     650
     651        TESTS::
     652
     653            sage: R.<t> = QQ[]
     654            sage: sigma = R.hom([t+1])
     655            sage: S.<x> = R['x',sigma]
     656            sage: a = x^2 + t; a
     657            x^2 + t
     658            sage: b = x^2 + (t + 1)*x; b
     659            x^2 + (t + 1)*x
     660            sage: a * b
     661            x^4 + (t + 3)*x^3 + t*x^2 + (t^2 + t)*x
     662            sage: a * b == b * a
     663            False
     664        """
     665        cdef list x = (<SkewPolynomial>self)._list_c()
     666        cdef list y = (<SkewPolynomial>right)._list_c()
     667        cdef Py_ssize_t i, k, start, end
     668        cdef Py_ssize_t dx = len(x)-1, dy = len(y)-1
     669        parent = self._parent
     670        if dx == -1:
     671            return self
     672        elif dy == -1:
     673            return right
     674        elif dx == 0:
     675            c = x[0]
     676            return self._new_c([c*a for a in y], parent, 0)
     677        cdef list coeffs = []
     678        for k from 0 <= k <= dx+dy:
     679            start = 0 if k <= dy else k-dy # max(0, k-dy)
     680            end = k if k <= dx else dx # min(k, dx)
     681            sum = x[start] * parent.twist_map(start)(y[k-start])
     682            for i from start < i <= end:
     683                sum += x[i] * parent.twist_map(i)(y[k-i])
     684            coeffs.append(sum)
     685        return self._new_c(coeffs, parent, 0)
     686
     687
     688    def square(self):
     689        """
     690        Return the square of self
     691
     692        TESTS::
     693
     694            sage: R.<t> = QQ[]
     695            sage: sigma = R.hom([t+1])
     696            sage: S.<x> = R['x',sigma]
     697            sage: a = x + t; a
     698            x + t
     699            sage: a.square()
     700            x^2 + (2*t + 1)*x + t^2
     701            sage: a.square() == a*a
     702            True
     703        """
     704        return self * self
     705
     706
     707    def conjugate(self,n):
     708        """
     709        INPUT:
     710
     711        -  `n` -- an integer
     712
     713        OUTPUT:
     714
     715        -  this skew polynomial conjugated by x^n (where x is
     716           the variable)
     717
     718        .. NOTE::
     719
     720            This conjugate is obtained from the skew polynomial by
     721            applying the `n`-th iterate of the twist map to each of
     722            its coefficients.
     723
     724        EXAMPLES::
     725
     726            sage: R.<t> = QQ[]
     727            sage: sigma = R.hom([t+1])
     728            sage: S.<x> = R['x',sigma]
     729            sage: a = t*x^3 + (t^2 + 1)*x^2 + 2*t
     730            sage: b = a.conjugate(2); b
     731            (t + 2)*x^3 + (t^2 + 4*t + 5)*x^2 + 2*t + 4
     732            sage: x^2*a == b*x^2
     733            True
     734
     735        In principle, negative values for `n` are allowed... but Sage
     736        needs to be able to invert the twist map::
     737
     738            sage: b = a.conjugate(-1)
     739            Traceback (most recent call last):
     740            ...
     741            NotImplementedError
     742
     743        Here is a working example::
     744
     745            sage: k.<t> = GF(5^3)
     746            sage: Frob = k.frobenius_endomorphism()
     747            sage: T.<y> = k['y',Frob]
     748            sage: u = T.random_element(); u
     749            (2*t^2 + 3)*y^2 + (4*t^2 + t + 4)*y + 2*t^2 + 2
     750            sage: v = u.conjugate(-1); v
     751            (3*t^2 + t)*y^2 + (4*t^2 + 2*t + 4)*y + 3*t^2 + t + 4
     752            sage: u*y == y*v
     753            True
     754        """
     755        return self._new_c([ self._parent.twist_map(n)(x) for x in (<SkewPolynomial>self)._list_c() ], self._parent, 0)
     756
     757
     758    # Other useful mathematical functions
     759    # -----------------------------------
     760
     761    def constant_coefficient(self):
     762        """
     763        Return the constant coefficient (i.e. the coefficient of degree
     764        `0`) of this skew polynomial.
     765
     766        OUTPUT: element of base ring
     767
     768        EXAMPLES::
     769
     770            sage: R.<t> = ZZ[]
     771            sage: sigma = R.hom([t+1])
     772            sage: S.<x> = R['x',sigma]
     773            sage: a = x + t^2 + 2
     774            sage: a.constant_coefficient()
     775            t^2 + 2
     776        """
     777        cdef x = (<SkewPolynomial>self)._list_c()
     778        if len(x) == 0:
     779            return self.base_ring()(0)
     780        else:
     781            return x[0]
     782
     783
     784    def leading_coefficient(self):
     785        """
     786        Return the leading coefficient of this skew polynomial.
     787
     788        OUTPUT: element of the base ring
     789
     790        EXAMPLES::
     791
     792            sage: R.<t> = ZZ[]
     793            sage: sigma = R.hom([t+1])
     794            sage: S.<x> = R['x',sigma]
     795            sage: a = x + (t+1)*x^5 + t^2*x^3 - x^5
     796            sage: a.leading_coefficient()
     797            t
     798        """
     799        cdef x = (<SkewPolynomial>self)._list_c()
     800        if len(x) == 0:
     801            raise ValueError("")
     802        return x[-1]
     803
     804
     805    def __nonzero__(self):
     806        return self.degree() >= 0
     807
     808
     809    def is_unit(self):
     810        """
     811        Return True if this skew polynomial is a unit.
     812
     813        Not yet implemented.
     814        """
     815        raise NotImplementedError
     816
     817
     818    def is_monic(self):
     819        """
     820        Returns True if this skew polynomial is monic. The zero polynomial
     821        is by definition not monic.
     822
     823        EXAMPLES::
     824
     825            sage: R.<t> = ZZ[]
     826            sage: sigma = R.hom([t+1])
     827            sage: S.<x> = R['x',sigma]
     828            sage: a = x + t
     829            sage: a.is_monic()
     830            True
     831            sage: a = 0*x
     832            sage: a.is_monic()
     833            False
     834            sage: a = t*x^3 + x^4 + (t+1)*x^2
     835            sage: a.is_monic()
     836            True
     837            sage: a = (t^2 + 2*t)*x^2 + x^3 + t^10*x^5
     838            sage: a.is_monic()
     839            False
     840        """
     841        return not self.is_zero() and self[self.degree()] == 1
     842
     843
     844    def lmonic(self):
     845        """
     846        Return the unique monic skew polynomial `a` of the same
     847        degree which divides this skew polynomial on the left
     848
     849        .. Note::
     850
     851            This skew polynomial is self dividing on the
     852            *right* by the `n`-th iterative (`n` is the degree of
     853            self) of the inverse of the twist map applied to the
     854            leading coefficient.
     855
     856        EXAMPLES::
     857
     858            sage: k.<t> = GF(5^3)
     859            sage: Frob = k.frobenius_endomorphism()
     860            sage: S.<x> = k['x',Frob]
     861            sage: a = (3*t^2 + 3*t + 2)*x^3 + (2*t^2 + 3)*x^2 + (4*t^2 + t + 4)*x + 2*t^2 + 2
     862            sage: b = a.lmonic(); b
     863            x^3 + (4*t^2 + 3*t)*x^2 + (4*t + 2)*x + 2*t^2 + 4*t + 3
     864
     865        Check list::
     866
     867            sage: b.degree() == a.degree()
     868            True
     869            sage: b.is_divisible_by(a,side=Left)
     870            True
     871            sage: twist = S.twist_map(-a.degree())
     872            sage: a == b * twist(a.leading_coefficient())
     873            True
     874
     875        Note that `b` does not divise `a` on the right::
     876
     877            sage: b.is_divisible_by(a,side=Right)
     878            False
     879
     880        This function does not work if the leading coefficient is not a
     881        unit::
     882
     883            sage: R.<t> = QQ[]
     884            sage: sigma = R.hom([t+1])
     885            sage: S.<x> = R['x',sigma]
     886            sage: a = t*x
     887            sage: a.lmonic()
     888            Traceback (most recent call last):
     889            ...
     890            NotImplementedError: the leading coefficient is not a unit
     891        """
     892        try:
     893            a = self.base_ring()(~self.leading_coefficient())
     894        except (ZeroDivisionError, TypeError):
     895            raise NotImplementedError("the leading coefficient is not a unit")
     896        return self*self._parent.twist_map(-self.degree())(a)
     897
     898
     899    def rmonic(self):
     900        """
     901        Return the unique monic skew polynomial `a` of the same
     902        degree which divides this skew polynomial on the right
     903
     904        .. Note::
     905
     906            This skew polynomial is self dividing on the *left*
     907            by its leading coefficient.
     908
     909        EXAMPLES::
     910
     911            sage: k.<t> = GF(5^3)
     912            sage: Frob = k.frobenius_endomorphism()
     913            sage: S.<x> = k['x',Frob]
     914            sage: a = (3*t^2 + 3*t + 2)*x^3 + (2*t^2 + 3)*x^2 + (4*t^2 + t + 4)*x + 2*t^2 + 2
     915            sage: b = a.rmonic(); b
     916            x^3 + (2*t^2 + 3*t + 4)*x^2 + (3*t^2 + 4*t + 1)*x + 2*t^2 + 4*t + 3
     917
     918        Check list::
     919
     920            sage: b.degree() == a.degree()
     921            True
     922            sage: b.is_divisible_by(a,side=Right)
     923            True
     924            sage: a == a.leading_coefficient() * b
     925            True
     926
     927        Note that `b` does not divise `a` on the right::
     928
     929            sage: b.is_divisible_by(a,side=Left)
     930            False
     931
     932        This function does not work if the leading coefficient is not a
     933        unit::
     934
     935            sage: R.<t> = QQ[]
     936            sage: sigma = R.hom([t+1])
     937            sage: S.<x> = R['x',sigma]
     938            sage: a = t*x
     939            sage: a.rmonic()
     940            Traceback (most recent call last):
     941            ...
     942            NotImplementedError: the leading coefficient is not a unit
     943        """
     944        try:
     945            a = self.base_ring()(~self.leading_coefficient())
     946        except (ZeroDivisionError, TypeError):
     947            raise NotImplementedError("the leading coefficient is not a unit")
     948        return a*self
     949
     950
     951    def monic(self, side=Right):
     952        """
     953        INPUT:
     954
     955        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     956
     957        OUTPUT:
     958
     959        -  the unique monic skew polynomial `a` of the same
     960           degree which divides this skew polynomial on ``side``
     961
     962        .. Note::
     963
     964            if ``side`` is ``Right``, this skew polynomial is self
     965            dividing on the *left* by its leading coefficient; if
     966            ``side`` is ``Left``, it is self dividing on the
     967            *right* by the `n`-th iterative (`n` is the degree of
     968            self) of the inverse of the twist map applied to the
     969            leading coefficient.
     970
     971        EXAMPLES::
     972
     973            sage: k.<t> = GF(5^3)
     974            sage: Frob = k.frobenius_endomorphism()
     975            sage: S.<x> = k['x',Frob]
     976            sage: a = (3*t^2 + 3*t + 2)*x^3 + (2*t^2 + 3)*x^2 + (4*t^2 + t + 4)*x + 2*t^2 + 2
     977            sage: b = a.monic(); b
     978            x^3 + (2*t^2 + 3*t + 4)*x^2 + (3*t^2 + 4*t + 1)*x + 2*t^2 + 4*t + 3
     979
     980        Check list::
     981
     982            sage: b.degree() == a.degree()
     983            True
     984            sage: b.is_divisible_by(a,side=Right)
     985            True
     986            sage: a == a.leading_coefficient() * b
     987            True
     988
     989        Note that `b` does not divise `a` on the left::
     990
     991            sage: b.is_divisible_by(a,side=Left)
     992            False
     993
     994        The same on the left::
     995
     996            sage: b = a.monic(side=Left); b
     997            x^3 + (4*t^2 + 3*t)*x^2 + (4*t + 2)*x + 2*t^2 + 4*t + 3
     998            sage: b.degree() == a.degree()
     999            True
     1000            sage: b.is_divisible_by(a,side=Left)
     1001            True
     1002            sage: twist = S.twist_map(-a.degree())
     1003            sage: a == b * twist(a.leading_coefficient())
     1004            True
     1005
     1006        This function does not work if the leading coefficient is not a
     1007        unit::
     1008
     1009            sage: R.<t> = QQ[]
     1010            sage: sigma = R.hom([t+1])
     1011            sage: S.<x> = R['x',sigma]
     1012            sage: a = t*x
     1013            sage: a.monic()
     1014            Traceback (most recent call last):
     1015            ...
     1016            NotImplementedError: the leading coefficient is not a unit
     1017        """
     1018        if self.is_monic():
     1019            return self
     1020        if side == Right:
     1021            return self.rmonic()
     1022        else:
     1023            return self.lmonic()
     1024
     1025
     1026    # Divisibility
     1027    # ------------
     1028
     1029
     1030    def lquo_rem(self, other):
     1031        """
     1032        INPUT:
     1033
     1034        -  ``other`` -- a skew polynomial ring over the same
     1035           base ring
     1036
     1037        OUTPUT:
     1038
     1039        -  the quotient and the remainder of the left euclidean
     1040           division of this skew polynomial by ``other``
     1041
     1042        .. NOTE::
     1043
     1044            Doesn't work if the leading coefficient of the divisor
     1045            is not a unit or if Sage can't invert the twist map.
     1046
     1047        EXAMPLES::
     1048
     1049            sage: k.<t> = GF(5^3)
     1050            sage: Frob = k.frobenius_endomorphism()
     1051            sage: S.<x> = k['x',Frob]
     1052            sage: a = (3*t^2 + 3*t + 2)*x^3 + (2*t^2 + 3)*x^2 + (4*t^2 + t + 4)*x + 2*t^2 + 2
     1053            sage: b = (3*t^2 + 4*t + 2)*x^2 + (2*t^2 + 4*t + 3)*x + 2*t^2 + t + 1
     1054            sage: q,r = a.lquo_rem(b)
     1055            sage: a == b*q + r
     1056            True
     1057
     1058        In the following example, Sage does not know the inverse
     1059        of the twist map::
     1060
     1061            sage: R.<t> = ZZ[]
     1062            sage: sigma = R.hom([t+1])
     1063            sage: S.<x> = R['x',sigma]
     1064            sage: a = (-2*t^2 - t + 1)*x^3 + (-t^2 + t)*x^2 + (-12*t - 2)*x - t^2 - 95*t + 1
     1065            sage: b = x^2 + (5*t - 6)*x - 4*t^2 + 4*t - 1
     1066            sage: a.lquo_rem(b)
     1067            Traceback (most recent call last):
     1068            ...
     1069            NotImplementedError
     1070        """
     1071        cdef list a = list((<SkewPolynomial>self)._list_c())
     1072        cdef list b = (<SkewPolynomial>other)._list_c()
     1073        cdef Py_ssize_t i, j
     1074        cdef Py_ssize_t da = self.degree(), db = other.degree()
     1075        if db < 0:
     1076            raise ZeroDivisionError
     1077        if da < db:
     1078            return self._new_c([],self._parent), self
     1079        try:
     1080            inv = self.base_ring()(~b[db])
     1081        except (ZeroDivisionError, TypeError):
     1082            raise NotImplementedError("the leading coefficient of the divisor is not invertible")
     1083        cdef list q = [ ]
     1084        parent = self._parent
     1085        for i from da-db >= i >= 0:
     1086            c = parent.twist_map(-db)(inv*a[i+db])
     1087            for j from 0 <= j < db:
     1088                a[i+j] -= b[j] * parent.twist_map(j)(c)
     1089            q.append(c)
     1090        q.reverse()
     1091        return self._new_c(q,parent), self._new_c(a[:db],parent,1)
     1092
     1093
     1094    def rquo_rem(self, other):
     1095        """
     1096        INPUT:
     1097
     1098        -  ``other`` -- a skew polynomial ring over the same
     1099           base ring
     1100
     1101        OUTPUT:
     1102
     1103        -  the quotient and the remainder of the left euclidean
     1104           division of this skew polynomial by ``other``
     1105
     1106        .. NOTE::
     1107
     1108            Doesn't work if the leading coefficient of the divisor
     1109            is not a unit.
     1110
     1111        EXAMPLES::
     1112
     1113            sage: R.<t> = ZZ[]
     1114            sage: sigma = R.hom([t+1])
     1115            sage: S.<x> = R['x',sigma]
     1116            sage: a = S.random_element(degree=4); a
     1117            t^2*x^4 + (-12*t^2 - 2*t - 1)*x^3 + (-95*t^2 + t + 2)*x^2 + (-t^2 + t)*x + 2*t - 8
     1118            sage: b = S.random_element(monic=True); b
     1119            x^2 + (4*t^2 - t - 2)*x - t^2 + t - 1
     1120            sage: q,r = a.rquo_rem(b)
     1121            sage: a == q*b + r
     1122            True
     1123
     1124        The leading coefficient of the divisor need to be invertible::
     1125
     1126            sage: c = S.random_element(); c
     1127            (-4*t^2 + t)*x^2 - 2*t^2*x + 5*t^2 - 6*t - 4
     1128            sage: a.rquo_rem(c)
     1129            Traceback (most recent call last):
     1130            ...
     1131            NotImplementedError: the leading coefficient of the divisor is not invertible
     1132        """
     1133        cdef list a = list((<SkewPolynomial>self)._list_c())
     1134        cdef list b = (<SkewPolynomial>other)._list_c()
     1135        cdef Py_ssize_t i, j
     1136        cdef Py_ssize_t da = self.degree(), db = other.degree()
     1137        parent = self._parent
     1138        if db < 0:
     1139            raise ZeroDivisionError
     1140        if da < db:
     1141            return self._new_c([],parent), self
     1142        try:
     1143            inv = self.base_ring()(~b[db])
     1144        except (ZeroDivisionError, TypeError):
     1145            raise NotImplementedError("the leading coefficient of the divisor is not invertible")
     1146        cdef list q = [ ]
     1147        parent = self._parent
     1148        for i from da-db >= i >= 0:
     1149            c = parent.twist_map(i)(inv) * a[i+db]
     1150            for j from 0 <= j < db:
     1151                a[i+j] -= c * parent.twist_map(i)(b[j])
     1152            q.append(c)
     1153        q.reverse()
     1154        return self._new_c(q,parent), self._new_c(a[:db],parent,1)
     1155
     1156
     1157    def quo_rem(self,other,side=Right):
     1158        r"""
     1159        INPUT:
     1160
     1161        -  ``other`` -- a skew polynomial ring over the same
     1162           base ring
     1163
     1164        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1165
     1166        OUTPUT:
     1167
     1168        -  the quotient and the remainder of the euclidean
     1169           division of this skew polynomial by ``other``
     1170
     1171        .. NOTE::
     1172
     1173            Doesn't work if the leading coefficient of the divisor
     1174            is not a unit.
     1175
     1176        EXAMPLES::
     1177
     1178            sage: R.<t> = ZZ[]
     1179            sage: sigma = R.hom([t+1])
     1180            sage: S.<x> = R['x',sigma]
     1181            sage: a = S.random_element(degree=4); a
     1182            t^2*x^4 + (-12*t^2 - 2*t - 1)*x^3 + (-95*t^2 + t + 2)*x^2 + (-t^2 + t)*x + 2*t - 8
     1183            sage: b = S.random_element(monic=True); b
     1184            x^2 + (4*t^2 - t - 2)*x - t^2 + t - 1
     1185            sage: q,r = a.quo_rem(b)  # right euclidean division
     1186            sage: a == q*b + r
     1187            True
     1188
     1189        The left euclidean division doesn't work over this `S` because
     1190        Sage cannot invert `\sigma`::
     1191
     1192            sage: q,r = a.quo_rem(b,side=Left)
     1193            Traceback (most recent call last):
     1194            ...
     1195            NotImplementedError
     1196
     1197        In any case, the leading coefficient of the divisor need to be
     1198        invertible::
     1199
     1200            sage: c = S.random_element(); c
     1201            (-4*t^2 + t)*x^2 - 2*t^2*x + 5*t^2 - 6*t - 4
     1202            sage: a.quo_rem(c)  # right euclidean division
     1203            Traceback (most recent call last):
     1204            ...
     1205            NotImplementedError: the leading coefficient of the divisor is not invertible
     1206
     1207        When the base ring is a finite field, everything works fine::
     1208
     1209            sage: k.<t> = GF(5^3)
     1210            sage: Frob = k.frobenius_endomorphism()
     1211            sage: S.<x> = k['x',Frob]
     1212            sage: a = 4*x^3 + (t^2 + 2*t + 4)*x^2 + (4*t^2 + 2*t + 4)*x + t^2 + 2
     1213            sage: b = (3*t + 2)*x^2 + 2*t*x + 3*t^2
     1214            sage: q,r = a.quo_rem(b)  # right euclidean division
     1215            sage: a == q*b + r
     1216            True
     1217            sage: q,r = a.quo_rem(b,side=Left)
     1218            sage: a == b*q + r
     1219            True
     1220        """
     1221        if side == Right:
     1222            return self.rquo_rem(other)
     1223        else:
     1224            return self.lquo_rem(other)
     1225
     1226
     1227    def rem(self,other,side=Right):
     1228        r"""
     1229        INPUT:
     1230
     1231        -  ``other`` -- a skew polynomial ring over the same
     1232           base ring
     1233
     1234        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1235
     1236        OUTPUT:
     1237
     1238        -  the remainder of the euclidean division of this
     1239           skew polynomial by ``other``
     1240
     1241        .. NOTE::
     1242
     1243            Doesn't work if the leading coefficient of the divisor
     1244            is not a unit.
     1245
     1246        EXAMPLES::
     1247
     1248            sage: R.<t> = ZZ[]
     1249            sage: sigma = R.hom([t+1])
     1250            sage: S.<x> = R['x',sigma]
     1251            sage: b = x^2 + 2*t*x + 2
     1252            sage: a = (x+t)*b + t*x + 1
     1253            sage: a.rem(b)   # right euclidean division
     1254            t*x + 1
     1255
     1256        The left euclidean division doesn't work over this `S` because
     1257        Sage cannot invert `\sigma`::
     1258
     1259            sage: a.rem(b,side=Left)
     1260            Traceback (most recent call last):
     1261            ...
     1262            NotImplementedError
     1263
     1264        In any case, the leading coefficient of the divisor need to be
     1265        invertible::
     1266
     1267            sage: (a*t).rem(b*t)
     1268            Traceback (most recent call last):
     1269            ...
     1270            NotImplementedError: the leading coefficient of the divisor is not invertible
     1271
     1272        When the base ring is a finite field, everything works fine::
     1273
     1274            sage: k.<t> = GF(5^3)
     1275            sage: Frob = k.frobenius_endomorphism()
     1276            sage: S.<x> = k['x',Frob]
     1277            sage: b = x^2 + 2*t*x + 2
     1278            sage: a = (x+t)*b + t*x + 1
     1279            sage: a.rem(b)  # right euclidean division
     1280            t*x + 1
     1281            sage: a.rem(b,side=Left)
     1282            (3*t^2 + 2)*x + 2*t^2
     1283        """
     1284        _,r = self.quo_rem(other,side=side)
     1285        return r
     1286
     1287
     1288    def __mod__(self,other):
     1289        """
     1290        Return the remainder in the *right* euclidean division of
     1291        this skew polynomial by ``other``
     1292
     1293        TESTS::
     1294
     1295            sage: R.<t> = ZZ[]
     1296            sage: sigma = R.hom([t+1])
     1297            sage: S.<x> = R['x',sigma]
     1298            sage: b = x^2 + 2*t*x + 2
     1299            sage: a = (x+t)*b + t*x + 1
     1300            sage: a % b
     1301            t*x + 1
     1302
     1303            sage: (a*t).rem(b*t)
     1304            Traceback (most recent call last):
     1305            ...
     1306            NotImplementedError: the leading coefficient of the divisor is not invertible
     1307        """
     1308        _,r = self.rquo_rem(other)
     1309        return r
     1310
     1311
     1312    def __floordiv__(self,right):
     1313        """
     1314        Return the quotient of the right euclidean division of self by right.
     1315
     1316        The algorithm fails if the leading coefficient of the divisor (right)
     1317        is not invertible.
     1318
     1319        .. SEEALSO::
     1320
     1321
     1322        TESTS::
     1323
     1324            sage: R.<t> = QQ[]
     1325            sage: sigma = R.hom([t+1])
     1326            sage: S.<x> = R['x',sigma]
     1327            sage: b = x^2 + t
     1328            sage: a = (x^2 + t*x + 1)*b + t^3*x
     1329            sage: a // b
     1330            x^2 + t*x + 1
     1331
     1332            sage: (t*a) // (t*b)
     1333            Traceback (most recent call last):
     1334            ...
     1335            NotImplementedError: the leading coefficient of the divisor is not invertible
     1336
     1337        """
     1338        q,_ = self.rquo_rem(right)
     1339        return q
     1340
     1341
     1342    cpdef RingElement _div_(self, RingElement right):
     1343        """
     1344        Not Implemented (since localization of Ore rings is
     1345        not yet implemented, see trac #13215).
     1346
     1347        Use the operator `//` even for exact division.
     1348
     1349        EXAMPLES::
     1350
     1351            sage: R.<t> = QQ[]
     1352            sage: sigma = R.hom([t+1])
     1353            sage: S.<x> = R['x',sigma]
     1354            sage: a = x^5 + (t + 2)*x^2 + t^2
     1355            sage: b = x^3 + 4*t
     1356            sage: c = a*b
     1357
     1358            sage: c / b
     1359            Traceback (most recent call last):
     1360            ...
     1361            NotImplementedError: Please use `//` (even for exact division)
     1362
     1363            sage: c // b == a
     1364            True
     1365        """
     1366        raise NotImplementedError("Please use `//` (even for exact division)")
     1367
     1368
     1369    def is_divisible_by(self,other,side=Right):
     1370        """
     1371        INPUT:
     1372
     1373        -  ``other`` -- a skew polynomial over the same base
     1374
     1375        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1376
     1377        OUTPUT:
     1378
     1379        Return True iff self is divisible by other on ``side``
     1380
     1381        EXAMPLES::
     1382
     1383            sage: k.<t> = GF(5^3)
     1384            sage: Frob = k.frobenius_endomorphism()
     1385            sage: S.<x> = k['x',Frob]
     1386            sage: a = x^2 + t*x + t^2 + 3
     1387            sage: b = x^3 + (t + 1)*x^2 + 1
     1388            sage: c = a*b
     1389            sage: c.is_divisible_by(a)   # on the right
     1390            False
     1391            sage: c.is_divisible_by(b)
     1392            True
     1393            sage: c.is_divisible_by(a,side=Left)
     1394            True
     1395            sage: c.is_divisible_by(b,side=Left)
     1396            False
     1397
     1398        Divisibility by 0 has no sense::
     1399
     1400            sage: c.is_divisible_by(S(0))
     1401            Traceback (most recent call last):
     1402            ...
     1403            ZeroDivisionError
     1404
     1405        This function does not work if the leading coefficient of the divisor
     1406        is not a unit::
     1407
     1408            sage: R.<t> = QQ[]
     1409            sage: sigma = R.hom([t+1])
     1410            sage: S.<x> = R['x',sigma]
     1411            sage: a = x^2 + 2*x + t
     1412            sage: b = (t+1)*x + t^2
     1413            sage: c = a*b
     1414            sage: c.is_divisible_by(b)
     1415            Traceback (most recent call last):
     1416            ...
     1417            NotImplementedError: the leading coefficient of the divisor is not invertible
     1418        """
     1419        return self.rem(other,side=side).is_zero()
     1420
     1421
     1422    def divides(self,other,side=Right):
     1423        """
     1424        INPUT:
     1425
     1426        -  ``other`` -- a skew polynomial over the same base
     1427
     1428        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1429
     1430        OUTPUT:
     1431
     1432        Return True iff self divides other on ``side``
     1433
     1434        EXAMPLES::
     1435
     1436            sage: k.<t> = GF(5^3)
     1437            sage: Frob = k.frobenius_endomorphism()
     1438            sage: S.<x> = k['x',Frob]
     1439            sage: a = x^2 + t*x + t^2 + 3
     1440            sage: b = x^3 + (t + 1)*x^2 + 1
     1441            sage: c = a*b
     1442            sage: a.divides(c)   # on the right
     1443            False
     1444            sage: b.divides(c)
     1445            True
     1446            sage: a.divides(c,side=Left)
     1447            True
     1448            sage: b.divides(c,side=Left)
     1449            False
     1450
     1451        Divisibility by 0 has no sense::
     1452
     1453            sage: S(0).divides(c)
     1454            Traceback (most recent call last):
     1455            ...
     1456            ZeroDivisionError
     1457
     1458        This function does not work if the leading coefficient of the divisor
     1459        is not a unit::
     1460
     1461            sage: R.<t> = QQ[]
     1462            sage: sigma = R.hom([t+1])
     1463            sage: S.<x> = R['x',sigma]
     1464            sage: a = x^2 + 2*x + t
     1465            sage: b = (t+1)*x + t^2
     1466            sage: c = a*b
     1467            sage: b.divides(c)
     1468            Traceback (most recent call last):
     1469            ...
     1470            NotImplementedError: the leading coefficient of the divisor is not invertible
     1471        """
     1472        return other.rem(self,side=side).is_zero()
     1473
     1474
     1475    # greastest commun divisor
     1476    # ------------------------
     1477
     1478    def lxgcd(self,other,monic=True):
     1479        """
     1480        INPUT:
     1481
     1482        -  ``other`` -- an other skew polynomial over the same
     1483           base
     1484
     1485        -  ``monic`` -- boolean (default: True)
     1486
     1487        OUTPUT:
     1488
     1489        - The left gcd of self and other, that is a skew polynomial
     1490        `g` with the following property: any skew polynomial is
     1491        divisible on the left by `g` iff it is divisible on the left
     1492        by both ``self`` and ``other``.
     1493        If monic is ``True``, `g` is in addition monic. (With this
     1494        extra condition, it is uniquely determined.)
     1495
     1496        - Two skew polynomials `u` and `v` such that:
     1497
     1498        .. MATH::
     1499
     1500            g = self * u + other * v
     1501
     1502        .. NOTE::
     1503
     1504            Works only if two following conditions are fulfilled
     1505            (otherwise left gcd do not exist in general):
     1506            1) the base ring is a field and 2) the twist map on
     1507            this field is bijective.
     1508
     1509        EXAMPLES::
     1510
     1511            sage: k.<t> = GF(5^3)
     1512            sage: Frob = k.frobenius_endomorphism()
     1513            sage: S.<x> = k['x',Frob]
     1514            sage: a = (x + t) * (x^2 + t*x + 1)
     1515            sage: b = 2 * (x + t) * (x^3 + (t+1)*x^2 + t^2)
     1516            sage: g,u,v = a.lxgcd(b); g
     1517            x + t
     1518            sage: a*u + b*v == g
     1519            True
     1520
     1521        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     1522
     1523            sage: g,u,v = a.lxgcd(b,monic=False); g
     1524            2*t*x + 4*t + 2
     1525            sage: a*u + b*v == g
     1526            True
     1527
     1528        The base ring need to be a field::
     1529
     1530            sage: R.<t> = QQ[]
     1531            sage: sigma = R.hom([t+1])
     1532            sage: S.<x> = R['x',sigma]
     1533            sage: a = (x + t) * (x^2 + t*x + 1)
     1534            sage: b = 2 * (x + t) * (x^3 + (t+1)*x^2 + t^2)
     1535            sage: a.lxgcd(b)
     1536            Traceback (most recent call last):
     1537            ...
     1538            TypeError: the base ring must be a field
     1539
     1540        And the twist map need to be bijective::
     1541
     1542            sage: FR = R.fraction_field()
     1543            sage: f = FR.hom([FR(t)^2])
     1544            sage: S.<x> = FR['x',f]
     1545            sage: a = (x + t) * (x^2 + t*x + 1)
     1546            sage: b = 2 * (x + t) * (x^3 + (t+1)*x^2 + t^2)
     1547            sage: a.lxgcd(b)
     1548            Traceback (most recent call last):
     1549            ...
     1550            NotImplementedError
     1551        """
     1552        if not isinstance(self.base_ring(),Field):
     1553            raise TypeError("the base ring must be a field")
     1554        G = self
     1555        U = self._parent(1)
     1556        if other.is_zero():
     1557            V = self._parent(0)
     1558        else:
     1559            V1 = self._parent(0)
     1560            V3 = other
     1561            while not V3.is_zero():
     1562                Q,R = G.lquo_rem(V3)
     1563                T = U - V1*Q
     1564                U = V1
     1565                G = V3
     1566                V1 = T
     1567                V3 = R
     1568            V,_ = (G-self*U).lquo_rem(other)
     1569        if monic:
     1570            lc = ~G.leading_coefficient()
     1571            lc = self._parent.twist_map(-G.degree())(lc)
     1572            G = G*lc
     1573            U = U*lc
     1574            V = V*lc
     1575        return G,U,V
     1576
     1577
     1578    def rxgcd(self,other,monic=True):
     1579        """
     1580        INPUT:
     1581
     1582        -  ``other`` -- an other skew polynomial over the same
     1583           base
     1584
     1585        -  ``monic`` -- boolean (default: True)
     1586
     1587        OUTPUT:
     1588
     1589        - The right gcd of self and other, that is a skew polynomial
     1590        `g` with the following property: any skew polynomial is
     1591        divisible on the right by `g` iff it is divisible on the right
     1592        by both ``self`` and ``other``.
     1593        If monic is ``True``, `g` is in addition monic. (With this
     1594        extra condition, it is uniquely determined.)
     1595
     1596        - Two skew polynomials `u` and `v` such that:
     1597
     1598        .. MATH::
     1599
     1600            g = u * self + v * other
     1601
     1602        .. NOTE::
     1603
     1604            Works only if the base ring is a field (otherwise right
     1605            gcd do not exist in general).
     1606
     1607        EXAMPLES::
     1608
     1609            sage: k.<t> = GF(5^3)
     1610            sage: Frob = k.frobenius_endomorphism()
     1611            sage: S.<x> = k['x',Frob]
     1612            sage: a = (x^2 + t*x + 1) * (x + t)
     1613            sage: b = 2 * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1614            sage: g,u,v = a.rxgcd(b); g
     1615            x + t
     1616            sage: u*a + v*b == g
     1617            True
     1618
     1619        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     1620
     1621            sage: g,u,v = a.rxgcd(b,monic=False); g
     1622            (4*t^2 + 4*t + 1)*x + 4*t^2 + 4*t + 3
     1623            sage: u*a + v*b == g
     1624            True
     1625
     1626        The base ring need to be a field::
     1627
     1628            sage: R.<t> = QQ[]
     1629            sage: sigma = R.hom([t+1])
     1630            sage: S.<x> = R['x',sigma]
     1631            sage: a = (x^2 + t*x + 1) * (x + t)
     1632            sage: b = 2 * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1633            sage: a.rxgcd(b)
     1634            Traceback (most recent call last):
     1635            ...
     1636            TypeError: the base ring must be a field
     1637        """
     1638        if not isinstance(self.base_ring(),Field):
     1639            raise TypeError("the base ring must be a field")
     1640        G = self
     1641        U = self._parent(1)
     1642        if other.is_zero():
     1643            V = self._parent(0)
     1644        else:
     1645            V1 = self._parent(0)
     1646            V3 = other
     1647            while not V3.is_zero():
     1648                Q, R = G.rquo_rem(V3)
     1649                T = U - Q*V1
     1650                U = V1
     1651                G = V3
     1652                V1 = T
     1653                V3 = R
     1654            V,_ = (G-U*self).rquo_rem(other)
     1655        if monic:
     1656            lc = ~G.leading_coefficient()
     1657            G = lc*G
     1658            U = lc*U
     1659            V = lc*V
     1660        return G,U,V
     1661
     1662
     1663    def xgcd(self,other,side=Right,monic=True):
     1664        """
     1665        INPUT:
     1666
     1667        -  ``other`` -- an other skew polynomial over the same
     1668           base
     1669
     1670        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1671
     1672        -  ``monic`` -- boolean (default: True)
     1673
     1674        OUTPUT:
     1675
     1676        - The left (resp. right) gcd of self and other, that is a
     1677        skew polynomial `g` with the following property: any skew
     1678        polynomial is divisible on the left (resp. right) by `g`
     1679        iff it is divisible on the left (resp. right) by both
     1680        ``self`` and ``other``.
     1681        If monic is ``True``, `g` is in addition monic. (With this
     1682        extra condition, it is uniquely determined.)
     1683
     1684        - Two skew polynomials `u` and `v` such that:
     1685
     1686            - (if side=Left) `g = self * u + other * v`
     1687
     1688            - (if side=Right) `g = u * self + v * other`
     1689
     1690        .. NOTE::
     1691
     1692            Works only if the base ring is a field and, when
     1693            ``side`` is ``Left`` if the twist map on this field
     1694            is bijective (otherwise gcd do not exist in general).
     1695
     1696        EXAMPLES::
     1697
     1698            sage: k.<t> = GF(5^3)
     1699            sage: Frob = k.frobenius_endomorphism()
     1700            sage: S.<x> = k['x',Frob]
     1701            sage: a = (x + 2*t) * (x^2 + t*x + 1) * (x + t)
     1702            sage: b = 2 * (x + 2*t) * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1703            sage: g,u,v = a.xgcd(b); g    # side = Right
     1704            x + t
     1705            sage: u*a + v*b == g
     1706            True
     1707            sage: g,u,v = a.xgcd(b,side=Left); g
     1708            x + 2*t
     1709            sage: a*u + b*v == g
     1710            True
     1711
     1712        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     1713
     1714            sage: g,u,v = a.xgcd(b,monic=False); g
     1715            (2*t^2 + 2*t + 4)*x + 2*t^2 + 3*t + 4
     1716            sage: u*a + v*b == g
     1717            True
     1718            sage: g,u,v = a.xgcd(b,side=Left,monic=False); g
     1719            (2*t^2 + 3*t + 4)*x + 2*t^2 + t + 3
     1720            sage: a*u + b*v == g
     1721            True
     1722
     1723        The base ring need to be a field::
     1724
     1725            sage: R.<t> = QQ[]
     1726            sage: sigma = R.hom([t+1])
     1727            sage: S.<x> = R['x',sigma]
     1728            sage: a = (x + 2*t) * (x^2 + t*x + 1) * (x + t)
     1729            sage: b = 2 * (x + 2*t) * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1730            sage: a.xgcd(b)
     1731            Traceback (most recent call last):
     1732            ...
     1733            TypeError: the base ring must be a field
     1734            sage: a.xgcd(b,side=Left)
     1735            Traceback (most recent call last):
     1736            ...
     1737            TypeError: the base ring must be a field
     1738
     1739        And the twist map need to be bijective::
     1740
     1741            sage: FR = R.fraction_field()
     1742            sage: f = FR.hom([FR(t)^2])
     1743            sage: S.<x> = FR['x',f]
     1744            sage: a = (x + 2*t) * (x^2 + t*x + 1) * (x + t)
     1745            sage: b = 2 * (x + 2*t) * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1746            sage: g,u,v = a.xgcd(b); g
     1747            x + t
     1748            sage: g,u,v = a.xgcd(b,side=Left); g
     1749            Traceback (most recent call last):
     1750            ...
     1751            NotImplementedError
     1752        """
     1753        if side == Right:
     1754            return self.rxgcd(other,monic=monic)
     1755        else:
     1756            return self.lxgcd(other,monic=monic)
     1757
     1758
     1759    def rgcd(self,other,monic=True):
     1760        """
     1761        INPUT:
     1762
     1763        -  ``other`` -- an other skew polynomial over the same
     1764           base
     1765
     1766        -  ``monic`` -- boolean (default: True)
     1767
     1768        OUTPUT:
     1769
     1770        - The right gcd of self and other, that is a skew polynomial
     1771        `g` with the following property: any skew polynomial is
     1772        divisible on the right by `g` iff it is divisible on the right
     1773        by both ``self`` and ``other``.
     1774        If monic is ``True``, `g` is in addition monic. (With this
     1775        extra condition, it is uniquely determined.)
     1776
     1777        .. NOTE::
     1778
     1779            Works only if the base ring is a field (otherwise right
     1780            gcd do not exist in general).
     1781
     1782        EXAMPLES::
     1783
     1784            sage: k.<t> = GF(5^3)
     1785            sage: Frob = k.frobenius_endomorphism()
     1786            sage: S.<x> = k['x',Frob]
     1787            sage: a = (x^2 + t*x + 1) * (x + t)
     1788            sage: b = 2 * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1789            sage: a.rgcd(b)
     1790            x + t
     1791
     1792        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     1793
     1794            sage: a.rgcd(b,monic=False)
     1795            (4*t^2 + 4*t + 1)*x + 4*t^2 + 4*t + 3
     1796
     1797        The base ring need to be a field::
     1798
     1799            sage: R.<t> = QQ[]
     1800            sage: sigma = R.hom([t+1])
     1801            sage: S.<x> = R['x',sigma]
     1802            sage: a = (x^2 + t*x + 1) * (x + t)
     1803            sage: b = 2 * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1804            sage: a.rgcd(b)
     1805            Traceback (most recent call last):
     1806            ...
     1807            TypeError: the base ring must be a field
     1808        """
     1809        if not isinstance(self.base_ring(),Field):
     1810            raise TypeError("the base ring must be a field")
     1811        if other.is_zero():
     1812            return self
     1813        A = self
     1814        B = other
     1815        while not B.is_zero():
     1816            A,B = B, A % B
     1817        if monic:
     1818            A = A.rmonic()
     1819        return A
     1820
     1821
     1822    def lgcd(self,other,monic=True):
     1823        """
     1824        INPUT:
     1825
     1826        -  ``other`` -- an other skew polynomial over the same
     1827           base
     1828
     1829        -  ``monic`` -- boolean (default: True)
     1830
     1831        OUTPUT:
     1832
     1833        - The left gcd of self and other, that is a skew polynomial
     1834        `g` with the following property: any skew polynomial is
     1835        divisible on the left by `g` iff it is divisible on the left
     1836        by both ``self`` and ``other``.
     1837        If monic is ``True``, `g` is in addition monic. (With this
     1838        extra condition, it is uniquely determined.)
     1839
     1840        .. NOTE::
     1841
     1842            Works only if two following conditions are fulfilled
     1843            (otherwise left gcd do not exist in general):
     1844            1) the base ring is a field and 2) the twist map on
     1845            this field is bijective.
     1846
     1847        EXAMPLES::
     1848
     1849            sage: k.<t> = GF(5^3)
     1850            sage: Frob = k.frobenius_endomorphism()
     1851            sage: S.<x> = k['x',Frob]
     1852            sage: a = (x + t) * (x^2 + t*x + 1)
     1853            sage: b = 2 * (x + t) * (x^3 + (t+1)*x^2 + t^2)
     1854            sage: a.lgcd(b)
     1855            x + t
     1856
     1857        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     1858
     1859            sage: a.lgcd(b,monic=False)
     1860            2*t*x + 4*t + 2
     1861
     1862        The base ring need to be a field::
     1863
     1864            sage: R.<t> = QQ[]
     1865            sage: sigma = R.hom([t+1])
     1866            sage: S.<x> = R['x',sigma]
     1867            sage: a = (x + t) * (x^2 + t*x + 1)
     1868            sage: b = 2 * (x + t) * (x^3 + (t+1)*x^2 + t^2)
     1869            sage: a.lgcd(b)
     1870            Traceback (most recent call last):
     1871            ...
     1872            TypeError: the base ring must be a field
     1873
     1874        And the twist map need to be bijective::
     1875
     1876            sage: FR = R.fraction_field()
     1877            sage: f = FR.hom([FR(t)^2])
     1878            sage: S.<x> = FR['x',f]
     1879            sage: a = (x + t) * (x^2 + t*x + 1)
     1880            sage: b = 2 * (x + t) * (x^3 + (t+1)*x^2 + t^2)
     1881            sage: a.lgcd(b)
     1882            Traceback (most recent call last):
     1883            ...
     1884            NotImplementedError
     1885        """
     1886        if not isinstance(self.base_ring(),Field):
     1887            raise TypeError("the base ring must be a field")
     1888        if other.is_zero():
     1889            return self
     1890        A = self
     1891        B = other
     1892        while not B.is_zero():
     1893            A,B = B, A.rem(B,side=Left)
     1894        if monic:
     1895            A = A.lmonic()
     1896        return A
     1897
     1898
     1899    def gcd(self,other,side=Right,monic=True):
     1900        """
     1901        INPUT:
     1902
     1903        -  ``other`` -- an other skew polynomial over the same
     1904           base
     1905
     1906        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1907
     1908        -  ``monic`` -- boolean (default: True)
     1909
     1910        OUTPUT:
     1911
     1912        - The left (resp. right) gcd of self and other, that is a
     1913        skew polynomial `g` with the following property: any skew
     1914        polynomial is divisible on the left (resp. right) by `g`
     1915        iff it is divisible on the left (resp. right) by both
     1916        ``self`` and ``other``.
     1917        If monic is ``True``, `g` is in addition monic. (With this
     1918        extra condition, it is uniquely determined.)
     1919
     1920        .. NOTE::
     1921
     1922            Works only if the base ring is a field and, when
     1923            ``side=Left`` if, in addition, the twist map on this
     1924            field is bijective.
     1925
     1926        EXAMPLES::
     1927
     1928            sage: k.<t> = GF(5^3)
     1929            sage: Frob = k.frobenius_endomorphism()
     1930            sage: S.<x> = k['x',Frob]
     1931            sage: a = (x + 2*t) * (x^2 + t*x + 1) * (x + t)
     1932            sage: b = 2 * (x + 2*t) * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1933            sage: a.gcd(b);    # side = Right
     1934            x + t
     1935            sage: a.gcd(b,side=Left)
     1936            x + 2*t
     1937
     1938        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     1939
     1940            sage: a.gcd(b,monic=False)
     1941            (2*t^2 + 2*t + 4)*x + 2*t^2 + 3*t + 4
     1942            sage: a.gcd(b,side=Left,monic=False)
     1943            (2*t^2 + 3*t + 4)*x + 2*t^2 + t + 3
     1944
     1945        The base ring need to be a field::
     1946
     1947            sage: R.<t> = QQ[]
     1948            sage: sigma = R.hom([t+1])
     1949            sage: S.<x> = R['x',sigma]
     1950            sage: a = (x + 2*t) * (x^2 + t*x + 1) * (x + t)
     1951            sage: b = 2 * (x + 2*t) * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1952            sage: a.gcd(b)
     1953            Traceback (most recent call last):
     1954            ...
     1955            TypeError: the base ring must be a field
     1956            sage: a.gcd(b,side=Left)
     1957            Traceback (most recent call last):
     1958            ...
     1959            TypeError: the base ring must be a field
     1960
     1961        And, when ``side=Left``, the twist map need to be bijective::
     1962
     1963            sage: FR = R.fraction_field()
     1964            sage: f = FR.hom([FR(t)^2])
     1965            sage: S.<x> = FR['x',f]
     1966            sage: a = (x + 2*t) * (x^2 + t*x + 1) * (x + t)
     1967            sage: b = 2 * (x + 2*t) * (x^3 + (t+1)*x^2 + t^2) * (x + t)
     1968            sage: a.gcd(b)
     1969            x + t
     1970            sage: a.gcd(b,side=Left)
     1971            Traceback (most recent call last):
     1972            ...
     1973            NotImplementedError
     1974        """
     1975        if side == Right:
     1976            return self.rgcd(other,monic=monic)
     1977        else:
     1978            return self.lgcd(other,monic=monic)
     1979
     1980
     1981    # Lowest common multiple
     1982    # ----------------------
     1983
     1984    def llcm(self,other,monic=True):
     1985        """
     1986        INPUT:
     1987
     1988        -  ``other`` -- an other skew polynomial over the same
     1989           base
     1990
     1991        -  ``monic`` -- boolean (default: True)
     1992
     1993        OUTPUT:
     1994
     1995        - The left lcm of self and other, that is a skew polynomial
     1996        `g` with the following property: any skew polynomial divides
     1997        `g` on the *right* iff it divides both ``self`` and ``other``
     1998        on the *right*.
     1999        If monic is ``True``, `g` is in addition monic. (With this
     2000        extra condition, it is uniquely determined.)
     2001
     2002        .. NOTE::
     2003
     2004            Works only if the base ring is a field (otherwise left
     2005            lcm do not exist in general).
     2006
     2007        EXAMPLES::
     2008
     2009            sage: k.<t> = GF(5^3)
     2010            sage: Frob = k.frobenius_endomorphism()
     2011            sage: S.<x> = k['x',Frob]
     2012            sage: a = (x + t^2) * (x + t)
     2013            sage: b = 2 * (x^2 + t + 1) * (x * t)
     2014            sage: c = a.llcm(b); c
     2015            x^5 + (2*t^2 + t + 4)*x^4 + (3*t^2 + 4)*x^3 + (3*t^2 + 3*t + 2)*x^2 + (t^2 + t + 2)*x
     2016            sage: c.is_divisible_by(a,side=Right)
     2017            True
     2018            sage: c.is_divisible_by(b,side=Right)
     2019            True
     2020            sage: a.degree() + b.degree() == c.degree() + a.rgcd(b).degree()
     2021            True
     2022
     2023        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     2024
     2025            sage: a.llcm(b,monic=False)
     2026            (t^2 + t)*x^5 + (4*t^2 + 4*t + 1)*x^4 + (t + 1)*x^3 + (t^2 + 2)*x^2 + (3*t + 4)*x
     2027
     2028        The base ring need to be a field::
     2029
     2030            sage: R.<t> = QQ[]
     2031            sage: sigma = R.hom([t+1])
     2032            sage: S.<x> = R['x',sigma]
     2033            sage: a = (x + t^2) * (x + t)
     2034            sage: b = 2 * (x^2 + t + 1) * (x * t)
     2035            sage: a.llcm(b)
     2036            Traceback (most recent call last):
     2037            ...
     2038            TypeError: the base ring must be a field
     2039        """
     2040        if not isinstance(self.base_ring(),Field):
     2041            raise TypeError("the base ring must be a field")
     2042        if self.is_zero() or other.is_zero():
     2043            raise ZeroDivisionError
     2044        U = self._parent(1)
     2045        G = self
     2046        V1 = self._parent(0)
     2047        V3 = other
     2048        while not V3.is_zero():
     2049            Q, R = G.rquo_rem(V3)
     2050            T = U - Q*V1
     2051            U = V1
     2052            G = V3
     2053            V1 = T
     2054            V3 = R
     2055        V1 = V1*self
     2056        if monic:
     2057            V1 = V1.rmonic()
     2058        return V1
     2059
     2060
     2061    def rlcm(self,other,monic=True):
     2062        """
     2063        INPUT:
     2064
     2065        -  ``other`` -- an other skew polynomial over the same
     2066           base
     2067
     2068        -  ``monic`` -- boolean (default: True)
     2069
     2070        OUTPUT:
     2071
     2072        - The right lcm of self and other, that is a skew polynomial
     2073        `g` with the following property: any skew polynomial divides
     2074        `g` on the *left* iff it divides both ``self`` and ``other``
     2075        on the *left*.
     2076        If monic is ``True``, `g` is in addition monic. (With this
     2077        extra condition, it is uniquely determined.)
     2078
     2079        .. NOTE::
     2080
     2081            Works only if two following conditions are fulfilled
     2082            (otherwise right lcm do not exist in general):
     2083            1) the base ring is a field and 2) the twist map on
     2084            this field is bijective.
     2085
     2086        EXAMPLES::
     2087
     2088            sage: k.<t> = GF(5^3)
     2089            sage: Frob = k.frobenius_endomorphism()
     2090            sage: S.<x> = k['x',Frob]
     2091            sage: a = (x + t) * (x + t^2)
     2092            sage: b = 2 * (x + t) * (x^2 + t + 1)
     2093            sage: c = a.rlcm(b); c
     2094            x^4 + (2*t^2 + t + 2)*x^3 + (3*t^2 + 4*t + 1)*x^2 + (3*t^2 + 4*t + 1)*x + t^2 + 4
     2095            sage: c.is_divisible_by(a,side=Left)
     2096            True
     2097            sage: c.is_divisible_by(b,side=Left)
     2098            True
     2099            sage: a.degree() + b.degree() == c.degree() + a.lgcd(b).degree()
     2100            True
     2101
     2102        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     2103
     2104            sage: a.rlcm(b,monic=False)
     2105            2*t*x^4 + (3*t + 1)*x^3 + (4*t^2 + 4*t + 3)*x^2 + (3*t^2 + 4*t + 2)*x + 3*t^2 + 2*t + 3
     2106
     2107        The base ring need to be a field::
     2108
     2109            sage: R.<t> = QQ[]
     2110            sage: sigma = R.hom([t+1])
     2111            sage: S.<x> = R['x',sigma]
     2112            sage: a = (x + t) * (x + t^2)
     2113            sage: b = 2 * (x + t) * (x^2 + t + 1)
     2114            sage: a.rlcm(b)
     2115            Traceback (most recent call last):
     2116            ...
     2117            TypeError: the base ring must be a field
     2118
     2119        And the twist map need to be bijective::
     2120
     2121            sage: FR = R.fraction_field()
     2122            sage: f = FR.hom([FR(t)^2])
     2123            sage: S.<x> = FR['x',f]
     2124            sage: a = (x + t) * (x + t^2)
     2125            sage: b = 2 * (x + t) * (x^2 + t + 1)
     2126            sage: a.rlcm(b)
     2127            Traceback (most recent call last):
     2128            ...
     2129            NotImplementedError
     2130        """
     2131        if not isinstance(self.base_ring(),Field):
     2132            raise TypeError("the base ring must be a field")
     2133        if self.is_zero() or other.is_zero():
     2134            raise ZeroDivisionError
     2135        R = self.parent()
     2136        U = R.one_element()
     2137        G = self
     2138        V1 = R.zero_element()
     2139        V3 = other
     2140        while not V3.is_zero():
     2141            Q, R = G.lquo_rem(V3)
     2142            T = U - V1*Q
     2143            U = V1
     2144            G = V3
     2145            V1 = T
     2146            V3 = R
     2147        V1 = self*V1
     2148        if monic:
     2149            V1 = V1.lmonic()
     2150        return V1
     2151
     2152
     2153    def lcm(self,other,side=Left,monic=True):
     2154        """
     2155        INPUT:
     2156
     2157        -  ``other`` -- an other skew polynomial over the same
     2158           base
     2159
     2160        -  ``side`` -- ``Left`` or ``Right`` (default: Left)
     2161
     2162        -  ``monic`` -- boolean (default: True)
     2163
     2164        OUTPUT:
     2165
     2166        - The left (resp. right) lcm of self and other, that is a
     2167        skew polynomial `g` with the following property: any skew
     2168        polynomial divides `g` on the right (resp. left) iff it
     2169        divides both ``self`` and ``other`` on the rignt (resp.
     2170        left)
     2171        If monic is ``True``, `g` is in addition monic. (With this
     2172        extra condition, it is uniquely determined.)
     2173
     2174        .. NOTE::
     2175
     2176            Works only if the base ring is a field and, when
     2177            ``side=Right`` if, in addition, the twist map on this
     2178            field is bijective.
     2179
     2180        EXAMPLES::
     2181
     2182            sage: k.<t> = GF(5^3)
     2183            sage: Frob = k.frobenius_endomorphism()
     2184            sage: S.<x> = k['x',Frob]
     2185            sage: a = (x + 2*t) * (x + t^2) * (x + t)
     2186            sage: b = 2 * (x + 2*t) * (x + t + 1) * (x + t)
     2187
     2188            sage: c = a.lcm(b); c      # side = Left
     2189            x^5 + (2*t^2 + 4*t)*x^4 + (2*t^2 + 3*t + 4)*x^3 + (t^2 + t + 3)*x^2 + (2*t^2 + 4*t + 3)*x + t^2 + t + 1
     2190            sage: c.is_divisible_by(a,side=Right)
     2191            True
     2192            sage: c.is_divisible_by(b,side=Right)
     2193            True
     2194
     2195            sage: c = a.lcm(b,side=Right); c
     2196            x^5 + (3*t^2 + 2*t + 1)*x^4 + (4*t^2 + t + 2)*x^3 + (4*t^2 + 2*t + 3)*x^2 + 4*t^2*x + t^2 + 2
     2197            sage: c.is_divisible_by(a,side=Left)
     2198            True
     2199            sage: c.is_divisible_by(b,side=Left)
     2200            True
     2201
     2202        Specifying ``monic=False``, we *can* get a nonmonic gcd::
     2203
     2204            sage: a.lcm(b,monic=False)
     2205            (3*t^2 + 3*t + 4)*x^5 + (2*t^2 + 4*t + 1)*x^4 + (t^2 + t + 1)*x^3 + (2*t^2 + t + 4)*x^2 + (t^2 + 3*t + 3)*x + t^2 + 1
     2206            sage: a.lcm(b,side=Right,monic=False)
     2207            (4*t + 4)*x^5 + (3*t^2 + t + 1)*x^4 + (t^2 + 4)*x^3 + (4*t^2 + 2*t + 4)*x^2 + (3*t^2 + t)*x + 2*t^2 + 2
     2208
     2209        The base ring need to be a field::
     2210
     2211            sage: R.<t> = QQ[]
     2212            sage: sigma = R.hom([t+1])
     2213            sage: S.<x> = R['x',sigma]
     2214            sage: a = (x + 2*t) * (x + t^2) * (x + t)
     2215            sage: b = 2 * (x + 2*t) * (x + t + 1) * (x + t)
     2216            sage: a.lcm(b)
     2217            Traceback (most recent call last):
     2218            ...
     2219            TypeError: the base ring must be a field
     2220            sage: a.lcm(b,side=Right)
     2221            Traceback (most recent call last):
     2222            ...
     2223            TypeError: the base ring must be a field
     2224
     2225        And the twist map need to be bijective::
     2226
     2227            sage: FR = R.fraction_field()
     2228            sage: f = FR.hom([FR(t)^2])
     2229            sage: S.<x> = FR['x',f]
     2230            sage: a = (x + 2*t) * (x + t^2) * (x + t)
     2231            sage: b = 2 * (x + 2*t) * (x + t + 1) * (x + t)
     2232            sage: a.lcm(b,side=Right)
     2233            Traceback (most recent call last):
     2234            ...
     2235            NotImplementedError
     2236        """
     2237        if side == Right:
     2238            return self.rlcm(other,monic=monic)
     2239        else:
     2240            return self.llcm(other,monic=monic)
     2241
     2242
     2243    # Printing
     2244    # --------
     2245
     2246    def _repr_(self):
     2247        """
     2248        Return string representation of this skew polynomial.
     2249
     2250        EXAMPLES::
     2251
     2252            sage: R.<t> = QQ[]
     2253            sage: sigma = R.hom([t+1])
     2254            sage: S.<x> = R['x',sigma]
     2255            sage: a = t^2 + 1/2*x*t;
     2256            sage: a._repr_()
     2257            '(1/2*t + 1/2)*x + t^2'
     2258        """
     2259        return self._repr()
     2260
     2261
     2262    def _repr(self,name=None):
     2263        """
     2264        Return string representation of this skew polynomial.
     2265
     2266        INPUT:
     2267
     2268        -  ``name`` -- the name of the variable (default: the
     2269           name given when the skew polynomial ring was created)
     2270
     2271        EXAMPLES::
     2272
     2273            sage: R.<t> = QQ[]
     2274            sage: sigma = R.hom([t+1])
     2275            sage: S.<x> = R['x',sigma]
     2276            sage: a = t^2 + 1/2*x*t;
     2277            sage: a._repr()
     2278            '(1/2*t + 1/2)*x + t^2'
     2279            sage: a._repr(name='y')
     2280            '(1/2*t + 1/2)*y + t^2'
     2281        """
     2282        s = " "
     2283        m = self.degree() + 1
     2284        if name is None:
     2285            name = self.parent().variable_name()
     2286        atomic_repr = self.parent().base_ring()._repr_option('element_is_atomic')
     2287        coeffs = self.list()
     2288        for n in reversed(xrange(m)):
     2289            x = coeffs[n]
     2290            if x:
     2291                if n != m-1:
     2292                    s += " + "
     2293                x = y = repr(x)
     2294                if y.find("-") == 0:
     2295                    y = y[1:]
     2296                if not atomic_repr and n > 0 and (y.find("+") != -1 or y.find("-") != -1):
     2297                    x = "(%s)"%x
     2298                if n > 1:
     2299                    var = "*%s^%s"%(name,n)
     2300                elif n==1:
     2301                    var = "*%s"%name
     2302                else:
     2303                    var = ""
     2304                s += "%s%s"%(x,var)
     2305        s = s.replace(" + -", " - ")
     2306        s = re.sub(r' 1(\.0+)?\*',' ', s)
     2307        s = re.sub(r' -1(\.0+)?\*',' -', s)
     2308        if s == " ":
     2309            return "0"
     2310        return s[1:]
     2311
     2312
     2313    def _latex_(self,name=None):
     2314        """
     2315        Return a latex representation of this skew polynomial.
     2316
     2317        INPUT:
     2318
     2319        -  ``name`` -- the name of the variable (default: the
     2320           name given when the skew polynomial ring was created)
     2321
     2322        EXAMPLES::
     2323
     2324            sage: R.<t> = QQ[]
     2325            sage: sigma = R.hom([t+1])
     2326            sage: S.<x> = R['x',sigma]
     2327            sage: a = t^2 + 1/2*x*t;
     2328            sage: a._latex_()
     2329            '\\left(\\frac{1}{2} t + \\frac{1}{2}\\right) x + t^{2}'
     2330            sage: a._latex_(name='y')
     2331            '\\left(\\frac{1}{2} t + \\frac{1}{2}\\right) y + t^{2}'
     2332        """
     2333        s = " "
     2334        coeffs = self.list()
     2335        m = len(coeffs)
     2336        if name is None:
     2337            name = self.parent().latex_variable_names()[0]
     2338        atomic_repr = self.parent().base_ring()._repr_option('element_is_atomic')
     2339        for n in reversed(xrange(m)):
     2340            x = coeffs[n]
     2341            x = y = latex(x)
     2342            if x != '0':
     2343                if n != m-1:
     2344                    s += " + "
     2345                if y.find("-") == 0:
     2346                    y = y[1:]
     2347                if not atomic_repr and n > 0 and (y.find("+") != -1 or y.find("-") != -1):
     2348                    x = "\\left(%s\\right)"%x
     2349                if n > 1:
     2350                    var = "|%s^{%s}"%(name,n)
     2351                elif n==1:
     2352                    var = "|%s"%name
     2353                else:
     2354                    var = ""
     2355                s += "%s %s"%(x,var)
     2356        s = s.replace(" + -", " - ")
     2357        s = re.sub(" 1(\.0+)? \|"," ", s)
     2358        s = re.sub(" -1(\.0+)? \|", " -", s)
     2359        s = s.replace("|","")
     2360        if s==" ":
     2361            return "0"
     2362        return s[1:].lstrip().rstrip()
     2363
     2364
     2365    # Misc
     2366    # ----
     2367
     2368    def _is_atomic(self):
     2369        return (self.degree() == self.valuation() and
     2370                self.leading_coefficient()._is_atomic())
     2371
     2372
     2373    def base_ring(self):
     2374        """
     2375        Return the base ring of this skew polynomial.
     2376
     2377        EXAMPLES::
     2378
     2379            sage: R.<t> = ZZ[]
     2380            sage: sigma = R.hom([t+1])
     2381            sage: S.<x> = R['x',sigma]
     2382            sage: a = S.random_element()
     2383            sage: a.base_ring()
     2384            Univariate Polynomial Ring in t over Integer Ring
     2385            sage: a.base_ring() is R
     2386            True
     2387        """
     2388        return self.parent().base_ring()
     2389
     2390
     2391    def shift(self, n):
     2392        """
     2393        Return this skew polynomial multiplied on the right by the
     2394        power `x^n`. If `n` is negative, terms below `x^n` will be
     2395        discarded.
     2396
     2397        EXAMPLES::
     2398
     2399            sage: R.<t> = QQ[]
     2400            sage: sigma = R.hom([t+1])
     2401            sage: S.<x> = R['x',sigma]
     2402            sage: a = x^5 + t^4*x^4 + t^2*x^2 + t^10
     2403            sage: a.shift(0)
     2404            x^5 + t^4*x^4 + t^2*x^2 + t^10
     2405            sage: a.shift(-1)
     2406            x^4 + t^4*x^3 + t^2*x
     2407            sage: a.shift(-5)
     2408            1
     2409            sage: a.shift(2)
     2410            x^7 + t^4*x^6 + t^2*x^4 + t^10*x^2
     2411
     2412        One can also use the infix shift operator::
     2413
     2414            sage: a >> 2
     2415            x^3 + t^4*x^2 + t^2
     2416            sage: a << 2
     2417            x^7 + t^4*x^6 + t^2*x^4 + t^10*x^2
     2418        """
     2419        if n == 0 or self.degree() < 0:
     2420            return self
     2421        if n > 0:
     2422            return self._parent(n*[self.base_ring().zero_element()] + self.list(), check=False)
     2423        if n < 0:
     2424            if n > self.degree():
     2425                return self._parent([])
     2426            else:
     2427                return self._parent(self.list()[-n:], check=False)
     2428
     2429
     2430    def __lshift__(self, k):
     2431        return self.shift(k)
     2432
     2433
     2434    def __rshift__(self, k):
     2435        return self.shift(-k)
     2436
     2437
     2438    def change_variable_name(self, var):
     2439        """
     2440        Return a new polynomial over the same base ring but in a different
     2441        variable.
     2442
     2443        INPUT:
     2444
     2445        -  ``var`` -- the name of the new variable
     2446
     2447        EXAMPLES::
     2448
     2449            sage: R.<t> = ZZ[]
     2450            sage: sigma = R.hom([t+1])
     2451            sage: S.<x> = R['x',sigma]
     2452            sage: a = x^3 + (2*t + 1)*x  + t^2 + 3*t + 5
     2453            sage: b = a.change_variable_name('y'); b
     2454            y^3 + (2*t + 1)*y  + t^2 + 3*t + 5
     2455
     2456        Remark that a new parent is created at the same time::
     2457
     2458            sage: b.parent()
     2459            Skew Polynomial Ring in y over Univariate Polynomial Ring in t over Integer Ring twisted by t |--> t + 1
     2460        """
     2461        parent = self._parent
     2462        R = parent.base_ring()[var,parent.twist_map()]
     2463        return R(self.list())
     2464
     2465
     2466    def __copy__(self):
     2467        """
     2468        Return a copy of this skew polynomial.
     2469        """
     2470        return self
     2471
     2472
     2473    def dict(self):
     2474        """
     2475        Return a sparse dictionary representation of this skew
     2476        polynomial.
     2477
     2478        EXAMPLES::
     2479
     2480            sage: R.<t> = ZZ[]
     2481            sage: sigma = R.hom([t+1])
     2482            sage: S.<x> = R['x',sigma]
     2483            sage: a = x^2012 + t*x^1006 + t^3 + 2*t
     2484            sage: a.dict()
     2485            {0: t^3 + 2*t, 2012: 1, 1006: t}
     2486        """
     2487        X = {}
     2488        Y = self.list()
     2489        for i in xrange(len(Y)):
     2490            c = Y[i]
     2491            if c:
     2492                X[i] = c
     2493        return X
     2494
     2495
     2496    def is_monomial(self):
     2497        """
     2498        Returns True if self is a monomial, i.e., a power of the generator.
     2499
     2500        EXAMPLES::
     2501
     2502            sage: R.<t> = ZZ[]
     2503            sage: sigma = R.hom([t+1])
     2504            sage: S.<x> = R['x',sigma]
     2505            sage: x.is_monomial()
     2506            True
     2507            sage: (x+1).is_monomial()
     2508            False
     2509            sage: (x^2).is_monomial()
     2510            True
     2511            sage: S(1).is_monomial()
     2512            True
     2513
     2514        The coefficient must be 1::
     2515
     2516            sage: (2*x^5).is_monomial()
     2517            False
     2518            sage: S(t).is_monomial()
     2519            False
     2520
     2521        To allow a non-1 leading coefficient, use is_term()::
     2522
     2523            sage: (2*x^5).is_term()
     2524            True
     2525            sage: S(t).is_term()
     2526            True
     2527        """
     2528        return len(self.exponents()) == 1 and self.leading_coefficient() == 1
     2529
     2530
     2531    def is_term(self):
     2532        """
     2533        Return True if self is an element of the base ring times a
     2534        power of the generator.
     2535
     2536        EXAMPLES::
     2537
     2538            sage: R.<t> = ZZ[]
     2539            sage: sigma = R.hom([t+1])
     2540            sage: S.<x> = R['x',sigma]
     2541            sage: x.is_term()
     2542            True
     2543            sage: R(1).is_term()
     2544            True
     2545            sage: (3*x^5).is_term()
     2546            True
     2547            sage: (1+3*x^5).is_term()
     2548            False
     2549
     2550        To require that the coefficient is 1, use is_monomial() instead::
     2551
     2552            sage: (3*x^5).is_monomial()
     2553            False
     2554        """
     2555        return len(self.exponents()) == 1
     2556
     2557
     2558    def is_gen(self):
     2559        return self._is_gen
     2560
     2561
     2562    def coefficients(self):
     2563        """
     2564        Return the nonzero coefficients of the monomials
     2565        appearing in self.
     2566
     2567        EXAMPLES::
     2568
     2569            sage: R.<t> = QQ[]
     2570            sage: sigma = R.hom([t+1])
     2571            sage: S.<x> = R['x',sigma]
     2572            sage: a = 1 + x^4 + (t+1)*x^2 + t^2
     2573            sage: a.coefficients()
     2574            [t^2 + 1, t + 1, 1]
     2575        """
     2576        return [c for c in self.list() if not c.is_zero()]
     2577
     2578
     2579    def exponents(self):
     2580        """
     2581        Return the exponents of the monomials appearing in self.
     2582
     2583        EXAMPLES::
     2584
     2585            sage: R.<t> = QQ[]
     2586            sage: sigma = R.hom([t+1])
     2587            sage: S.<x> = R['x',sigma]
     2588            sage: a = 1 + x^4 + (t+1)*x^2 + t^2
     2589            sage: a.exponents()
     2590            [0, 2, 4]
     2591        """
     2592        l = self.list()
     2593        return [i for i in range(len(l)) if not l[i].is_zero()]
     2594
     2595
     2596    def prec(self):
     2597        """
     2598        Return the precision of this polynomial. This is always infinity,
     2599        since polynomials are of infinite precision by definition (there is
     2600        no big-oh).
     2601
     2602        EXAMPLES::
     2603
     2604            sage: R.<t> = QQ[]
     2605            sage: sigma = R.hom([t+1])
     2606            sage: S.<x> = R['x',sigma]
     2607            sage: x.prec()
     2608            +Infinity
     2609        """
     2610        return infinity.infinity
     2611
     2612
     2613    def padded_list(self,n=None):
     2614        """
     2615        Return list of coefficients of self up to (but not including)
     2616        `q^n`.
     2617
     2618        Includes 0's in the list on the right so that the list has length
     2619        `n`.
     2620
     2621        INPUT:
     2622
     2623
     2624        -  ``n`` - (default: None); if given, an integer that
     2625           is at least 0
     2626
     2627
     2628        EXAMPLES::
     2629
     2630            sage: R.<t> = QQ[]
     2631            sage: sigma = R.hom([t+1])
     2632            sage: S.<x> = R['x',sigma]
     2633            sage: a = 1 + t*x^3 + t^2*x^5
     2634            sage: a.padded_list()
     2635            [1, 0, 0, t, 0, t^2]
     2636            sage: a.padded_list(10)
     2637            [1, 0, 0, t, 0, t^2, 0, 0, 0, 0]
     2638            sage: len(a.padded_list(10))
     2639            10
     2640            sage: a.padded_list(3)
     2641            [1, 0, 0]
     2642            sage: a.padded_list(0)
     2643            []
     2644            sage: a.padded_list(-1)
     2645            Traceback (most recent call last):
     2646            ...
     2647            ValueError: n must be at least 0
     2648        """
     2649        v = self.list()
     2650        if n is None:
     2651            return v
     2652        if n < 0:
     2653            raise ValueError("n must be at least 0")
     2654        if len(v) < n:
     2655            z = self._parent.base_ring().zero_element()
     2656            return v + [z]*(n - len(v))
     2657        else:
     2658            return v[:int(n)]
     2659
     2660
     2661    def coeffs(self):
     2662        r"""
     2663        Returns ``self.list()``.
     2664
     2665        (It is potentially slightly faster to use
     2666        ``self.list()`` directly.)
     2667
     2668        EXAMPLES::
     2669
     2670            sage: R.<t> = QQ[]
     2671            sage: sigma = R.hom([t+1])
     2672            sage: S.<x> = R['x',sigma]
     2673            sage: a = 1 + x^4 + (t+1)*x^2 + t^2
     2674            sage: a.coeffs()
     2675            [t^2 + 1, 0, t + 1, 0, 1]
     2676        """
     2677        return self.list()
     2678
     2679
     2680    def variable_name(self):
     2681        """
     2682        Return name of variable used in this polynomial as a string.
     2683
     2684        EXAMPLES::
     2685
     2686            sage: R.<t> = QQ[]
     2687            sage: sigma = R.hom([t+1])
     2688            sage: S.<x> = R['x',sigma]
     2689            sage: a = x + t
     2690            sage: a.variable_name()
     2691            'x'
     2692        """
     2693        return self.parent().variable_name()
     2694
     2695
     2696cdef class SkewPolynomial_generic_dense(SkewPolynomial):
     2697    def __init__(self, parent, x=None, int check=1, is_gen=False, int construct=0, **kwds):
     2698        """
     2699        TESTS::
     2700
     2701            sage: R.<t> = QQ[]
     2702            sage: sigma = R.hom([t+1])
     2703            sage: S.<x> = R['x',sigma]
     2704
     2705        We create a skew polynomial from a list::
     2706
     2707            sage: S([t,1])
     2708            x + t
     2709
     2710        from another skew polynomial::
     2711
     2712            sage: S(x^2 + t)
     2713            x^2 + t
     2714
     2715        from a constant::
     2716
     2717            sage: x = S(t^2 + 1); x
     2718            t^2 + 1
     2719            sage: x.parent() is S
     2720            True
     2721        """
     2722        SkewPolynomial.__init__(self, parent, is_gen=is_gen)
     2723        if x is None:
     2724            self.__coeffs = []
     2725            return
     2726
     2727        R = parent.base_ring()
     2728        if PY_TYPE_CHECK(x, list):
     2729            if check:
     2730                self.__coeffs = [R(t) for t in x]
     2731                self.__normalize()
     2732            else:
     2733                self.__coeffs = x
     2734            return
     2735
     2736        if PY_TYPE_CHECK(x, SkewPolynomial):
     2737            if (<Element>x)._parent is self._parent:
     2738                x = list(x.list())
     2739            elif R.has_coerce_map_from((<Element>x)._parent):# is R or (<Element>x)._parent == R:
     2740                try:
     2741                    if x.is_zero():
     2742                        self.__coeffs = []
     2743                        return
     2744                except (AttributeError, TypeError):
     2745                    pass
     2746                x = [x]
     2747            else:
     2748                self.__coeffs = [R(a, **kwds) for a in x.list()]
     2749                if check:
     2750                    self.__normalize()
     2751                return
     2752
     2753        elif PY_TYPE_CHECK(x, int) and x == 0:
     2754            self.__coeffs = []
     2755            return
     2756
     2757        elif isinstance(x, dict):
     2758            x = self._dict_to_list(x, R.zero_element())
     2759
     2760        elif not isinstance(x, list):
     2761            # We trust that the element constructors do not send x=0
     2762#            if x:
     2763            x = [x] # constant polynomials
     2764#            else:
     2765#                x = [] # zero polynomial
     2766        if check:
     2767            self.__coeffs = [R(z, **kwds) for z in x]
     2768            self.__normalize()
     2769        else:
     2770            self.__coeffs = x
     2771
     2772
     2773    cdef list _list_c(self):
     2774        """
     2775        Return the list of the underlying elements of self.
     2776
     2777        .. WARNING::
     2778
     2779            It is a priori not a copy; do not modify this list!
     2780        """
     2781        return self.__coeffs
     2782
     2783
     2784    cdef SkewPolynomial _new_c(self, list coeffs, Parent P, char check=0):
     2785        """
     2786        Fast creation of a new skew polynomial
     2787        """
     2788        cdef SkewPolynomial_generic_dense f = <SkewPolynomial_generic_dense>PY_NEW_SAME_TYPE(self)
     2789        f._parent = P
     2790        f.__coeffs = coeffs
     2791        if check:
     2792            f.__normalize()
     2793        return f
     2794
     2795
     2796    cdef void __normalize(self):
     2797        x = self.__coeffs
     2798        cdef Py_ssize_t n = len(x) - 1
     2799        while n >= 0 and not x[n]:
     2800            del x[n]
     2801            n -= 1
     2802
     2803
     2804    # Basic operations in place
     2805    # -------------------------
     2806
     2807    cdef void _inplace_add(self, SkewPolynomial_generic_dense right):
     2808        """
     2809        Replace self by self+right (only for internal use).
     2810        """
     2811        cdef Py_ssize_t i, min
     2812        x = (<SkewPolynomial_generic_dense>self).__coeffs
     2813        y = (<SkewPolynomial_generic_dense>right).__coeffs
     2814        if len(x) > len(y):
     2815            for i from 0 <= i < len(y):
     2816                x[i] += y[i]
     2817        else:
     2818            x += y[len(x):]
     2819            for i from 0 <= i < len(x):
     2820                x[i] += y[i]
     2821            x += y[len(x):]
     2822        if len(x) == len(y):
     2823            self.__normalize()
     2824
     2825    cdef void _inplace_sub(self, SkewPolynomial_generic_dense right):
     2826        """
     2827        Replace self by self-right (only for internal use).
     2828        """
     2829        cdef Py_ssize_t i, min
     2830        cdef list x = (<SkewPolynomial_generic_dense>self).__coeffs
     2831        cdef list y = (<SkewPolynomial_generic_dense>right).__coeffs
     2832        if len(x) >= len(y):
     2833            for i from 0 <= i < len(y):
     2834                x[i] -= y[i]
     2835        else:
     2836            for i from 0 <= i < len(x):
     2837                x[i] -= y[i]
     2838            x += [-c for c in y[len(x):]]
     2839        if len(x) == len(y):
     2840            self.__normalize()
     2841
     2842    cdef void _inplace_rmul(self, SkewPolynomial_generic_dense right):
     2843        """
     2844        Replace self by self*right (only for internal use).
     2845        """
     2846        cdef list x = (<SkewPolynomial_generic_dense>self).__coeffs
     2847        cdef list y = (<SkewPolynomial_generic_dense>right).__coeffs
     2848        cdef Py_ssize_t i, k, start, end
     2849        cdef Py_ssize_t d1 = len(x)-1, d2 = len(y)-1
     2850        parent = self._parent
     2851        if d2 == -1:
     2852            (<SkewPolynomial_generic_dense>self).__coeffs = [ ]
     2853        elif d1 >= 0:
     2854            for k from d1 < k <= d1+d2:
     2855                start = 0 if k <= d2 else k-d2 # max(0, k-d2)
     2856                sum = x[start] * parent.twist_map(start)(y[k-start])
     2857                for i from start < i <= d1:
     2858                    sum += x[i] * parent.twist_map(i)(y[k-i])
     2859                x.append(sum)
     2860            for k from d1 >= k >= 0:
     2861                start = 0 if k <= d2 else k-d2 # max(0, k-d2)
     2862                end = k if k <= d1 else d1 # min(k, d1)
     2863                sum = x[start] * parent.twist_map(start)(y[k-start])
     2864                for i from start < i <= end:
     2865                    sum += x[i] * parent.twist_map(i)(y[k-i])
     2866                x[k] = sum
     2867
     2868    cdef void _inplace_lmul(self, SkewPolynomial_generic_dense left):
     2869        """
     2870        Replace self by left*self (only for internal use).
     2871        """
     2872        cdef list x = (<SkewPolynomial_generic_dense>self).__coeffs
     2873        cdef list y = (<SkewPolynomial_generic_dense>left).__coeffs
     2874        cdef Py_ssize_t i, k, start, end
     2875        cdef Py_ssize_t d1 = len(x)-1, d2 = len(y)-1
     2876        parent = self._parent
     2877        if d2 == -1:
     2878            (<SkewPolynomial_generic_dense>self).__coeffs = [ ]
     2879        elif d1 >= 0:
     2880            for k from d1 < k <= d1+d2:
     2881                start = 0 if k <= d2 else k-d2 # max(0, k-d2)
     2882                sum = parent.twist_map(k-start)(x[start]) * y[k-start]
     2883                for i from start < i <= d1:
     2884                    sum += parent.twist_map(k-i)(x[i]) * y[k-i]
     2885                x.append(sum)
     2886            for k from d1 >= k >= 0:
     2887                start = 0 if k <= d2 else k-d2 # max(0, k-d2)
     2888                end = k if k <= d1 else d1 # min(k, d1)
     2889                sum = parent.twist_map(k-start)(x[start]) * y[k-start]
     2890                for i from start < i <= end:
     2891                    sum += parent.twist_map(k-i)(x[i]) * y[k-i]
     2892                x[k] = sum
     2893
     2894
     2895    # Fast exponentiation
     2896    # -------------------
     2897
     2898    cdef void _inplace_pow(self, Py_ssize_t n):
     2899        """
     2900        Replace self by self**n.
     2901        """
     2902        while n & 1 == 0:
     2903            self._inplace_rmul(self)
     2904            n = n >> 1
     2905        cdef SkewPolynomial_generic_dense selfpow = <SkewPolynomial_generic_dense>self._new_c(list(self.__coeffs),self._parent)
     2906        n = n >> 1
     2907        while n != 0:
     2908            selfpow._inplace_rmul(selfpow)
     2909            if n&1 == 1:
     2910                self._inplace_rmul(selfpow)
     2911            n = n >> 1
     2912
     2913
     2914    cpdef _pow_(self,exp,modulus=None,side=Right):
     2915        """
     2916        INPUT:
     2917
     2918        -  ``exp`` -- an Integer
     2919
     2920        -  ``modulus`` -- a skew polynomial over the same ring (default: None)
     2921
     2922        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     2923
     2924        OUTPUT:
     2925
     2926        If ``modulus`` is None, return ``self**exp``.
     2927
     2928        Otherwise, return the remainder of self**exp in the ``side``
     2929        euclidean division by ``modulus``.
     2930
     2931        REMARK:
     2932
     2933        The quotient of the underlying skew polynomial ring by the
     2934        principal ideal generated by ``modulus`` is in general *not*
     2935        a ring.
     2936
     2937        As a consequence, Sage first computes exactly ``self**exp``
     2938        and then reduce it modulo ``modulus``.
     2939
     2940        However, if the base ring is a finite field, Sage uses the
     2941        following optimized algorithm:
     2942
     2943        #. One first compute a central skew polynomial `N` which is
     2944           divisible by ``modulus``. (Since `N` lies in center, the
     2945           quotient `K[X,\sigma]/N` inherits a ring structure.)
     2946
     2947        #. One compute ``self**exp`` in the quotient ring `K[X,\sigma]/N`
     2948
     2949        #. One reduce modulo ``modulus`` the result computed in the
     2950           previous step
     2951
     2952        EXAMPLES::
     2953
     2954            sage: k.<t> = GF(5^3)
     2955            sage: Frob = k.frobenius_endomorphism()
     2956            sage: S.<x> = k['x',Frob]
     2957            sage: a = x + t
     2958            sage: b = a^10  # short form for ``a._pow_(10)``
     2959            sage: b == a*a*a*a*a*a*a*a*a*a
     2960            True
     2961
     2962            sage: modulus = x^3 + t*x^2 + (t+3)*x - 2
     2963            sage: br = a._pow_(10,modulus); br
     2964            (t^2 + t)*x^2 + (3*t^2 + 1)*x + t^2 + t
     2965            sage: br == b.rem(modulus)
     2966            True
     2967            sage: bl = a._pow_(10,modulus,side=Left); bl
     2968            (4*t^2 + 2*t + 3)*x^2 + (3*t^2 + 1)*x + 2*t + 3
     2969            sage: bl == b.rem(modulus,side=Left)
     2970            True
     2971
     2972            sage: a._pow_(10^100,modulus)  # quite fast
     2973            (3*t^2 + 3)*x^2 + (t^2 + 2*t + 4)*x + 4*t^2 + 2*t + 1
     2974        """
     2975        cdef SkewPolynomial_generic_dense r
     2976
     2977        if not PY_TYPE_CHECK_EXACT(exp, Integer) or \
     2978                PY_TYPE_CHECK_EXACT(exp, int):
     2979                    try:
     2980                        exp = Integer(exp)
     2981                    except TypeError:
     2982                        raise TypeError("non-integral exponents not supported")
     2983
     2984        if self.degree() <= 0:
     2985            return self.parent()(self[0]**exp)
     2986        if exp == 0:
     2987            return self.parent()(1)
     2988        if exp < 0:
     2989            return (~self).pow(-exp,modulus,side=side)
     2990
     2991        if self == self.parent().gen(): # special case x**n should be faster!
     2992            P = self.parent()
     2993            R = P.base_ring()
     2994            v = [R.zero_element()]*exp + [R.one_element()]
     2995            r = <SkewPolynomial_generic_dense>self._parent(v)
     2996        else:
     2997            r = <SkewPolynomial_generic_dense>self._new_c(list(self.__coeffs),self._parent)
     2998            sig_on()
     2999            r._inplace_pow(exp)
     3000            sig_off()
     3001
     3002        if modulus:
     3003            sig_on()
     3004            if side == Right:
     3005                r._inplace_rrem(modulus)
     3006            else:
     3007                r._inplace_lrem(modulus)
     3008            sig_off()
     3009
     3010        return r
     3011
     3012
     3013    def __pow__(self,exp,modulus):
     3014        """
     3015        INPUT:
     3016
     3017        -  ``exp`` -- an Integer
     3018
     3019        -  ``modulus`` -- a skew polynomial over the same ring
     3020           (default: None)
     3021
     3022        OUTPUT:
     3023
     3024        If ``modulus`` is None, return ``self**exp``.
     3025
     3026        Otherwise, return the remainder of self**exp in the right
     3027        euclidean division by ``modulus``.
     3028
     3029        .. SEEALSO::
     3030
     3031            :meth:`~sage.rings.polynomial.skew_polynomial_element._pow_`
     3032
     3033        EXAMPLES::
     3034
     3035            sage: k.<t> = GF(5^3)
     3036            sage: Frob = k.frobenius_endomorphism()
     3037            sage: S.<x> = k['x',Frob]
     3038            sage: a = x + t
     3039            sage: b = a^10
     3040            sage: b == a*a*a*a*a*a*a*a*a*a
     3041            True
     3042
     3043            sage: modulus = x^3 + t*x^2 + (t+3)*x - 2
     3044            sage: bmod = a._pow_(10,modulus); bmod
     3045            (t^2 + t)*x^2 + (3*t^2 + 1)*x + t^2 + t
     3046            sage: bmod == b.rem(modulus)
     3047            True
     3048        """
     3049        return self._pow_(exp,modulus)
     3050
     3051
     3052def make_generic_skew_polynomial(parent, coeffs):
     3053    return parent(coeffs)
     3054
     3055
     3056cdef class ConstantSkewPolynomialSection(Map):
     3057    cpdef Element _call_(self, x):
     3058        if x.degree() <= 0:
     3059            try:
     3060                return <Element>(x.constant_coefficient())
     3061            except AttributeError:
     3062                return <Element>((<SkewPolynomial>x).constant_coefficient())
     3063        else:
     3064            raise TypeError("not a constant polynomial")
     3065
     3066
     3067cdef class SkewPolynomialBaseringInjection(RingHomomorphism):
     3068    def __init__(self, domain, codomain):
     3069        assert codomain.base_ring() is domain, "domain must be basering"
     3070        RingHomomorphism.__init__(self, Hom(domain,codomain))
     3071        self._an_element = codomain.gen()
     3072        self._repr_type_str = "Polynomial base injection"
     3073        self._new_constant_poly_ = self._an_element._new_constant_poly
     3074
     3075    cpdef Element _call_(self, x):
     3076        return self._new_constant_poly_(x, self._codomain)
     3077
     3078    cpdef Element _call_with_args(self, x, args=(), kwds={}):
     3079        try:
     3080            return self._codomain._element_constructor_(x, *args, **kwds)
     3081        except AttributeError:
     3082            # if there is no element constructor, there is a custom call method.
     3083            return self._codomain(x, *args, **kwds)
     3084
     3085    def section(self):
     3086        return ConstantSkewPolynomialSection(self._codomain, self._domain)
  • new file sage/rings/polynomial/skew_polynomial_finite_field.pxd

    diff --git a/sage/rings/polynomial/skew_polynomial_finite_field.pxd b/sage/rings/polynomial/skew_polynomial_finite_field.pxd
    new file mode 100644
    - +  
     1from sage.rings.integer cimport Integer
     2
     3from skew_polynomial_element cimport SkewPolynomial_generic_dense
     4from sage.matrix.matrix_dense cimport Matrix_dense
     5from skew_polynomial_element cimport CenterSkewPolynomial_generic_dense
     6
     7from polynomial_element cimport Polynomial
     8from sage.structure.element cimport RingElement
     9
     10
     11cdef class SkewPolynomial_finite_field_dense (SkewPolynomial_generic_dense):
     12    # cache
     13    cdef list _conjugates
     14    cdef Polynomial _norm
     15    cdef _norm_factor
     16    cdef _optbound
     17    cdef dict _rdivisors
     18    cdef dict _types
     19    cdef _factorization
     20
     21    cdef inline void _init_cache(self)
     22
     23    # Karatsuba
     24    #cpdef RingElement _mul_karatsuba(self, RingElement other, cutoff=*)
     25    cpdef SkewPolynomial_finite_field_dense _mul_central(self, SkewPolynomial_finite_field_dense right)
     26    cpdef RingElement _mul_(self, RingElement right)
     27    cpdef rquo_rem_karatsuba(self, RingElement other, cutoff=*)
     28
     29    cdef SkewPolynomial_finite_field_dense _rgcd(self,SkewPolynomial_finite_field_dense other)
     30    cdef SkewPolynomial_finite_field_dense _coeff_llcm(self, SkewPolynomial_finite_field_dense other)
     31
     32    # Inplace functions
     33    cdef void _inplace_lrem(self, SkewPolynomial_finite_field_dense other)
     34    cdef void _inplace_rrem(self, SkewPolynomial_finite_field_dense other)
     35    cdef void _inplace_lfloordiv(self, SkewPolynomial_finite_field_dense other)
     36    cdef void _inplace_rfloordiv(self, SkewPolynomial_finite_field_dense other)
     37    cdef void _inplace_pow_mod(self, Integer n, SkewPolynomial_finite_field_dense mod)
     38    cdef void _inplace_lmonic(self)
     39    cdef void _inplace_rmonic(self)
     40    cdef void _inplace_rgcd(self,SkewPolynomial_finite_field_dense other)
     41    cdef Py_ssize_t _val_inplace_unit(self)
     42    cdef SkewPolynomial_finite_field_dense _rquo_inplace_rem(self, SkewPolynomial_finite_field_dense other)
     43
     44    # Specific functions
     45    cdef Matrix_dense _matphir_c(self)
     46    cdef Matrix_dense _matmul_c(self)
     47
     48    # Finding divisors
     49    cdef SkewPolynomial_finite_field_dense _rdivisor_c(P, CenterSkewPolynomial_generic_dense N)
     50
     51    # Finding factorizations
     52    cdef _factor_c(self)
     53    cdef _factor_uniform_c(self)
     54
     55
     56cdef class SkewPolynomial_finite_field_karatsuba:
     57    cdef _parent
     58    cdef Py_ssize_t _order
     59    cdef Py_ssize_t _cutoff
     60    cdef RingElement _zero
     61    cdef _twist
     62    cdef char _algo_matrix
     63    cdef RingElement _t
     64    cdef Matrix_dense _T, _Tinv
     65
     66    cdef list mul_step (self, list x, list y)
     67    cdef list mul_step_matrix(self, list x, list y)
     68    cdef list mul_iter(self, list x, list y, char flag)
     69    cdef list _twinv
     70    cdef list div_step(self, list a, Py_ssize_t ia, Py_ssize_t da, list b, Py_ssize_t ib, Py_ssize_t db)
     71    cdef list div_iter(self, list a, Py_ssize_t ia, Py_ssize_t da, list b, Py_ssize_t ib, Py_ssize_t db)
  • new file sage/rings/polynomial/skew_polynomial_finite_field.pyx

    diff --git a/sage/rings/polynomial/skew_polynomial_finite_field.pyx b/sage/rings/polynomial/skew_polynomial_finite_field.pyx
    new file mode 100644
    - +  
     1r"""
     2This module implements skew polynomials over finite fields.
     3
     4Let `k` be a finite field and `\sigma` be a ring automorphism
     5of `k` (i.e. a power of the Frobenius endomorphism). Let
     6Put `S = k[X,\sigma]`: as an addtive group, it is the usual ring
     7of polynomials with coefficients in `k` and the multiplication
     8on `S` is defined by the rule `X * a = \sigma(a) * X`.
     9
     10We recall that:
     11
     12#. `S` is a left (resp. right) euclidean noncommutative ring
     13
     14#. in particular, every left (resp. right) ideal is principal
     15
     16Since `k` is a finite field, we have the additional following
     17properties:
     18
     19#. the center of `S`, denoted by `Z`, is the univariate
     20   polynomial ring over `k^\sigma` (`k` fixed by `\sigma`)
     21   in the variable `X^r` where `r` is the order of `\sigma`
     22
     23#. `S[1/X]` is an Azumaya algebra over `Z[1/X^r]` (i.e. etale
     24   locally `S[1/X]` is a matrix algebra over `Z[1/X^r]`)
     25
     26#. in particular, we have a reduced norm map `N` from `S[1/X]`
     27   to `Z[1/X^r]` (etale locally, it is the determinant); one
     28   can prove that it maps `S` to `Z`
     29
     30#. `N` has very good properties regarding to factorizations; in
     31   particular:
     32
     33   #. if `a` is a skew polynomial, `a` always divides `N(a)`
     34
     35   #. if `a` is a skew polynomial, any factorization of `N(a)`
     36      (in any order) lifts to a factorization of `a` (and we
     37      have a precise control on the number of such lifts); as
     38      a consequence there is an explicit (but complicated)
     39      formula counting the number of factorizations of a skew
     40      polynomial.
     41
     42
     43.. TODO::
     44
     45    Try to replace as possible ``finite field`` by ``field
     46    endowed with a finite order twist morphism``. It may cause
     47    new phenomena due to the non trivality of the Brauer group.
     48
     49EXAMPLES::
     50
     51We illustrate some properties listed above::
     52
     53    sage: k.<t> = GF(5^3)
     54    sage: Frob = k.frobenius_endomorphism()
     55    sage: S.<x> = k['x',Frob]; S
     56    Skew Polynomial Ring in x over Finite Field in t of size 5^3 twisted by t |--> t^5
     57    sage: Z = S.center(); Z
     58    Center of Skew Polynomial Ring in x over Finite Field in t of size 5^3 twisted by t |--> t^5:
     59    Univariate Polynomial Ring in (x^3) over Finite Field of size 5
     60    sage: a = x^5 + (2*t^2 + t + 1)*x^4 + (3*t^2 + 3*t + 2)*x^3 + (2*t^2 + 3)*x^2 + (4*t^2 + t + 4)*x + 2*t^2 + 2
     61
     62We compute the reduced norm of `a`::
     63
     64    sage: N = a.reduced_norm(); N
     65    (x^3)^5 + 3*(x^3)^4 + 2*(x^3)^3 + 3*(x^3) + 4
     66
     67Note that the parent of `N` is the center `Z` (and not `S`)::
     68
     69    sage: N.parent()
     70    Center of Skew Polynomial Ring in x over Finite Field in t of size 5^3 twisted by t |--> t^5:
     71    Univariate Polynomial Ring in (x^3) over Finite Field of size 5
     72    sage: N.parent() == Z
     73    True
     74
     75    sage: S(N)  # coercion of N into S
     76    x^15 + 3*x^12 + 2*x^9 + 3*x^3 + 4
     77
     78`N` is a mutiple of `a`::
     79
     80    sage: S(N).is_divisible_by(a)
     81    True
     82    sage: S(N).is_divisible_by(a,side=Left)
     83    True
     84
     85.. NOTE::
     86
     87    We really need to coerce first `N` into `S`. Otherwise an
     88    error occurs::
     89
     90        sage: N.is_divisible_by(a)
     91        Traceback (most recent call last):
     92        ...
     93        AttributeError: 'sage.rings.polynomial.skew_polynomial_element.CenterSkewPolynomial_generic_dense' object has no attribute 'is_divisible_by'
     94
     95As a polynomial in `x^r` (here `r = 3`), ``N`` factors as a produit
     96of two irreducible polynomials::
     97
     98    sage: N.factor()
     99    ((x^3)^2 + (x^3) + 1) * ((x^3)^3 + 2*(x^3)^2 + 4*(x^3) + 4)
     100
     101And so does `a`::
     102
     103    sage: F = a.factor(); F
     104    (x^3 + (t^2 + 2*t + 1)*x^2 + (4*t^2 + t + 2)*x + 3*t^2 + t + 4) * (x^2 + (t^2 + 4*t)*x + 2*t)
     105
     106We can check that each of these two factors of `a` corresponds
     107to a factor of `N`::
     108
     109    sage: F[0][0].reduced_norm()
     110    (x^3)^3 + 2*(x^3)^2 + 4*(x^3) + 4
     111    sage: F[1][0].reduced_norm()
     112    (x^3)^2 + (x^3) + 1
     113
     114Actually, `a` has exactly two factorizations corresponding to the
     115two possible orderings of the irreducible factors of `N`::
     116
     117    sage: a.count_factorizations()
     118    2
     119    sage: for F in a.factorizations(): print F
     120    (x^3 + (t^2 + 2*t + 1)*x^2 + (4*t^2 + t + 2)*x + 3*t^2 + t + 4) * (x^2 + (t^2 + 4*t)*x + 2*t)
     121    (x^2 + (2*t + 2)*x + t + 1) * (x^3 + (3*t^2 + 1)*x^2 + (2*t + 4)*x + 4*t^2 + 3*t + 4)
     122
     123AUTHOR::
     124
     125- Xavier Caruso (2012-06-29)
     126"""
     127
     128#############################################################################
     129#    Copyright (C) 2012 Xavier Caruso <xavier.caruso@normalesup.org>
     130#
     131#  Distributed under the terms of the GNU General Public License (GPL)
     132#
     133#                  http://www.gnu.org/licenses/
     134#****************************************************************************
     135
     136
     137include "../../ext/stdsage.pxi"
     138include "../../ext/interrupt.pxi"
     139
     140import operator, copy, re
     141
     142from sage.rings.integer cimport Integer
     143from sage.rings.integer_ring import ZZ
     144from sage.functions.other import ceil
     145
     146from sage.structure.parent cimport Parent
     147
     148from skew_polynomial_element cimport SkewPolynomial_generic_dense
     149from sage.matrix.constructor import matrix
     150from sage.matrix.constructor import zero_matrix
     151from sage.matrix.matrix_dense cimport Matrix_dense
     152from sage.matrix.matrix_space import MatrixSpace
     153from sage.structure.factorization import Factorization
     154
     155from skew_polynomial_element cimport SkewPolynomial
     156from polynomial_ring import PolynomialRing_general
     157from polynomial_ring_constructor import PolynomialRing
     158
     159from sage.rings.ring cimport Ring
     160from sage.structure.element cimport RingElement
     161
     162from sage.structure.side import Left, Right
     163from sage.rings.infinity import Infinity
     164
     165from sage.combinat.q_analogues import q_jordan
     166from sage.functions.other import factorial
     167from sage.combinat.permutation import Permutations
     168from sage.combinat.partition import Partition
     169from sage.misc.mrange import xmrange_iter
     170
     171
     172cdef class SkewPolynomial_finite_field_dense (SkewPolynomial_generic_dense):
     173    def __init__(self, parent, x=None, int check=1, is_gen=False, int construct=0, **kwds):
     174        SkewPolynomial_generic_dense.__init__ (self, parent, x, check, is_gen, construct, **kwds)
     175        self._init_cache()
     176
     177
     178    cdef inline void _init_cache(self):
     179        """
     180        Initialize cached variables (set them to None).
     181        """
     182        self._conjugates = [ self.__coeffs ]
     183        self._norm = None
     184        self._norm_factor = None
     185        self._optbound = None
     186        self._rdivisors = None
     187        self._types = None
     188        self._factorization = None
     189
     190
     191    cdef SkewPolynomial _new_c(self, list coeffs, Parent P, char check=0):
     192        """
     193        Fast creation of a new skew polynomial
     194        """
     195        cdef SkewPolynomial_finite_field_dense f = <SkewPolynomial_finite_field_dense>PY_NEW_SAME_TYPE(self)
     196        f._parent = P
     197        f.__coeffs = coeffs
     198        f._init_cache()
     199        if check:
     200            f.__normalize()
     201        return f
     202
     203
     204    # Skew multiplication
     205    # -------------------
     206
     207    def _mul_karatsuba(self,right,cutoff=None):
     208        """
     209        Karatsuba multiplication
     210
     211        INPUT:
     212
     213        - ``right`` -- an other skew polynomial in the same ring
     214
     215        - ``cutoff`` -- ``None``, an integer or Infinity (default: None)
     216
     217        .. WARNING::
     218
     219            ``cutoff`` need to be greater than or equal to the order of the
     220            twist map acting on the base ring of the underlying skew polynomial
     221            ring.
     222
     223        OUTPUT:
     224
     225        The result of the product self*right (computed by a variant of
     226        Karatsuba`s algorithm)
     227
     228        .. NOTE::
     229
     230            if ``cutoff`` is None, use the default cutoff which is the
     231            maximum between 150 and the order of the twist map.
     232
     233        EXAMPLES::
     234
     235            sage: k.<t> = GF(5^3)
     236            sage: Frob = k.frobenius_endomorphism()
     237            sage: S.<x> = k['x',Frob]
     238            sage: a = S.random_element(5000)
     239            sage: b = S.random_element(5000)
     240            sage: timeit("c = a._mul_karatsuba(b)")  # random, long time
     241            5 loops, best of 3: 659 ms per loop
     242            sage: timeit("c = a._mul_classical(b)")  # random, long time
     243            5 loops, best of 3: 1.9 s per loop
     244            sage: a._mul_karatsuba(b) == a._mul_classical(b)
     245            True
     246
     247        The operator ``*`` performs Karatsuba multiplication::
     248
     249            sage: timeit("c = a*b")  # random, long time
     250            5 loops, best of 3: 653 ms per loop
     251        """
     252        karatsuba_class = self._parent._karatsuba_class
     253        if cutoff != None:
     254            save_cutoff = karatsuba_class.get_cutoff()
     255            karatsuba_class.set_cutoff(cutoff)
     256        res = karatsuba_class.mul(self,right)
     257        if cutoff != None:
     258            karatsuba_class.set_cutoff(save_cutoff)
     259        return res
     260
     261
     262    def _mul_karatsuba_matrix(self,right):
     263        """
     264        Karatsuba multiplication with multiplication step based
     265        on an isomorphism between a quotient of the underlying
     266        skew polynomial ring and a ring of matrices.
     267
     268        INPUT:
     269
     270        - ``right`` -- an other skew polynomial in the same ring
     271
     272        OUTPUT:
     273
     274        The result of the product self*right
     275
     276        EXAMPLES::
     277
     278            sage: k.<t> = GF(5^3)
     279            sage: Frob = k.frobenius_endomorphism()
     280            sage: S.<x> = k['x',Frob]
     281            sage: a = S.random_element(degree=20)
     282            sage: b = S.random_element(degree=20)
     283            sage: a._mul_karatsuba_matrix(b) == a*b
     284            True
     285
     286        This routine is only efficient when the twisting map (here
     287        ``Frob``) has a large order `r` and the degrees of ``self``
     288        and ``other`` have a very special shape (just below a power
     289        of `2` times `r * floor(r/2)`)::
     290
     291            sage: k.<t> = GF(5^40)
     292            sage: Frob = k.frobenius_endomorphism()
     293            sage: S.<x> = k['x',Frob]
     294            sage: a = S.random_element(degree=799)
     295            sage: b = S.random_element(degree=799)
     296            sage: timeit("c = a*b")  # random, long time
     297            5 loops, best of 3: 23.1 s per loop
     298            sage: timeit("c = a._mul_karatsuba_matrix(b)")  # random, long time
     299            5 loops, best of 3: 12.2 s per loop
     300        """
     301        karatsuba_class = self._parent._karatsuba_class
     302        return karatsuba_class.mul_matrix(self,right)
     303
     304
     305    cpdef SkewPolynomial_finite_field_dense _mul_central(self, SkewPolynomial_finite_field_dense right):
     306        r"""
     307        Return self * right
     308
     309        .. WARNING::
     310
     311            Do you use this function! It is very slow due to a quite
     312            slow interface with ``polynomial_zz_pex``.
     313
     314        ALGORITHM::
     315
     316        Notations::
     317
     318        -  `S` is the underlyling skew polynomial ring
     319
     320        -  `x` is the variable on `S`
     321
     322        -  `k` is the base ring of `S` (it is a finite field)
     323
     324        -  `\sigma` is the twisting automorphism acting on `k`
     325
     326        -  `r` is the order of `\sigma`
     327
     328        -  `t` is a generator of `k` over `k^\sigma`
     329
     330        #. We decompose the polynomial ``right`` as follows::
     331
     332           .. MATH::
     333
     334               right = \sum_{i=0}^{r-1} \sum_{j=0}^{r-1} y_{i,j} t^j x^i
     335
     336           where `y_{i,j}` are polynomials in the center `k^\sigma[x^r]`.
     337
     338        #. We compute all products `z_{i,j} = left * y_{i,j}`; since
     339           all `y_{i,j}` lie in the center, we can compute all these
     340           products as if `left` was a commutative polynomial (and we
     341           can therefore use fast algorithms like FFT and/or fast
     342           implementations)
     343
     344        #. We compute and return the sum
     345
     346           .. MATH::
     347
     348               \sum_{i=0}^{r-1} \sum_{j=0}^{r-1} z_{i,j} t^j x^i
     349
     350        EXAMPLES::
     351
     352            sage: k.<t> = GF(5^3)
     353            sage: Frob = k.frobenius_endomorphism()
     354            sage: S.<x> = k['x',Frob]
     355            sage: a = S.random_element(degree=10)
     356            sage: b = S.random_element(degree=10)
     357            sage: a._mul_central(b) == a*b
     358            True
     359
     360        TESTS::
     361
     362        Here is an example where `k^\sigma` is not a prime field::
     363
     364            sage: k.<t> = GF(5^6)
     365            sage: Frob = k.frobenius_endomorphism(2)
     366            sage: S.<x> = k['x',Frob]
     367            sage: a = S.random_element(degree=10)
     368            sage: b = S.random_element(degree=10)
     369            sage: a._mul_central(b) == a*b
     370            True
     371        """
     372        skew_ring = self._parent
     373        base_ring = skew_ring.base_ring()
     374        commutative_ring = PolynomialRing(skew_ring.base_ring(),name='x')
     375        cdef RingElement c
     376        cdef RingElement zero = base_ring(0)
     377        cdef Py_ssize_t i, j, k
     378        cdef Py_ssize_t order = skew_ring._order
     379        cdef Py_ssize_t degree = base_ring.degree()
     380
     381        left = commutative_ring(self.__coeffs)
     382        cdef list y = [ c.polynomial() for c in right.__coeffs ]
     383        cdef Py_ssize_t leny = len(y)
     384        cdef list yc = leny * [zero]
     385        cdef list res = (leny + len(self.__coeffs) - 1) * [zero]
     386        cdef list term
     387        cdef list twist = [ base_ring.gen() ]
     388        for i from 0 <= i < order-1:
     389            twist.append(skew_ring.twist_map(1)(twist[i]))
     390        for i from 0 <= i < order:
     391            for j from 0 <= j < degree:
     392                for k from i <= k < leny by order:
     393                    yc[k] = y[k][j]
     394                term = (left * commutative_ring(yc)).list()
     395                for k from i <= k < len(term):
     396                    res[k] += term[k] * twist[(k-i)%order]**j
     397            for k from i <= k < leny by order:
     398                yc[k] = zero
     399        return self._new_c(res,skew_ring,1)
     400
     401
     402    cpdef RingElement _mul_(self, RingElement right):
     403        """
     404        Compute self * right (in this order)
     405
     406        .. NOTE::
     407
     408            Use skew Karatsuba's algorithm for skew
     409            polynomials of large degrees.
     410
     411        INPUT:
     412
     413        -  right -- a skew polynomial in the same ring
     414
     415        TESTS::
     416
     417            sage: k.<t> = GF(5^3)
     418            sage: Frob = k.frobenius_endomorphism()
     419            sage: S.<x> = k['x',Frob]
     420            sage: a = x^2 + t; a
     421            x^2 + t
     422            sage: b = x^2 + (t + 1)*x; b
     423            x^2 + (t + 1)*x
     424            sage: a * b
     425            x^4 + (3*t^2 + 2)*x^3 + t*x^2 + (t^2 + t)*x
     426            sage: a * b == b * a
     427            False
     428        """
     429        return self._parent._karatsuba_class.mul(self,right)
     430
     431
     432    def _mul_classical(self,right):
     433        """
     434        Compute self * right (in this order) using the
     435        skew SchoolBook algorithm.
     436
     437        INPUT:
     438
     439        -  ``right`` -- a skew polynomial in the same ring
     440
     441        TESTS::
     442
     443            sage: k.<t> = GF(5^3)
     444            sage: Frob = k.frobenius_endomorphism()
     445            sage: S.<x> = k['x',Frob]
     446            sage: a = x^2 + t; a
     447            x^2 + t
     448            sage: b = x^2 + (t + 1)*x; b
     449            x^2 + (t + 1)*x
     450            sage: a._mul_classical(b)
     451            x^4 + (3*t^2 + 2)*x^3 + t*x^2 + (t^2 + t)*x
     452            sage: a * b == b * a
     453            False
     454        """
     455        karatsuba_class = self._parent._karatsuba_class
     456        save_cutoff = karatsuba_class.get_cutoff()
     457        karatsuba_class.set_cutoff(Infinity)
     458        res = karatsuba_class.mul(self,right)
     459        karatsuba_class.set_cutoff(save_cutoff)
     460        return res
     461
     462
     463    cpdef rquo_rem_karatsuba(self, RingElement other, cutoff=None):
     464        """
     465        Right euclidean division based on Karatsuba's algorithm.
     466
     467        DO NOT USE THIS! It is not efficient for usual degrees!
     468
     469        .. TODO::
     470
     471            Try to understand why...
     472
     473        TESTS::
     474
     475            sage: k.<t> = GF(5^3)
     476            sage: Frob = k.frobenius_endomorphism()
     477            sage: S.<x> = k['x',Frob]
     478
     479            sage: a = S.random_element(2000)
     480            sage: b = S.random_element(1000)
     481            sage: timeit("q,r = a.rquo_rem_karatsuba(b)")  # random, long time
     482            5 loops, best of 3: 104 ms per loop
     483            sage: timeit("q,r = a.rquo_rem(b)")  # random, long time
     484            5 loops, best of 3: 79.6 ms per loop
     485            sage: a.rquo_rem(b) == a.rquo_rem_karatsuba(b)
     486            True
     487
     488            sage: a = S.random_element(10000)
     489            sage: b = S.random_element(5000)
     490            sage: timeit("q,r = a.rquo_rem_karatsuba(b)")  # random, long time
     491            5 loops, best of 3: 1.79 s per loop
     492            sage: timeit("q,r = a.rquo_rem(b)")  # random, long time
     493            5 loops, best of 3: 1.93 s per loop
     494            sage: a.rquo_rem(b) == a.rquo_rem_karatsuba(b)
     495            True
     496        """
     497        karatsuba_class = self._parent._karatsuba_class
     498        if cutoff != None:
     499            save_cutoff = karatsuba_class.get_cutoff()
     500            karatsuba_class.set_cutoff(cutoff)
     501        res = karatsuba_class.div(self,other)
     502        if cutoff != None:
     503            karatsuba_class.set_cutoff(save_cutoff)
     504        return res
     505
     506
     507    # We improve some functions
     508    # -------------------------
     509
     510    def rquo_rem(self,other):
     511        """
     512        DEFINITION:
     513
     514        Let `a` and `b` be two skew polynomials over the same
     515        ring. The *right euclidean division* of `a` by `b` is
     516        a couple `(q,r)` such that
     517
     518        -  `a = q*b + r`
     519
     520        -  the degree of `r` is less than the degree of `b`
     521
     522        `q` (resp. `r`) is called the *quotient* (resp. the
     523        remainder) of this euclidean division.
     524
     525        If the leading coefficient of `b` is a unit (e.g. if
     526        `b` is monic) then `q` and `r` exist and are unique.
     527
     528        INPUT:
     529
     530        -  ``other`` -- a skew polynomial ring over the same
     531           base ring
     532
     533        OUTPUT:
     534
     535        -  the quotient and the remainder of the left euclidean
     536           division of this skew polynomial by ``other``
     537
     538        .. NOTE::
     539
     540            Doesn't work if the leading coefficient of the divisor
     541            is not a unit.
     542
     543        EXAMPLES::
     544
     545            sage: R.<t> = ZZ[]
     546            sage: sigma = R.hom([t+1])
     547            sage: S.<x> = R['x',sigma]
     548            sage: a = S.random_element(degree=4); a
     549            t^2*x^4 + (-12*t^2 - 2*t - 1)*x^3 + (-95*t^2 + t + 2)*x^2 + (-t^2 + t)*x + 2*t - 8
     550            sage: b = S.random_element(monic=True); b
     551            x^2 + (4*t^2 - t - 2)*x - t^2 + t - 1
     552            sage: q,r = a.rquo_rem(b)
     553            sage: a == q*b + r
     554            True
     555
     556        The leading coefficient of the divisor need to be invertible::
     557
     558            sage: c = S.random_element(); c
     559            (-4*t^2 + t)*x^2 - 2*t^2*x + 5*t^2 - 6*t - 4
     560            sage: a.rquo_rem(c)
     561            Traceback (most recent call last):
     562            ...
     563            NotImplementedError: the leading coefficient of the divisor is not invertible
     564        """
     565        cdef list a = self.list()
     566        cdef list b = other.list()
     567        cdef Py_ssize_t i, j
     568        cdef Py_ssize_t da = len(a)-1, db = len(b)-1
     569        parent = self._parent
     570        if db < 0:
     571            raise ZeroDivisionError
     572        if da < db:
     573            return parent(0), self
     574        cdef RingElement inv = ~b[db]
     575        cdef list q = [ ]
     576        cdef Py_ssize_t order = parent._order
     577        cdef list twinv = [ inv ], twb = [ b ]
     578        cdef RingElement c, x
     579        for i from 0 <= i < min(da-db,order-1):
     580            twinv.append(parent.twist_map()(twinv[i]))
     581            twb.append([ parent.twist_map()(x) for x in twb[i] ])
     582        for i from da-db >= i >= 0:
     583            c = twinv[i%order] * a[i+db]
     584            for j from 0 <= j < db:
     585                a[i+j] -= c * twb[i%order][j]
     586            q.append(c)
     587        q.reverse()
     588        return parent(q), parent(a[:db])
     589
     590
     591    cdef SkewPolynomial_finite_field_dense _rgcd(self,SkewPolynomial_finite_field_dense other):
     592        """
     593        Fast right gcd.
     594        """
     595        cdef SkewPolynomial_finite_field_dense A = self
     596        cdef SkewPolynomial_finite_field_dense B = other
     597        cdef SkewPolynomial_finite_field_dense swap
     598        if len(B.__coeffs):
     599            A = <SkewPolynomial_finite_field_dense>self._new_c(A.__coeffs[:],A._parent)
     600            B = <SkewPolynomial_finite_field_dense>B._new_c(B.__coeffs[:],B._parent)
     601            while len(B.__coeffs):
     602                A._inplace_rrem(B)
     603                swap = A; A = B; B = swap
     604            return A
     605        else:
     606            return self
     607
     608
     609    cdef void _inplace_rmul(self, SkewPolynomial_generic_dense right):
     610        """
     611        Replace self by self*right
     612        """
     613        self.__coeffs = self._parent._karatsuba_class.mul_list(self.__coeffs,right.__coeffs)
     614        self._init_cache()
     615
     616
     617    cdef void _inplace_lmul(self, SkewPolynomial_generic_dense left):
     618        """
     619        Replace self by left*self
     620        """
     621        self.__coeffs = self._parent._karatsuba_class.mul_list(left.__coeffs,self.__coeffs)
     622        self._init_cache()
     623
     624
     625    cpdef _pow_(self,right,modulus=None,side=Right):
     626        """
     627        INPUT:
     628
     629        -  ``exp`` -- an Integer
     630
     631        -  ``modulus`` -- a skew polynomial over the same ring (default: None)
     632
     633        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     634
     635        OUTPUT:
     636
     637        If ``modulus`` is None, return ``self**exp``.
     638
     639        Otherwise, return the remainder of self**exp in the ``side``
     640        euclidean division by ``modulus``.
     641
     642        REMARK:
     643
     644        The quotient of the underlying skew polynomial ring by the
     645        principal ideal generated by ``modulus`` is in general *not*
     646        a ring.
     647
     648        As a consequence, Sage first computes exactly ``self**exp``
     649        and then reduce it modulo ``modulus``.
     650
     651        However, if the base ring is a finite field, Sage uses the
     652        following optimized algorithm:
     653
     654        #. One first compute a central skew polynomial `N` which is
     655           divisible by ``modulus``. (Since `N` lies in center, the
     656           quotient `K[X,\sigma]/N` inherits a ring structure.)
     657
     658        #. One compute ``self**exp`` in the quotient ring `K[X,\sigma]/N`
     659
     660        #. One reduce modulo ``modulus`` the result computed in the
     661           previous step
     662
     663        EXAMPLES::
     664
     665            sage: k.<t> = GF(5^3)
     666            sage: Frob = k.frobenius_endomorphism()
     667            sage: S.<x> = k['x',Frob]
     668            sage: a = x + t
     669            sage: b = a^10  # short form for ``a._pow_(10)``
     670            sage: b == a*a*a*a*a*a*a*a*a*a
     671            True
     672
     673            sage: modulus = x^3 + t*x^2 + (t+3)*x - 2
     674            sage: br = a._pow_(10,modulus); br
     675            (t^2 + t)*x^2 + (3*t^2 + 1)*x + t^2 + t
     676            sage: br == b.rem(modulus)
     677            True
     678            sage: bl = a._pow_(10,modulus,side=Left); bl
     679            (4*t^2 + 2*t + 3)*x^2 + (3*t^2 + 1)*x + 2*t + 3
     680            sage: bl == b.rem(modulus,side=Left)
     681            True
     682
     683            sage: a._pow_(10^100,modulus)  # rather fast
     684            (3*t^2 + 3)*x^2 + (t^2 + 2*t + 4)*x + 4*t^2 + 2*t + 1
     685        """
     686        cdef SkewPolynomial_finite_field_dense r
     687
     688        if not isinstance(right, Integer) or isinstance(right, int):
     689            try:
     690                right = Integer(right)
     691            except TypeError:
     692                raise TypeError("non-integral exponents not supported")
     693
     694        if self.degree() <= 0:
     695            return self.parent()(self[0]**right)
     696        if right == 0:
     697            return self.parent()(1)
     698        if right < 0:
     699            return (~self).pow(-right,modulus,side=side)
     700
     701        if self == self.parent().gen(): # special case x**n should be faster!
     702            P = self.parent()
     703            R = P.base_ring()
     704            v = [R.zero_element()]*right + [R.one_element()]
     705            r = <SkewPolynomial_generic_dense>self._parent(v)
     706            sig_on()
     707            if modulus:
     708                if side is Right:
     709                    r._inplace_rrem(modulus)
     710                else:
     711                    r._inplace_lrem(modulus)
     712                r._init_cache()
     713            sig_off()
     714            return r
     715
     716        mod = modulus
     717        if not modulus is None:
     718            try:
     719                mod = self.parent()(mod.bound())
     720            except NotImplementedError:
     721                mod = None
     722        r = <SkewPolynomial_generic_dense>self._new_c(self.__coeffs,self._parent)
     723        sig_on()
     724        if mod:
     725            r._inplace_pow_mod(right,mod)
     726        else:
     727            r._inplace_pow(right)
     728        if (not modulus is None) and modulus != mod:
     729            if side is Right:
     730                r._inplace_rrem(modulus)
     731            else:
     732                r._inplace_lrem(modulus)
     733            r._init_cache()
     734        sig_off()
     735        return r
     736
     737
     738    # Inplace functions
     739    # -----------------
     740
     741    #cdef void _inplace_conjugate(self,n):
     742    #    cdef Py_ssize_t i
     743    #    cdef Morphism twist = <Morphism>self._parent.twist_map(n)
     744    #    cdef RingElement x
     745    #    for i from 0 <= i < len(self.__coeffs):
     746    #        x = twist(self.__coeffs[i])
     747    #        self.__coeffs[i] = x
     748
     749
     750    cdef void _inplace_lrem(self, SkewPolynomial_finite_field_dense other):
     751        """
     752        Replace self by the remainder in the left euclidean division
     753        of self by other (only for internal use).
     754        """
     755        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     756        cdef list b = (<SkewPolynomial_finite_field_dense>other).__coeffs
     757        cdef Py_ssize_t da = len(a)-1, db = len(b)-1
     758        cdef Py_ssize_t i, j
     759        cdef RingElement c, inv
     760        parent = self._parent
     761        if db < 0:
     762            raise ZeroDivisionError
     763        if da >= db:
     764            inv = ~b[db]
     765            for i from da-db >= i >= 0:
     766                c = parent.twist_map(-db)(inv*a[i+db])
     767                for j from 0 <= j < db:
     768                    a[i+j] -= b[j] * parent.twist_map(j)(c)
     769            del a[db:]
     770            self.__normalize()
     771        self._init_cache()
     772
     773
     774    cdef void _inplace_rrem(self, SkewPolynomial_finite_field_dense other):
     775        """
     776        Replace self by the remainder in the right euclidean division
     777        of self by other (only for internal use).
     778        """
     779        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     780        cdef list b = (<SkewPolynomial_finite_field_dense>other).__coeffs
     781        cdef Py_ssize_t da = len(a)-1, db = len(b)-1
     782        cdef Py_ssize_t i, j, order
     783        cdef RingElement c, x, inv
     784        cdef list twinv, twb
     785        parent = self._parent
     786        if db < 0:
     787            raise ZeroDivisionError
     788        if da >= db:
     789            order = parent._order
     790            inv = ~b[db]
     791            twinv = [ inv ]
     792            for i from 0 <= i < min(da-db,order-1):
     793                twinv.append(parent.twist_map()(twinv[i]))
     794            twb = (<SkewPolynomial_finite_field_dense>other)._conjugates
     795            for i from len(twb)-1 <= i < min(da-db,order-1):
     796                twb.append([ parent.twist_map()(x) for x in twb[i] ])
     797            for i from da-db >= i >= 0:
     798                c = twinv[i%order] * a[i+db]
     799                for j from 0 <= j < db:
     800                    a[i+j] -= c * twb[i%order][j]
     801            del a[db:]
     802            self.__normalize()
     803        self._init_cache()
     804
     805
     806    cdef void _inplace_lfloordiv(self, SkewPolynomial_finite_field_dense other):
     807        """
     808        Replace self by the quotient in the left euclidean division
     809        of self by other (only for internal use).
     810        """
     811        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     812        cdef list b = (<SkewPolynomial_finite_field_dense>other).__coeffs
     813        cdef Py_ssize_t da = len(a)-1, db = len(b)-1
     814        cdef Py_ssize_t i, j, deb
     815        cdef RingElement c, inv
     816        parent = self._parent
     817        if db < 0:
     818            raise ZeroDivisionError
     819        if da < db:
     820            (<SkewPolynomial_finite_field_dense>self).__coeffs = [ ]
     821        else:
     822            inv = ~b[db]
     823            for i from da-db >= i >= 0:
     824                c = a[i+db] = parent.twist_map(-db)(inv*a[i+db])
     825                if i < db: deb = db
     826                else: deb = i
     827                for j from deb <= j < db+i:
     828                    a[j] -= b[j-i] * parent.twist_map(j-i)(c)
     829            del a[:db]
     830            self.__normalize()
     831        self._init_cache()
     832
     833
     834    cdef void _inplace_rfloordiv(self, SkewPolynomial_finite_field_dense other):
     835        """
     836        Replace self by the quotient in the right euclidean division
     837        of self by other (only for internal use).
     838        """
     839        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     840        cdef list b = (<SkewPolynomial_finite_field_dense>other).__coeffs
     841        cdef Py_ssize_t da = len(a)-1, db = len(b)-1
     842        cdef Py_ssize_t i, j, deb, order
     843        cdef RingElement c, x, inv
     844        parent = self._parent
     845        if db < 0:
     846            raise ZeroDivisionError
     847        if da < db:
     848            (<SkewPolynomial_finite_field_dense>self).__coeffs = [ ]
     849        else:
     850            order = parent._order
     851            inv = ~b[db]
     852            twinv = [ inv ]
     853            for i from 0 <= i < min(da-db,order-1):
     854                twinv.append(parent.twist_map()(twinv[i]))
     855            twb = (<SkewPolynomial_finite_field_dense>other)._conjugates
     856            for i from len(twb)-1 <= i < min(da-db,order-1):
     857                twb.append([ parent.twist_map()(x) for x in twb[i] ])
     858            for i from da-db >= i >= 0:
     859                c = a[i+db] = twinv[i%order] * a[i+db]
     860                if i < db: deb = db
     861                else: deb = i
     862                for j from deb <= j < db+i:
     863                    a[j] -= c * twb[i%order][j-i]
     864            del a[:db]
     865            self.__normalize()
     866        self._init_cache()
     867
     868
     869    cdef void _inplace_pow_mod(self, Integer n, SkewPolynomial_finite_field_dense mod):
     870        """
     871        Replace self by the remainder in the euclidean
     872        division of self**n by mod (only for internal use).
     873
     874        .. WARNING::
     875
     876            Assume that mod is central
     877        """
     878        while n&1 == 0:
     879            self._inplace_rmul(self)
     880            self._inplace_rrem(mod)
     881            n = n >> 1
     882        cdef SkewPolynomial_finite_field_dense selfpow = <SkewPolynomial_finite_field_dense>self._new_c(self.__coeffs[:], self._parent)
     883        n = n >> 1
     884        while n != 0:
     885            selfpow._inplace_rmul(selfpow)
     886            selfpow._inplace_rrem(mod)
     887            if n&1 == 1:
     888                self._inplace_rmul(selfpow)
     889                self._inplace_rrem(mod)
     890            n = n >> 1
     891
     892
     893    cdef void _inplace_lmonic(self):
     894        """
     895        Replace self by ``self.lmonic()`` (only for internal use).
     896        """
     897        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     898        cdef Py_ssize_t da = len(a)-1, i
     899        cdef RingElement inv = ~a[da]
     900        parent = self._parent
     901        a[da] = parent.base_ring()(1)
     902        for i from 0 <= i < da:
     903            a[i] *= parent.twist_map(i-da)(inv)
     904        self._init_cache()
     905
     906
     907    cdef void _inplace_rmonic(self):
     908        """
     909        Replace self by ``self.rmonic()`` (only for internal use).
     910        """
     911        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     912        cdef Py_ssize_t da = len(a)-1, i
     913        cdef RingElement inv = ~a[da]
     914        a[da] = self._parent.base_ring()(1)
     915        for i from 0 <= i < da:
     916            a[i] *= inv
     917        self._init_cache()
     918
     919
     920    cdef void _inplace_rgcd(self,SkewPolynomial_finite_field_dense other):
     921        """
     922        Replace self by its right gcd with other (only for internal use).
     923        """
     924        cdef SkewPolynomial_finite_field_dense B
     925        cdef list swap
     926        if len(other.__coeffs):
     927            B = <SkewPolynomial_finite_field_dense>self._new_c(other.__coeffs[:],other._parent)
     928            while len(B.__coeffs):
     929                B._conjugates = [ B.__coeffs ]
     930                self._inplace_rrem(B)
     931                swap = self.__coeffs
     932                self.__coeffs = B.__coeffs
     933                B.__coeffs = swap
     934        self._init_cache()
     935
     936
     937    cdef SkewPolynomial_finite_field_dense _rquo_inplace_rem(self, SkewPolynomial_finite_field_dense other):
     938        """
     939        Replace self by the remainder in the right euclidean division
     940        of self by other and return the quotient (only for internal use).
     941        """
     942        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     943        cdef list b = (<SkewPolynomial_finite_field_dense>other).__coeffs
     944        cdef Py_ssize_t da = len(a)-1, db = len(b)-1
     945        cdef Py_ssize_t i, j
     946        cdef RingElement c, inv
     947        cdef list q
     948        parent = self._parent
     949        if db < 0:
     950            raise ZeroDivisionError
     951        if da < db:
     952            return self._new_c([],self._parent)
     953        inv = ~b[db]
     954        q = [ ]
     955        for i from da-db >= i >= 0:
     956            c = parent.twist_map(i)(inv) * a[i+db]
     957            q.append(c)
     958            for j from 0 <= j < db:
     959                a[i+j] -= c * parent.twist_map(i)(b[j])
     960        del a[db:]
     961        self.__normalize()
     962        self._init_cache()
     963        q.reverse()
     964        return self._new_c(q,self._parent)
     965
     966
     967    cdef Py_ssize_t _val_inplace_unit(self):
     968        """
     969        Return `v` the valuation of self and replace self by
     970        self >> v (only for internal use).
     971        """
     972        cdef list a = (<SkewPolynomial_finite_field_dense>self).__coeffs
     973        cdef Py_ssize_t val = 0
     974        if len(a) < 0:
     975            return -1
     976        while a[0].is_zero():
     977            del a[0]
     978            val += 1
     979        self._init_cache()
     980        return val
     981
     982
     983    # lowest common divisor: other functions
     984    # --------------------------------------
     985
     986    cdef SkewPolynomial_finite_field_dense _coeff_llcm(self, SkewPolynomial_finite_field_dense other):
     987        """
     988        (only for internal use)
     989        """
     990        cdef SkewPolynomial_finite_field_dense A = <SkewPolynomial_finite_field_dense>self
     991        cdef SkewPolynomial_finite_field_dense B = <SkewPolynomial_finite_field_dense>other
     992        if not len(A.__coeffs) and not len(B.__coeffs):
     993            raise ZeroDivisionError
     994        A = <SkewPolynomial_finite_field_dense>A._new_c(A.__coeffs[:],A._parent)
     995        B = <SkewPolynomial_finite_field_dense>A._new_c(B.__coeffs[:],B._parent)
     996        cdef parent = self._parent
     997        cdef RingElement one = self.base_ring().one_element()
     998        cdef SkewPolynomial_finite_field_dense U = self._new_c([one],parent)
     999        cdef SkewPolynomial_finite_field_dense V = self._new_c([],parent)
     1000        cdef SkewPolynomial_finite_field_dense Q, R
     1001        cdef list swap
     1002        while len(B.__coeffs):
     1003            Q = A._rquo_inplace_rem(B)
     1004            Q._inplace_rmul(V)
     1005            U._inplace_sub(Q)
     1006            swap = A.__coeffs
     1007            A.__coeffs = B.__coeffs
     1008            B.__coeffs = swap
     1009            swap = U.__coeffs
     1010            U.__coeffs = V.__coeffs
     1011            V.__coeffs = swap
     1012        return V
     1013
     1014
     1015    # Specific functions
     1016    # ------------------
     1017
     1018    cdef Matrix_dense _matphir_c(self):
     1019        r"""
     1020        Return the matrix of the multiplication by `X^r` on
     1021        the quotient `K[X,\sigma] / K[X,\sigma]*self`.
     1022        """
     1023        cdef Py_ssize_t i, j, col, exp, n
     1024        cdef Py_ssize_t d = self.degree()
     1025        cdef Py_ssize_t r = self.parent()._order
     1026        parent = self._parent
     1027        cdef k = parent.base_ring()
     1028        cdef Matrix_dense phir = <Matrix_dense?>zero_matrix(k,d,d)
     1029        cdef RingElement zero = k(0)
     1030        cdef RingElement one = k(1)
     1031        if r < d:
     1032            for i from 0 <= i < d-r:
     1033                phir.set_unsafe(r+i,i,one)
     1034            col = d-r
     1035            exp = d
     1036        else:
     1037            col = 0
     1038            exp = r
     1039        cdef SkewPolynomial_finite_field_dense powx = <SkewPolynomial_finite_field_dense>self._new_c([zero,one], parent)
     1040        cdef SkewPolynomial_finite_field_dense v
     1041        if (exp % 2 == 1):
     1042            v = <SkewPolynomial_finite_field_dense>self._new_c([zero,one], parent)
     1043            if self.degree() == 1:
     1044                v._inplace_rrem(self)
     1045        else:
     1046            v = <SkewPolynomial_finite_field_dense>self._new_c([one], parent)
     1047        exp = exp >> 1
     1048        n = 1
     1049        while exp != 0:
     1050            powx._inplace_lmul(powx.conjugate(n))
     1051            powx._inplace_rrem(self)
     1052            n = n << 1
     1053            if (exp % 2 == 1):
     1054                v = v.conjugate(n)
     1055                v._inplace_rmul(powx)
     1056                v._inplace_rrem(self)
     1057            exp = exp >> 1
     1058        l = v.list()
     1059        for i from 0 <= i < len(l):
     1060            phir.set_unsafe(i,col,l[i])
     1061        for j from col+1 <= j < d:
     1062            v <<= 1
     1063            v = v.conjugate(1)
     1064            v._inplace_rrem(self)
     1065            for i from 0 <= i <= v.degree():
     1066                phir.set_unsafe(i,j,v.__coeffs[i])
     1067        return phir
     1068
     1069    def smurf(self):
     1070        return self._matphir_c()
     1071
     1072
     1073    cdef Matrix_dense _matmul_c(self):
     1074        r"""
     1075        Return the matrix of the multiplication by self on
     1076        `K[X,\sigma]` considered as a free module over `K[X^r]`
     1077        (here `r` is the order of `\sigma`).
     1078
     1079        .. WARNING::
     1080
     1081            Does not work if self is not monic.
     1082        """
     1083        cdef Py_ssize_t i, j, deb, k, r = self.parent()._order
     1084        cdef Py_ssize_t d = self.degree ()
     1085        cdef Ring base_ring = <Ring?>self.parent().base_ring()
     1086        cdef RingElement minusone = <RingElement?>base_ring(-1)
     1087        cdef RingElement zero = <RingElement?>base_ring(0)
     1088        cdef Polk = PolynomialRing (base_ring, 'xr')
     1089        cdef Matrix_dense M = <Matrix_dense?>zero_matrix(Polk,r,r)
     1090        cdef list l = self.list()
     1091        for j from 0 <= j < r:
     1092            for i from 0 <= i < r:
     1093                if i < j:
     1094                    pol = [ zero ]
     1095                    deb = i-j+r
     1096                else:
     1097                    pol = [ ]
     1098                    deb = i-j
     1099                for k from deb <= k <= d by r:
     1100                    pol.append(l[k])
     1101                M.set_unsafe(i,j,Polk(pol))
     1102            for i from 0 <= i <= d:
     1103                l[i] = self._parent.twist_map()(l[i])
     1104        return M
     1105
     1106
     1107    def reduced_norm(self):
     1108        r"""
     1109        Return the reduced norm of this skew polynomial.
     1110
     1111        .. NOTE::
     1112
     1113            The result is cached.
     1114
     1115        ALGORITHM:
     1116
     1117        If `r` (= the order of the twist map) is small compared
     1118        to `d` (= the degree of this skew polynomial), the reduced
     1119        norm is computed as the determinant of the multiplication
     1120        by `P` (= this skew polynomial) acting on `K[X,\sigma]`
     1121        (= the underlying skew ring) viewed as a free module of
     1122        rank `r` over `K[X^r]`.
     1123
     1124        Otherwise, the reduced norm is computed as the characteristic
     1125        polynomial (considered as a polynomial of the variable `X^r`)
     1126        of the left multiplication by `X` on the quotient
     1127        `K[X,\sigma] / K[X,\sigma]*P` (which is a `K`-vector space
     1128        of dimension `d`).
     1129
     1130        EXAMPLES::
     1131
     1132            sage: k.<t> = GF(5^3)
     1133            sage: Frob = k.frobenius_endomorphism()
     1134            sage: S.<x> = k['x',Frob]
     1135            sage: a = S.random_element(degree=3,monic=True); a
     1136            x^3 + (2*t^2 + 3)*x^2 + (4*t^2 + t + 4)*x + 2*t^2 + 2
     1137            sage: N = a.reduced_norm(); N
     1138            (x^3)^3 + 4*(x^3)^2 + 4
     1139
     1140        Note that the parent of `N` is the center of the `S`
     1141        (and not `S` itself)::
     1142
     1143            sage: N.parent()
     1144            Center of Skew Polynomial Ring in x over Finite Field in t of size 5^3 twisted by t |--> t^5:
     1145            Univariate Polynomial Ring in (x^3) over Finite Field of size 5
     1146            sage: N.parent() == S.center()
     1147            True
     1148
     1149        In any case, coercion works fine::
     1150
     1151            sage: S(N)
     1152            x^9 + 4*x^6 + 4
     1153            sage: N + a
     1154            x^9 + 4*x^6 + x^3 + (2*t^2 + 3)*x^2 + (4*t^2 + t + 4)*x + 2*t^2 + 1
     1155
     1156        We check that `N` is a multiple of `a`::
     1157
     1158            sage: S(N).is_divisible_by(a)
     1159            True
     1160            sage: S(N).is_divisible_by(a,side=Left)
     1161            True
     1162
     1163        .. NOTE::
     1164
     1165            We really need to coerce first `N` into `S`. Otherwise an
     1166            error occurs::
     1167
     1168                 sage: N.is_divisible_by(a)
     1169                 Traceback (most recent call last):
     1170                 ...
     1171                 AttributeError: 'sage.rings.polynomial.skew_polynomial_element.CenterSkewPolynomial_generic_dense' object has no attribute 'is_divisible_by'
     1172
     1173        We check that the reduced norm is a multiplicative map::
     1174
     1175            sage: a = S.random_element(degree=5)
     1176            sage: b = S.random_element(degree=7)
     1177            sage: a.reduced_norm() * b.reduced_norm() == (a*b).reduced_norm()
     1178            True
     1179        """
     1180        if self._norm is None:
     1181            center = self.parent().center()
     1182            if self.is_zero():
     1183                self._norm = center(0)
     1184            else:
     1185                section = center._embed_basering.section()
     1186                exp = (self.parent().base_ring().cardinality() - 1) / (center.base_ring().cardinality() - 1)
     1187                order = self.parent()._order
     1188                lc = section(self.leading_coefficient()**exp)
     1189                if order < self.degree():
     1190                    M = self._matmul_c()
     1191                    self._norm = center([ lc*section(x) for x in M.determinant().monic().list() ])
     1192                else:
     1193                    charpoly = self._matphir_c().characteristic_polynomial()
     1194                    self._norm = center([ lc*section(x) for x in charpoly.list() ])
     1195        return self._norm
     1196
     1197
     1198    def reduced_norm_factor(self):
     1199        """
     1200        Return the reduced norm of this polynomial
     1201        factorized in the centre.
     1202
     1203        EXAMPLES:
     1204
     1205            sage: k.<t> = GF(5^3)
     1206            sage: Frob = k.frobenius_endomorphism()
     1207            sage: S.<x> = k['x',Frob]
     1208
     1209            sage: a = (x^2 + 1) * (x+3)
     1210            sage: a.reduced_norm_factor()
     1211            ((x^3) + 3) * ((x^3) + 2)^2
     1212        """
     1213        if self._norm_factor is None:
     1214            self._norm_factor = self.reduced_norm().factor()
     1215        return self._norm_factor
     1216
     1217
     1218    def is_central(self):
     1219        """
     1220        Return True if this skew polynomial lies in the center.
     1221
     1222        EXAMPLES::
     1223
     1224            sage: k.<t> = GF(5^3)
     1225            sage: Frob = k.frobenius_endomorphism()
     1226            sage: S.<x> = k['x',Frob]
     1227
     1228            sage: x.is_central()
     1229            False
     1230            sage: (t*x^3).is_central()
     1231            False
     1232            sage: (x^6 + x^3).is_central()
     1233            True
     1234        """
     1235        center = self.parent().center()
     1236        try:
     1237            center(self)
     1238            return True
     1239        except ValueError:
     1240            return False
     1241
     1242
     1243    def bound(self):
     1244        """
     1245        Return a bound of this skew polynomial (i.e. a multiple
     1246        of this skew polynomial lying in the center).
     1247
     1248        .. NOTE::
     1249
     1250            Since `b` is central, it divides a skew polynomial
     1251            on the left iff it divides it on the right
     1252
     1253        ALGORITHM:
     1254
     1255        #. Sage first checks whether ``self`` is itself in the
     1256           center. It if is, it returns ``self``
     1257
     1258        #. If an optimal bound was previously computed and
     1259           cached, Sage returns it
     1260
     1261        #. Otherwise, Sage returns the reduced norm of ``self``
     1262
     1263        As a consequence, the output of this function may depend
     1264        on previous computations (an example is given below).
     1265
     1266        EXAMPLES::
     1267
     1268            sage: k.<t> = GF(5^3)
     1269            sage: Frob = k.frobenius_endomorphism()
     1270            sage: S.<x> = k['x',Frob]
     1271            sage: Z = S.center()
     1272
     1273            sage: a = x^2 + (4*t + 2)*x + 4*t^2 + 3
     1274            sage: b = a.bound(); b
     1275            (x^3)^2 + (x^3) + 4
     1276
     1277        Note that the parent of `b` is the center of the `S`
     1278        (and not `S` itself)::
     1279
     1280            sage: b.parent()
     1281            Center of Skew Polynomial Ring in x over Finite Field in t of size 5^3 twisted by t |--> t^5:
     1282            Univariate Polynomial Ring in (x^3) over Finite Field of size 5
     1283            sage: b.parent() == Z
     1284            True
     1285
     1286        We check that `b` is divisible by `a`::
     1287
     1288            sage: S(b).is_divisible_by(a)
     1289            True
     1290            sage: S(b).is_divisible_by(a,side=Left)
     1291            True
     1292
     1293        Actually, `b` is the reduced norm of `a`::
     1294
     1295            sage: b == a.reduced_norm()
     1296            True
     1297
     1298        Now, we compute the optimal bound of `a` and see that
     1299        it affects the behaviour of ``bound()``::
     1300
     1301            sage: a.optimal_bound()
     1302            (x^3) + 3
     1303            sage: a.bound()
     1304            (x^3) + 3
     1305
     1306        We finally check that if `a` is a central skew polynomial,
     1307        then ``a.bound()`` returns simply `a`::
     1308
     1309            sage: a = S(Z.random_element(degree=4)); a
     1310            2*x^12 + x^9 + 2*x^3
     1311            sage: b = a.bound(); b
     1312            2*(x^3)^4 + (x^3)^3 + 2*(x^3)
     1313            sage: a == b
     1314            True
     1315        """
     1316        center = self.parent().center()
     1317        try:
     1318            return center(self)
     1319        except ValueError:
     1320            pass
     1321        if not self._optbound is None:
     1322            return center(self._optbound)
     1323        return self.reduced_norm()
     1324
     1325
     1326    def optimal_bound(self):
     1327        """
     1328        Return the optimal bound of this skew polynomial (i.e.
     1329        the monic multiple of this skew polynomial of minimal
     1330        degree lying in the center).
     1331
     1332        .. NOTE::
     1333
     1334            The result is cached.
     1335
     1336        EXAMPLES::
     1337
     1338            sage: k.<t> = GF(5^3)
     1339            sage: Frob = k.frobenius_endomorphism()
     1340            sage: S.<x> = k['x',Frob]
     1341            sage: Z = S.center()
     1342
     1343            sage: a = x^2 + (4*t + 2)*x + 4*t^2 + 3
     1344            sage: b = a.optimal_bound(); b
     1345            (x^3) + 3
     1346
     1347        Note that the parent of `b` is the center of the `S`
     1348        (and not `S` itself)::
     1349
     1350            sage: b.parent()
     1351            Center of Skew Polynomial Ring in x over Finite Field in t of size 5^3 twisted by t |--> t^5:
     1352            Univariate Polynomial Ring in (x^3) over Finite Field of size 5
     1353            sage: b.parent() == Z
     1354            True
     1355
     1356        We check that `b` is divisible by `a`::
     1357
     1358            sage: S(b).is_divisible_by(a)
     1359            True
     1360            sage: S(b).is_divisible_by(a,side=Left)
     1361            True
     1362        """
     1363        center = self.parent().center()
     1364        if self._optbound is None:
     1365            try:
     1366                self._optbound = center(self).monic()
     1367            except ValueError:
     1368                bound = self._matphir_c().minimal_polynomial()
     1369                section = center._embed_basering.section()
     1370                self._optbound = [ section(x) for x in bound.list() ]
     1371        return center(self._optbound)
     1372
     1373
     1374    def is_irreducible(self):
     1375        """
     1376        Return true if this skew polynomial is irreducible.
     1377
     1378        EXAMPLES::
     1379
     1380            sage: k.<t> = GF(5^3)
     1381            sage: Frob = k.frobenius_endomorphism()
     1382            sage: S.<x> = k['x',Frob]
     1383
     1384            sage: a = x^2 + t*x + 1
     1385            sage: a.is_irreducible()
     1386            False
     1387            sage: a.factor()
     1388            (x + 4*t^2 + 4*t + 1) * (x + 3*t + 2)
     1389
     1390            sage: a = x^2 + t*x + t + 1
     1391            sage: a.is_irreducible()
     1392            True
     1393            sage: a.factor()
     1394            x^2 + t*x + t + 1
     1395
     1396        Skew polynomials of degree `1` are of course irreducible::
     1397
     1398            sage: a = x + t
     1399            sage: a.is_irreducible()
     1400            True
     1401
     1402        A random irreducible skew polynomial is irreducible::
     1403
     1404            sage: a = S.random_irreducible(degree=4,monic=True); a   # random
     1405            x^4 + (t + 1)*x^3 + (3*t^2 + 2*t + 3)*x^2 + 3*t*x + 3*t
     1406            sage: a.is_irreducible()
     1407            True
     1408
     1409        By convention, constant skew polynomials are not irreducible::
     1410
     1411            sage: S(1).is_irreducible()
     1412            False
     1413            sage: S(0).is_irreducible()
     1414            False
     1415        """
     1416        return self.reduced_norm().is_irreducible()
     1417
     1418
     1419    def type(self,N):
     1420        """
     1421        INPUT:
     1422
     1423        -  ``N`` -- an irreducible polynomial in the
     1424           center of the underlying skew polynomial ring
     1425
     1426        OUTPUT:
     1427
     1428        The `N`-type of this skew polynomial
     1429
     1430        .. NOTE::
     1431
     1432            The result is cached.
     1433
     1434        DEFINITION:
     1435
     1436        The `N`-type of a skew polynomial `a` is the Partition
     1437        `(t_0, t_1, t_2, ...)` defined by
     1438
     1439        .. MATH::
     1440
     1441            t_0 + \cdots + t_i = \frac{\deg gcd(a,N^i)}{\deg N}
     1442
     1443        where `\deg N` is the degree of `N` considered as an
     1444        element in the center.
     1445
     1446        This notion has an important mathematic interest because
     1447        it corresponds to the Jordan type of the `N`-typical part
     1448        of the associated Galois representation.
     1449
     1450        EXAMPLES::
     1451
     1452            sage: k.<t> = GF(5^3)
     1453            sage: Frob = k.frobenius_endomorphism()
     1454            sage: S.<x> = k['x',Frob]
     1455            sage: Z = S.center(); x3 = Z.gen()
     1456
     1457            sage: a = x^4 + x^3 + (4*t^2 + 4)*x^2 + (t^2 + 2)*x + 4*t^2
     1458            sage: N = x3^2 + x3 + 1
     1459            sage: a.type(N)
     1460            [1]
     1461            sage: N = x3 + 1
     1462            sage: a.type(N)
     1463            [2]
     1464
     1465            sage: a = x^3 + (3*t^2 + 1)*x^2 + (3*t^2 + t + 1)*x + t + 1
     1466            sage: N = x3 + 1
     1467            sage: a.type(N)
     1468            [2, 1]
     1469
     1470        If `N` does not divide the reduced map of `a`, the type
     1471        is empty::
     1472
     1473            sage: N = x3 + 2
     1474            sage: a.type(N)
     1475            []
     1476
     1477        If `a = N`, the type is just `[r]` where `r` is the order
     1478        of the twist map ``Frob``::
     1479
     1480            sage: N = x3^2 + x3 + 1
     1481            sage: S(N).type(N)
     1482            [3]
     1483
     1484        `N` must be irreducible::
     1485
     1486            sage: N = (x3 + 1) * (x3 + 2)
     1487            sage: a.type(N)
     1488            Traceback (most recent call last):
     1489            ...
     1490            ValueError: N is not irreducible
     1491        """
     1492        try:
     1493            return self._types[N]
     1494        except (KeyError, TypeError):
     1495            if not N.is_irreducible():
     1496                raise ValueError("N is not irreducible")
     1497            skew_ring = self._parent
     1498            if self._norm_factor is None:
     1499                m = -1
     1500            else:
     1501                i = [ n for n,_ in self._norm_factor ].index(N)
     1502                m = self._norm_factor[i][1]
     1503            NS = skew_ring(N)
     1504            type = [ ]
     1505            degN = N.degree()
     1506            while True:
     1507                d = self.gcd(NS)
     1508                deg = d.degree()/degN
     1509                if deg == 0:
     1510                    break
     1511                if m >= 0:
     1512                    if deg == 1:
     1513                        type += m * [1]
     1514                        break
     1515                    m -= deg
     1516                self = self // d
     1517                type.append(deg)
     1518            type = Partition(type)
     1519            if self._types is None:
     1520                self._types = { N: type }
     1521            else:
     1522                self._types[N] = type
     1523            return type
     1524
     1525
     1526    # Finding divisors
     1527    # ----------------
     1528
     1529    cdef SkewPolynomial_finite_field_dense _rdivisor_c(P, CenterSkewPolynomial_generic_dense N):
     1530        """
     1531        cython procedure computing an irreducible monic right divisor
     1532        of `P` whose reduced norm is `N`
     1533
     1534        .. WARNING::
     1535
     1536            `N` needs to be an irreducible factor of the
     1537            reduced norm of `P`. This function does not check
     1538            this (and his behaviour is not defined if the
     1539            require property doesn't hold).
     1540        """
     1541        cdef skew_ring = P._parent
     1542        cdef Py_ssize_t d = N.degree()
     1543        cdef Py_ssize_t e = P.degree()/d
     1544        cdef SkewPolynomial_finite_field_dense D
     1545        if e == 1:
     1546            D = <SkewPolynomial_finite_field_dense>P._new_c(list(P.__coeffs),skew_ring)
     1547            D._inplace_rmonic()
     1548            return D
     1549
     1550        E = N.parent().base_ring().extension(N,name='xr')
     1551        PE = PolynomialRing(E,name='T')
     1552        cdef Integer exp
     1553        if skew_ring.characteristic() != 2:
     1554            exp = Integer((E.cardinality()-1)/2)
     1555        cdef SkewPolynomial_finite_field_dense NS = <SkewPolynomial_finite_field_dense>skew_ring(N)
     1556        cdef SkewPolynomial_finite_field_dense Q = <SkewPolynomial_finite_field_dense>(NS // P)
     1557        cdef SkewPolynomial_finite_field_dense R, X
     1558        cdef Matrix_dense M = <Matrix_dense?>MatrixSpace(E,e,e)(0)
     1559        cdef Matrix_dense V = <Matrix_dense?>MatrixSpace(E,e,1)(0)
     1560        cdef Matrix_dense W
     1561        cdef Py_ssize_t i, j, t, r = skew_ring._order
     1562        cdef Polynomial dd, xx, yy, zz
     1563
     1564        while 1:
     1565            R = <SkewPolynomial_finite_field_dense>skew_ring.random_element((e*r-1,e*r-1))
     1566            R._inplace_lmul(Q)
     1567            X = <SkewPolynomial_finite_field_dense>Q._new_c(Q.__coeffs[:],Q._parent)
     1568            for j from 0 <= j < e:
     1569                for i from 0 <= i < e:
     1570                    M.set_unsafe(i, j, E([skew_ring._retraction(X[t*r+i]) for t in range(d)]))
     1571                X._inplace_lmul(R)
     1572                X._inplace_rrem(NS)
     1573            for i from 0 <= i < e:
     1574                V.set_unsafe(i, 0, E([skew_ring._retraction(X[t*r+i]) for t in range(d)]))
     1575            W = M._solve_right_nonsingular_square(V)
     1576            if M*W != V:
     1577                skew_ring._new_retraction_map()
     1578                continue
     1579            xx = PE(W.list()+[E(-1)])
     1580            if skew_ring.characteristic() == 2:
     1581                yy = PE.gen()
     1582                zz = PE.gen()
     1583                for i from 1 <= i < d:
     1584                    zz = (zz*zz) % xx
     1585                    yy += zz
     1586                dd = xx.gcd(yy)
     1587                if dd.degree() != 1: continue
     1588            else:
     1589                yy = PE.gen().__pow__(exp,xx) - 1
     1590                dd = xx.gcd(yy)
     1591                if dd.degree() != 1:
     1592                    yy += 2
     1593                    dd = xx.gcd(yy)
     1594                    if dd.degree() != 1: continue
     1595            D = P._rgcd(R + skew_ring.center()((dd[0]/dd[1]).list()))
     1596            if D.degree() == 0:
     1597                continue
     1598            D._inplace_rmonic()
     1599            D._init_cache()
     1600            return D
     1601
     1602
     1603    def irreducible_divisor(self,side=Right,distribution=None):
     1604        """
     1605        INPUT:
     1606
     1607        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1608
     1609        -  ``distribution`` -- None (default) or ``uniform``
     1610
     1611           - None: no particular specification
     1612
     1613           - ``uniform``: the returned irreducible divisor is
     1614             uniformly distributed
     1615
     1616        .. NOTE::
     1617
     1618            ``uniform`` is a little bit slower.
     1619
     1620        OUTPUT:
     1621
     1622        -  an irreducible monic ``side`` divisor of ``self``
     1623
     1624        EXAMPLES::
     1625
     1626            sage: k.<t> = GF(5^3)
     1627            sage: Frob = k.frobenius_endomorphism()
     1628            sage: S.<x> = k['x',Frob]
     1629            sage: a = x^6 + 3*t*x^5 + (3*t + 1)*x^3 + (4*t^2 + 3*t + 4)*x^2 + (t^2 + 2)*x + 4*t^2 + 3*t + 3
     1630
     1631            sage: dr = a.irreducible_divisor(); dr  # random
     1632            x^3 + (2*t^2 + t + 4)*x^2 + (4*t + 1)*x + 4*t^2 + t + 1
     1633            sage: a.is_divisible_by(dr)
     1634            True
     1635
     1636            sage: dl = a.irreducible_divisor(side=Left); dl  # random
     1637            x^3 + (2*t^2 + t + 1)*x^2 + (4*t^2 + 3*t + 3)*x + 4*t^2 + 2*t + 1
     1638            sage: a.is_divisible_by(dl,side=Left)
     1639            True
     1640
     1641        Right divisors are cached. Hence, if we ask again for a
     1642        right divisor, we will get the same answer::
     1643
     1644            sage: a.irreducible_divisor()  # random
     1645            x^3 + (2*t^2 + t + 4)*x^2 + (4*t + 1)*x + 4*t^2 + t + 1
     1646
     1647        However the algorithm is probabilistic. Hence, if we first
     1648        reinitialiaze `a`, we may get a different answer::
     1649
     1650            sage: a = x^6 + 3*t*x^5 + (3*t + 1)*x^3 + (4*t^2 + 3*t + 4)*x^2 + (t^2 + 2)*x + 4*t^2 + 3*t + 3
     1651            sage: a.irreducible_divisor()  # random
     1652            x^3 + (t^2 + 3*t + 4)*x^2 + (t + 2)*x + 4*t^2 + t + 1
     1653
     1654        We can also generate uniformly distributed irreducible monic
     1655        divisors as follows::
     1656
     1657            sage: a.irreducible_divisor(distribution="uniform")  # random
     1658            x^3 + (4*t + 2)*x^2 + (2*t^2 + 2*t + 2)*x + 2*t^2 + 2
     1659            sage: a.irreducible_divisor(distribution="uniform")  # random
     1660            x^3 + (t^2 + 2)*x^2 + (3*t^2 + 1)*x + 4*t^2 + 2*t
     1661            sage: a.irreducible_divisor(distribution="uniform")  # random
     1662            x^3 + x^2 + (4*t^2 + 2*t + 4)*x + t^2 + 3
     1663
     1664        By convention, the zero skew polynomial has no irreducible
     1665        divisor:
     1666
     1667            sage: S(0).irreducible_divisor()
     1668            Traceback (most recent call last):
     1669            ...
     1670            ValueError: 0 has no irreducible divisor
     1671        """
     1672        if self.is_zero():
     1673            raise ValueError("0 has no irreducible divisor")
     1674        if not (distribution is None or distribution == "uniform"):
     1675            raise ValueError("distribution must be None or 'uniform'")
     1676        if distribution == "uniform":
     1677            skew_ring = self._parent
     1678            center = skew_ring.center()
     1679            cardcenter = center.base_ring().cardinality()
     1680            gencenter = center.gen()
     1681            count = [ ]
     1682            total = 0
     1683            F = self.reduced_norm_factor()
     1684            for n,_ in F:
     1685                if n == gencenter:
     1686                    total += 1
     1687                else:
     1688                    degn = n.degree()
     1689                    P = self.gcd(skew_ring(n))
     1690                    m = P.degree()/degn
     1691                    cardL = cardcenter**degn
     1692                    total += (cardL**m - 1)/(cardL - 1)
     1693                count.append(total)
     1694            if total == 0:
     1695                raise ValueError("No irreducible divisor having given reduced norm")
     1696            random = ZZ.random_element(total)
     1697            for i in range(len(F)):
     1698                if random < count[i]:
     1699                    N = F[i][0]
     1700                    break
     1701        else:
     1702            N = self.reduced_norm_factor()[0][0]
     1703        return self.irreducible_divisor_with_norm(N,side=side,distribution=distribution)
     1704
     1705
     1706    def irreducible_divisor_with_norm(self,N,side=Right,distribution=None): # Ajouter side
     1707        """
     1708        INPUT:
     1709
     1710        -  ``N`` -- an irreducible polynomial in the center
     1711           of the underlying skew polynomial ring
     1712
     1713        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1714
     1715        -  ``distribution`` -- None (default) or ``uniform``
     1716
     1717           - None: no particular specification
     1718
     1719           - ``uniform``: the returned irreducible divisor is
     1720             uniformly distributed
     1721
     1722        .. NOTE::
     1723
     1724            ``uniform`` is a little bit slower.
     1725
     1726        OUTPUT:
     1727
     1728        -  an irreducible monic ``side`` divisor of ``self``
     1729           whose reduced norm is similar to `N` (i.e. `N` times
     1730           a unit).
     1731
     1732        EXAMPLES::
     1733
     1734            sage: k.<t> = GF(5^3)
     1735            sage: Frob = k.frobenius_endomorphism()
     1736            sage: S.<x> = k['x',Frob]
     1737            sage: Z = S.center(); x3 = Z.gen()
     1738            sage: a = x^6 + 3*x^3 + 2
     1739
     1740            sage: d1 = a.irreducible_divisor_with_norm(x3+1); d1   # random
     1741            x + t^2 + 3*t
     1742            sage: a.is_divisible_by(d1)
     1743            True
     1744            sage: d1.reduced_norm()
     1745            (x^3) + 1
     1746
     1747            sage: d2 = a.irreducible_divisor_with_norm(x3+2); d2   # random
     1748            x + 2*t^2 + 3*t + 2
     1749            sage: a.is_divisible_by(d2)
     1750            True
     1751            sage: d2.reduced_norm()
     1752            (x^3) + 2
     1753
     1754            sage: d3 = a.irreducible_divisor_with_norm(x3+3)
     1755            Traceback (most recent call last):
     1756            ...
     1757            ValueError: No irreducible divisor having given reduced norm
     1758
     1759        We can also generate uniformly distributed irreducible monic
     1760        divisors as follows::
     1761
     1762            sage: a.irreducible_divisor_with_norm(x3+1,distribution="uniform")   # random
     1763            x + 3*t^2 + 3*t + 1
     1764            sage: a.irreducible_divisor_with_norm(x3+1,distribution="uniform")   # random
     1765            x + 1
     1766            sage: a.irreducible_divisor_with_norm(x3+1,distribution="uniform")   # random
     1767            x + 2*t^2 + 4*t
     1768        """
     1769        cdef SkewPolynomial_finite_field_dense cP1
     1770        cdef CenterSkewPolynomial_generic_dense cN
     1771        if self.is_zero():
     1772            raise "No irreducible divisor having given reduced norm"
     1773        skew_ring = self._parent
     1774        center = skew_ring.center()
     1775        try:
     1776            N = center(N)
     1777        except TypeError:
     1778            raise TypeError("N must be a polynomial in the center")
     1779        cardcenter = center.base_ring().cardinality()
     1780        gencenter = center.gen()
     1781
     1782        if N == gencenter:
     1783            if self[0] == 0:
     1784                return skew_ring.gen()
     1785            else:
     1786                raise ValueError("No irreducible divisor having given reduced norm")
     1787
     1788        D = None
     1789        try:
     1790            D = self._rdivisors[N]
     1791        except (KeyError, TypeError):
     1792            if N.is_irreducible():
     1793                cP1 = <SkewPolynomial_finite_field_dense>self._rgcd(self._parent(N))
     1794                cN = <CenterSkewPolynomial_generic_dense>N
     1795                if cP1.degree() > 0:
     1796                    D = cP1._rdivisor_c(cN)
     1797            if self._rdivisors is None:
     1798                self._rdivisors = { N: D }
     1799            else:
     1800                self._rdivisors[N] = D
     1801            distribution = ""
     1802        if D is None:
     1803            raise ValueError("No irreducible divisor having given reduced norm")
     1804
     1805        NS = self._parent(N)
     1806        degN = N.degree()
     1807        if side is Right:
     1808            if distribution == "uniform":
     1809                P1 = self._rgcd(NS)
     1810                if P1.degree() != degN:
     1811                    Q1 = NS // P1
     1812                    deg = P1.degree()-1
     1813                    while True:
     1814                        R = Q1*skew_ring.random_element((deg,deg))
     1815                        if P1.gcd(R) == 1:
     1816                            break
     1817                    D = P1.gcd(D*R)
     1818            return D
     1819        else:
     1820            deg = NS.degree()-1
     1821            P1 = self.lgcd(NS)
     1822            while True:
     1823                if distribution == "uniform":
     1824                    while True:
     1825                        R = skew_ring.random_element((deg,deg))
     1826                        if NS.gcd(R) == 1:
     1827                            break
     1828                    D = NS.gcd(D*R)
     1829                Dp = NS // D
     1830                LDp = P1.gcd(Dp)
     1831                LD = P1 // LDp
     1832                if LD.degree() == degN:
     1833                    return LD
     1834                distribution = "uniform"
     1835
     1836
     1837    def irreducible_divisors(self,side=Right):
     1838        """
     1839        INPUT:
     1840
     1841        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1842
     1843        OUTPUT:
     1844
     1845        An iterator over all irreducible monic ``side`` divisors
     1846        of this skew polynomial
     1847
     1848        EXAMPLES:
     1849
     1850            sage: k.<t> = GF(5^3)
     1851            sage: Frob = k.frobenius_endomorphism()
     1852            sage: S.<x> = k['x',Frob]
     1853            sage: a = x^4 + 2*t*x^3 + 3*t^2*x^2 + (t^2 + t + 1)*x + 4*t + 3
     1854            sage: iter = a.irreducible_divisors(); iter
     1855            <generator object at 0x...>
     1856            sage: iter.next()   # random
     1857            x + 2*t^2 + 4*t + 4
     1858            sage: iter.next()   # random
     1859            x + 3*t^2 + 4*t + 1
     1860
     1861        We can use this function to build the list of all monic
     1862        irreducible divisors of `a`::
     1863
     1864            sage: rightdiv = [ d for d in a.irreducible_divisors() ]
     1865            sage: leftdiv = [ d for d in a.irreducible_divisors(side=Left) ]
     1866
     1867        We do some checks::
     1868
     1869            sage: len(rightdiv) == a.count_irreducible_divisors()
     1870            True
     1871            sage: len(rightdiv) == len(Set(rightdiv))  # check no duplicates
     1872            True
     1873            sage: for d in rightdiv:
     1874            ...       if not a.is_divisible_by(d):
     1875            ...           print "Found %s which is not a right divisor" % d
     1876            ...       elif not d.is_irreducible():
     1877            ...           print "Found %s which is not irreducible" % d
     1878
     1879            sage: len(leftdiv) == a.count_irreducible_divisors(side=Left)
     1880            True
     1881            sage: len(leftdiv) == len(Set(leftdiv))  # check no duplicates
     1882            True
     1883            sage: for d in leftdiv:
     1884            ...       if not a.is_divisible_by(d,side=Left):
     1885            ...           print "Found %s which is not a left divisor" % d
     1886            ...       elif not d.is_irreducible():
     1887            ...           print "Found %s which is not irreducible" % d
     1888
     1889        Note that left divisors and right divisors differ::
     1890
     1891            sage: Set(rightdiv) == Set(leftdiv)
     1892            False
     1893
     1894        Note that the algorithm is probabilistic. As a consequence, if we
     1895        build again the list of right monic irreducible divisors of `a`, we
     1896        may get a different ordering::
     1897
     1898            sage: rightdiv2 = [ d for d in a.irreducible_divisors() ]
     1899            sage: rightdiv == rightdiv2
     1900            False
     1901            sage: Set(rightdiv) == Set(rightdiv2)
     1902            True
     1903        """
     1904        return self._irreducible_divisors(side)
     1905
     1906
     1907    def _irreducible_divisors(self,side): # prendre side en compte
     1908        """
     1909        Return an iterator over all irreducible monic
     1910        divisors of this skew polynomial.
     1911
     1912        Do not use this function. Use instead
     1913        ``self.irreducible_divisors()``.
     1914        """
     1915        if self.is_zero():
     1916            return
     1917        skew_ring = self._parent
     1918        center = skew_ring.center()
     1919        kfixed = center.base_ring()
     1920        F = self.reduced_norm_factor()
     1921        oppside = side.opposite()
     1922        for N,_ in F:
     1923            if N == center.gen():
     1924                yield skew_ring.gen()
     1925                continue
     1926            degN = N.degree()
     1927            NS = skew_ring(N)
     1928            P = self.gcd(NS,side=side)
     1929            m = P.degree()/degN
     1930            if m == 1:
     1931                yield P
     1932                continue
     1933            degrandom = P.degree() - 1
     1934            Q,_ = NS.quo_rem(P,side=side)
     1935            P1 = self.irreducible_divisor_with_norm(N,side=side)
     1936            Q1,_ = P.quo_rem(P1,side=side)
     1937            while True:
     1938                R = skew_ring.random_element((degrandom,degrandom))
     1939                if side is Right:
     1940                    g = (R*Q).rem(P,side=Left)
     1941                else:
     1942                    g = (Q*R).rem(P)
     1943                if g.gcd(P,side=oppside) != 1: continue
     1944                L = Q1
     1945                V = L
     1946                for i in range(1,m):
     1947                    if side is Right:
     1948                        L = (g*L).gcd(P,side=Left)
     1949                    else:
     1950                        L = (L*g).gcd(P)
     1951                    V = V.gcd(L,side=oppside)
     1952                if V == 1: break
     1953            rng = xmrange_iter([kfixed]*degN,center)
     1954            for i in range(m):
     1955                for pol in xmrange_iter([rng]*i):
     1956                    f = skew_ring(1)
     1957                    for j in range(i):
     1958                        coeff = pol.pop()
     1959                        f = (g*f+coeff).rem(P,side=oppside)
     1960                    if side is Right:
     1961                        d = (f*Q1).gcd(P,side=Left)
     1962                    else:
     1963                        d = (Q1*f).gcd(P)
     1964                    d,_ = P.quo_rem(d,side=oppside)
     1965                    yield d
     1966
     1967
     1968    def count_irreducible_divisors(self,side=Right):
     1969        """
     1970        INPUT:
     1971
     1972        -  ``side`` -- ``Left`` or ``Right`` (default: Right)
     1973
     1974        OUTPUT:
     1975
     1976        The number of irreducible monic ``side`` divisors of
     1977        this skew polynomial.
     1978
     1979        .. NOTE::
     1980
     1981            Actually, one can prove that there are always as
     1982            many left irreducible monic divisors as right
     1983            irreducible monic divisors.
     1984
     1985        EXAMPLES::
     1986
     1987            sage: k.<t> = GF(5^3)
     1988            sage: Frob = k.frobenius_endomorphism()
     1989            sage: S.<x> = k['x',Frob]
     1990
     1991        We illustrate that a skew polynomial may have a number of irreducible
     1992        divisors greater than its degree.
     1993
     1994            sage: a = x^4 + (4*t + 3)*x^3 + t^2*x^2 + (4*t^2 + 3*t)*x + 3*t
     1995            sage: a.count_irreducible_divisors()
     1996            12
     1997            sage: a.count_irreducible_divisors(side=Left)
     1998            12
     1999
     2000        We illustrate that an irreducible polynomial in the center have
     2001        in general a lot of irreducible divisors in the skew polynomial
     2002        ring::
     2003
     2004            sage: Z = S.center(); x3 = Z.gen()
     2005            sage: N = x3^5 + 4*x3^4 + 4*x3^2 + 4*x3 + 3; N
     2006            (x^3)^5 + 4*(x^3)^4 + 4*(x^3)^2 + 4*(x^3) + 3
     2007            sage: N.is_irreducible()
     2008            True
     2009            sage: S(N).count_irreducible_divisors()
     2010            9768751
     2011        """
     2012        if self.is_zero():
     2013            return 0
     2014        skew_ring = self.parent()
     2015        cardcenter = skew_ring.center().base_ring().cardinality()
     2016        gencenter = skew_ring.center().gen()
     2017        F = self.reduced_norm_factor()
     2018        val = self.valuation()
     2019        self >>= val
     2020        count = 0
     2021        if val > 0:
     2022            count = 1
     2023        for N,_ in F:
     2024            if N == gencenter:
     2025                continue
     2026            degN = N.degree()
     2027            P = self.gcd(skew_ring(N), side=side)
     2028            m = P.degree()/degN
     2029            cardL = cardcenter**degN
     2030            count += (cardL**m - 1)/(cardL - 1)
     2031        return count
     2032
     2033
     2034    # Finding factorizations
     2035    # ----------------------
     2036
     2037    cdef _factor_c(self):
     2038        """
     2039        Compute a factorization of ``self``
     2040        """
     2041        cdef skew_ring = self._parent
     2042        cdef Py_ssize_t degQ, degrandom, m, mP, i
     2043        cdef CenterSkewPolynomial_generic_dense N
     2044        cdef SkewPolynomial_finite_field_dense poly = <SkewPolynomial_finite_field_dense>self.rmonic()
     2045        cdef val = poly._val_inplace_unit()
     2046        if val == -1:
     2047            return Factorization([], sort=False, unit=skew_ring.zero_element())
     2048        cdef list factors = [ (skew_ring.gen(), val) ]
     2049        cdef SkewPolynomial_finite_field_dense P, Q, P1, NS, g, right, Pn
     2050        cdef SkewPolynomial_finite_field_dense right2 = skew_ring(1) << val
     2051        cdef RingElement unit = <RingElement>self.leading_coefficient()
     2052        cdef Polynomial gencenter = skew_ring.center().gen()
     2053        cdef Py_ssize_t p = skew_ring.characteristic()
     2054        cdef F = self.reduced_norm_factor()
     2055
     2056        for N,m in F:
     2057            if N == gencenter:
     2058                continue
     2059            degN = N.degree()
     2060            if poly.degree() == degN:
     2061                factors.append((poly,1))
     2062                break
     2063            NS = <SkewPolynomial_finite_field_dense>skew_ring(N)
     2064            P1 = None
     2065            while 1:
     2066                P = <SkewPolynomial_finite_field_dense>poly._rgcd(NS)
     2067                P._inplace_rmonic()
     2068                mP = P.degree() / degN
     2069                if mP == 0: break
     2070                if mP == 1:
     2071                    factors.append((P,1))
     2072                    poly._inplace_rfloordiv(P)
     2073                    for i from 1 <= i < m:
     2074                        if poly.degree() == degN:
     2075                            factors.append((poly,1))
     2076                            break
     2077                        P = poly._rgcd(NS)
     2078                        P._inplace_rmonic()
     2079                        factors.append((P,1))
     2080                        poly._inplace_rfloordiv(P)
     2081                    break
     2082                if P1 is None:
     2083                    P1 = P._rdivisor_c(N)
     2084                Q = <SkewPolynomial_finite_field_dense>NS._new_c(NS.__coeffs[:], NS._parent)
     2085                Q._inplace_rfloordiv(P)
     2086                Q._inplace_lmul(P1)
     2087                factors.append((P1,1))
     2088                right = <SkewPolynomial_finite_field_dense>P1._new_c(P1.__coeffs[:], P1._parent)
     2089                m -= (mP-1)
     2090                degrandom = P.degree()
     2091                while mP > 2:
     2092                    while 1:
     2093                        g = <SkewPolynomial_finite_field_dense>skew_ring.random_element((degrandom,degrandom))
     2094                        g._inplace_lmul(Q)
     2095                        g._inplace_rgcd(P)
     2096                        Pn = right._coeff_llcm(g)
     2097                        if len(Pn.__coeffs)-1 == degN: break
     2098                    Pn._inplace_rmonic()
     2099                    factors.append((Pn,1))
     2100                    right._inplace_lmul(Pn)
     2101                    degrandom -= degN
     2102                    mP -= 1
     2103                poly._inplace_rfloordiv(right)
     2104                P1,_ = P.rquo_rem(right)
     2105        factors.reverse()
     2106        return Factorization(factors, sort=False, unit=unit)
     2107
     2108
     2109    cdef _factor_uniform_c(self):
     2110        """
     2111        Compute a uniformly distrbuted factorization of ``self``
     2112        """
     2113        skew_ring = self._parent
     2114        cdef Integer cardE, cardcenter = skew_ring.center().base_ring().cardinality()
     2115        cdef CenterSkewPolynomial_generic_dense gencenter = <CenterSkewPolynomial_generic_dense>skew_ring.center().gen()
     2116        cdef SkewPolynomial_finite_field_dense gen = <SkewPolynomial_finite_field_dense>skew_ring.gen()
     2117
     2118        cdef list factorsN = [ ]
     2119        cdef dict dict_divisor = { }
     2120        cdef dict dict_type = { }
     2121        cdef dict dict_right = { }
     2122        cdef CenterSkewPolynomial_generic_dense N
     2123        cdef Py_ssize_t m
     2124        cdef list type
     2125
     2126        for N,m in self.reduced_norm_factor():
     2127            factorsN += m * [N]
     2128            if N == gencenter: continue
     2129            type = list(self.type(N))
     2130            dict_type[N] = type
     2131            if type[0] > 1:
     2132                dict_divisor[N] = self.irreducible_divisor_with_norm(N)
     2133                dict_right[N] = skew_ring(1)
     2134        cdef list indices = list(Permutations(len(factorsN)).random_element())
     2135
     2136        cdef RingElement unit = self.leading_coefficient()
     2137        cdef SkewPolynomial_finite_field_dense left = self._new_c(self.__coeffs[:],skew_ring)
     2138        left._inplace_rmonic()
     2139        cdef SkewPolynomial_finite_field_dense right = <SkewPolynomial_finite_field_dense>skew_ring(1)
     2140        cdef SkewPolynomial_finite_field_dense L, R
     2141        cdef SkewPolynomial_finite_field_dense NS, P, Q, D, D1, D2, d
     2142        cdef list factors = [ ]
     2143        cdef list maxtype
     2144        cdef Py_ssize_t i, j, degN, deg
     2145        cdef count, maxcount
     2146
     2147        for i in indices:
     2148            N = factorsN[i-1]
     2149            if N == gencenter:
     2150                D1 = gen
     2151            else:
     2152                type = dict_type[N]
     2153                NS = skew_ring(N)
     2154                P = left.gcd(NS)
     2155                if type[0] == 1:
     2156                    D1 = P
     2157                else:
     2158                    R = right._new_c(right.__coeffs[:],skew_ring)
     2159                    R._inplace_rfloordiv(dict_right[N])
     2160                    D = R._coeff_llcm(dict_divisor[N])
     2161                    maxtype = list(type)
     2162                    maxtype[-1] -= 1
     2163                    degN = N.degree()
     2164                    cardE = cardcenter ** degN
     2165                    maxcount = q_jordan(Partition(maxtype),cardE)
     2166                    Q = NS // P
     2167                    deg = P.degree()-1
     2168                    while 1:
     2169                        while 1:
     2170                            R = <SkewPolynomial_finite_field_dense>skew_ring.random_element((deg,deg))
     2171                            R._inplace_lmul(Q)
     2172                            if P._rgcd(R).degree() == 0:
     2173                                break
     2174                        D1 = P._rgcd(D*R)
     2175                        D1._inplace_rmonic()
     2176
     2177                        L = left._new_c(list(left.__coeffs),skew_ring)
     2178                        L._inplace_rfloordiv(D1)
     2179                        degN = N.degree()
     2180                        for j in range(len(type)):
     2181                            if type[j] == 1:
     2182                                newtype = type[:-1]
     2183                                break
     2184                            d = L._rgcd(NS)
     2185                            d._inplace_rmonic()
     2186                            deg = d.degree() / degN
     2187                            if deg < type[j]:
     2188                                newtype = type[:]
     2189                                newtype[j] = deg
     2190                                break
     2191                            L._inplace_rfloordiv(d)
     2192                        count = q_jordan(Partition(newtype),cardE)
     2193                        if ZZ.random_element(maxcount) < count:
     2194                            break
     2195                    dict_type[N] = newtype
     2196
     2197                    D2 = D._new_c(list(D.__coeffs),skew_ring)
     2198                    D2._inplace_rmonic()
     2199                    while D2 == D1:
     2200                        while 1:
     2201                            R = <SkewPolynomial_finite_field_dense>skew_ring.random_element((deg,deg))
     2202                            R._inplace_lmul(Q)
     2203                            if P._rgcd(R).degree() == 0:
     2204                                break
     2205                        D2 = P._rgcd(D*R)
     2206                        D2._inplace_rmonic()
     2207                    dict_divisor[N] = D1._coeff_llcm(D2)
     2208            factors.append((D1,1))
     2209            left._inplace_rfloordiv(D1)
     2210            right._inplace_lmul(D1)
     2211            dict_right[N] = right._new_c(list(right.__coeffs),skew_ring)
     2212
     2213        factors.reverse()
     2214        return Factorization(factors,sort=False,unit=unit)
     2215
     2216
     2217    def factor(self,distribution=None):
     2218        """
     2219        Return a factorization of this skew polynomial.
     2220
     2221        INPUT:
     2222
     2223        -  ``distribution`` -- None (default) or ``uniform``
     2224
     2225           - None: no particular specification
     2226
     2227           - ``uniform``: the returned factorization is uniformly
     2228             distributed among all possible factorizations
     2229
     2230        .. NOTE::
     2231
     2232            ``uniform`` is a little bit slower.
     2233
     2234        OUTPUT:
     2235
     2236        -  ``Factorization`` -- a factorization of self as a
     2237           product of a unit and a product of irreducible monic
     2238           factors
     2239
     2240        EXAMPLES::
     2241
     2242            sage: k.<t> = GF(5^3)
     2243            sage: Frob = k.frobenius_endomorphism()
     2244            sage: S.<x> = k['x',Frob]
     2245            sage: a = x^3 + (t^2 + 4*t + 2)*x^2 + (3*t + 3)*x + t^2 + 1
     2246            sage: F = a.factor(); F  # random
     2247            (x + t^2 + 4) * (x + t + 3) * (x + t)
     2248            sage: F.value() == a
     2249            True
     2250
     2251        The result of the factorization is cached. Hence, if we try
     2252        again to factor `a`, we will get the same answer::
     2253
     2254            sage: a.factor()  # random
     2255            (x + t^2 + 4) * (x + t + 3) * (x + t)
     2256
     2257        However, the algorithm is probabilistic. Hence if we first
     2258        reinitialiaze `a`, we may get a different answer::
     2259
     2260            sage: a = x^3 + (t^2 + 4*t + 2)*x^2 + (3*t + 3)*x + t^2 + 1
     2261            sage: F = a.factor(); F   # random
     2262            (x + t^2 + t + 2) * (x + 2*t^2 + t + 4) * (x + t)
     2263            sage: F.value() == a
     2264            True
     2265
     2266        There is no guarantee on the distribution of the factorizations
     2267        we get that way. (For this particular `a` for example, we get the
     2268        uniform distribution on the subset of all factorizations ending
     2269        by the factor `x + t`.)
     2270
     2271        If we rather want uniform distribution among all factorizations,
     2272        we need to specify it as follows::
     2273
     2274            sage: a.factor(distribution="uniform")   # random
     2275            (x + t^2 + 4) * (x + t) * (x + t + 3)
     2276            sage: a.factor(distribution="uniform")   # random
     2277            (x + 2*t^2) * (x + t^2 + t + 1) * (x + t^2 + t + 2)
     2278            sage: a.factor(distribution="uniform")   # random
     2279            (x + 2*t^2 + 3*t) * (x + 4*t + 2) * (x + 2*t + 2)
     2280
     2281        By convention, the zero skew polynomial has no factorization:
     2282
     2283            sage: S(0).factor()
     2284            Traceback (most recent call last):
     2285            ...
     2286            ValueError: factorization of 0 not defined
     2287        """
     2288        if not (distribution is None or distribution == "uniform"):
     2289            raise ValueError("distribution must be None or 'uniform'")
     2290        if self.is_zero():
     2291            raise ValueError("factorization of 0 not defined")
     2292        sig_on()
     2293        if distribution is None:
     2294            if self._factorization is None:
     2295                self._factorization = self._factor_c()
     2296            F = self._factorization
     2297        else:
     2298            F = self._factor_uniform_c()
     2299            if self._factorization is None:
     2300                self._factorization = F
     2301        sig_off()
     2302        return F
     2303
     2304
     2305    def count_factorizations(self):
     2306        """
     2307        Return the number of factorizations (as a product of a
     2308        unit and a product of irreducible monic factors) of this
     2309        skew polynomial.
     2310
     2311        EXAMPLES::
     2312
     2313            sage: k.<t> = GF(5^3)
     2314            sage: Frob = k.frobenius_endomorphism()
     2315            sage: S.<x> = k['x',Frob]
     2316            sage: a = x^4 + (4*t + 3)*x^3 + t^2*x^2 + (4*t^2 + 3*t)*x + 3*t
     2317            sage: a.count_factorizations()
     2318            216
     2319
     2320        We illustrate that an irreducible polynomial in the center have
     2321        in general a lot of distinct factorizations in the skew polynomial
     2322        ring::
     2323
     2324            sage: Z = S.center(); x3 = Z.gen()
     2325            sage: N = x3^5 + 4*x3^4 + 4*x3^2 + 4*x3 + 3; N
     2326            (x^3)^5 + 4*(x^3)^4 + 4*(x^3)^2 + 4*(x^3) + 3
     2327            sage: N.is_irreducible()
     2328            True
     2329            sage: S(N).count_factorizations()
     2330            30537115626
     2331        """
     2332        if self.is_zero():
     2333            raise ValueError("factorization of 0 not defined")
     2334        cardcenter = self._parent.center().base_ring().cardinality()
     2335        gencenter = self._parent.center().gen()
     2336        F = self.reduced_norm_factor()
     2337        summ = 0
     2338        count = 1
     2339        for N,m in F:
     2340            summ += m
     2341            if m == 1: continue
     2342            if N != gencenter:
     2343                count *= q_jordan(self.type(N),cardcenter**N.degree())
     2344            count /= factorial(m)
     2345        return count * factorial(summ)
     2346
     2347    def count_factorisations(self):
     2348        """
     2349        Return the number of factorisations (as a product of a
     2350        unit and a product of irreducible monic factors) of this
     2351        skew polynomial.
     2352
     2353        EXAMPLES::
     2354
     2355            sage: k.<t> = GF(5^3)
     2356            sage: Frob = k.frobenius_endomorphism()
     2357            sage: S.<x> = k['x',Frob]
     2358            sage: a = x^4 + (4*t + 3)*x^3 + t^2*x^2 + (4*t^2 + 3*t)*x + 3*t
     2359            sage: a.count_factorisations()
     2360            216
     2361
     2362        We illustrate that an irreducible polynomial in the center have
     2363        in general a lot of distinct factorisations in the skew polynomial
     2364        ring::
     2365
     2366            sage: Z = S.center(); x3 = Z.gen()
     2367            sage: N = x3^5 + 4*x3^4 + 4*x3^2 + 4*x3 + 3; N
     2368            (x^3)^5 + 4*(x^3)^4 + 4*(x^3)^2 + 4*(x^3) + 3
     2369            sage: N.is_irreducible()
     2370            True
     2371            sage: S(N).count_factorisations()
     2372            30537115626
     2373        """
     2374        return self.count_factorizations()
     2375
     2376
     2377    # Not optimized (many calls to reduced_norm, reduced_norm_factor,_rdivisor_c, which are slow)
     2378    def factorizations(self):
     2379        """
     2380        Return an iterator over all factorizations (as a product
     2381        of a unit and a product of irreducible monic factors) of
     2382        this skew polynomial.
     2383
     2384        .. NOTE::
     2385
     2386            The algorithm is probabilistic. As a consequence, if
     2387            we execute two times with the same input we can get
     2388            the list of all factorizations in two differents orders.
     2389
     2390        EXAMPLES::
     2391
     2392            sage: k.<t> = GF(5^3)
     2393            sage: Frob = k.frobenius_endomorphism()
     2394            sage: S.<x> = k['x',Frob]
     2395            sage: a = x^3 + (t^2 + 1)*x^2 + (2*t + 3)*x + t^2 + t + 2
     2396            sage: iter = a.factorizations(); iter
     2397            <generator object at 0x...>
     2398            sage: iter.next()   # random
     2399            (x + 3*t^2 + 4*t) * (x + 2*t^2) * (x + 4*t^2 + 4*t + 2)
     2400            sage: iter.next()   # random
     2401            (x + 3*t^2 + 4*t) * (x + 3*t^2 + 2*t + 2) * (x + 4*t^2 + t + 2)
     2402
     2403        We can use this function to build the list of factorizations
     2404        of `a`::
     2405
     2406            sage: factorizations = [ F for F in a.factorizations() ]
     2407
     2408        We do some checks::
     2409
     2410            sage: len(factorizations) == a.count_factorizations()
     2411            True
     2412            sage: len(factorizations) == len(Set(factorizations))  # check no duplicates
     2413            True
     2414            sage: for F in factorizations:
     2415            ...       if F.value() != a:
     2416            ...           print "Found %s which is not a correct factorization" % d
     2417            ...           continue
     2418            ...       for d,_ in F:
     2419            ...           if not d.is_irreducible():
     2420            ...               print "Found %s which is not a correct factorization" % d
     2421        """
     2422        if self.is_zero():
     2423            raise ValueError("factorization of 0 not defined")
     2424        unit = self.leading_coefficient()
     2425        poly = self.rmonic()
     2426        for factors in self._factorizations_rec():
     2427            yield Factorization(factors,sort=False,unit=unit)
     2428    def _factorizations_rec(self):
     2429        if self.is_irreducible():
     2430            yield [ (self,1) ]
     2431        else:
     2432            for div in self._irreducible_divisors(Right):
     2433                poly = self // div
     2434                # Here, we should update poly._norm, poly._norm_factor, poly._rdivisors
     2435                for factors in poly._factorizations_rec():
     2436                    factors.append((div,1))
     2437                    yield factors
     2438
     2439
     2440    def factorisations(self):
     2441        """
     2442        Return an iterator over all factorisations (as a product
     2443        of a unit and a product of irreducible monic factors) of
     2444        this skew polynomial.
     2445
     2446        .. NOTE::
     2447
     2448            The algorithm is probabilistic. As a consequence, if
     2449            we execute two times with the same input we can get
     2450            the list of all factorizations in two differents orders.
     2451
     2452        EXAMPLES::
     2453
     2454            sage: k.<t> = GF(5^3)
     2455            sage: Frob = k.frobenius_endomorphism()
     2456            sage: S.<x> = k['x',Frob]
     2457            sage: a = x^3 + (t^2 + 1)*x^2 + (2*t + 3)*x + t^2 + t + 2
     2458            sage: iter = a.factorisations(); iter
     2459            <generator object at 0x...>
     2460            sage: iter.next()   # random
     2461            (x + 3*t^2 + 4*t) * (x + 2*t^2) * (x + 4*t^2 + 4*t + 2)
     2462            sage: iter.next()   # random
     2463            (x + 3*t^2 + 4*t) * (x + 3*t^2 + 2*t + 2) * (x + 4*t^2 + t + 2)
     2464
     2465        We can use this function to build the list of factorizations
     2466        of `a`::
     2467
     2468            sage: factorisations = [ F for F in a.factorisations() ]
     2469
     2470        We do some checks::
     2471
     2472            sage: len(factorisations) == a.count_factorisations()
     2473            True
     2474            sage: len(factorisations) == len(Set(factorisations))  # check no duplicates
     2475            True
     2476            sage: for F in factorisations:
     2477            ...       if F.value() != a:
     2478            ...           print "Found %s which is not a correct factorization" % d
     2479            ...           continue
     2480            ...       for d,_ in F:
     2481            ...           if not d.is_irreducible():
     2482            ...               print "Found %s which is not a correct factorization" % d
     2483        """
     2484        return self.factorizations()
     2485
     2486
     2487# Karatsuba class
     2488#################
     2489
     2490cdef class SkewPolynomial_finite_field_karatsuba:
     2491    """
     2492    A special class implementing Karatsuba multiplication and
     2493    euclidean division on skew polynomial rings over finite fields.
     2494
     2495    Only for internal use!
     2496
     2497    TESTS::
     2498
     2499        sage: k.<t> = GF(5^3)
     2500        sage: Frob = k.frobenius_endomorphism()
     2501        sage: S.<x> = k['x',Frob]
     2502
     2503    An instance of this class is automatically created when the skew
     2504    ring is created. It is accessible via the key work ``_karatsuba_class``::
     2505
     2506        sage: KarClass = S._karatsuba_class; KarClass
     2507        Special class for Karatsuba multiplication and euclidean division on
     2508        Skew Polynomial Ring in x over Finite Field in t of size 5^3 twisted by t |--> t^5
     2509
     2510    Each instance of this class is equipped with a ``cutoff`` variable. The
     2511    default value is the maximum between 150 and the order of the twist map
     2512    acting on the finite field::
     2513
     2514        sage: KarClass.get_cutoff()
     2515        150
     2516
     2517    We can set a new value to this variable as follows::
     2518
     2519        sage: KarClass.set_cutoff(27)
     2520        sage: KarClass.get_cutoff()
     2521        27
     2522    """
     2523    def __init__(self,parent,cutoff=0):
     2524        """
     2525        Initialize a new instance of this class.
     2526        """
     2527        self._parent = parent
     2528        self._order = parent._order
     2529        self._zero = parent.base_ring()(0)
     2530        self.set_cutoff(cutoff)
     2531        self._algo_matrix = 0
     2532        self._t = None
     2533        self._T = None
     2534        self._Tinv = None
     2535
     2536
     2537    def __repr__(self):
     2538        """
     2539        Return a representation of ``self``
     2540        """
     2541        return "Special class for Karatsuba multiplication and euclidean division on\n%s" % self._parent
     2542
     2543
     2544    def set_cutoff(self,cutoff=0):
     2545        """
     2546        Set a new cutoff for all Karatsuba multiplications
     2547        and euclidean divisions performed by this class.
     2548
     2549        INPUT:
     2550
     2551        -  ``cutoff`` -- a nonnegative integer or +Infinity
     2552           (default: 0)
     2553
     2554        .. NOTE::
     2555
     2556            If ``cutoff`` is `0`, the new cutoff is set to the
     2557            default value which is the maximum between 150 and
     2558            the order of the twist map acting on the underlying
     2559            finite field.
     2560
     2561        EXAMPLES::
     2562
     2563            sage: k.<t> = GF(5^3)
     2564            sage: Frob = k.frobenius_endomorphism()
     2565            sage: S.<x> = k['x',Frob]
     2566            sage: KarClass = S._karatsuba_class
     2567
     2568            sage: KarClass.set_cutoff(27)
     2569            sage: KarClass.get_cutoff()
     2570            27
     2571
     2572            sage: KarClass.set_cutoff(Infinity)
     2573            sage: KarClass.get_cutoff()
     2574            +Infinity
     2575
     2576            sage: KarClass.set_cutoff(0)  # set the default value
     2577            sage: KarClass.get_cutoff()
     2578            150
     2579
     2580        The cutoff can't be less than the order of the twist map::
     2581
     2582            sage: KarClass.set_cutoff(2)
     2583            Traceback (most recent call last):
     2584            ...
     2585            ValueError: cutoff must be 0 or >= 3
     2586        """
     2587        if cutoff == 0:
     2588            self._cutoff = max(150,self._order)
     2589        else:
     2590            if cutoff < self._order:
     2591                raise ValueError("cutoff must be 0 or >= %s" % self._order)
     2592            if cutoff == Infinity:
     2593                self._cutoff = -1
     2594            else:
     2595                self._cutoff = cutoff
     2596
     2597
     2598    def get_cutoff(self):
     2599        """
     2600        Return the current cutoff
     2601
     2602        EXAMPLES::
     2603
     2604            sage: k.<t> = GF(5^3)
     2605            sage: Frob = k.frobenius_endomorphism()
     2606            sage: S.<x> = k['x',Frob]
     2607            sage: KarClass = S._karatsuba_class
     2608            sage: KarClass.get_cutoff()
     2609            150
     2610            sage: KarClass.set_cutoff(27)
     2611            sage: KarClass.get_cutoff()
     2612            27
     2613        """
     2614        if self._cutoff < 0:
     2615            return Infinity
     2616        else:
     2617            return self._cutoff
     2618
     2619
     2620    def mul(self,left,right):
     2621        """
     2622        INPUT:
     2623
     2624        -  ``left`` -- a skew polynomial in the skew polynomial
     2625           right attached to the class
     2626
     2627        -  ``right`` -- a skew polynomial in the skew polynomial
     2628           right attached to the class
     2629
     2630        OUTPUT:
     2631
     2632        The product left * right (computed by Karatsuba's algorithm)
     2633
     2634        EXAMPLES::
     2635
     2636            sage: k.<t> = GF(5^3)
     2637            sage: Frob = k.frobenius_endomorphism()
     2638            sage: S.<x> = k['x',Frob]
     2639            sage: KarClass = S._karatsuba_class
     2640
     2641            sage: a = S.random_element(degree=500)
     2642            sage: b = S.random_element(degree=500)
     2643            sage: c = KarClass.mul(a,b)
     2644            sage: c == a*b
     2645            True
     2646
     2647        .. NOTE::
     2648
     2649            Behind the scene, the operator ``*`` calls actually
     2650            the function ``KarClass.mul()``.
     2651        """
     2652        cdef Py_ssize_t dx = left.degree(), dy = right.degree()
     2653        cdef list x = left.list()
     2654        cdef list y = right.list()
     2655        cdef list res
     2656        if self._cutoff < 0:
     2657            res = self.mul_step(x,y)
     2658        else:
     2659            if dx < dy:
     2660                res = self.mul_iter(y,x,1)
     2661            else:
     2662                res = self.mul_iter(x,y,0)
     2663        return self._parent(res)
     2664
     2665
     2666    def mul_matrix(self,left,right):
     2667        """
     2668        INPUT:
     2669
     2670        -  ``left`` -- an element in the skew polynomial ring
     2671           attached to the class
     2672
     2673        -  ``right`` -- an element in the skew polynomial ring
     2674           attached to the class
     2675
     2676        OUTPUT:
     2677
     2678        The product left * right (computed by Karatsuba-matrix
     2679        algorithm)
     2680
     2681        EXAMPLES::
     2682
     2683            sage: k.<t> = GF(5^3)
     2684            sage: Frob = k.frobenius_endomorphism()
     2685            sage: S.<x> = k['x',Frob]
     2686            sage: KarClass = S._karatsuba_class
     2687
     2688            sage: a = S.random_element(degree=500)
     2689            sage: b = S.random_element(degree=500)
     2690            sage: c = KarClass.mul_matrix(a,b)
     2691            sage: c == a*b
     2692            True
     2693        """
     2694        cdef Py_ssize_t dx = left.degree(), dy = right.degree()
     2695        cdef list x = left.list()
     2696        cdef list y = right.list()
     2697        cdef list res
     2698        cdef Py_ssize_t save_cutoff = self._cutoff
     2699        cdef Py_ssize_t i, j
     2700        self._cutoff = int(self._order * self._order / 2)
     2701        self._algo_matrix = 1
     2702        if self._t is None:
     2703            self._t = self._parent.base_ring().gen()
     2704        if self._T is None:
     2705            self._T = <Matrix_dense>zero_matrix(self._parent.base_ring(),self._order)
     2706            for i from 0 <= i < self._order:
     2707                map = self._parent.twist_map(-i)
     2708                for j from 0 <= j < self._order:
     2709                    self._T.set_unsafe(i,j,map(self._t**j))
     2710        if self._Tinv is None:
     2711            self._Tinv = <Matrix_dense>self._T.inverse()
     2712        if dx < dy:
     2713            res = self.mul_iter(y,x,1)
     2714        else:
     2715            res = self.mul_iter(x,y,0)
     2716        self._cutoff = save_cutoff
     2717        self._algo_matrix = 0
     2718        return self._parent(res)
     2719
     2720
     2721    def mul_list(self,x,y):
     2722        """
     2723        INPUT:
     2724
     2725        -  ``x`` -- a list of coefficients
     2726
     2727        -  ``y`` -- a list of coefficients
     2728
     2729        OUTPUT:
     2730
     2731        The list of coefficients of the product `a * b`
     2732        where `a` (resp. `b`) is the skew polynomial whose
     2733        coefficients are given by `x` (resp. `y`).
     2734
     2735        EXAMPLES::
     2736
     2737            sage: k.<t> = GF(5^3)
     2738            sage: Frob = k.frobenius_endomorphism()
     2739            sage: S.<u> = k['u',Frob]
     2740            sage: KarClass = S._karatsuba_class
     2741
     2742            sage: x = [ k.random_element() for _ in range(500) ]
     2743            sage: y = [ k.random_element() for _ in range(500) ]
     2744            sage: z = KarClass.mul_list(x,y)
     2745            sage: S(z) == S(x)*S(y)
     2746            True
     2747        """
     2748        if self._cutoff < 0:
     2749            return self.mul_step(x,y)
     2750        cdef Py_ssize_t dx = len(x)-1, dy = len(y)-1
     2751        if dx < dy:
     2752            return self.mul_iter(y,x,1)
     2753        else:
     2754            return self.mul_iter(x,y,0)
     2755
     2756
     2757    cdef list mul_step(self, list x, list y):
     2758        """
     2759        Multiplication step in Karatsuba algorithm
     2760        """
     2761        cdef Py_ssize_t dx = len(x)-1, dy = len(y) - 1
     2762        if dx < 0 or dy < 0: return [ ]
     2763        cdef Py_ssize_t i, j, start, end = dx if dx < self._order else self._order-1
     2764        cdef list twists = [ y ]
     2765        for i from 0 <= i < end:
     2766            twists.append([ self._parent.twist_map()(a) for a in twists[i] ])
     2767        cdef list res = [ ]
     2768        for j from 0 <= j <= dx+dy:
     2769            start = 0 if j <= dy else j-dy
     2770            end = j if j <= dx else dx
     2771            sum = x[start] * twists[start % self._order][j-start]
     2772            for i from start < i <= end:
     2773                sum += x[i] * twists[i % self._order][j-i]
     2774            res.append(sum)
     2775        return res
     2776
     2777
     2778    cdef list mul_step_matrix(self, list x, list y):
     2779        r"""
     2780        Multiplication step in Karatsuba-matrix algorithm. It
     2781        is based on the ring isomorphism:
     2782
     2783        .. MATH::
     2784
     2785            S / NS \simeq M_r(k^\sigma)
     2786
     2787        with the standard notations::
     2788
     2789        -  `S` is the underlying skew polynomial ring
     2790
     2791        -  `k` is the base ring of `S` (it's a finite field)
     2792
     2793        -  `\sigma` is the twisting automorphism on `k`
     2794
     2795        -  `r` is the order of `\sigma`
     2796
     2797        -  `N` is a polynomial of definition of the extension
     2798           `k / k^\sigma`
     2799
     2800        .. WARNING::
     2801
     2802            The polynomials `x` and `y` must have degree
     2803            at most `r^2/2`
     2804        """
     2805        cdef Py_ssize_t dx = len(x)-1, dy = len(y) - 1
     2806        cdef Py_ssize_t i, j
     2807        cdef Py_ssize_t r = self._order
     2808        cdef Mx, My, M
     2809        cdef list row
     2810        cdef RingElement c
     2811        k = self._parent.base_ring()
     2812        Mx = self._T * matrix(k, r, r, x + [self._zero] * (r*r-len(x)))
     2813        My = self._T * matrix(k, r, r, y + [self._zero] * (r*r-len(y)))
     2814        for i from 0 <= i < r:
     2815            frob = self._parent.twist_map(i)
     2816            row = Mx.row(i).list()
     2817            row = map(frob,row)
     2818            row = [ self._t*c for c in row[r-i:] ] + row[:r-i]
     2819            Mx.set_row(i,row)
     2820            row = My.row(i).list()
     2821            row = map(frob,row)
     2822            row = [ self._t*c for c in row[r-i:] ] + row[:r-i]
     2823            My.set_row(i,row)
     2824        M = Mx * My
     2825        for i from 0 <= i < r:
     2826            frob = self._parent.twist_map(-i)
     2827            row = M.row(i).list()
     2828            row = row[i:] + [ c/self._t for c in row[:i] ]
     2829            row = map(frob,row)
     2830            M.set_row(i,row)
     2831        M = self._Tinv * M
     2832        return M.list()[:dx+dy+1]
     2833
     2834
     2835    cdef list mul_iter(self, list x, list y, char flag): # Assume dx >= dy
     2836        """
     2837        Karatsuba's recursive iteration
     2838
     2839        .. WARNING::
     2840
     2841            This function assumes that len(x) >= len(y).
     2842        """
     2843        cdef Py_ssize_t i, j, k
     2844        cdef Py_ssize_t dx = len(x)-1, dy = len(y)-1
     2845        if self._algo_matrix:
     2846            if dx < self._cutoff:
     2847                if flag:
     2848                    return self.mul_step_matrix(y,x)
     2849                else:
     2850                    return self.mul_step_matrix(x,y)
     2851        else:
     2852            if dy < self._cutoff:
     2853                if flag:
     2854                    return self.mul_step(y,x)
     2855                else:
     2856                    return self.mul_step(x,y)
     2857        cdef Py_ssize_t dp = self._order * (1 + ceil(dx/self._order/2))
     2858        cdef list x1 = x[:dp], x2 = x[dp:]
     2859        cdef list y1 = y[:dp], y2 = y[dp:]
     2860        cdef list p1, p2, res
     2861        res = self.mul_iter(x1,y1,flag)
     2862        if dy >= dp:
     2863            res.append(self._zero)
     2864            p1 = self.mul_iter(x2,y2,flag)
     2865            res.extend(p1)
     2866            for i from 0 <= i <= dx-dp: x1[i] += x2[i]
     2867            for i from 0 <= i <= dy-dp: y1[i] += y2[i]
     2868            p2 = self.mul_iter(x1,y1,flag)
     2869            j = dx + dy - 2*dp
     2870            k = min(2*dp-2, len(res)-dp-1)
     2871            for i from k >= i > j:
     2872                res[i+dp] += p2[i] - res[i]
     2873            for i from j >= i >= 0:
     2874                res[i+dp] += p2[i] - p1[i] - res[i]
     2875        else:
     2876            if dx-dp < dy:
     2877                p2 = self.mul_iter(y1,x2,not flag)
     2878            else:
     2879                p2 = self.mul_iter(x2,y1,flag)
     2880            for i from 0 <= i < dy:
     2881                res[i+dp] += p2[i]
     2882            res.extend(p2[dy:])
     2883        return res
     2884
     2885
     2886    def div(self,left,right):
     2887        """
     2888        INPUT:
     2889
     2890        -  ``left`` -- a skew polynomial in the skew polynomial
     2891           right attached to the class
     2892
     2893        -  ``right`` -- a skew polynomial in the skew polynomial
     2894           right attached to the class
     2895
     2896        OUTPUT:
     2897
     2898        The quotient and the remainder of the right euclidien division
     2899        of left by right (computed by Karatsuba's algorithm).
     2900
     2901        .. WARNING::
     2902
     2903            This algorithm is only efficient for very very large
     2904            degrees. Do not use it!
     2905
     2906        .. TODO::
     2907
     2908            Try to understand why...
     2909
     2910        EXAMPLES::
     2911
     2912            sage: k.<t> = GF(5^3)
     2913            sage: Frob = k.frobenius_endomorphism()
     2914            sage: S.<x> = k['x',Frob]
     2915            sage: KarClass = S._karatsuba_class
     2916
     2917            sage: a = S.random_element(degree=1000)
     2918            sage: b = S.random_element(degree=500)
     2919            sage: q,r = KarClass.div(a,b)
     2920            sage: q2,r2 = a.quo_rem(b)
     2921            sage: q == q2
     2922            True
     2923            sage: r == r2
     2924            True
     2925        """
     2926        cdef list a = left.list()
     2927        cdef list b = right.list()
     2928        cdef Py_ssize_t da = len(a)-1, db = len(b)-1, i
     2929        cdef list q
     2930        self._twinv = [ ~b[db] ]
     2931        for i from 0 <= i < min(da-db,self._order-1):
     2932            self._twinv.append(self._parent.twist_map()(self._twinv[i]))
     2933        if self._cutoff < 0:
     2934            q = self.div_step(a,0,da,b,0,db)
     2935        else:
     2936            q = self.div_iter(a,0,da,b,0,db)
     2937        return self._parent(q), self._parent(a[:db])
     2938
     2939
     2940    cdef list div_step(self, list a, Py_ssize_t ia, Py_ssize_t da, list b, Py_ssize_t ib, Py_ssize_t db):
     2941        """
     2942        Division step in Karatsuba's algorithm
     2943        """
     2944        cdef Py_ssize_t i, j
     2945        cdef list q = [ ]
     2946        cdef list twb = [ b[ib:ib+db] ]
     2947        for i from 0 <= i < min(da-db,self._order-1):
     2948            twb.append([ self._parent.twist_map()(x) for x in twb[i] ])
     2949        for i from da-db >= i >= 0:
     2950            c = self._twinv[i%self._order] * a[ia+i+db]
     2951            for j from 0 <= j < db:
     2952                a[ia+i+j] -= c * twb[i%self._order][j]
     2953            q.append(c)
     2954        q.reverse()
     2955        return q
     2956
     2957
     2958    cdef list div_iter(self, list a, Py_ssize_t ia, Py_ssize_t da, list b, Py_ssize_t ib, Py_ssize_t db):
     2959        """
     2960        Karatsuba recursive iteration
     2961        """
     2962        cdef Py_ssize_t delta = da - db
     2963        if delta < self._cutoff:
     2964            return self.div_step(a,ia,da,b,ib,db)
     2965        cdef Py_ssize_t i
     2966        cdef Py_ssize_t dp = self._order * (1 + int(delta/self._order/2))
     2967        cdef list q = self.div_iter(a,ia+db,delta,b,ib+db-dp,dp), pr
     2968        if db < delta:
     2969            pr = self.mul_iter(q,b[ib:ib+db-dp],0)
     2970        else:
     2971            pr = self.mul_iter(b[ib:ib+db-dp],q,1)
     2972        for i from 0 <= i < len(pr):
     2973            a[i+ia+dp] -= pr[i]
     2974        cdef list qq = self.div_iter(a,ia,db+dp-1,b,ib,db)
     2975        qq.extend(q)
     2976        return qq
  • new file sage/rings/polynomial/skew_polynomial_ring.py

    diff --git a/sage/rings/polynomial/skew_polynomial_ring.py b/sage/rings/polynomial/skew_polynomial_ring.py
    new file mode 100644
    - +  
     1r"""
     2Skew Univariate Polynomial Rings
     3
     4Sage implements dense skew univariate polynomials over commutative rings.
     5
     6DEFINITION:
     7
     8Given a ring `R` and a ring endomorphism `\sigma` of `R`, the ring of
     9skew polynomials `R[x,\sigma]` is the usual abelian group polynomial
     10`R[x]` equipped with the modification multiplication deduced from the
     11rule `X a = \sigma(a) X`.
     12
     13.. TODO::
     14
     15    Add derivations.
     16
     17EXAMPLES::
     18
     19    sage: R.<t> = ZZ[]
     20    sage: sigma = R.hom([t+1])
     21    sage: S.<x> = SkewPolynomialRing(R,sigma); S
     22    Skew Polynomial Ring in x over Univariate Polynomial Ring in t over Integer Ring twisted by t |--> t + 1
     23
     24One can also use a shorter syntax::
     25
     26    sage: S.<x> = R['x',sigma]; S
     27    Skew Polynomial Ring in x over Univariate Polynomial Ring in t over Integer Ring twisted by t |--> t + 1
     28
     29Be careful, with the latter syntax, one cannot omit the name of the
     30variable neither in LHS nor in RHS. If we omit it in LHS, the variable
     31is not created::
     32
     33    sage: Sy = R['y',sigma]; Sy
     34    Skew Polynomial Ring in y over Univariate Polynomial Ring in t over Integer Ring twisted by t |--> t + 1
     35    sage: y.parent()
     36    Traceback (most recent call last):
     37    ...
     38    NameError: name 'y' is not defined
     39
     40If we omit it in RHS, sage tries to create a polynomial ring and fails::
     41
     42    sage: Sz.<z> = R[sigma]
     43    Traceback (most recent call last):
     44    ...
     45    ValueError: variable names must be alphanumeric, but one is 'Ring endomorphism of Univariate Polynomial Ring in t over Integer Ring
     46      Defn: t |--> t + 1' which is not.
     47
     48As for polynomials, skew polynomial rings with different variable names
     49are not equal::
     50
     51    sage: R['x',sigma] == R['y',sigma]
     52    False
     53
     54Of course, skew polynomial rings with different twist maps are not
     55equal as well
     56
     57    sage: R['x',sigma] == R['x',sigma^2]
     58    False
     59
     60Saving and loading of polynomial rings works::
     61
     62    sage: loads(dumps(R['x',sigma])) == R['x',sigma]
     63    True
     64
     65There is a coercion map from the base ring of the skew polynomial rings::
     66
     67    sage: S.has_coerce_map_from(R)
     68    True
     69    sage: x.parent()
     70    Skew Polynomial Ring in x over Univariate Polynomial Ring in t over Integer Ring twisted by t |--> t + 1
     71    sage: t.parent()
     72    Univariate Polynomial Ring in t over Integer Ring
     73    sage: y = x+t; y
     74    x + t
     75    sage: y.parent() is S
     76    True
     77
     78AUTHOR:
     79
     80- Xavier Caruso (2012-06-29)
     81"""
     82
     83#############################################################################
     84#    Copyright (C) 2012 Xavier Caruso <xavier.caruso@normalesup.org>
     85#
     86#  Distributed under the terms of the GNU General Public License (GPL)
     87#
     88#                  http://www.gnu.org/licenses/
     89#****************************************************************************
     90
     91from sage.structure.unique_representation import UniqueRepresentation
     92
     93from sage.rings.infinity import Infinity
     94
     95from sage.structure.element import Element
     96import sage.algebras.algebra
     97import sage.categories.basic as categories
     98import sage.rings.ring as ring
     99import sage.rings.ring_element as ring_element
     100import sage.rings.integral_domain as integral_domain
     101import sage.rings.principal_ideal_domain as principal_ideal_domain
     102import sage.rings.polynomial.polynomial_element_generic as polynomial_element_generic
     103import sage.rings.rational_field as rational_field
     104from sage.rings.integer_ring import is_IntegerRing, IntegerRing
     105from sage.rings.integer import Integer
     106from sage.libs.pari.all import pari_gen
     107
     108from sage.structure.parent_gens import normalize_names
     109from sage.misc.prandom import randint
     110
     111from sage.categories.morphism import Morphism
     112from sage.categories.morphism import IdentityMorphism
     113from sage.categories.homset import Hom
     114from sage.categories.map import Section
     115from sage.rings.morphism import RingHomomorphism
     116
     117from sage.matrix.matrix_space import MatrixSpace
     118
     119import sage.misc.latex as latex
     120
     121from sage.rings.polynomial.polynomial_ring import PolynomialRing_general
     122from sage.rings.polynomial.polynomial_element import Polynomial_generic_dense
     123from sage.rings.polynomial.polynomial_element import PolynomialBaseringInjection
     124import copy, re
     125
     126
     127def is_SkewPolynomialRing(S):
     128    """
     129    Return True if S is a skew polynomial ring.
     130
     131    EXAMPLES::
     132
     133        sage: from sage.rings.polynomial.skew_polynomial_ring import is_SkewPolynomialRing
     134        sage: is_SkewPolynomialRing(QQ['x'])
     135        False
     136
     137        sage: k.<t> = GF(5^3)
     138        sage: Frob = k.frobenius_endomorphism()
     139        sage: is_SkewPolynomialRing(k['x',Frob])
     140        True
     141    """
     142    return isinstance(S, SkewPolynomialRing_general)
     143
     144
     145#########################################################################################
     146
     147class SectionSkewPolynomialCenterInjection(Section):
     148    def _call_ (self, x):
     149        order = self.inverse()._order
     150        section = self.inverse()._embed.section()
     151        lx = x.list()
     152        l = [ ]
     153        mod = 0
     154        for c in lx:
     155            if mod == 0:
     156                l.append(section(c))
     157            else:
     158                if not c.is_zero():
     159                    raise ValueError("%s is not in the center" % x)
     160            mod += 1
     161            if mod == order:
     162                mod = 0
     163        return self.codomain()(l)
     164
     165
     166class SkewPolynomialCenterInjection(RingHomomorphism):
     167    def __init__(self,domain,codomain,embed,order):
     168        RingHomomorphism.__init__(self,Hom(domain,codomain))
     169        self._embed = embed
     170        self._order = order
     171        self._codomain = codomain
     172        self._section = SectionSkewPolynomialCenterInjection(self)
     173
     174    def _repr_(self):
     175        return "Embedding of the center of %s into this ring" % self._codomain
     176
     177    def _call_(self,x):
     178        k = self._codomain.base_ring ()
     179        l = [ ]
     180        lz = [ k(0) ] * (self._order-1)
     181        for c in x.list():
     182            l += [ self._embed(c) ] + lz
     183        return self._codomain (l)
     184
     185    def section(self):
     186        return self._section
     187
     188
     189class CenterSkewPolynomialRing(PolynomialRing_general):
     190    """
     191    A specific class for the center of a skew polynomial ring.
     192    """
     193
     194    def __init__ (self, skew_ring, names=None, sparse=False, element_class=None):
     195        if not isinstance (skew_ring, SkewPolynomialRing_general):
     196            raise TypeError("%s is not a Skew Polynomial Ring" % skew_ring)
     197        self._skew_ring = skew_ring
     198        base_ring = skew_ring.base_ring()
     199        kfixed, embed = skew_ring._map.fixed_points()
     200        self._embed_basering = embed
     201        order = skew_ring._map.order()
     202        if order == Infinity:
     203            raise NotImplementedError
     204        self.__is_sparse = sparse
     205        self._PolynomialRing_general__is_sparse = sparse
     206        if element_class:
     207            self._polynomial_class = element_class
     208        else:
     209            if sparse:
     210                raise NotImplementedError("sparse skew polynomials are not implemented")
     211            else:
     212                self._polynomial_class = sage.rings.polynomial.skew_polynomial_element.CenterSkewPolynomial_generic_dense
     213
     214        # Algebra.__init__ also calls __init_extra__ of Algebras(...).parent_class, which
     215        # tries to provide a conversion from the base ring, if it does not exist.
     216        # This is for algebras that only do the generic stuff in their initialisation.
     217        # But here, we want to use PolynomialBaseringInjection. Hence, we need to
     218        # wipe the memory and construct the conversion from scratch.
     219        sage.algebras.algebra.Algebra.__init__(self, kfixed, names=names, normalize=True, category=None)
     220
     221        if names is None:
     222            if order == 1:
     223                self._variable_name = skew_ring.variable_name()
     224                self._latex_variable_name = skew_ring.latex_variable_names()[0]
     225                self._parenthesis = False
     226            else:
     227                self._variable_name = skew_ring.variable_name () + "^" + str(order)
     228                self._latex_variable_name = skew_ring.latex_variable_names()[0] + "^{" + str (order) + "}"
     229                self._parenthesis = True
     230        else:
     231            self._variable_name = sage.algebras.algebra.Algebra.variable_name(self)
     232            self._latex_variable_name = sage.algebras.algebra.Algebra.latex_variable_names(self)[0]
     233            self._parenthesis = False
     234        self._names = [ self._variable_name ]
     235        self.__generator = self._polynomial_class (self, [0,1], is_gen=True)
     236        base_inject = PolynomialBaseringInjection(kfixed,self)
     237        center_inject = SkewPolynomialCenterInjection (self, skew_ring, embed, order)
     238        self._unset_coercions_used()
     239        self._populate_coercion_lists_(
     240            coerce_list = [base_inject],
     241            convert_list = [list, base_inject],
     242            embedding = center_inject)
     243
     244    def _repr_ (self):
     245        """
     246        Return a string representation of this ring.
     247        """
     248        s = "Center of %s:\n" % self._skew_ring
     249        s += PolynomialRing_general._repr_(self)
     250        return s
     251
     252    def _latex_ (self):
     253        """
     254        Return a latex representation of this ring.
     255        """
     256        return "%s[%s]"%(latex.latex(self.base_ring()), self._latex_variable_name)
     257
     258    def gen (self,n=0):
     259        """
     260        Return the generator of this ring.
     261
     262        EXAMPLES::
     263
     264            sage: k.<t> = GF(2^10)
     265      &n