Ticket #4965: zmod_poly.patch

File zmod_poly.patch, 100.3 KB (added by Martin Albrecht, 14 years ago)

going down to 230

  • module_list.py

    # HG changeset patch
    # User Martin Albrecht <malb@informatik.uni-bremen.de>
    # Date 1232431397 28800
    # Node ID c253f3eb84d063f2ec2bd86bffc921f7d08ae2af
    # Parent  4ddd93e16f5e39475c5d792857cbeb49c51b07f8
    univariate polynomials over Z/nZ where n <= sys.maxint via FLINT's zmod_poly
    
    This patch also improves polynomial_template.pxi (which was used to provide the functionality
    of this patch) and provides a simple technology-preview level interface to David Harvey's zn_poly
    library. Finally, some pxi files became pxd files in the process.
    
    This patch contains/folds in Robert Bradshaw's patches:
     * 4965-shift-fix.patch
     * 4965-pxi-to-pxd.patch
    
    diff -r 4ddd93e16f5e -r c253f3eb84d0 module_list.py
    a b ext_modules = [ 
    10491049              language = 'c++',
    10501050              include_dirs = ['sage/libs/ntl/']),
    10511051
     1052    Extension('sage.rings.polynomial.polynomial_zmod_flint',
     1053              sources = ['sage/rings/polynomial/polynomial_zmod_flint.pyx'],
     1054              libraries = ["csage", "flint", "gmp", "gmpxx", "ntl", "zn_poly"],
     1055              extra_compile_args=["-std=c99", "-D_XPG6"],
     1056              include_dirs = [SAGE_ROOT+'/local/include/FLINT/'],
     1057              depends = [SAGE_ROOT + "/local/include/FLINT/flint.h"]),
     1058
    10521059    Extension('sage.rings.polynomial.polynomial_integer_dense_flint',
    10531060              sources = ['sage/rings/polynomial/polynomial_integer_dense_flint.pyx'],
    10541061              language = 'c++',
  • new file sage/libs/flint/flint.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/flint.pxd
    - +  
     1cdef extern from "flint.h":
     2
     3    cdef long FLINT_BITS
     4    cdef long FLINT_D_BITS
     5
     6    cdef unsigned long FLINT_BIT_COUNT(unsigned long)
  • deleted file sage/libs/flint/flint.pxi

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/flint.pxi
    + -  
    1 cdef extern from "flint.h":
    2     pass
  • sage/libs/flint/fmpz_poly.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/fmpz_poly.pxd
    a b include "../../ext/stdsage.pxi" 
    11include "../../ext/stdsage.pxi"
    22include "../../ext/cdefs.pxi"
    33
    4 include "flint.pxi"
     4from flint cimport *
    55include "fmpz_poly.pxi"
    66
    77from sage.structure.sage_object cimport SageObject
  • new file sage/libs/flint/long_extras.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/long_extras.pxd
    - +  
     1include "../../ext/stdsage.pxi"
     2include "../../ext/cdefs.pxi"
     3
     4from sage.libs.flint.flint cimport *
     5
     6
     7cdef extern from "FLINT/long_extras.h":
     8
     9    cdef unsigned long z_randint(unsigned long limit)
     10
     11    cdef unsigned long z_randbits(unsigned long bits)
     12
     13    cdef unsigned long z_randprime(unsigned long bits)
     14
     15    cdef double z_precompute_inverse(unsigned long n)
     16
     17    cdef double z_precompute_inverse2(unsigned long n)
     18
     19    cdef double z_ll_precompute_inverse(unsigned long n)
     20
     21    cdef unsigned long z_addmod(unsigned long a, unsigned long b, unsigned long p)
     22
     23    cdef unsigned long z_submod(unsigned long a, unsigned long b, unsigned long p)
     24
     25    cdef unsigned long z_negmod(unsigned long a, unsigned long p)
     26
     27    cdef unsigned long z_mod_precomp(unsigned long a, unsigned long n, double ninv)
     28
     29    cdef unsigned long z_div2_precomp(unsigned long a, unsigned long n, double ninv)
     30
     31    cdef unsigned long z_mod2_precomp(unsigned long a, unsigned long n, double ninv)
     32
     33    cdef unsigned long z_ll_mod_precomp(unsigned long a_hi, unsigned long a_lo, unsigned long n, double ninv)
     34
     35    cdef unsigned long z_mulmod_precomp(unsigned long a, unsigned long b, unsigned long n, double ninv)
     36
     37    cdef unsigned long z_mulmod2_precomp(unsigned long a, unsigned long b, unsigned long n, double ninv)
     38
     39    cdef unsigned long z_powmod(unsigned long a, long exp, unsigned long n)
     40
     41    cdef unsigned long z_powmod2(unsigned long a, long exp, unsigned long n)
     42
     43    cdef unsigned long z_powmod_precomp(unsigned long a, long exp, unsigned long n, double ninv)
     44
     45    cdef unsigned long z_powmod2_precomp(unsigned long a, long exp, unsigned long n, double ninv)
     46
     47    cdef int z_legendre_precomp(unsigned long a, unsigned long p, double pinv)
     48
     49    cdef int z_jacobi(long x, unsigned long y)
     50
     51    cdef int z_ispseudoprime_fermat(unsigned long n, unsigned long i)
     52
     53    cdef int z_isprime(unsigned long n)
     54
     55    cdef int z_isprime_precomp(unsigned long n, double ninv)
     56
     57    cdef unsigned long z_nextprime(unsigned long n)
     58
     59    cdef int z_isprime_pocklington(unsigned long n, unsigned long iterations)
     60
     61    cdef int z_ispseudoprime_lucas_ab(unsigned long n, int a, int b)
     62
     63    cdef int z_ispseudoprime_lucas(unsigned long n)
     64
     65    cdef unsigned long z_pow(unsigned long a, unsigned long exp)
     66
     67    cdef unsigned long z_sqrtmod(unsigned long a, unsigned long p)
     68
     69    cdef unsigned long z_cuberootmod(unsigned long * cuberoot1, unsigned long a, unsigned long p)
     70
     71    cdef unsigned long z_invert(unsigned long a, unsigned long p)
     72   
     73    cdef long z_gcd_invert(long* a, long x, long y)
     74
     75    cdef long z_extgcd(long* a, long* b, long x, long y)
     76
     77    cdef unsigned long z_gcd(long x, long y)
     78
     79    cdef unsigned long z_intsqrt(unsigned long r)
     80
     81    cdef int z_issquare(long x)
     82
     83    cdef unsigned long z_CRT(unsigned long x1, unsigned long n1, unsigned long x2, unsigned long n2)
     84
     85    cdef int z_issquarefree(unsigned long n)
     86
     87    cdef int z_remove_precomp(unsigned long * n, unsigned long p, double pinv)
     88
     89    cdef int z_remove(unsigned long * n, unsigned long p)
     90
     91    cdef unsigned long z_factor_SQUFOF(unsigned long n)
     92
     93    cdef unsigned long z_primitive_root(unsigned long p)
     94
     95    cdef unsigned long z_primitive_root_precomp(unsigned long p, double p_inv)
  • new file sage/libs/flint/ntl_interface.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/ntl_interface.pxd
    - +  
     1include "../../ext/cdefs.pxi"
     2include "fmpz.pxi"
     3
     4from sage.libs.flint.fmpz_poly cimport fmpz_poly_t
     5from sage.libs.ntl.ntl_ZZ_decl cimport ZZ_c
     6from sage.libs.ntl.ntl_ZZX_decl cimport ZZX_c
     7
     8cdef extern from "FLINT/NTL-interface.h":
     9    unsigned long ZZ_limbs(ZZ_c z)
     10
     11    void fmpz_poly_to_ZZX(ZZX_c output, fmpz_poly_t poly)
     12    void ZZX_to_fmpz_poly(fmpz_poly_t output, ZZX_c poly)
     13
     14    void fmpz_to_mpz(mpz_t res, fmpz_t f)
     15    void ZZ_to_fmpz(fmpz_t output, ZZ_c z)
     16   
  • deleted file sage/libs/flint/ntl_interface.pxi

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/ntl_interface.pxi
    + -  
    1 cdef extern from "FLINT/NTL-interface.h":
    2     unsigned long ZZ_limbs(ZZ_c z)
    3 
    4     void fmpz_poly_to_ZZX(ZZX_c output, fmpz_poly_t poly)
    5     void ZZX_to_fmpz_poly(fmpz_poly_t output, ZZX_c poly)
    6 
    7     void fmpz_to_mpz(mpz_t res, fmpz_t f)
    8     void ZZ_to_fmpz(fmpz_t output, ZZ_c z)
    9    
  • new file sage/libs/flint/zmod_poly.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/zmod_poly.pxd
    - +  
     1include "../../ext/stdsage.pxi"
     2include "../../ext/cdefs.pxi"
     3
     4from sage.libs.flint.flint cimport *
     5
     6from flint import *
     7
     8cdef extern from "FLINT/zmod_poly.h":
     9    ctypedef struct zmod_poly_struct:
     10        unsigned long *coeffs
     11        unsigned long alloc
     12        unsigned long length
     13        unsigned long p
     14        double p_inv
     15
     16    ctypedef zmod_poly_struct* zmod_poly_t
     17
     18    cdef void zmod_poly_init(zmod_poly_t poly, unsigned long p)
     19    cdef void zmod_poly_init_precomp(zmod_poly_t poly, unsigned long p, double p_inv)
     20    cdef void zmod_poly_init2(zmod_poly_t poly, unsigned long p, unsigned long alloc)
     21    cdef void zmod_poly_init2_precomp(zmod_poly_t poly, unsigned long p, double p_inv, unsigned long alloc)
     22    cdef void zmod_poly_clear(zmod_poly_t poly)
     23
     24    cdef void zmod_poly_realloc(zmod_poly_t poly, unsigned long alloc)
     25    # _bits_ only applies to newly allocated coefficients, not existing ones...
     26
     27    # this non-inlined version REQUIRES that alloc > poly->alloc
     28    void __zmod_poly_fit_length(zmod_poly_t poly, unsigned long alloc)
     29
     30    # this is arranged so that the initial comparison (very frequent) is inlined,
     31    # but the actual allocation (infrequent) is not
     32    cdef void zmod_poly_fit_length(zmod_poly_t poly, unsigned long alloc)
     33
     34    # ------------------------------------------------------
     35    # Setting/retrieving coefficients
     36
     37    cdef unsigned long zmod_poly_get_coeff_ui(zmod_poly_t poly, unsigned long n)
     38
     39    cdef unsigned long _zmod_poly_get_coeff_ui(zmod_poly_t poly, unsigned long n)
     40
     41    cdef void zmod_poly_set_coeff_ui(zmod_poly_t poly, unsigned long n, unsigned long c)
     42
     43    cdef void _zmod_poly_set_coeff_ui(zmod_poly_t poly, unsigned long n, unsigned long c)
     44
     45    # ------------------------------------------------------
     46    # String conversions and I/O
     47
     48    cdef int zmod_poly_from_string(zmod_poly_t poly, char* s)
     49    cdef char* zmod_poly_to_string(zmod_poly_t poly)
     50    cdef void zmod_poly_print(zmod_poly_t poly)
     51    cdef void zmod_poly_fprint(zmod_poly_t poly, FILE* f)
     52    cdef int zmod_poly_read(zmod_poly_t poly)
     53    cdef int zmod_poly_fread(zmod_poly_t poly, FILE* f)
     54
     55    # ------------------------------------------------------
     56    # Length and degree
     57
     58    cdef void __zmod_poly_normalise(zmod_poly_t poly)
     59    cdef int __zmod_poly_normalised(zmod_poly_t poly)
     60    cdef void zmod_poly_truncate(zmod_poly_t poly, unsigned long length)
     61
     62    cdef unsigned long zmod_poly_length(zmod_poly_t poly)
     63
     64    cdef long zmod_poly_degree(zmod_poly_t poly)
     65
     66    cdef unsigned long zmod_poly_modulus(zmod_poly_t poly)
     67
     68    cdef double zmod_poly_precomputed_inverse(zmod_poly_t poly)
     69
     70    # ------------------------------------------------------
     71    # Assignment
     72
     73    cdef void _zmod_poly_set(zmod_poly_t res, zmod_poly_t poly)
     74    cdef void zmod_poly_set(zmod_poly_t res, zmod_poly_t poly)
     75
     76    cdef void zmod_poly_zero(zmod_poly_t poly)
     77
     78    cdef void zmod_poly_swap(zmod_poly_t poly1, zmod_poly_t poly2)
     79
     80    #
     81    # Subpolynomials
     82    #
     83
     84    cdef void _zmod_poly_attach(zmod_poly_t output, zmod_poly_t input)
     85
     86    cdef void zmod_poly_attach(zmod_poly_t output, zmod_poly_t input)
     87
     88    #
     89    # Attach input shifted right by n to output
     90    #
     91
     92    cdef void _zmod_poly_attach_shift(zmod_poly_t output, zmod_poly_t input, unsigned long n)
     93
     94    cdef void zmod_poly_attach_shift(zmod_poly_t output, zmod_poly_t input, unsigned long n)
     95
     96    #
     97    # Attach input to first n coefficients of input
     98    #
     99
     100    cdef void _zmod_poly_attach_truncate(zmod_poly_t output,  zmod_poly_t input, unsigned long n)
     101
     102    cdef void zmod_poly_attach_truncate(zmod_poly_t output, zmod_poly_t input, unsigned long n)
     103
     104    #
     105    # Comparison functions
     106    #
     107
     108    cdef int zmod_poly_equal(zmod_poly_t poly1, zmod_poly_t poly2)
     109
     110    cdef int zmod_poly_is_one(zmod_poly_t poly1)
     111
     112    #
     113    # Reversal
     114    #
     115
     116    cdef void _zmod_poly_reverse(zmod_poly_t output, zmod_poly_t input, unsigned long length)
     117    cdef void zmod_poly_reverse(zmod_poly_t output, zmod_poly_t input, unsigned long length)
     118
     119    #
     120    # Monic polys
     121    #
     122
     123    cdef void zmod_poly_make_monic(zmod_poly_t output, zmod_poly_t pol)
     124
     125    #
     126    # Addition and subtraction
     127    #
     128
     129    cdef void zmod_poly_add(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     130    cdef void zmod_poly_add_without_mod(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     131    cdef void zmod_poly_sub(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     132    cdef void _zmod_poly_sub(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     133    cdef void zmod_poly_neg(zmod_poly_t res, zmod_poly_t poly)
     134
     135    #
     136    # Shifting functions
     137    #
     138
     139    cdef void zmod_poly_left_shift(zmod_poly_t res, zmod_poly_t poly, unsigned long k)
     140    cdef void zmod_poly_right_shift(zmod_poly_t res, zmod_poly_t poly, unsigned long k)
     141
     142    #
     143    # Polynomial multiplication
     144    #
     145    # All multiplication functions require that the modulus be no more than FLINT_BITS-1 bits
     146    #
     147
     148    cdef void zmod_poly_mul(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     149    cdef void zmod_poly_sqr(zmod_poly_t res, zmod_poly_t poly)
     150
     151    # Requires that poly1 bits + poly2 bits + log_length is not greater than 2*FLINT_BITS
     152
     153    cdef void zmod_poly_mul_KS(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits_input)
     154    cdef void _zmod_poly_mul_KS(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits_input)
     155
     156    cdef void zmod_poly_mul_KS_trunc(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits_input, unsigned long trunc)
     157    cdef void _zmod_poly_mul_KS_trunc(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits_input, unsigned long trunc)
     158
     159    cdef void _zmod_poly_mul_classical(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     160    cdef void __zmod_poly_mul_classical_mod_last(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits)
     161    cdef void __zmod_poly_mul_classical_mod_throughout(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits)
     162    cdef void zmod_poly_mul_classical(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     163    cdef void _zmod_poly_sqr_classical(zmod_poly_t res, zmod_poly_t poly)
     164    cdef void zmod_poly_sqr_classical(zmod_poly_t res, zmod_poly_t poly)
     165
     166    cdef void _zmod_poly_mul_classical_trunc(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long trunc)
     167    cdef void __zmod_poly_mul_classical_trunc_mod_last(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits, unsigned long trunc)
     168    cdef void __zmod_poly_mul_classical_trunc_mod_throughout(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits, unsigned long trunc)
     169    cdef void zmod_poly_mul_classical_trunc(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long trunc)
     170
     171    cdef void _zmod_poly_mul_classical_trunc_left(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long trunc)
     172    cdef void __zmod_poly_mul_classical_trunc_left_mod_last(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits, unsigned long trunc)
     173    cdef void __zmod_poly_mul_classical_trunc_left_mod_throughout(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long bits, unsigned long trunc)
     174    cdef void zmod_poly_mul_classical_trunc_left(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long trunc)
     175
     176    cdef void zmod_poly_mul_trunc_n(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long trunc)
     177    cdef void zmod_poly_mul_trunc_left_n(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2, unsigned long trunc)
     178
     179    #
     180    # Bit packing functions
     181    #
     182
     183    cdef unsigned long zmod_poly_bits(zmod_poly_t poly)
     184    cdef void _zmod_poly_bit_pack_mpn(mp_limb_t * res, zmod_poly_t poly, unsigned long bits, unsigned long length)
     185    cdef void _zmod_poly_bit_unpack_mpn(zmod_poly_t poly, mp_limb_t *mpn, unsigned long length, unsigned long bits)
     186
     187    cdef void print_binary(unsigned long n, unsigned long len)
     188    cdef void print_binary2(unsigned long n, unsigned long len, unsigned long space_bit)
     189
     190    #
     191    # Scalar multiplication
     192    #
     193
     194    cdef void _zmod_poly_scalar_mul(zmod_poly_t res, zmod_poly_t poly, unsigned long scalar)
     195    cdef void zmod_poly_scalar_mul(zmod_poly_t res, zmod_poly_t poly, unsigned long scalar)
     196    cdef void __zmod_poly_scalar_mul_without_mod(zmod_poly_t res, zmod_poly_t poly, unsigned long scalar)
     197
     198    #
     199    # Division
     200    #
     201
     202    cdef void zmod_poly_divrem_classical(zmod_poly_t Q, zmod_poly_t R, zmod_poly_t A, zmod_poly_t B)
     203    cdef void __zmod_poly_divrem_classical_mod_last(zmod_poly_t Q, zmod_poly_t R, zmod_poly_t A, zmod_poly_t B)
     204    cdef void zmod_poly_div_classical(zmod_poly_t Q, zmod_poly_t A, zmod_poly_t B)
     205    cdef void __zmod_poly_div_classical_mod_last(zmod_poly_t Q, zmod_poly_t A, zmod_poly_t B)
     206    cdef void zmod_poly_div_divconquer_recursive(zmod_poly_t Q, zmod_poly_t BQ, zmod_poly_t A, zmod_poly_t B)
     207    cdef void zmod_poly_divrem_divconquer(zmod_poly_t Q, zmod_poly_t R, zmod_poly_t A, zmod_poly_t B)
     208    cdef void zmod_poly_div_divconquer(zmod_poly_t Q, zmod_poly_t A, zmod_poly_t B)
     209
     210    #
     211    # Newton Inversion
     212    #
     213
     214    cdef void zmod_poly_newton_invert_basecase(zmod_poly_t Q_inv, zmod_poly_t Q, unsigned long n)
     215    cdef void zmod_poly_newton_invert(zmod_poly_t Q_inv, zmod_poly_t Q, unsigned long n)
     216
     217    #
     218    # Newton Division
     219    #
     220
     221    cdef void zmod_poly_div_series(zmod_poly_t Q, zmod_poly_t A, zmod_poly_t B, unsigned long n)
     222    cdef void zmod_poly_div_newton(zmod_poly_t Q, zmod_poly_t A, zmod_poly_t B)
     223    cdef void zmod_poly_divrem_newton(zmod_poly_t Q, zmod_poly_t R, zmod_poly_t A, zmod_poly_t B)
     224
     225    cdef void zmod_poly_divrem(zmod_poly_t Q, zmod_poly_t R, zmod_poly_t A, zmod_poly_t B)
     226
     227    cdef void zmod_poly_div(zmod_poly_t Q, zmod_poly_t A, zmod_poly_t B)
     228
     229    #
     230    # Resultant
     231    #
     232
     233    cdef unsigned long zmod_poly_resultant_euclidean(zmod_poly_t a, zmod_poly_t b)
     234
     235    cdef unsigned long zmod_poly_resultant(zmod_poly_t a, zmod_poly_t b)
     236
     237    #
     238    # GCD
     239    #
     240
     241    cdef void zmod_poly_gcd(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     242    cdef int zmod_poly_gcd_invert(zmod_poly_t res, zmod_poly_t poly1, zmod_poly_t poly2)
     243    cdef void zmod_poly_xgcd(zmod_poly_t res, zmod_poly_t s, zmod_poly_t t, zmod_poly_t poly1, zmod_poly_t poly2)
     244
     245
     246    # FLINT 1.1 will have:
     247    # cdef void zmod_poly_powmod(zmod_poly_t res, zmod_poly_t pol, long exp, zmod_poly_t f)
  • new file sage/libs/flint/zmod_poly_linkage.pxi

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/flint/zmod_poly_linkage.pxi
    - +  
     1r"""
     2Linkage for arithmetic with FLINT's zmod_poly_t elements.
     3
     4This file provides the backend for \class{Polynomial_zmod_flint} via
     5templating.
     6
     7AUTHOR:
     8    -- Martin Albrecht (2009-01) another initial implementation
     9    -- Burcin Erocal (2008-11) initial implementation
     10"""
     11#*****************************************************************************
     12#       Copyright (C) 2008 Burcin Erocal <burcin@erocal.org>
     13#       Copyright (C) 2009 Martin Albrecht <M.R.Albrecht@rhul.ac.uk>
     14#
     15#  Distributed under the terms of the GNU General Public License (GPL)
     16#                  http://www.gnu.org/licenses/
     17#*****************************************************************************
     18
     19from sage.libs.flint.zmod_poly cimport *, zmod_poly_t
     20from sage.libs.flint.long_extras cimport *
     21
     22include "../../ext/stdsage.pxi"
     23
     24cdef inline celement *celement_new(unsigned long n):
     25    cdef celement *g = <celement *>sage_malloc(sizeof(zmod_poly_t))
     26    zmod_poly_init(g, n)
     27    return g
     28
     29cdef inline int celement_delete(zmod_poly_t e, unsigned long n):
     30    zmod_poly_clear(e)
     31    sage_free(e)
     32
     33cdef inline int celement_construct(zmod_poly_t e, unsigned long n):
     34    """
     35    EXAMPLE:
     36        sage: P.<x> = GF(32003)[]
     37
     38        sage: Q.<x> = GF(7)[]
     39    """
     40    zmod_poly_init(e, n)
     41
     42cdef inline int celement_destruct(zmod_poly_t e, unsigned long n):
     43    """
     44    EXAMPLE:
     45        sage: P.<x> = GF(32003)[]
     46        sage: del x
     47
     48        sage: Q.<x> = GF(7)[]
     49        sage: del x
     50    """
     51    zmod_poly_clear(e)
     52
     53cdef inline int celement_gen(zmod_poly_t e, long i, unsigned long n) except -2:
     54    """
     55    EXAMPLE:
     56        sage: P.<x> = GF(32003)[]
     57
     58        sage: Q.<x> = GF(7)[]
     59    """
     60    zmod_poly_zero(e)
     61    zmod_poly_set_coeff_ui(e, 1, 1)
     62
     63cdef object celement_repr(zmod_poly_t e, unsigned long n):
     64    raise NotImplementedError
     65
     66cdef inline int celement_set(zmod_poly_t res, zmod_poly_t a, unsigned long n) except -2:
     67    """
     68    EXAMPLE:
     69        sage: P.<x> = GF(32003)[]
     70        sage: y = copy(x)
     71        sage: y is x
     72        False
     73        sage: y == x
     74        True
     75
     76        sage: Q.<x> = GF(7)[]
     77        sage: y = copy(x)
     78        sage: y is x
     79        False
     80        sage: y == x
     81        True
     82    """
     83    zmod_poly_set(res, a)
     84
     85cdef inline int celement_set_si(zmod_poly_t res, long i, unsigned long n) except -2:
     86    """
     87    EXAMPLE:
     88        sage: P.<x> = GF(32003)[]
     89        sage: P(32003)
     90        0
     91        sage: P(1)
     92        1
     93        sage: P(32004)
     94        1
     95
     96        sage: Q.<x> = GF(7)[]
     97        sage: Q(7)
     98        0
     99        sage: Q(1)
     100        1
     101        sage: Q(8)
     102        1
     103    """
     104    while i < 0:
     105        i += n
     106    zmod_poly_zero(res)
     107    if i:
     108        zmod_poly_set_coeff_ui(res, 0, <unsigned long>i)
     109
     110cdef inline long celement_get_si(zmod_poly_t res, unsigned long n) except -2:
     111    raise NotImplementedError
     112   
     113cdef inline bint celement_is_zero(zmod_poly_t a, unsigned long n) except -2:
     114    """
     115    EXAMPLE:
     116        sage: P.<x> = GF(32003)[]
     117        sage: P(1).is_zero()
     118        False
     119        sage: P(0).is_zero()
     120        True
     121
     122        sage: Q.<x> = GF(7)[]
     123        sage: Q(1).is_zero()
     124        False
     125        sage: Q(0).is_zero()
     126        True
     127    """
     128    # is_zero doesn't exist
     129    return zmod_poly_degree(a) == -1
     130
     131cdef inline bint celement_is_one(zmod_poly_t a, unsigned long n) except -2:
     132    """
     133    EXAMPLE:
     134        sage: P.<x> = GF(32003)[]
     135        sage: P(1).is_one()
     136        True
     137        sage: P(0).is_one()
     138        False
     139
     140        sage: Q.<x> = GF(7)[]
     141        sage: Q(1).is_one()
     142        True
     143        sage: Q(0).is_one()
     144        False
     145    """
     146
     147    return zmod_poly_is_one(a)
     148
     149cdef inline bint celement_equal(zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     150    """
     151    EXAMPLE:
     152        sage: P.<x> = GF(32003)[]
     153        sage: (3*2)*x == 3*(2*x)
     154        True
     155        sage: (3*2)*x + 2 == 3*(2*x) + 1 + 1
     156        True
     157        sage: (3*2)*x + 7 == 3*(2*x) + 1 + 1
     158        False
     159
     160        sage: Q.<x> = GF(7)[]
     161        sage: (3*2)*x == 3*(2*x)
     162        True
     163        sage: (3*2)*x + 2 == 3*(2*x) + 1 + 1
     164        True
     165        sage: (3*2)*x + 7 == 3*(2*x) + 1 + 1
     166        False
     167    """
     168    return zmod_poly_equal(a, b)
     169
     170cdef inline int celement_cmp(zmod_poly_t l, zmod_poly_t r, unsigned long n) except -2:
     171    """
     172    EXAMPLE:
     173        sage: P.<x> = GF(32003)[]
     174        sage: x > x
     175        False
     176        sage: x^2 > x
     177        True
     178        sage: 3*x > x
     179        True
     180
     181        sage: Q.<x> = GF(7)[]
     182        sage: x > x
     183        False
     184        sage: x^2 > x
     185        True
     186        sage: 3*x > x
     187        True
     188
     189        sage: f = x^64 + x^20 + 1
     190        sage: g = x^63 + x^20 + 1
     191        sage: f > g
     192        True
     193
     194        sage: f = x^64 + x^10 + 1
     195        sage: g = x^64 + x^20 + 1
     196        sage: f < g
     197        True
     198
     199        sage: f = x^64 + x^20
     200        sage: g = x^64 + x^20 + 1
     201        sage: f < g
     202        True
     203    """
     204    cdef int deg_right = zmod_poly_degree(r)
     205    cdef int degdiff = deg_right - zmod_poly_degree(l)
     206    cdef int i
     207    cdef unsigned long rcoeff, lcoeff
     208    if degdiff > 0:
     209        return -1
     210    elif degdiff < 0:
     211        return 1
     212    else:
     213        if zmod_poly_equal(l, r):
     214            return 0
     215        i = deg_right
     216        rcoeff = zmod_poly_get_coeff_ui(r, i)
     217        lcoeff = zmod_poly_get_coeff_ui(l, i)
     218        while rcoeff == lcoeff and i > 0:
     219            i -= 1
     220            rcoeff = zmod_poly_get_coeff_ui(r, i)
     221            lcoeff = zmod_poly_get_coeff_ui(l, i)
     222        return cmp(lcoeff, rcoeff)
     223
     224cdef long celement_len(zmod_poly_t a, unsigned long n) except -2:
     225    """
     226    EXAMPLE:
     227        sage: P.<x> = GF(32003)[]
     228        sage: (x + 1).degree()
     229        1
     230        sage: (x).degree()
     231        1
     232        sage: P(0).degree()
     233        -1
     234
     235        sage: Q.<x> = GF(7)[]
     236        sage: (x + 1).degree()
     237        1
     238        sage: (x).degree()
     239        1
     240        sage: P(0).degree()
     241        -1
     242    """
     243    return <long>zmod_poly_length(a)
     244
     245cdef inline int celement_add(zmod_poly_t res, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     246    """
     247    EXAMPLE:
     248        sage: P.<x> = GF(32003)[]
     249        sage: x + 1
     250        x + 1
     251
     252        sage: Q.<x> = GF(7)[]
     253        sage: x + 1
     254        x + 1
     255    """
     256    zmod_poly_add(res, a, b)
     257
     258cdef inline int celement_sub(zmod_poly_t res, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     259    """
     260    EXAMPLE:
     261        sage: P.<x> = GF(32003)[]
     262        sage: x - 1
     263        x + 32002
     264
     265        sage: Q.<x> = GF(7)[]
     266        sage: x - 1
     267        x + 6
     268    """
     269    zmod_poly_sub(res, a, b)
     270
     271cdef inline int celement_neg(zmod_poly_t res, zmod_poly_t a, unsigned long n) except -2:
     272    """
     273    EXAMPLE:
     274        sage: P.<x> = GF(32003)[]
     275        sage: -(x + 2)
     276        32002*x + 32001
     277
     278        sage: Q.<x> = GF(7)[]
     279        sage: -(x + 2)
     280        6*x + 5       
     281    """
     282    zmod_poly_neg(res, a)
     283
     284cdef inline int celement_mul(zmod_poly_t res, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     285    """
     286    EXAMPLE:
     287        sage: P.<x> = GF(32003)[]
     288        sage: (x + 1) * (x + 2)
     289        x^2 + 3*x + 2
     290
     291        sage: Q.<x> = GF(7)[]
     292        sage: (x + 1) * (x + 2)
     293        x^2 + 3*x + 2
     294    """
     295    zmod_poly_mul(res, a, b)
     296
     297cdef inline int celement_div(zmod_poly_t res, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     298    raise NotImplementedError
     299
     300cdef inline int celement_floordiv(zmod_poly_t res, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     301    """
     302    EXAMPLE:
     303        sage: P.<x> = GF(32003)[]
     304        sage: (x + 1) // (x + 2)
     305        1
     306        sage: (3*x^2 + 1) // (x + 2)
     307        3*x + 31997
     308        sage: (x^2 + 3*x + 2)//(x + 1)
     309        x + 2
     310        sage: (x^2 + 3*x + 2)//(x + 2)
     311        x + 1
     312
     313        sage: Q.<x> = GF(7)[]
     314        sage: (x + 1) // (x + 2)
     315        1
     316        sage: (3*x^2 + 1) // (x + 2)
     317        3*x + 1
     318        sage: (x^2 + 3*x + 2)//(x + 1)
     319        x + 2
     320        sage: (x^2 + 3*x + 2)//(x + 2)
     321        x + 1
     322    """
     323    zmod_poly_div(res, a, b)
     324
     325cdef inline int celement_mod(zmod_poly_t res, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     326    """
     327    EXAMPLE:
     328        sage: P.<x> = GF(32003)[]
     329        sage: f = 24998*x^2 + 29761*x + 2252
     330        sage: g = 20778*x^2 + 15346*x + 12697
     331
     332        sage: f % g
     333        5815*x + 10280
     334
     335        sage: f^5 % g
     336        7231*x + 17274
     337
     338        sage: R.<x> = Integers(81)[]
     339        sage: f = x^7 + x + 1; g = x^3
     340        sage: r = f % g; r
     341        x + 1
     342        sage: g * x^4 + r
     343        x^7 + x + 1
     344
     345        sage: f = x^3 + 1
     346        sage: f % 3
     347        Traceback (most recent call last):
     348        ...
     349        ValueError: Leading coefficient of a must be invertible.
     350
     351        sage: f % 0
     352        Traceback (most recent call last):
     353        ...
     354        ZeroDivisionError
     355    """
     356    cdef zmod_poly_t q
     357    cdef unsigned long leadcoeff, modulus
     358   
     359    zmod_poly_init(q, n)
     360    leadcoeff = zmod_poly_get_coeff_ui(b, zmod_poly_degree(b))
     361    modulus = zmod_poly_modulus(b)
     362    if (leadcoeff > 1 and z_gcd(modulus,leadcoeff) != 1):
     363        raise ValueError("Leading coefficient of a must be invertible.")
     364
     365    zmod_poly_divrem(q, res, a, b)
     366    zmod_poly_clear(q)
     367
     368cdef inline int celement_quorem(zmod_poly_t q, zmod_poly_t r, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     369    """
     370    EXAMPLES:
     371        sage: R.<x> = Integers(125)[]
     372        sage: f = x^5+1; g = (x+1)^2
     373        sage: q, r = f.quo_rem(g)
     374        sage: q
     375        x^3 + 123*x^2 + 3*x + 121
     376        sage: r
     377        5*x + 5
     378        sage: q*g + r
     379        x^5 + 1
     380
     381        sage: x.quo_rem(5*x)
     382        Traceback (most recent call last):
     383        ...
     384        ValueError: Leading coefficient of a must be invertible.
     385
     386        sage: x.quo_rem(0)
     387        Traceback (most recent call last):
     388        ...
     389        ZeroDivisionError
     390    """
     391    cdef unsigned long leadcoeff, modulus
     392
     393    leadcoeff = zmod_poly_get_coeff_ui(b, zmod_poly_degree(b))
     394    modulus = zmod_poly_modulus(b)
     395    if (leadcoeff > 1 and z_gcd(modulus,leadcoeff) != 1):
     396        raise ValueError("Leading coefficient of a must be invertible.")
     397
     398    zmod_poly_divrem(q, r, a, b)
     399
     400cdef inline int celement_inv(zmod_poly_t res, zmod_poly_t a, unsigned long n) except -2:
     401    raise NotImplementedError
     402
     403cdef inline int celement_pow(zmod_poly_t res, zmod_poly_t x, long e, zmod_poly_t modulus, unsigned long n) except -2:
     404    """
     405    EXAMPLE:
     406        sage: P.<x> = GF(32003)[]
     407        sage: f = 24998*x^2 + 29761*x + 2252
     408
     409        sage: f*f
     410        9426*x^4 + 15477*x^3 + 6531*x^2 + 14980*x + 15030
     411        sage: f^2
     412        9426*x^4 + 15477*x^3 + 6531*x^2 + 14980*x + 15030
     413
     414        sage: f*f*f
     415        25062*x^6 + 30670*x^5 + 15816*x^4 + 20746*x^3 + 9142*x^2 + 5697*x + 20389
     416        sage: f^3
     417        25062*x^6 + 30670*x^5 + 15816*x^4 + 20746*x^3 + 9142*x^2 + 5697*x + 20389
     418
     419        sage: f*f*f*f*f
     420        20269*x^10 + 20535*x^9 + 7313*x^8 + 7311*x^7 + 16853*x^6 + 142*x^5 + 23853*x^4 + 12065*x^3 + 516*x^2 + 8473*x + 17945
     421        sage: f^5
     422        20269*x^10 + 20535*x^9 + 7313*x^8 + 7311*x^7 + 16853*x^6 + 142*x^5 + 23853*x^4 + 12065*x^3 + 516*x^2 + 8473*x + 17945
     423
     424        sage: f^0
     425        1
     426
     427        sage: f^1
     428        24998*x^2 + 29761*x + 2252
     429
     430        sage: f^-1
     431        1/(24998*x^2 + 29761*x + 2252)
     432
     433        sage: f^-5
     434        1/(20269*x^10 + 20535*x^9 + 7313*x^8 + 7311*x^7 + 16853*x^6 + 142*x^5 + 23853*x^4 + 12065*x^3 + 516*x^2 + 8473*x + 17945)
     435
     436     Testing the modulus:
     437
     438        sage: g = 20778*x^2 + 15346*x + 12697
     439
     440        sage: pow(f, 2, g)
     441        15328*x + 6968
     442        sage: f^2 % g
     443        15328*x + 6968
     444
     445        sage: pow(f, -2, g)
     446        1/(15328*x + 6968)
     447        sage: (f^2 % g)^-1
     448        1/(15328*x + 6968)
     449
     450        sage: pow(f, 5, g)
     451        7231*x + 17274
     452        sage: f^5 % g
     453        7231*x + 17274
     454    """
     455    cdef zmod_poly_t pow2
     456    cdef zmod_poly_t q
     457    cdef zmod_poly_t tmp
     458
     459    zmod_poly_init(q, n)
     460    zmod_poly_init(tmp, n)
     461
     462    if zmod_poly_degree(x) == 1 and zmod_poly_get_coeff_ui(x,0) == 0 and zmod_poly_get_coeff_ui(x,1) == 1:
     463        zmod_poly_zero(res)
     464        zmod_poly_set_coeff_ui(res,e,1)
     465    elif e == 0:
     466        zmod_poly_zero(res)
     467        zmod_poly_set_coeff_ui(res,0,1)
     468    elif e == 1:
     469        zmod_poly_set(res, x)
     470    elif e == 2:
     471        zmod_poly_sqr(res, x)
     472    else:
     473        if res == x:
     474            zmod_poly_set(tmp, x)
     475            x = tmp
     476        zmod_poly_init(pow2, n)
     477        zmod_poly_set(pow2, x)
     478        if e % 2:
     479            zmod_poly_set(res, x)
     480        else:
     481            zmod_poly_zero(res)
     482            zmod_poly_set_coeff_ui(res, 0, 1)
     483        e = e >> 1
     484        while(e != 0):
     485            zmod_poly_sqr(pow2, pow2)
     486            if e % 2:
     487                zmod_poly_mul(res, res, pow2)
     488            e = e >> 1
     489            if modulus != NULL:
     490                zmod_poly_divrem(q, res, res, modulus)
     491        zmod_poly_clear(pow2)
     492
     493    if modulus != NULL:
     494        zmod_poly_divrem(q, res, res, modulus)
     495    zmod_poly_clear(q)
     496    zmod_poly_clear(tmp)
     497
     498cdef inline int celement_gcd(zmod_poly_t res, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     499    """
     500    EXAMPLE:
     501        sage: P.<x> = GF(32003)[]
     502        sage: f = P.random_element(degree=4); f
     503        16660*x^4 + 10640*x^3 + 1430*x^2 + 16460*x + 3566
     504        sage: g = P.random_element(degree=3); g
     505        28452*x^3 + 2561*x^2 + 22429*x + 5847
     506        sage: h = P.random_element(degree=2); h
     507        24731*x^2 + 28238*x + 18622
     508        sage: F = f*g; F
     509        13887*x^7 + 19164*x^6 + 25146*x^5 + 25986*x^4 + 21143*x^3 + 14830*x^2 + 14916*x + 16449
     510        sage: G = f*h; G
     511        11838*x^6 + 10154*x^5 + 15609*x^4 + 26164*x^3 + 11353*x^2 + 8656*x + 31830
     512        sage: d = (F).gcd(G); d
     513        x^4 + 18557*x^3 + 22917*x^2 + 30813*x + 4914
     514        sage: (F//d)*d == F
     515        True
     516        sage: (G//d)*d == G
     517        True
     518
     519        sage: Q.<x> = GF(7)[]
     520        sage: f = Q.random_element(degree=4); f
     521        5*x^4 + 3*x^3 + 6*x^2 + 6*x + 1
     522        sage: g = Q.random_element(degree=3); g
     523        2*x^3 + 5*x^2 + 2*x + 3
     524        sage: h = Q.random_element(degree=2); h
     525        4*x^2 + 4*x + 6
     526        sage: F = f*g; F
     527        3*x^7 + 3*x^6 + 2*x^5 + 4*x^3 + 6*x + 3
     528        sage: G = f*h; G
     529        6*x^6 + 4*x^5 + 3*x^4 + 3*x^3 + x^2 + 5*x + 6
     530        sage: d = (F).gcd(G); d
     531        x^4 + 2*x^3 + 4*x^2 + 4*x + 3
     532        sage: (F//d)*d == F
     533        True
     534        sage: (G//d)*d == G
     535        True
     536    """
     537    if celement_is_zero(b, n):
     538        zmod_poly_set(res, a)
     539        return 0
     540
     541    zmod_poly_gcd(res, a, b)
     542    cdef unsigned long leadcoeff = zmod_poly_get_coeff_ui(res, zmod_poly_degree(res))
     543    cdef unsigned long modulus = zmod_poly_modulus(res)
     544    if z_gcd(modulus,leadcoeff) == 1:
     545        zmod_poly_make_monic(res, res)
     546
     547cdef inline int celement_xgcd(zmod_poly_t res, zmod_poly_t s, zmod_poly_t t, zmod_poly_t a, zmod_poly_t b, unsigned long n) except -2:
     548    """
     549    EXAMPLE:
     550        sage: P.<x> = GF(32003)[]
     551        sage: f = P.random_element(degree=4); f
     552        16660*x^4 + 10640*x^3 + 1430*x^2 + 16460*x + 3566
     553        sage: g = P.random_element(degree=3); g
     554        28452*x^3 + 2561*x^2 + 22429*x + 5847
     555        sage: h = P.random_element(degree=2); h
     556        24731*x^2 + 28238*x + 18622
     557        sage: F = f*g; F
     558        13887*x^7 + 19164*x^6 + 25146*x^5 + 25986*x^4 + 21143*x^3 + 14830*x^2 + 14916*x + 16449
     559        sage: G = f*h; G
     560        11838*x^6 + 10154*x^5 + 15609*x^4 + 26164*x^3 + 11353*x^2 + 8656*x + 31830
     561        sage: d,s,t = (F).xgcd(G); d
     562        x^4 + 18557*x^3 + 22917*x^2 + 30813*x + 4914
     563        sage: (F//d)*d == F
     564        True
     565        sage: (G//d)*d == G
     566        True
     567
     568        sage: Q.<x> = GF(7)[]
     569        sage: f = Q.random_element(degree=4); f
     570        5*x^4 + 3*x^3 + 6*x^2 + 6*x + 1
     571        sage: g = Q.random_element(degree=3); g
     572        2*x^3 + 5*x^2 + 2*x + 3
     573        sage: h = Q.random_element(degree=2); h
     574        4*x^2 + 4*x + 6
     575        sage: F = f*g; F
     576        3*x^7 + 3*x^6 + 2*x^5 + 4*x^3 + 6*x + 3
     577        sage: G = f*h; G
     578        6*x^6 + 4*x^5 + 3*x^4 + 3*x^3 + x^2 + 5*x + 6
     579        sage: d,s,t = (F).xgcd(G); d
     580        x^4 + 2*x^3 + 4*x^2 + 4*x + 3
     581        sage: (F//d)*d == F
     582        True
     583        sage: (G//d)*d == G
     584        True
     585    """
     586    zmod_poly_xgcd(res, s, t, a, b)
     587
     588
     589cdef inline unsigned long zmod_poly_evaluate_horner(zmod_poly_t _poly, unsigned long c):
     590    """
     591    Evaluate _poly at c using Horner's rule.
     592
     593    AUTHOR:
     594        -- Bill Hart (2009-01)
     595
     596    EXAMPLE:
     597        sage: P.<x> = PolynomialRing(GF(7))
     598        sage: f= x^2 + 1
     599        sage: f(0)
     600        1
     601        sage: f(2)
     602        5
     603        sage: f(3)
     604        3
     605
     606    NOTE:
     607        This function should be removed as soon as we update the the
     608        new FLINT version which has this function natively.
     609    """
     610    cdef zmod_poly_struct *poly = <zmod_poly_struct *>_poly
     611    cdef long i
     612    cdef unsigned long bits
     613    if poly.length == 0:
     614        return 0
     615
     616    cdef unsigned long p = poly.p
     617    cdef double pinv = poly.p_inv
     618
     619    cdef unsigned long value = poly.coeffs[poly.length - 1]
     620
     621    if FLINT_BITS == 64:
     622        bits = FLINT_BIT_COUNT(p)
     623        if bits < FLINT_D_BITS:
     624            for i in xrange(poly.length - 2, -1, -1):
     625                value = z_addmod(z_mulmod_precomp(value, c, p, pinv), poly.coeffs[i], p)
     626        else:
     627            for i in xrange(poly.length - 2, -1, -1):
     628                value = z_addmod(z_mulmod2_precomp(value, c, p, pinv), poly.coeffs[i], p)
     629    else:
     630        for i in xrange(poly.length - 2, -1, -1):
     631            value = z_addmod(z_mulmod2_precomp(value, c, p, pinv), poly.coeffs[i], p)
     632    return value
  • sage/libs/ntl/ntl_GF2X_linkage.pxi

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/libs/ntl/ntl_GF2X_linkage.pxi
    a b AUTHOR: 
    1616
    1717from sage.libs.ntl.ntl_GF2_decl cimport *, GF2_c
    1818from sage.libs.ntl.ntl_GF2X_decl cimport *, GF2X_c, GF2XModulus_c
     19
     20
     21cdef GF2X_c *celement_new(long parent):
     22    """
     23    EXAMPLE:
     24        sage: P.<x> = GF(2)[]
     25    """
     26    cdef GF2X_c *e = GF2X_new()
     27    return e
     28
     29cdef int celement_delete(GF2X_c *e, long parent):
     30    """
     31    EXAMPLE:
     32        sage: P.<x> = GF(2)[]
     33        sage: del x
     34    """
     35    GF2X_delete(e)
    1936
    2037cdef int celement_construct(GF2X_c *e, long parent):
    2138    """
    cdef inline int celement_cmp(GF2X_c *a,  
    122139        False
    123140        sage: x > 1
    124141        True
     142
     143        sage: f = x^64 + x^20 + 1
     144        sage: g = x^63 + x^20 + 1
     145        sage: f > g
     146        True
     147
     148        sage: f = x^64 + x^10 + 1
     149        sage: g = x^64 + x^20 + 1
     150        sage: f < g
     151        True
     152
     153        sage: f = x^64 + x^20
     154        sage: g = x^64 + x^20 + 1
     155        sage: f < g
     156        True
    125157    """
    126158    cdef bint t
    127159    cdef long diff
     160    cdef long ca, cb
    128161    diff = GF2X_NumBits(a[0]) - GF2X_NumBits(b[0])
    129162    if diff > 0:
    130163        return 1
    131     elif diff == 0:
    132         t = GF2X_equal(a[0], b[0])
    133         return not t
     164    elif diff < 0:
     165        return -1
    134166    else:
    135         return -1
    136 
    137 cdef inline long celement_hash(GF2X_c *a, long parent) except -2:
    138     """
    139     EXAMPLE:
    140         sage: P.<x> = GF(2)[]
    141         sage: {x:1}
    142         {x: 1}
    143     """
    144     cdef long _hex = GF2XHexOutput_c[0]
    145     GF2XHexOutput_c[0] = 1
    146     s = GF2X_to_PyString(a)
    147     GF2XHexOutput_c[0] = _hex
    148     return hash(s)
     167        for i in xrange(GF2X_NumBits(a[0])-1, -1, -1):
     168            ca = GF2_conv_to_long(GF2X_coeff(a[0], i))
     169            cb = GF2_conv_to_long(GF2X_coeff(b[0], i))
     170            if ca < cb:
     171                return -1
     172            elif ca > cb:
     173                return 1
     174    return 0
    149175
    150176cdef long celement_len(GF2X_c *a, long parent) except -2:
    151177    """
    152178    EXAMPLE:
    153179        sage: P.<x> = GF(2)[]
    154         sage: len(x)
    155         2
    156         sage: len(x+1)
    157         2
     180        sage: x.degree()
     181        1
     182        sage: (x+1).degree()
     183        1
    158184    """
    159185    return int(GF2X_NumBits(a[0]))
    160186
    cdef inline int celement_gcd(GF2X_c* res 
    286312        x
    287313    """
    288314    GF2X_GCD(res[0], a[0], b[0])   
     315
     316cdef inline int celement_xgcd(GF2X_c* res, GF2X_c* s, GF2X_c *t, GF2X_c* a, GF2X_c *b, long parent) except -2:
     317    raise NotImplementedError
  • sage/modular/modsym/apply.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/modular/modsym/apply.pxd
    a b include "../../ext/stdsage.pxi" 
    11include "../../ext/stdsage.pxi"
    22include "../../ext/cdefs.pxi"
    3 include "../../libs/flint/flint.pxi"
     3
     4from sage.libs.flint.flint cimport *
    45include "../../libs/flint/fmpz_poly.pxi"
    56
    67
  • sage/modular/modsym/heilbronn.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/modular/modsym/heilbronn.pyx
    a b include '../../ext/cdefs.pxi' 
    2424include '../../ext/cdefs.pxi'
    2525include '../../ext/interrupt.pxi'
    2626include '../../ext/stdsage.pxi'
    27 include "../../libs/flint/flint.pxi"
     27from sage.libs.flint.flint cimport *
    2828include "../../libs/flint/fmpz_poly.pxi"
    2929
    3030cimport p1list
  • sage/rings/finite_field.py

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/finite_field.py
    a b class FiniteFieldFactory(UniqueFactory): 
    177177        sage: f = K.modulus(); f
    178178        x^5 + 4*x + 1
    179179        sage: type(f)
    180         <type 'sage.rings.polynomial.polynomial_modn_dense_ntl.Polynomial_dense_mod_p'>
     180         <type 'sage.rings.polynomial.polynomial_zmod_flint.Polynomial_zmod_flint'>
    181181
    182182    The modulus must be irreducible:
    183183        sage: K.<a> = GF(5**5, name='a', modulus=x^5 - x )
  • sage/rings/finite_field_ext_pari.py

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/finite_field_ext_pari.py
    a b class FiniteField_ext_pari(FiniteField_g 
    615615        Return the hash of this field.
    616616
    617617        EXAMPLES:
    618             sage: hash(GF(9,'a'))
    619             -417021630                     # 32-bit
    620             1006006598732398914            # 64-bit
    621             sage: hash(GF(9,'b'))
    622             995034271                      # 32-bit
    623             8600900932991911071            # 64-bit
     618            sage: {GF(9,'a'): 1} # indirect doctest
     619            {Finite Field in a of size 3^2: 1}
     620            sage: {GF(9,'b'): 1} # indirect doctest
     621            {Finite Field in b of size 3^2: 1}
    624622        """
    625623        return hash((self.__order, self.variable_name(), self.__modulus))
    626624
  • sage/rings/finite_field_givaro.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/finite_field_givaro.pyx
    a b cdef class FiniteField_givaro(FiniteFiel 
    876876        characterstic polynomial and the string 'givaro'
    877877
    878878        EXAMPLES:
    879             sage: hash(GF(3^4, 'a'))
    880             -417021630                # 32-bit
    881             1006006598732398914       # 64-bit
     879            sage: {GF(3^4, 'a'):1} #indirect doctest
     880            {Finite Field in a of size 3^4: 1}
    882881        """
    883882        if self._hash is None:
    884883            pass
    cdef class FiniteField_givaroElement(Fin 
    17791778            sage: f = (b^2+1).polynomial(); f
    17801779            b + 4
    17811780            sage: type(f)
    1782             <type 'sage.rings.polynomial.polynomial_modn_dense_ntl.Polynomial_dense_mod_p'>
     1781            <type 'sage.rings.polynomial.polynomial_zmod_flint.Polynomial_zmod_flint'>
    17831782            sage: parent(f)
    17841783            Univariate Polynomial Ring in b over Finite Field of size 5       
    17851784        """
  • sage/rings/integer_mod.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/integer_mod.pyx
    a b cdef class IntegerMod_abstract(sage.stru 
    503503            sage: a.polynomial()
    504504            1
    505505            sage: type(a.polynomial())
    506             <type 'sage.rings.polynomial.polynomial_modn_dense_ntl.Polynomial_dense_mod_p'>
     506            <type 'sage.rings.polynomial.polynomial_zmod_flint.Polynomial_zmod_flint'>
    507507        """
    508508        R = self.parent()[var]
    509509        return R(self)
  • sage/rings/polynomial/polynomial_element.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_element.pyx
    a b cdef class Polynomial(CommutativeAlgebra 
    22642264            if unit is None:
    22652265                c *= f.leading_coefficient()
    22662266            F.append((f,e))
    2267            
    22682267        if unit is None:
    2269            
    22702268            unit = R.base_ring()(self.leading_coefficient()/c)
    2271            
    22722269        if not unit.is_unit():
    22732270           
    22742271            F.append((R(unit), ZZ(1)))
  • sage/rings/polynomial/polynomial_gf2x.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_gf2x.pyx
    a b cdef class Polynomial_GF2X(Polynomial_te 
    6060            sage: f[1]
    6161            0
    6262        """
    63         cdef long c
     63        cdef long c = 0
    6464        if 0 <= i < GF2X_NumBits(self.x):
    6565            c = GF2_conv_to_long(GF2X_coeff(self.x, i))
    6666        return self._parent.base_ring()(c)
  • sage/rings/polynomial/polynomial_integer_dense_flint.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_integer_dense_flint.pxd
    a b include "../../ext/cdefs.pxi" 
    11include "../../ext/cdefs.pxi"
    22include "../../libs/flint/fmpz.pxi"
    33include "../../libs/flint/fmpz_poly.pxi"
    4 include "../../libs/flint/ntl_interface.pxi"
     4
     5from sage.libs.flint.ntl_interface cimport *
    56
    67from sage.rings.polynomial.polynomial_element cimport Polynomial
    78from sage.rings.integer cimport Integer
  • sage/rings/polynomial/polynomial_modn_dense_ntl.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_modn_dense_ntl.pyx
    a b cdef class Polynomial_dense_mod_n(Polyno 
    316316
    317317        EXAMPLES:
    318318            sage: R.<x> = PolynomialRing(Integers(100))
    319             sage: R([1,-2,3])
     319            sage: from sage.rings.polynomial.polynomial_modn_dense_ntl import Polynomial_dense_mod_n as poly_modn_dense
     320            sage: poly_modn_dense(R, ([1,-2,3]))
    320321            3*x^2 + 98*x + 1
    321             sage: f = R(0)
     322            sage: f = poly_modn_dense(R, 0)
    322323            sage: f.ntl_set_directly([1,-2,3])
    323324            sage: f
    324325            3*x^2 + 98*x + 1
    cdef class Polynomial_dense_mod_n(Polyno 
    357358    def resultant(self, other, variable=None):
    358359        return polynomial_singular_interface.resultant_func(self, other, variable)
    359360
     361    def small_roots(self, *args, **kwds):
     362        r"""
     363        See \code{sage.rings.polynomial.polynomial_modn_dense_ntl.small_roots}
     364        for the documentation of this function.
    360365
    361     def small_roots(self, X=None, beta=1.0, epsilon=None, **kwds):
    362         r"""
    363         Let $N$ be the characteristic of the base ring this polynomial
    364         is defined over: \code{N = self.base_ring().characteristic()}.
    365         This method returns small roots of this polynomial modulo some
    366         factor $b$ of $N$ with the constraint that $b >= N^\beta$.
    367         Small in this context means that if $x$ is a root of $f$
    368         modulo $b$ then $|x| < X$. This $X$ is either provided by the
    369         user or the maximum $X$ is chosen such that this algorithm
    370         terminates in polynomial time. If $X$ is chosen automatically
    371         it is $X = ceil(1/2 N^{\beta^2/\delta - \epsilon})$.
    372         The algorithm may also return some roots which are larger than $X$.
    373         `This algorithm' in this context means Coppersmith's algorithm for
    374         finding small roots using the LLL algorithm. The
    375         implementation of this algorithm follows Alexander May's PhD
    376         thesis referenced below.
    377 
    378         INPUT:
    379           X -- an absolute bound for the root (default: see above)
    380           beta -- compute a root mod $b$ where $b$ is a factor of $N$
    381                    and $b >= N^\beta$.  (default: 1.0 => $b = N$)
    382           epsilon -- the parameter $\epsilon$ described above
    383                      (default: beta/8)
    384           **kwds -- passed through to LLL method of
    385                     \code{Matrix_integer_dense}.
    386 
    387         EXAMPLES:
    388 
    389         First consider a small example:
    390 
     366        EXAMPLE:
    391367            sage: N = 10001
    392368            sage: K = Zmod(10001)
    393369            sage: P.<x> = PolynomialRing(K)
    394370            sage: f = x^3 + 10*x^2 + 5000*x - 222
    395 
    396         This polynomial has no roots without modular reduction
    397         (i.e. over $\ZZ$):
    398 
    399             sage: f.change_ring(ZZ).roots()
    400             []
    401 
    402         To compute its roots we need to factor the modulus $N$ and use
    403         the chinese remainder theorem:
    404 
    405             sage: p,q = N.prime_divisors()
    406             sage: f.change_ring(GF(p)).roots()
    407             [(4, 1)]
    408             sage: f.change_ring(GF(q)).roots()
    409             [(4, 1)]
    410 
    411             sage: crt(4, 4, p, q)
    412             4
    413 
    414         This root is quite small compared to $N$ so we can attempt to
    415         recover it without factoring $N$ using Coppersmith's small root
    416         method:
    417 
    418371            sage: f.small_roots()
    419372            [4]
     373        """
     374        return small_roots(self, *args, **kwds)
    420375
    421         An application of this method is to consider RSA. We are using
    422         512-bit RSA with public exponent $e=3$ to encrypt a 56-bit DES
    423         key. Because it would be easy to attack this setting if no
    424         padding was used we pad the key $K$ with 1s to get a large
    425         number.
     376def small_roots(self, X=None, beta=1.0, epsilon=None, **kwds):
     377    r"""
     378    Let $N$ be the characteristic of the base ring this polynomial
     379    is defined over: \code{N = self.base_ring().characteristic()}.
     380    This method returns small roots of this polynomial modulo some
     381    factor $b$ of $N$ with the constraint that $b >= N^\beta$.
     382    Small in this context means that if $x$ is a root of $f$
     383    modulo $b$ then $|x| < X$. This $X$ is either provided by the
     384    user or the maximum $X$ is chosen such that this algorithm
     385    terminates in polynomial time. If $X$ is chosen automatically
     386    it is $X = ceil(1/2 N^{\beta^2/\delta - \epsilon})$.
     387    The algorithm may also return some roots which are larger than $X$.
     388    `This algorithm' in this context means Coppersmith's algorithm for
     389    finding small roots using the LLL algorithm. The
     390    implementation of this algorithm follows Alexander May's PhD
     391    thesis referenced below.
    426392
    427             sage: Nbits, Kbits = 512, 56
    428             sage: e = 3
     393    INPUT:
     394      X -- an absolute bound for the root (default: see above)
     395      beta -- compute a root mod $b$ where $b$ is a factor of $N$
     396               and $b >= N^\beta$.  (default: 1.0 => $b = N$)
     397      epsilon -- the parameter $\epsilon$ described above
     398                 (default: beta/8)
     399      **kwds -- passed through to LLL method of
     400                \code{Matrix_integer_dense}.
    429401
    430         We choose two primes of size 256-bit each.
     402    EXAMPLES:
    431403
    432             sage: p = 2^256 + 2^8 + 2^5 + 2^3 + 1
    433             sage: q = 2^256 + 2^8 + 2^5 + 2^3 + 2^2 + 1
    434             sage: N = p*q
    435             sage: ZmodN = Zmod( N )
     404    First consider a small example:
    436405
    437         We choose a random key
     406        sage: N = 10001
     407        sage: K = Zmod(10001)
     408        sage: P.<x> = PolynomialRing(K)
     409        sage: f = x^3 + 10*x^2 + 5000*x - 222
    438410
    439             sage: K = ZZ.random_element(0, 2^Kbits)
     411    This polynomial has no roots without modular reduction
     412    (i.e. over $\ZZ$):
    440413
    441         and pad it with 512-56=456 1s
     414        sage: f.change_ring(ZZ).roots()
     415        []
    442416
    443             sage: Kdigits = K.digits(2)
    444             sage: M = [0]*Kbits + [1]*(Nbits-Kbits)
    445             sage: for i in range(len(Kdigits)): M[i] = Kdigits[i]
     417    To compute its roots we need to factor the modulus $N$ and use
     418    the chinese remainder theorem:
    446419
    447             sage: M = ZZ(M, 2)
     420        sage: p,q = N.prime_divisors()
     421        sage: f.change_ring(GF(p)).roots()
     422        [(4, 1)]
     423        sage: f.change_ring(GF(q)).roots()
     424        [(4, 1)]
    448425
    449         Now we encrypt the resulting message:
     426        sage: crt(4, 4, p, q)
     427        4
    450428
    451             sage: C = ZmodN(M)^e
     429    This root is quite small compared to $N$ so we can attempt to
     430    recover it without factoring $N$ using Coppersmith's small root
     431    method:
    452432
    453         To recover $K$ we consider the following polynomial modulo $N$:
     433        sage: f.small_roots()
     434        [4]
    454435
    455             sage: P.<x> = PolynomialRing(ZmodN)
    456             sage: f = (2^Nbits - 2^Kbits + x)^e - C
     436    An application of this method is to consider RSA. We are using
     437    512-bit RSA with public exponent $e=3$ to encrypt a 56-bit DES
     438    key. Because it would be easy to attack this setting if no
     439    padding was used we pad the key $K$ with 1s to get a large
     440    number.
    457441
    458         and recover its small roots:
     442        sage: Nbits, Kbits = 512, 56
     443        sage: e = 3
    459444
    460             sage: Kbar = f.small_roots()[0]
    461             sage: K == Kbar
    462             True
     445    We choose two primes of size 256-bit each.
    463446
    464         The same algorithm can be used to factor $N = pq$ if partial
    465         knowledge about $q$ is available. This example is from the
    466         Magma handbook:
     447        sage: p = 2^256 + 2^8 + 2^5 + 2^3 + 1
     448        sage: q = 2^256 + 2^8 + 2^5 + 2^3 + 2^2 + 1
     449        sage: N = p*q
     450        sage: ZmodN = Zmod( N )
    467451
    468         First, we set up $p$,$q$ and $N$.
     452    We choose a random key
    469453
    470             sage: length = 512
    471             sage: hidden = 110
    472             sage: p = next_prime(2^int(round(length/2)))
    473             sage: q = next_prime( round(pi.n()*p) )
    474             sage: N = p*q
     454        sage: K = ZZ.random_element(0, 2^Kbits)
    475455
    476         Now we disturb the low 110 bits of $q$
     456    and pad it with 512-56=456 1s
    477457
    478             sage: qbar = q + ZZ.random_element(0,2^hidden-1)
     458        sage: Kdigits = K.digits(2)
     459        sage: M = [0]*Kbits + [1]*(Nbits-Kbits)
     460        sage: for i in range(len(Kdigits)): M[i] = Kdigits[i]
    479461
    480         And try to recover $q$ from it:
     462        sage: M = ZZ(M, 2)
    481463
    482             sage: F.<x> = PolynomialRing(Zmod(N))
    483             sage: f = x - qbar
     464    Now we encrypt the resulting message:
    484465
    485         We know that the error is $<= 2^{hidden}-1$ and that the modulus
    486         we are looking for is $>= sqrt(N)$.
     466        sage: C = ZmodN(M)^e
    487467
    488             sage: set_verbose(2)
    489             sage: d = f.small_roots(X=2^hidden-1, beta=0.5)[0] # time random
    490             verbose 2 (<module>) m = 4
    491             verbose 2 (<module>) t = 4
    492             verbose 2 (<module>) X = 1298074214633706907132624082305023
    493             verbose 1 (<module>) LLL of 8x8 matrix (algorithm fpLLL:wrapper)
    494             verbose 1 (<module>) LLL finished (time = 0.006998)
    495             sage: q == qbar - d
    496             True
     468    To recover $K$ we consider the following polynomial modulo $N$:
    497469
    498         REFERENCES:
    499             Don Coppersmith. Finding a small root of a univariate
    500                 modular equation.  In Advances in Cryptology,
    501                 EuroCrypt 1996, volume 1070 of Lecture Notes in
    502                 Computer Science, p. 155--165. Springer, 1996.
    503                 http://cr.yp.to/bib/2001/coppersmith.pdf
     470        sage: P.<x> = PolynomialRing(ZmodN)
     471        sage: f = (2^Nbits - 2^Kbits + x)^e - C
    504472
    505             Alexander May. New RSA Vulnerabilities Using Lattice
    506                 Reduction Methods.  PhD thesis, University of
    507                 Paderborn, 2003
    508                 http://www.informatik.tu-darmstadt.de/KP/publications/03/bp.ps
    509         """
    510         from sage.misc.misc import verbose
    511         from sage.matrix.constructor import Matrix
    512         from sage.rings.all import RR
    513      
    514         N = self.parent().characteristic()
     473    and recover its small roots:
    515474
    516         if not self.is_monic():
    517             raise ArithmeticError, "Polynomial must be monic."
     475        sage: Kbar = f.small_roots()[0]
     476        sage: K == Kbar
     477        True
    518478
    519         beta = RR(beta)
    520         if beta <= 0.0 or beta > 1.0:
    521             raise ValueError, "0.0 < beta <= 1.0 not satisfied."
     479    The same algorithm can be used to factor $N = pq$ if partial
     480    knowledge about $q$ is available. This example is from the
     481    Magma handbook:
    522482
    523         f = self.change_ring(ZZ)
     483    First, we set up $p$,$q$ and $N$.
    524484
    525         P,(x,) = f.parent().objgens()
     485        sage: length = 512
     486        sage: hidden = 110
     487        sage: p = next_prime(2^int(round(length/2)))
     488        sage: q = next_prime( round(pi.n()*p) )
     489        sage: N = p*q
    526490
    527         delta = f.degree()
     491    Now we disturb the low 110 bits of $q$
    528492
    529         if epsilon is None:
    530             epsilon = beta/8
    531         verbose("epsilon = %d"%epsilon, level=2)
     493        sage: qbar = q + ZZ.random_element(0,2^hidden-1)
    532494
    533         m = max(beta**2/(delta * epsilon), 7*beta/delta).ceil()
    534         verbose("m = %d"%m, level=2)
     495    And try to recover $q$ from it:
    535496
    536         t = int( ( delta*m*(1/beta -1) ).floor() )
    537         verbose("t = %d"%t, level=2)
     497        sage: F.<x> = PolynomialRing(Zmod(N))
     498        sage: f = x - qbar
    538499
    539         if X is None:
    540             X = (0.5 * N**(beta**2/delta - epsilon)).ceil()
    541         verbose("X = %s"%X, level=2)
     500    We know that the error is $<= 2^{hidden}-1$ and that the modulus
     501    we are looking for is $>= sqrt(N)$.
    542502
    543         # we could do this much faster, but this is a cheap step
    544         # compared to LLL
    545         g  = [x**j * N**(m-i) * f**i for i in range(m) for j in range(delta) ]
    546         g.extend([x**i * f**m for i in range(t)]) # h
     503        sage: set_verbose(2)
     504        sage: d = f.small_roots(X=2^hidden-1, beta=0.5)[0] # time random
     505        verbose 2 (<module>) m = 4
     506        verbose 2 (<module>) t = 4
     507        verbose 2 (<module>) X = 1298074214633706907132624082305023
     508        verbose 1 (<module>) LLL of 8x8 matrix (algorithm fpLLL:wrapper)
     509        verbose 1 (<module>) LLL finished (time = 0.006998)
     510        sage: q == qbar - d
     511        True
    547512
    548         B = Matrix(ZZ, len(g), delta*m + max(delta,t) )
    549         for i in range(B.nrows()):
    550             for j in range( g[i].degree()+1 ):
    551                 B[i,j] = g[i][j]*X**j
     513    REFERENCES:
     514        Don Coppersmith. Finding a small root of a univariate
     515            modular equation.  In Advances in Cryptology,
     516            EuroCrypt 1996, volume 1070 of Lecture Notes in
     517            Computer Science, p. 155--165. Springer, 1996.
     518            http://cr.yp.to/bib/2001/coppersmith.pdf
    552519
    553         B =  B.LLL(**kwds)
     520        Alexander May. New RSA Vulnerabilities Using Lattice
     521            Reduction Methods.  PhD thesis, University of
     522            Paderborn, 2003
     523            http://www.informatik.tu-darmstadt.de/KP/publications/03/bp.ps
     524    """
     525    from sage.misc.misc import verbose
     526    from sage.matrix.constructor import Matrix
     527    from sage.rings.all import RR
    554528
    555         f = sum([ZZ(B[0,i]//X**i)*x**i for i in range(B.ncols())])
    556         R = f.roots()
     529    N = self.parent().characteristic()
    557530
    558         ZmodN = self.base_ring()
    559         roots = set([ZmodN(r) for r,m in R if abs(r) <= X])
    560         Nbeta = N**beta
    561         return [root for root in roots if N.gcd(ZZ(self(root))) >= Nbeta]
     531    if not self.is_monic():
     532        raise ArithmeticError, "Polynomial must be monic."
     533
     534    beta = RR(beta)
     535    if beta <= 0.0 or beta > 1.0:
     536        raise ValueError, "0.0 < beta <= 1.0 not satisfied."
     537
     538    f = self.change_ring(ZZ)
     539
     540    P,(x,) = f.parent().objgens()
     541
     542    delta = f.degree()
     543
     544    if epsilon is None:
     545        epsilon = beta/8
     546    verbose("epsilon = %d"%epsilon, level=2)
     547
     548    m = max(beta**2/(delta * epsilon), 7*beta/delta).ceil()
     549    verbose("m = %d"%m, level=2)
     550
     551    t = int( ( delta*m*(1/beta -1) ).floor() )
     552    verbose("t = %d"%t, level=2)
     553
     554    if X is None:
     555        X = (0.5 * N**(beta**2/delta - epsilon)).ceil()
     556    verbose("X = %s"%X, level=2)
     557
     558    # we could do this much faster, but this is a cheap step
     559    # compared to LLL
     560    g  = [x**j * N**(m-i) * f**i for i in range(m) for j in range(delta) ]
     561    g.extend([x**i * f**m for i in range(t)]) # h
     562
     563    B = Matrix(ZZ, len(g), delta*m + max(delta,t) )
     564    for i in range(B.nrows()):
     565        for j in range( g[i].degree()+1 ):
     566            B[i,j] = g[i][j]*X**j
     567
     568    B =  B.LLL(**kwds)
     569
     570    f = sum([ZZ(B[0,i]//X**i)*x**i for i in range(B.ncols())])
     571    R = f.roots()
     572
     573    ZmodN = self.base_ring()
     574    roots = set([ZmodN(r) for r,m in R if abs(r) <= X])
     575    Nbeta = N**beta
     576    return [root for root in roots if N.gcd(ZZ(self(root))) >= Nbeta]
    562577
    563578cdef class Polynomial_dense_modn_ntl_zz(Polynomial_dense_mod_n):
    564579   
    cdef class Polynomial_dense_modn_ntl_zz( 
    614629       
    615630        EXAMPLES:
    616631            sage: R.<x> = Integers(100)[]
    617             sage: f = x^3 + 5
     632            sage: from sage.rings.polynomial.polynomial_modn_dense_ntl import Polynomial_dense_mod_n as poly_modn_dense
     633            sage: f = poly_modn_dense(R,[5,0,0,1])
    618634            sage: f.int_list()
    619635            [5, 0, 0, 1]
    620636            sage: [type(a) for a in f.int_list()]
    cdef class Polynomial_dense_modn_ntl_zz( 
    810826            5*x + 5
    811827            sage: q*g + r
    812828            x^5 + 1
    813 
    814             sage: x.quo_rem(5*x)
    815             Traceback (most recent call last):
    816             ...
    817             RuntimeError: InvMod: inverse undefined
    818829        """
    819830        if PY_TYPE(self) != PY_TYPE(right) or self._parent is not (<Element>right)._parent:
    820831            self, right = canonical_coercion(self, right)
    cdef class Polynomial_dense_modn_ntl_zz( 
    863874            x + 1
    864875            sage: g * x^4 + r
    865876            x^7 + x + 1
    866 
    867             sage: f = x^3 + 1
    868             sage: f % 3
    869             Traceback (most recent call last):
    870             ...
    871             RuntimeError: InvMod: inverse undefined
    872877        """
    873878        if PY_TYPE(self) != PY_TYPE(right) or (<Element>self)._parent is not (<Element>right)._parent:
    874879            self, right = canonical_coercion(self, right)
    cdef class Polynomial_dense_modn_ntl_zz( 
    959964        cdef Polynomial_dense_modn_ntl_zz r = self._new()
    960965        zz_pX_diff(r.x, self.x)
    961966        return r
    962        
    963        
     967               
    964968    def reverse(self):
    965969        """
    966970        Reverses the coefficients of self. The reverse of f(x) is x^n f(1/x).
    cdef class Polynomial_dense_modn_ntl_ZZ( 
    11111115            self.c.restore_c()
    11121116        ZZ_pX_destruct(&self.x)
    11131117
    1114     def ntl_set_directly(self, v):
    1115         # TODO: Get rid of this
    1116         Polynomial_dense_mod_n.ntl_set_directly(self, v)
    1117         # verbatim from __init__
    1118         cdef ntl_ZZ_pX ntl = self.__poly
    1119         self.__poly = None # this will eventually go away
    1120         self.x = ntl.x
    1121         self.c = ntl.c
     1118#     def ntl_set_directly(self, v):
     1119#         # TODO: Get rid of this
     1120#         Polynomial_dense_mod_n.ntl_set_directly(self, v)
     1121#         # verbatim from __init__
     1122#         cdef ntl_ZZ_pX ntl = self.__poly
     1123#         self.__poly = None # this will eventually go away
     1124#         self.x = ntl.x
     1125#         self.c = ntl.c
    11221126
    11231127    cdef Polynomial_dense_modn_ntl_ZZ _new(self):
    11241128        cdef Polynomial_dense_modn_ntl_ZZ y = <Polynomial_dense_modn_ntl_ZZ>PY_NEW(Polynomial_dense_modn_ntl_ZZ)
  • sage/rings/polynomial/polynomial_ring.py

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_ring.py
    a b TESTS: 
    6666#
    6767#                  http://www.gnu.org/licenses/
    6868#*****************************************************************************
     69
     70import sys # for maxint
    6971
    7072from sage.structure.element import is_Element
    7173import sage.algebras.algebra
    class PolynomialRing_general(sage.algebr 
    490492            sage: type(RR['x'].0)
    491493            <type 'sage.rings.polynomial.polynomial_real_mpfr_dense.PolynomialRealDense'>
    492494            sage: type(Integers(4)['x'].0)
    493             <type 'sage.rings.polynomial.polynomial_modn_dense_ntl.Polynomial_dense_modn_ntl_zz'>
     495            <type 'sage.rings.polynomial.polynomial_zmod_flint.Polynomial_zmod_flint'>
    494496            sage: type(Integers(5*2^100)['x'].0)
    495497            <type 'sage.rings.polynomial.polynomial_modn_dense_ntl.Polynomial_dense_modn_ntl_ZZ'>
    496498            sage: type(CC['x'].0)
    class PolynomialRing_dense_mod_n(Polynom 
    11421144        """
    11431145        if hasattr(x, '_polynomial_'):
    11441146            return x._polynomial_(self)
    1145         if self.__modulus < ZZ_sage(polynomial_modn_dense_ntl.zz_p_max):
     1147        elif self.__modulus <= sys.maxint:
     1148            import sage.rings.polynomial.polynomial_zmod_flint as polynomial_zmod_flint
     1149            return polynomial_zmod_flint.Polynomial_zmod_flint(self, x, check, is_gen,construct=construct)
     1150        elif self.__modulus < ZZ_sage(polynomial_modn_dense_ntl.zz_p_max):
    11461151            return polynomial_modn_dense_ntl.Polynomial_dense_modn_ntl_zz(self, x, check, is_gen, construct=construct)
    11471152        else:
    11481153            return polynomial_modn_dense_ntl.Polynomial_dense_modn_ntl_ZZ(self, x, check, is_gen, construct=construct)
    11491154
    1150 class PolynomialRing_dense_mod_p(PolynomialRing_dense_mod_n,
    1151                                  PolynomialRing_singular_repr,
    1152                                  principal_ideal_domain.PrincipalIdealDomain):
     1155class PolynomialRing_dense_mod_p(PolynomialRing_field,
     1156                                 PolynomialRing_dense_mod_n,
     1157                                 PolynomialRing_singular_repr):
    11531158    def __init__(self, base_ring, name="x"):
    11541159        from sage.rings.polynomial.polynomial_singular_interface import can_convert_to_singular
    11551160        self.__modulus = base_ring.order()
    class PolynomialRing_dense_mod_p(Polynom 
    11711176                raise TypeError,"Unable to coerce string"
    11721177        elif hasattr(x, '_polynomial_'):
    11731178            return x._polynomial_(self)
    1174         if self.modulus() == 2:
     1179        modulus = self.modulus()
     1180        if modulus == 2:
    11751181            import sage.rings.polynomial.polynomial_gf2x as polynomial_gf2x
    11761182            return polynomial_gf2x.Polynomial_GF2X(self, x, check, is_gen,construct=construct)
     1183        elif modulus <= sys.maxint:
     1184            import sage.rings.polynomial.polynomial_zmod_flint as polynomial_zmod_flint
     1185            return polynomial_zmod_flint.Polynomial_zmod_flint(self, x, check, is_gen,construct=construct)
    11771186        else:
    11781187            return polynomial_modn_dense_ntl.Polynomial_dense_mod_p(self, x, check, is_gen,construct=construct)
    11791188
  • sage/rings/polynomial/polynomial_template.pxi

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_template.pxi
    a b include "../../ext/stdsage.pxi" 
    3333
    3434from sage.rings.polynomial.polynomial_element cimport Polynomial
    3535from sage.structure.element cimport ModuleElement, Element, RingElement
     36from sage.rings.fraction_field_element import FractionFieldElement
    3637from sage.rings.integer cimport Integer
    3738from sage.libs.all import pari_gen
     39
     40from sage.interfaces.all import singular as singular_default
    3841
    3942def make_element(parent, args):
    4043    return parent(*args)
    cdef class Polynomial_template(Polynomia 
    5558            sage: P(map(GF(2),[1,0,1]))
    5659            x^2 + 1
    5760        """
    58         cdef celement gen, monomial, *coeff
     61        cdef celement *gen, *monomial
    5962        cdef Py_ssize_t deg
    6063        cdef cparent _parent
    6164
    6265        Polynomial.__init__(self, parent, is_gen=is_gen)
    6366
    64         celement_construct(&self.x, get_cparent(parent))
    65 
    6667        if is_gen:
     68            celement_construct(&self.x, get_cparent(parent))
    6769            celement_gen(&self.x, 0, get_cparent(parent))
    6870
    6971        elif PY_TYPE_CHECK(x, Polynomial_template):
    7072            try:
     73                celement_construct(&self.x, get_cparent(parent))
    7174                celement_set(&self.x, &(<Polynomial_template>x).x, get_cparent(parent))
    7275            except NotImplementedError:
    7376                raise TypeError("%s not understood."%x)
    7477
    7578        elif PY_TYPE_CHECK(x, int) or PY_TYPE_CHECK(x, Integer):
    7679            try:
     80                celement_construct(&self.x, get_cparent(parent))
    7781                celement_set_si(&self.x, int(x), get_cparent(parent))
    7882            except NotImplementedError:
    7983                raise TypeError("%s not understood."%x)
    cdef class Polynomial_template(Polynomia 
    8286            parent = (<Polynomial_template>self)._parent
    8387            _parent = get_cparent(parent)
    8488
     89            celement_construct(&self.x, get_cparent(parent))
     90            gen = celement_new(_parent)
     91            monomial = celement_new(_parent)
     92
    8593            celement_set_si(&self.x, 0, _parent)
    86             celement_gen(&gen, 0, _parent)
     94            celement_gen(gen, 0, _parent)
    8795
    8896            deg = 0
    8997            for e in x:
    9098                # r += parent(e)*power
    91                 celement_pow(&monomial, &gen, deg, NULL, _parent)
    92                 coeff = &(<Polynomial_template>parent(e)).x
    93                 celement_mul(&monomial, coeff, &monomial, _parent)
    94                 celement_add(&self.x, &self.x, &monomial, _parent)
     99                celement_pow(monomial, gen, deg, NULL, _parent)
     100                celement_mul(monomial, &(<Polynomial_template>self.__class__(parent, e)).x, monomial, _parent)
     101                celement_add(&self.x, &self.x, monomial, _parent)
    95102                deg += 1
     103
     104            celement_delete(gen, _parent)
     105            celement_delete(monomial, _parent)
     106
     107        elif PY_TYPE_CHECK(x, dict):
     108            parent = (<Polynomial_template>self)._parent
     109            _parent = get_cparent(parent)
     110
     111            celement_construct(&self.x, get_cparent(parent))
     112
     113            gen = celement_new(_parent)
     114            monomial = celement_new(_parent)
     115
     116            celement_set_si(&self.x, 0, _parent)
     117            celement_gen(gen, 0, _parent)
     118
     119            for deg, coef in x.iteritems():
     120                celement_pow(monomial, gen, deg, NULL, _parent)
     121                celement_mul(monomial, &(<Polynomial_template>parent(coef)).x, monomial, _parent)
     122                celement_add(&self.x, &self.x, monomial, _parent)
     123
     124            celement_delete(gen, _parent)
     125            celement_delete(monomial, _parent)
    96126
    97127        elif PY_TYPE_CHECK(x, pari_gen):
    98128            k = (<Polynomial_template>self)._parent.base_ring()
    99129            x = [k(w) for w in x.Vecrev()]
    100             Polynomial_template.__init__(self, parent, x, check=True, is_gen=False, construct=construct)
     130            self.__class__.__init__(self, parent, x, check=True, is_gen=False, construct=construct)
    101131        elif PY_TYPE_CHECK(x, Polynomial):
    102132            k = (<Polynomial_template>self)._parent.base_ring()
    103133            x = [k(w) for w in list(x)]
    104134            Polynomial_template.__init__(self, parent, x, check=True, is_gen=False, construct=construct)
     135        elif PY_TYPE_CHECK(x, FractionFieldElement) and (x.parent().base() is parent or x.parent().base() == parent) and x.denominator() == 1:
     136            x = x.numerator()
     137            self.__class__.__init__(self, parent, x, check=check, is_gen=is_gen, construct=construct)
    105138        else:
    106             raise TypeError("Coercion from parent %s not supported."%x.parent())
     139            x = parent.base_ring()(x)
     140            self.__class__.__init__(self, parent, x, check=check, is_gen=is_gen, construct=construct)           
    107141
    108142    def __reduce__(self):
    109143        """
    cdef class Polynomial_template(Polynomia 
    127161        cdef cparent _parent = get_cparent((<Polynomial_template>self)._parent)
    128162        return [self[i] for i in range(celement_len(&self.x, _parent))]
    129163
    130 #     def __cinit__(self):
    131 #         """
    132 #         EXAMPLE:
    133 #             sage: P.<x> = GF(2)[]
    134 #         """
    135 #         self.x = <celement>NULL
    136 
    137164    def __dealloc__(self):
    138165        """
    139166        EXAMPLE:
    cdef class Polynomial_template(Polynomia 
    150177            x + 1
    151178        """
    152179        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     180        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    153181        r._parent = (<Polynomial_template>self)._parent
    154182        celement_add(&r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>right).x, get_cparent((<Polynomial_template>self)._parent))
     183        #assert(r._parent(pari(self) + pari(right)) == r)
    155184        return r
    156185
    157186    cpdef ModuleElement _sub_(self, ModuleElement right):
    cdef class Polynomial_template(Polynomia 
    162191            x + 1
    163192        """
    164193        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     194        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    165195        r._parent = (<Polynomial_template>self)._parent
    166         celement_add(&r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>right).x, get_cparent((<Polynomial_template>self)._parent))
     196        celement_sub(&r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>right).x, get_cparent((<Polynomial_template>self)._parent))
     197        #assert(r._parent(pari(self) - pari(right)) == r)
    167198        return r
    168199
    169200    def __neg__(self):
    cdef class Polynomial_template(Polynomia 
    174205            x
    175206        """
    176207        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     208        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    177209        r._parent = (<Polynomial_template>self)._parent
    178210        celement_neg(&r.x, &self.x, get_cparent((<Polynomial_template>self)._parent))
     211        #assert(r._parent(-pari(self)) == r)
    179212        return r
    180213
    181214    cpdef RingElement _mul_(self, RingElement right):
    cdef class Polynomial_template(Polynomia 
    186219            x^2 + x
    187220        """
    188221        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     222        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    189223        r._parent = (<Polynomial_template>self)._parent
    190224        celement_mul(&r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>right).x, get_cparent((<Polynomial_template>self)._parent))
     225        #assert(r._parent(pari(self) * pari(right)) == r)
    191226        return r
    192227
    193228    def gcd(self, Polynomial_template other):
    194229        """
     230        Return the greatest common divisor of self and other.
     231
    195232        EXAMPLE:
    196233            sage: P.<x> = GF(2)[]
    197234            sage: f = x*(x+1)
    cdef class Polynomial_template(Polynomia 
    200237            sage: f.gcd(x^2)
    201238            x
    202239        """
     240        if(celement_is_zero(&self.x, get_cparent((<Polynomial_template>self)._parent))):
     241            return other
     242        if(celement_is_zero(&other.x, get_cparent((<Polynomial_template>self)._parent))):
     243            return self
     244
    203245        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     246        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    204247        r._parent = (<Polynomial_template>self)._parent
    205248        celement_gcd(&r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>other).x, get_cparent((<Polynomial_template>self)._parent))
     249        #assert(r._parent(pari(self).gcd(pari(other))) == r)
    206250        return r
     251
     252    def xgcd(self, Polynomial_template other):
     253        """
     254        Computes extended gcd of self and other.
     255
     256        EXAMPLE:
     257            sage: P.<x> = GF(7)[]
     258            sage: f = x*(x+1)
     259            sage: f.xgcd(x+1)
     260            (x + 1, 0, 1)
     261            sage: f.xgcd(x^2)
     262            (x, 1, 6)
     263        """
     264        if(celement_is_zero(&self.x, get_cparent((<Polynomial_template>self)._parent))):
     265            return other, self._parent(0), self._parent(1)
     266        if(celement_is_zero(&other.x, get_cparent((<Polynomial_template>self)._parent))):
     267            return self, self._parent(1), self._parent(0)
     268
     269        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     270        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
     271        r._parent = (<Polynomial_template>self)._parent
     272
     273        cdef Polynomial_template s = <Polynomial_template>PY_NEW(self.__class__)
     274        celement_construct(&s.x, get_cparent((<Polynomial_template>self)._parent))
     275        s._parent = (<Polynomial_template>self)._parent
     276
     277        cdef Polynomial_template t = <Polynomial_template>PY_NEW(self.__class__)
     278        celement_construct(&t.x, get_cparent((<Polynomial_template>self)._parent))
     279        t._parent = (<Polynomial_template>self)._parent
     280
     281        celement_xgcd(&r.x, &s.x, &t.x, &(<Polynomial_template>self).x, &(<Polynomial_template>other).x, get_cparent((<Polynomial_template>self)._parent))
     282        #rp, sp, tp = pari(self).xgcd(pari(other))
     283        #assert(r._parent(rp) == r)
     284        #assert(s._parent(sp) == s)
     285        #assert(t._parent(tp) == t)
     286        return r,s,t
    207287
    208288    def __floordiv__(self, right):
    209289        """
    cdef class Polynomial_template(Polynomia 
    214294            sage: (x + 1)//x
    215295            1
    216296        """
    217         right = (<Polynomial_template>self)._parent._coerce_(right)
     297        cdef Polynomial_template _right = <Polynomial_template>(<Polynomial_template>self)._parent._coerce_(right)
     298        if celement_is_zero(&_right.x, get_cparent((<Polynomial_template>self)._parent)):
     299            raise ZeroDivisionError
    218300        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     301        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    219302        r._parent = (<Polynomial_template>self)._parent
     303        #assert(r._parent(pari(self) // pari(right)) == r)
    220304        celement_floordiv(&r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>right).x, get_cparent((<Polynomial_template>self)._parent))
    221305        return r
    222306
    cdef class Polynomial_template(Polynomia 
    227311            sage: (x^2 + 1) % x^2
    228312            1
    229313        """
    230         other = (<Polynomial_template>self)._parent._coerce_(other)
     314        cdef Polynomial_template _other = <Polynomial_template>(<Polynomial_template>self)._parent._coerce_(other)
     315        if celement_is_zero(&_other.x, get_cparent((<Polynomial_template>self)._parent)):
     316            raise ZeroDivisionError
     317
    231318        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     319        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    232320        r._parent = (<Polynomial_template>self)._parent
    233         celement_mod(&r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>other).x, get_cparent((<Polynomial_template>self)._parent))
     321        celement_mod(&r.x, &(<Polynomial_template>self).x, &_other.x, get_cparent((<Polynomial_template>self)._parent))
     322        #assert(r._parent(pari(self) % pari(other)) == r)
    234323        return r
    235324
    236325    def quo_rem(self, right):
    cdef class Polynomial_template(Polynomia 
    241330            sage: f.quo_rem(x + 1)
    242331            (x, 1)
    243332        """
    244         right = (<Polynomial_template>self)._parent._coerce_(right)
     333        cdef Polynomial_template _right = <Polynomial_template>(<Polynomial_template>self)._parent._coerce_(right)
     334
     335        if celement_is_zero(&_right.x, get_cparent((<Polynomial_template>self)._parent)):
     336            raise ZeroDivisionError
     337
    245338        cdef Polynomial_template q = <Polynomial_template>PY_NEW(self.__class__)
     339        celement_construct(&q.x, get_cparent((<Polynomial_template>self)._parent))
    246340        q._parent = (<Polynomial_template>self)._parent
     341
    247342        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     343        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    248344        r._parent = (<Polynomial_template>self)._parent
    249         celement_quorem(&q.x, &r.x, &(<Polynomial_template>self).x, &(<Polynomial_template>right).x, get_cparent((<Polynomial_template>self)._parent))
     345
     346        celement_quorem(&q.x, &r.x, &(<Polynomial_template>self).x, &_right.x, get_cparent((<Polynomial_template>self)._parent))
    250347        return q,r
    251348       
    252349    def __long__(self):
    cdef class Polynomial_template(Polynomia 
    303400            sage: {x:1}
    304401            {x: 1}
    305402        """
    306         return celement_hash(&self.x, get_cparent((<Polynomial_template>self)._parent))
     403        cdef long result = 0 # store it in a c-int and just let the overflowing additions wrap
     404        cdef long result_mon
     405        cdef long c_hash
     406        cdef long var_name_hash
     407        cdef int i
     408        for i from 0<= i <= self.degree():
     409            if i == 1:
     410                # we delay the hashing until now to not waste it one a constant poly
     411                var_name_hash = hash(self.variable_name())
     412            #  I'm assuming (incorrectly) that hashes of zero indicate that the element is 0.
     413            # This assumption is not true, but I think it is true enough for the purposes and it
     414            # it allows us to write fast code that omits terms with 0 coefficients.  This is
     415            # important if we want to maintain the '==' relationship with sparse polys.
     416            c_hash = hash(self[i])
     417            if c_hash != 0:
     418                if i == 0:
     419                    result += c_hash
     420                else:
     421                    # Hash (self[i], generator, i) as a tuple according to the algorithm.
     422                    result_mon = c_hash
     423                    result_mon = (1000003 * result_mon) ^ var_name_hash
     424                    result_mon = (1000003 * result_mon) ^ i
     425                    result += result_mon
     426        if result == -1:
     427            return -2
     428        return result
    307429
    308     def __len__(self):
    309         """
    310         EXAMPLE:
    311             sage: P.<x> = GF(2)[]
    312             sage: P.<x> = GF(2)[]
    313             sage: len(x)
    314             2
    315             sage: len(x+1)
    316             2
    317         """
    318         return celement_len(&self.x, get_cparent((<Polynomial_template>self)._parent))
    319430
    320431    def __pow__(self, ee, modulus):
    321432        """
    cdef class Polynomial_template(Polynomia 
    338449        if not PY_TYPE_CHECK(self, Polynomial_template):
    339450            raise NotImplementedError("%s^%s not defined."%(ee,self))
    340451        cdef bint recip = 0, do_sig
    341         cdef long e = ee
     452
     453        cdef long e
     454        try:
     455            e = ee
     456        except OverflowError:
     457            return Polynomial.__pow__(self, ee, modulus)
    342458        if e != ee:
    343459            raise TypeError("Only integral powers defined.")
    344460        elif e < 0:
    cdef class Polynomial_template(Polynomia 
    348464            if e == 0:
    349465                raise ArithmeticError, "0^0 is undefined."
    350466        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     467        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
    351468        r._parent = (<Polynomial_template>self)._parent
    352469
    353470        if modulus is None:
    cdef class Polynomial_template(Polynomia 
    355472        else:
    356473            modulus = (<Polynomial_template>self)._parent._coerce_(modulus)
    357474            celement_pow(&r.x, &(<Polynomial_template>self).x, e, &(<Polynomial_template>modulus).x, get_cparent((<Polynomial_template>self)._parent))
     475       
     476        #assert(r._parent(pari(self)**ee) == r)
    358477        if recip:
    359478            return ~r
    360479        else:
    361480            return r
    362481       
     482    def __copy__(self):
     483        """
     484        EXAMPLE:
     485            sage: P.<x> = GF(2)[]
     486            sage: copy(x) is x
     487            False
     488            sage: copy(x) == x
     489            True
     490        """
     491        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     492        celement_construct(&r.x, get_cparent((<Polynomial_template>self)._parent))
     493        r._parent = (<Polynomial_template>self)._parent
     494        celement_set(&r.x, &self.x, get_cparent((<Polynomial_template>self)._parent))
     495        return r
     496
     497    def __getslice__(self, i, j):
     498        """
     499        Returns the monomials of self of degree i <= n < j.
     500       
     501        EXAMPLES:
     502            sage: R.<x> = Integers(100)[]
     503            sage: f = (x+2)^7
     504            sage: f[3:6]
     505            84*x^5 + 80*x^4 + 60*x^3
     506            sage: f[-5:50] == f
     507            True
     508            sage: f[6:]
     509            x^7 + 14*x^6       
     510        """
     511        cdef cparent _parent = get_cparent((<Polynomial_template>self)._parent)
     512        if i < 0:
     513            i = 0
     514        if j > celement_len(&self.x,_parent):
     515            j = celement_len(&self.x, _parent)
     516        x = (<Polynomial_template>self)._parent.gen()
     517        v = [self[t] for t from i <= t < j]
     518
     519        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     520        Polynomial_template.__init__(r, (<Polynomial_template>self)._parent, v)
     521        return r << i
     522
    363523    def is_gen(self):
    364524        """
    365525        EXAMPLE:
    cdef class Polynomial_template(Polynomia 
    369529            sage: (x+1).is_gen()
    370530            False
    371531        """
    372         cdef celement gen
    373         celement_gen(&gen, 0, get_cparent((<Polynomial_template>self)._parent))
    374         return celement_equal(&self.x, &gen, get_cparent((<Polynomial_template>self)._parent))
     532        cdef celement *gen = celement_new(get_cparent((<Polynomial_template>self)._parent))
     533        celement_gen(gen, 0, get_cparent((<Polynomial_template>self)._parent))
     534        cdef bint r = celement_equal(&self.x, gen, get_cparent((<Polynomial_template>self)._parent))
     535        celement_delete(gen, get_cparent((<Polynomial_template>self)._parent))
     536        return r
    375537
    376538    def shift(self, int n):
    377539        """
    cdef class Polynomial_template(Polynomia 
    383545            sage: f.shift(-1)
    384546            x^2 + x
    385547        """
    386         cdef celement gen
    387548        cdef Polynomial_template r
    388549        if n == 0:
    389550            return self
    390551
    391552        parent = (<Polynomial_template>self)._parent
    392553        cdef cparent _parent = get_cparent(parent)
    393         celement_gen(&gen, 0, _parent)
    394         celement_pow(&gen, &gen, abs(n), NULL, _parent)
     554        cdef celement *gen = celement_new(_parent)
     555        celement_gen(gen, 0, _parent)
     556        celement_pow(gen, gen, abs(n), NULL, _parent)
    395557        r = <Polynomial_template>PY_NEW(self.__class__)
     558        celement_construct(&r.x, _parent)
    396559        r._parent = parent
    397560       
    398561        if n > 0:
    399             celement_mul(&r.x, &self.x, &gen, _parent)
     562            celement_mul(&r.x, &self.x, gen, _parent)
    400563        else:
    401             celement_floordiv(&r.x, &self.x, &gen, _parent)
     564            celement_floordiv(&r.x, &self.x, gen, _parent)
     565        celement_delete(gen, _parent)
    402566        return r
    403567
    404568    def __lshift__(self, int n):
    cdef class Polynomial_template(Polynomia 
    411575        """
    412576        if not PY_TYPE_CHECK(self, Polynomial_template):
    413577            raise TypeError("Cannot %s << %n."%(self, n))
    414         cdef celement gen
     578        cdef celement *gen, *tmp
    415579        cdef Polynomial_template r
    416580        if n == 0:
    417581            return self
    cdef class Polynomial_template(Polynomia 
    420584
    421585        parent = (<Polynomial_template>self)._parent
    422586        cdef cparent _parent = get_cparent(parent)
    423         celement_gen(&gen, 0, _parent)
    424         celement_pow(&gen, &gen, n, NULL, _parent)
     587        gen = celement_new(_parent)
     588        celement_gen(gen, 0, _parent)
     589        celement_pow(gen, gen, n, NULL, _parent)
    425590        r = <Polynomial_template>PY_NEW(self.__class__)
     591        celement_construct(&r.x, _parent)
    426592        r._parent = parent
    427         celement_mul(&r.x, &(<Polynomial_template>self).x, &gen, _parent)
     593        celement_mul(&r.x, &(<Polynomial_template>self).x, gen, _parent)
     594        celement_delete(gen, _parent)
    428595        return r
    429596
    430597    def __rshift__(self, int n):
    431598        """
    432599        EXAMPLE:
    433600            sage: P.<x> = GF(2)[]
     601            sage: x>>1
     602            1
     603            sage: (x^2 + x)>>1
     604            x + 1
     605            sage: (x^2 + x) >> -1
     606            x^3 + x^2
    434607        """
    435608        if not PY_TYPE_CHECK(self, Polynomial_template):
    436609            raise TypeError("Cannot %s >> %n."%(self, n))
    437         cdef celement gen
     610        cdef celement *gen
    438611        cdef Polynomial_template r
    439612        if n == 0:
    440613            return self
    441614        elif n < 0:
    442             return self >> -n
     615            return self << -n
    443616
    444617        parent = (<Polynomial_template>self)._parent
    445618        cdef cparent _parent = get_cparent(parent)
    446         celement_gen(&gen, 0, _parent)
    447         celement_pow(&gen, &gen, n, NULL, _parent)
     619        gen = celement_new(_parent)
     620        celement_gen(gen, 0, _parent)
     621        celement_pow(gen, gen, n, NULL, _parent)
    448622        r = <Polynomial_template>PY_NEW(self.__class__)
     623        celement_construct(&r.x, _parent)
    449624        r._parent = parent
    450625       
    451         celement_floordiv(&r.x, &(<Polynomial_template>self).x, &gen, _parent)
     626        celement_floordiv(&r.x, &(<Polynomial_template>self).x, gen, _parent)
     627        celement_delete(gen, _parent)
    452628        return r
    453629
    454630    def is_zero(self):
    455631        """
    456632        EXAMPLE:
    457633            sage: P.<x> = GF(2)[]
     634            sage: x.is_zero()
     635            False
    458636        """
    459637        return celement_is_zero(&self.x, get_cparent((<Polynomial_template>self)._parent))
    460638
    cdef class Polynomial_template(Polynomia 
    462640        """
    463641        EXAMPLE:
    464642            sage: P.<x> = GF(2)[]
     643            sage: P(1).is_one()
     644            True
    465645        """
    466646        return celement_is_one(&self.x, get_cparent((<Polynomial_template>self)._parent))
    467647
    cdef class Polynomial_template(Polynomia 
    469649        """
    470650        EXAMPLE:
    471651            sage: P.<x> = GF(2)[]
     652            sage: x.degree()
     653            1
     654            sage: P(1).degree()
     655            0
     656            sage: P(0).degree()
     657            -1
    472658        """
    473659        return Integer(celement_len(&self.x, get_cparent((<Polynomial_template>self)._parent))-1)
    474660
    cdef class Polynomial_template(Polynomia 
    483669            sage: f.truncate(6)
    484670            x^5 + x^4 + x^3 + x^2 + x + 1
    485671        """
    486         if n < 0:
    487             raise ValueError(" n must be >= 0.")
    488672        parent = (<Polynomial_template>self)._parent
    489673        cdef cparent _parent = get_cparent(parent)
    490         cdef celement gen
    491         celement_gen(&gen, 0, _parent)
    492         celement_pow(&gen, &gen, n, NULL, _parent)
    493674
    494675        cdef Polynomial_template r = <Polynomial_template>PY_NEW(self.__class__)
     676        celement_construct(&r.x, _parent)
    495677        r._parent = parent
    496         celement_mod(&r.x, &self.x, &gen, _parent)
     678
     679        if n <= 0:
     680            return r
     681
     682        cdef celement *gen = celement_new(_parent)
     683        celement_gen(gen, 0, _parent)
     684        celement_pow(gen, gen, n, NULL, _parent)
     685
     686        celement_mod(&r.x, &self.x, gen, _parent)
     687        celement_delete(gen, _parent)
    497688        return r
     689
     690    def _singular_(self, singular=singular_default, have_ring=False, force=False):
     691        r"""
     692        Return \Singular representation of this polynomial
     693
     694        INPUT:
     695            singular -- \Singular interpreter (default: default interpreter)
     696            have_ring -- set to True if the ring was already set in \Singular
     697            force -- ignored.
     698
     699        EXAMPLE:
     700            sage: P.<x> = PolynomialRing(GF(7))
     701            sage: f = 3*x^2 + 2*x + 5
     702            sage: singular(f)
     703            3*x^2+2*x-2
     704        """
     705        if not have_ring:
     706            self.parent()._singular_(singular,force=force).set_ring() #this is expensive
     707        return singular(self._singular_init_())
     708
     709    def _derivative(self, var=None):
     710        r"""
     711        Returns the formal derivative of self with respect to var.
     712
     713        var must be either the generator of the polynomial ring to which
     714        this polynomial belongs, or None (either way the behaviour is the
     715        same).
     716
     717        SEE ALSO:
     718            self.derivative()
     719
     720        EXAMPLES:
     721            sage: R.<x> = Integers(77)[]
     722            sage: f = x^4 - x - 1
     723            sage: f._derivative()
     724            4*x^3 + 76
     725            sage: f._derivative(None)
     726            4*x^3 + 76
     727
     728            sage: f._derivative(2*x)
     729            Traceback (most recent call last):
     730            ...
     731            ValueError: cannot differentiate with respect to 2*x
     732
     733            sage: y = var("y")
     734            sage: f._derivative(y)
     735            Traceback (most recent call last):
     736            ...
     737            ValueError: cannot differentiate with respect to y
     738        """
     739        if var is not None and var is not self._parent.gen():
     740            raise ValueError, "cannot differentiate with respect to %s" % var
     741
     742        P = self.parent()
     743        x = P.gen()
     744        res = P(0)
     745        for i,c in enumerate(self.list()[1:]):
     746            res += (i+1)*c*x**i
     747        return res
  • new file sage/rings/polynomial/polynomial_zmod_flint.pxd

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_zmod_flint.pxd
    - +  
     1from sage.libs.flint.zmod_poly cimport zmod_poly_t, zmod_poly_struct
     2
     3ctypedef zmod_poly_struct celement
     4
     5include "polynomial_template_header.pxi"
     6
     7cdef class Polynomial_zmod_flint(Polynomial_template):
     8    cdef _set_list(self, x)
     9   
  • new file sage/rings/polynomial/polynomial_zmod_flint.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/polynomial/polynomial_zmod_flint.pyx
    - +  
     1"""
     2Univariate Polynomials over Z/nZ for n <= sys.maxint via FLINT.
     3
     4AUTHOR:
     5    -- Martin Albrecht (2009-01) another initial implementation
     6    -- Burcin Erocal (2008-11) initial implementation
     7"""
     8
     9from sage.libs.ntl.ntl_lzz_pX import ntl_zz_pX
     10
     11# We need to define this stuff before including the templating stuff
     12# to make sure the function get_cparent is found since it is used in
     13# 'polynomial_template.pxi'.
     14
     15ctypedef unsigned long cparent
     16
     17cdef inline cparent get_cparent(parent):
     18    if parent is None:
     19        return 0
     20    return <unsigned long>(parent.modulus())
     21
     22# first we include the definitions
     23include "../../libs/flint/zmod_poly_linkage.pxi"
     24
     25# and then the interface
     26include "polynomial_template.pxi"
     27
     28from sage.libs.all import pari
     29
     30cdef extern from "zn_poly/zn_poly.h":
     31    ctypedef struct zn_mod_struct:
     32        pass
     33    cdef void zn_mod_init(zn_mod_struct *mod, unsigned long m)
     34    cdef void zn_mod_clear(zn_mod_struct *mod)
     35    cdef void zn_array_mul(unsigned long* res, unsigned long* op1, size_t n1, unsigned long* op2, size_t n2, zn_mod_struct *mod)
     36
     37cdef class Polynomial_zmod_flint(Polynomial_template):
     38    def __init__(self, parent, x=None, check=True, is_gen=False, construct=False):
     39        """
     40        EXAMPLE:
     41            sage: P.<x> = GF(32003)[]
     42            sage: f = 24998*x^2 + 29761*x + 2252
     43        """
     44        cdef long nlen
     45
     46        if PY_TYPE_CHECK(x, list) or PY_TYPE_CHECK(x, tuple):
     47            k = parent._base
     48            if check:
     49                lst = [k(i) for i in x]
     50            else:
     51                lst = x
     52            # remove trailing zeroes
     53            nlen = len(lst)
     54            while nlen and lst[nlen-1] == 0:
     55                nlen -= 1
     56            lst = lst[:nlen]
     57            Polynomial_template.__init__(self, parent, 0, check, is_gen, construct)
     58            self._set_list(lst)
     59            return
     60        else:
     61            if PY_TYPE_CHECK(x, ntl_zz_pX):
     62                x = x.list()
     63            try:
     64                if x.parent() is parent.base_ring() or x.parent() == parent.base_ring():
     65                    x = int(x) % parent.modulus()
     66            except AttributeError:
     67                pass
     68        Polynomial_template.__init__(self, parent, x, check, is_gen, construct)
     69
     70    cdef _set_list(self, x):
     71        """
     72        EXAMPLES:
     73            sage: P.<a>=GF(7)[]
     74            sage: P([2^60,0,1])
     75            a^2 + 1
     76            sage: P([])
     77            0
     78            sage: P(range(15))
     79            6*a^13 + 5*a^12 + 4*a^11 + 3*a^10 + 2*a^9 + a^8 + 6*a^6 + 5*a^5 + 4*a^4 + 3*a^3 + 2*a^2 + a
     80        """
     81        cdef list l_in = x
     82        cdef unsigned long length = len(l_in)
     83        cdef unsigned long modulus = zmod_poly_modulus(&self.x)
     84        cdef int i
     85        if length == 0:
     86            zmod_poly_zero(&self.x)
     87            return
     88
     89        # resize to length of list
     90        _sig_on
     91        zmod_poly_realloc(&self.x, length)
     92        _sig_off
     93
     94        for i from 0 <= i < length:
     95            _zmod_poly_set_coeff_ui(&self.x, i, l_in[i])
     96        # we need to set the length manually, we used _zmod_poly_set_coeff_ui
     97        self.x.length = length
     98
     99    def __getitem__(self, i):
     100        """
     101        EXAMPLE:
     102            sage: P.<x> = GF(32003)[]
     103            sage: f = 24998*x^2 + 29761*x + 2252
     104            sage: f[100]
     105            0
     106            sage: f[1]
     107            29761
     108            sage: f[0]
     109            2252
     110            sage: f[-1]
     111            0
     112        """
     113        cdef unsigned long c = 0
     114        if 0 <= i < zmod_poly_length(&self.x):
     115            c = zmod_poly_get_coeff_ui(&self.x, i)
     116        return self._parent.base_ring()(c)
     117
     118    def __call__(self, *x, **kwds):
     119        """
     120        Evaluate polynomial at x=a.
     121
     122        INPUT:
     123            a -- ring element a; need not be in the coefficient
     124                 ring of the polynomial.
     125          -- or --
     126            a dictionary for kwds:value pairs.  If the variable
     127            name of the polynomial is a kwds it is substituted in;
     128            otherwise this polynomial is returned unchanged.
     129                 
     130        EXAMPLE:
     131            sage: P.<x> = PolynomialRing(GF(7))
     132            sage: f= x^2 + 1
     133            sage: f(0)
     134            1
     135            sage: f(2)
     136            5
     137            sage: f(3)
     138            3
     139        """
     140        K = self._parent.base_ring()
     141        if len(kwds) == 0 and len(x) == 1:
     142            try:
     143                x = K._coerce_(x[0])
     144                return K(zmod_poly_evaluate_horner(&self.x, int(x)))
     145            except TypeError:
     146                pass
     147        return Polynomial.__call__(self, *x, **kwds)
     148
     149    def resultant(self, Polynomial_zmod_flint other):
     150        """
     151        Returns the resultant of self and other, which must lie in the same
     152        polynomial ring.
     153
     154        INPUT:
     155            other -- a polynomial
     156
     157        OUTPUT:
     158            an element of the base ring of the polynomial ring
     159
     160        EXAMPLES:
     161            sage: R.<x> = GF(19)['x']
     162            sage: f = x^3 + x + 1;  g = x^3 - x - 1
     163            sage: r = f.resultant(g); r
     164            11
     165            sage: r.parent() is GF(19)
     166            True
     167        """
     168        other = self.parent()._coerce_(other)
     169        res = zmod_poly_resultant(&(<Polynomial_template>self).x, &(<Polynomial_template>other).x)
     170        return self.parent().base_ring()(res)
     171
     172    def small_roots(self, *args, **kwds):
     173        r"""
     174        See \code{sage.rings.polynomial.polynomial_modn_dense_ntl.small_roots}
     175        for the documentation of this function.
     176
     177        EXAMPLE:
     178            sage: N = 10001
     179            sage: K = Zmod(10001)
     180            sage: P.<x> = PolynomialRing(K)
     181            sage: f = x^3 + 10*x^2 + 5000*x - 222
     182            sage: f.small_roots()
     183            [4]
     184        """
     185        from sage.rings.polynomial.polynomial_modn_dense_ntl import small_roots
     186        return small_roots(self, *args, **kwds)
     187
     188    def _unsafe_mutate(self, n, value):
     189        r"""
     190        Never use this unless you really know what you are doing.
     191
     192        INPUT:
     193            n -- degree
     194            value -- coefficient
     195
     196        WARNING: This could easily introduce subtle bugs, since \Sage
     197        assumes everywhere that polynomials are immutable.  It's OK to
     198        use this if you really know what you're doing.
     199
     200        EXAMPLES:
     201            sage: R.<x> = GF(7)[]
     202            sage: f = (1+2*x)^2; f
     203            4*x^2 + 4*x + 1
     204            sage: f._unsafe_mutate(1, -5)
     205            sage: f
     206            4*x^2 + 2*x + 1
     207        """
     208        cdef cparent _parent = get_cparent((<Polynomial_template>self)._parent)
     209
     210        n = int(n)
     211        value = self.base_ring()(value)
     212        if n >= 0:
     213            zmod_poly_set_coeff_ui(&self.x, n, int(value)%zmod_poly_modulus(&self.x))
     214        else:
     215            raise IndexError, "Polynomial coefficient index must be nonnegative."
     216       
     217    def _mul_zn_poly(self, other):
     218        r"""           
     219        Returns the product of two polynomials using the zn_poly library.
     220
     221        See \url{http://www.math.harvard.edu/~dmharvey/zn_poly/} for
     222        details on zn_poly.
     223
     224        INPUT:
     225           self: Polynomial
     226           right: Polynomial (over same base ring as self)
     227
     228        OUTPUT: Polynomial
     229           The product self*right.
     230
     231
     232        EXAMPLE:
     233            sage: P.<x> = PolynomialRing(GF(next_prime(2^30)))
     234            sage: f = P.random_element(1000)
     235            sage: g = P.random_element(1000)
     236            sage: f*g == f._mul_zn_poly(g)
     237            True
     238
     239            sage: P.<x> = PolynomialRing(Integers(100))
     240            sage: P
     241            Univariate Polynomial Ring in x over Ring of integers modulo 100
     242            sage: r = (10*x)._mul_zn_poly(10*x); r
     243            0
     244            sage: r.degree()
     245            -1
     246       
     247        ALGORITHM:
     248           uses David Harvey's zn_poly library.
     249
     250        NOTE: This function is a technology preview. It might
     251        disappear or be replaced without a deprecation warning.
     252        """
     253        cdef Polynomial_zmod_flint _other = <Polynomial_zmod_flint>self._parent._coerce_(other)
     254       
     255        cdef Polynomial_zmod_flint r = <Polynomial_zmod_flint>PY_NEW(self.__class__)
     256        r._parent = (<Polynomial_zmod_flint>self)._parent
     257
     258        cdef unsigned long p = self._parent.modulus()
     259        cdef unsigned long n1 = self.x.length
     260        cdef unsigned long n2 = _other.x.length
     261
     262        cdef zn_mod_struct zn_mod
     263
     264        zmod_poly_init2(&r.x, p, n1 + n2 -1 )
     265
     266        zn_mod_init(&zn_mod, p)
     267        zn_array_mul(r.x.coeffs, self.x.coeffs, n1, _other.x.coeffs, n2, &zn_mod)
     268        r.x.length = n1 + n2 -1
     269        __zmod_poly_normalise(&r.x)
     270        zn_mod_clear(&zn_mod)
     271        return r
     272
  • sage/rings/residue_field.pyx

    diff -r 4ddd93e16f5e -r c253f3eb84d0 sage/rings/residue_field.pyx
    a b class ReductionMap: 
    482482        # The reduction map is just x |--> F(to_vs(x) * (PB**(-1))) if
    483483        # either x is integral or the denominator of x is coprime to
    484484        # p; otherwise we work harder.
    485 
    486485        x = self.__K(x)
    487486        p = self.__F.p
    488487