Ticket #13215: trac_13642_modular_exp_polynomial.patch

File trac_13642_modular_exp_polynomial.patch, 3.7 KB (added by caruso, 8 years ago)
  • sage/rings/polynomial/polynomial_element.pxd

    # HG changeset patch
    # User Xavier Caruso <xavier.caruso@normalesup.org>
    # Date 1350912681 -7200
    # Node ID 67381b4251f22964ab68d7cc6e3cc69708f89b06
    # Parent  8ca0c60550327733689a85e812c6321ca13773f1
    trac #13642: modular exponentiation of polynomials
    
    diff --git a/sage/rings/polynomial/polynomial_element.pxd b/sage/rings/polynomial/polynomial_element.pxd
    a b  
    33include '../../ext/stdsage.pxi'
    44
    55
     6from sage.rings.integer cimport Integer
    67from sage.structure.element import Element, CommutativeAlgebraElement
    78from sage.structure.element cimport Element, CommutativeAlgebraElement, ModuleElement
    89from sage.structure.parent cimport Parent
     
    1718    cdef long _hash_c(self)
    1819    cpdef constant_coefficient(self)
    1920    cpdef Polynomial _new_constant_poly(self, a, Parent P)
     21    cdef Polynomial _pow_modulus(Polynomial self, Integer right, modulus)
    2022
    2123    # UNSAFE, only call from an inplace operator
    2224    # may return a new element if not possible to modify inplace
  • sage/rings/polynomial/polynomial_element.pyx

    diff --git a/sage/rings/polynomial/polynomial_element.pyx b/sage/rings/polynomial/polynomial_element.pyx
    a b  
    14071407            sage: f^3
    14081408            x^3 - 3*x^2 + 3*x - 1           
    14091409       
    1410             sage: R = PolynomialRing(GF(2), x)
    1411             sage: f = R(x^9 + x^7 + x^6 + x^5 + x^4 + x^2 + x)
    1412             sage: h = R(x^10 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + 1)
     1410            sage: R.<x> = GF(2)[]
     1411            sage: f = x^9 + x^7 + x^6 + x^5 + x^4 + x^2 + x
     1412            sage: h = x^10 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + 1
    14131413            sage: pow(f, 2, h)
    14141414            x^9 + x^8 + x^7 + x^5 + x^3
    14151415
     
    14351435                        raise TypeError("non-integral exponents not supported")
    14361436           
    14371437        if self.degree() <= 0:
    1438             r = self.parent()(self[0]**right)
     1438            return self.parent()(self[0]**right)
    14391439        elif right < 0:
    1440             r = (~self)**(-right)
     1440            return (~self)**(-right)
     1441        elif modulus:
     1442            return (<Polynomial>self)._pow_modulus(<Integer>right, modulus)
    14411443        elif (<Polynomial>self) == self.parent().gen():   # special case x**n should be faster!
    14421444            P = self.parent()
    14431445            R = P.base_ring()
     
    14451447                v = {right:R.one_element()}
    14461448            else:
    14471449                v = [R.zero_element()]*right + [R.one_element()]
    1448             r = self.parent()(v, check=False)
     1450            return self.parent()(v, check=False)
    14491451        else:
    1450             r = generic_power(self, right)
    1451         if modulus:
    1452             return r % modulus
    1453         else:
    1454             return r
     1452            return generic_power(self,right)
     1453
     1454    cdef Polynomial _pow_modulus(Polynomial self, Integer right, modulus):
     1455        """
     1456        TESTS::
     1457
     1458            sage: R = PolynomialRing(GF(2), x)
     1459            sage: f = R(x^9 + x^7 + x^6 + x^5 + x^4 + x^2 + x)
     1460            sage: h = R(x^10 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + 1)
     1461            sage: pow(f, 2, h)
     1462            x^9 + x^8 + x^7 + x^5 + x^3
     1463        """
     1464        cdef Polynomial pow = self
     1465        cdef Polynomial r
     1466        cdef Integer n = right
     1467        while n & 1 == 0:
     1468            pow = (pow*pow) % modulus
     1469            n = n >> 1
     1470        r = pow
     1471        n = n >> 1
     1472        while n != 0:
     1473            pow = (pow*pow) % modulus
     1474            if n & 1 != 0:
     1475                r = (r*pow) % modulus
     1476            n = n >> 1
     1477        return r
    14551478       
    14561479    def _pow(self, right):
    14571480        # TODO: fit __pow__ into the arithmetic structure