Ticket #7715: vector_mod2_dense.patch

File vector_mod2_dense.patch, 47.9 KB (added by malb, 12 years ago)
  • module_list.py

    # HG changeset patch
    # User Martin Albrecht <malb@informatik.uni-bremen.de>
    # Date 1260985109 0
    # Node ID b29fdccf32e16bfb81fd27d414f352ac6dae11b9
    # Parent  7e790bf432473ff2e87ca3b13ce0f5723ae3da0e
    implement dense vectors over GF(2) as m4ri's mzd_t with one row.
    
    diff -r 7e790bf43247 -r b29fdccf32e1 module_list.py
    a b  
    885885
    886886    Extension('sage.modules.vector_modn_dense',
    887887              sources = ['sage/modules/vector_modn_dense.pyx']),
     888
     889    Extension('sage.modules.vector_mod2_dense',
     890              sources = ['sage/modules/vector_mod2_dense.pyx'],
     891              libraries = ['gmp','m4ri', 'png12', 'gd'],
     892              depends = [SAGE_ROOT + "/local/include/png.h"]),
    888893   
    889894    Extension('sage.modules.vector_rational_dense',
    890895              sources = ['sage/modules/vector_rational_dense.pyx'],
  • new file sage/libs/m4ri.pxd

    diff -r 7e790bf43247 -r b29fdccf32e1 sage/libs/m4ri.pxd
    - +  
     1cdef extern from "m4ri/m4ri.h":
     2    ctypedef unsigned long long m4ri_word "word"
     3
     4    ctypedef struct mzd_t:
     5        int nrows
     6        int ncols
     7        int width
     8        m4ri_word **rows
     9
     10    ctypedef struct mzp_t:
     11        int *values
     12        int size
     13
     14    ctypedef int BIT
     15
     16    cdef int RADIX
     17
     18    ##############
     19    # Maintainance
     20    ##############
     21
     22    # builds all gray codes up to a certain size
     23    cdef void m4ri_build_all_codes()
     24    cdef void m4ri_destroy_all_codes()
     25
     26    ##############
     27    # Constructors
     28    ##############
     29 
     30    # create empty matrix
     31    cdef mzd_t *mzd_init(int,int)
     32
     33    # create the identity permutation
     34    cdef mzp_t *mzp_init(long)
     35   
     36    # free memory for the matrix
     37    cdef void mzd_free(mzd_t *)
     38
     39    # free memory for the permutation
     40    cdef void mzp_free(mzp_t *)
     41
     42    # filled uniformly random
     43    cdef void mzd_randomize(mzd_t *)
     44
     45    # identity matrix if i%2
     46    cdef void mzd_set_ui(mzd_t *, unsigned int i)
     47
     48    # [A],[B] -> [AB]
     49    cdef mzd_t *mzd_concat(mzd_t *C, mzd_t *A, mzd_t *B)
     50
     51    # [A],[B] -> | A |
     52    #            | B |
     53    cdef mzd_t *mzd_stack(mzd_t *C, mzd_t *A, mzd_t *B)
     54   
     55    # returns a submatrix from a
     56    cdef mzd_t *mzd_submatrix(mzd_t *S, mzd_t *A, int lowr, int lowc, int highr, int highc)
     57
     58    # return a matrix window to A
     59    cdef mzd_t *mzd_init_window( mzd_t *A, int lowr, int lowc, int highr, int highc)
     60
     61    cdef void mzd_free_window(mzd_t *)
     62
     63    # deep copy
     64    cdef mzd_t *mzd_copy(mzd_t *, mzd_t *)
     65
     66    ##############
     67    # Bit Level IO
     68    ##############
     69
     70    # set BIT
     71    cdef void mzd_write_bit( mzd_t *m, int row, int col, BIT value)
     72
     73    # get BIT
     74    cdef BIT mzd_read_bit( mzd_t *m, int row, int col )
     75
     76    # get BITs (n<=64)
     77    cdef m4ri_word mzd_read_bits( mzd_t *m, int row, int col, int n)
     78   
     79    #####################
     80    # Row/Column Based IO
     81    #####################
     82   
     83    cdef void mzd_row_swap(mzd_t *, int, int)
     84
     85    cdef void mzd_col_swap(mzd_t *, int, int)
     86
     87    cdef void mzd_row_clear_offset(mzd_t *m, int, int)
     88
     89    cdef void mzd_row_add_offset(mzd_t *m, int, int, int)
     90
     91    ############
     92    # Arithmetic
     93    ############
     94
     95    # matrix addition
     96    cdef mzd_t *mzd_add(mzd_t *, mzd_t *, mzd_t *)
     97
     98    # naive cubic matrix multiply
     99    cdef mzd_t *mzd_mul_naive(mzd_t *, mzd_t *, mzd_t *)
     100
     101    # naive cubic matrix multiply (b is pre-transposed)
     102    cdef mzd_t *_mzd_mul_naive(mzd_t *, mzd_t *, mzd_t *, int)
     103
     104    # matrix multiply using Gray codes
     105    cdef mzd_t *mzd_mul_m4rm(mzd_t *, mzd_t *, mzd_t *, int k)
     106
     107    # matrix multiply using Gray codes (transposed)
     108    cdef mzd_t *mzd_mul_m4rm_t(mzd_t *, mzd_t *, mzd_t *, int k)
     109
     110    # matrix multiply and addition using Gray codes: C = C + AB
     111    cdef mzd_t *mzd_addmul_m4rm(mzd_t *, mzd_t *, mzd_t *, int k)
     112
     113    # matrix multiplication via Strassen's formula
     114    cdef mzd_t *mzd_mul(mzd_t *, mzd_t *, mzd_t *, int cutoff)
     115
     116    # C = C + AB via Strassen's formula
     117    cdef mzd_t *mzd_addmul(mzd_t *, mzd_t *, mzd_t *, int cutoff)
     118
     119    # equality testing
     120    cdef int mzd_equal(mzd_t *, mzd_t *)
     121
     122    # returns -1,0,1
     123    cdef int mzd_cmp(mzd_t *, mzd_t *)
     124
     125    # transpose
     126    cdef mzd_t *mzd_transpose(mzd_t *, mzd_t *)
     127
     128    # density with given resolution
     129    cdef double mzd_density(mzd_t *, int resolution)
     130
     131    ########################
     132    # LAPACK Level Functions
     133    ########################
     134
     135    # cubic Gaussian elimination
     136    cdef int mzd_echelonize_naive(mzd_t *, int full)
     137
     138    # row echelon form using Gray codes
     139    cdef int mzd_echelonize_m4ri(mzd_t *m, int full, int k)
     140
     141    # reduced row echelon form from upper triangular form
     142    cdef void mzd_top_echelonize_m4ri(mzd_t *m, int k)
     143
     144    # matrix inversion using Gray codes
     145    cdef mzd_t *mzd_invert_m4ri(mzd_t *m, mzd_t *identity, int k)
     146
     147    # asymptotically fast PLUQ factorization
     148    cdef long mzd_pluq(mzd_t *A, mzp_t *P, mzp_t *Q, int cutoff)
     149
     150    # PLUQ factorization using Gray codes
     151    cdef long _mzd_pluq_mmpf(mzd_t *A, mzp_t *P, mzp_t *Q, int k)
     152
     153    # cubic PLUQ factorization
     154    cdef long _mzd_pluq_naive(mzd_t *A, mzp_t *P, mzp_t *Q)
     155
     156    # asymptotically fast LQUP factorization
     157    cdef long mzd_lqup(mzd_t *A, mzp_t *P, mzp_t *Q, int cutoff)
     158
     159    # LQUP factorization using Gray codes
     160    cdef long _mzd_lqup_mmpf(mzd_t *A, mzp_t *P, mzp_t *Q, int k)
     161
     162    # cubic LQUP factorization
     163    cdef long _mzd_lqup_naive(mzd_t *A, mzp_t *P, mzp_t *Q)
     164
     165    # reduced row echelon form using PLUQ factorization
     166    cdef long mzd_echelonize_pluq(mzd_t *A, int full)
     167
     168    # reduced row echelon form using PLUQ factorization
     169    cdef mzd_t *mzd_kernel_left_pluq(mzd_t *A, int cutoff)
     170
     171    ##################################
     172    # Internal Functions for Debugging
     173    ##################################
     174
     175    cdef void mzd_write_zeroed_bits(mzd_t *m, int x, int y, int n, m4ri_word values)
     176     
     177    cdef void mzd_clear_bits(mzd_t *m, int x, int y, int n)
     178
  • sage/matrix/matrix_mod2_dense.pxd

    diff -r 7e790bf43247 -r b29fdccf32e1 sage/matrix/matrix_mod2_dense.pxd
    a b  
    11# choose: dense or sparse
    22
    3 cdef extern from "m4ri/m4ri.h":
    4     ctypedef unsigned long long word "word"
     3cimport matrix_dense
    54
    6     ctypedef struct mzd_t:
    7         int nrows
    8         int ncols
    9         int width
    10         word **rows
    11 
    12     ctypedef struct mzp_t:
    13         int *values
    14         int size
    15 
    16     ctypedef int BIT
    17 
    18     cdef int RADIX
    19 
    20     ##############
    21     # Maintainance
    22     ##############
    23 
    24     # builds all gray codes up to a certain size
    25     cdef void m4ri_build_all_codes()
    26     cdef void m4ri_destroy_all_codes()
    27 
    28     ##############
    29     # Constructors
    30     ##############
    31  
    32     # create empty matrix
    33     cdef mzd_t *mzd_init(int,int)
    34 
    35     # create the identity permutation
    36     cdef mzp_t *mzp_init(long)
    37    
    38     # free memory for the matrix
    39     cdef void mzd_free(mzd_t *)
    40 
    41     # free memory for the permutation
    42     cdef void mzp_free(mzp_t *)
    43 
    44     # filled uniformly random
    45     cdef void mzd_randomize(mzd_t *)
    46 
    47     # identity matrix if i%2
    48     cdef void mzd_set_ui(mzd_t *, unsigned int i)
    49 
    50     # [A],[B] -> [AB]
    51     cdef mzd_t *mzd_concat(mzd_t *C, mzd_t *A, mzd_t *B)
    52 
    53     # [A],[B] -> | A |
    54     #            | B |
    55     cdef mzd_t *mzd_stack(mzd_t *C, mzd_t *A, mzd_t *B)
    56    
    57     # returns a submatrix from a
    58     cdef mzd_t *mzd_submatrix(mzd_t *S, mzd_t *A, int lowr, int lowc, int highr, int highc)
    59 
    60     # return a matrix window to A
    61     cdef mzd_t *mzd_init_window( mzd_t *A, int lowr, int lowc, int highr, int highc)
    62 
    63     cdef void mzd_free_window(mzd_t *)
    64 
    65     # deep copy
    66     cdef mzd_t *mzd_copy(mzd_t *, mzd_t *)
    67 
    68     ##############
    69     # Bit Level IO
    70     ##############
    71 
    72     # set BIT
    73     cdef void mzd_write_bit( mzd_t *m, int row, int col, BIT value)
    74 
    75     # get BIT
    76     cdef BIT mzd_read_bit( mzd_t *m, int row, int col )
    77 
    78     # get BITs (n<=64)
    79     cdef word mzd_read_bits( mzd_t *m, int row, int col, int n)
    80    
    81     #####################
    82     # Row/Column Based IO
    83     #####################
    84    
    85     cdef void mzd_row_swap(mzd_t *, int, int)
    86 
    87     cdef void mzd_col_swap(mzd_t *, int, int)
    88 
    89     cdef void mzd_row_clear_offset(mzd_t *m, int, int)
    90 
    91     cdef void mzd_row_add_offset(mzd_t *m, int, int, int)
    92 
    93     ############
    94     # Arithmetic
    95     ############
    96 
    97     # matrix addition
    98     cdef mzd_t *mzd_add(mzd_t *, mzd_t *, mzd_t *)
    99 
    100     # naive cubic matrix multiply
    101     cdef mzd_t *mzd_mul_naive(mzd_t *, mzd_t *, mzd_t *)
    102 
    103     # matrix multiply using Gray codes
    104     cdef mzd_t *mzd_mul_m4rm(mzd_t *, mzd_t *, mzd_t *, int k)
    105 
    106     # matrix multiply using Gray codes (transposed)
    107     cdef mzd_t *mzd_mul_m4rm_t(mzd_t *, mzd_t *, mzd_t *, int k)
    108 
    109     # matrix multiply and addition using Gray codes: C = C + AB
    110     cdef mzd_t *mzd_addmul_m4rm(mzd_t *, mzd_t *, mzd_t *, int k)
    111 
    112     # matrix multiplication via Strassen's formula
    113     cdef mzd_t *mzd_mul(mzd_t *, mzd_t *, mzd_t *, int cutoff)
    114 
    115     # C = C + AB via Strassen's formula
    116     cdef mzd_t *mzd_addmul(mzd_t *, mzd_t *, mzd_t *, int cutoff)
    117 
    118     # equality testing
    119     cdef int mzd_equal(mzd_t *, mzd_t *)
    120 
    121     # returns -1,0,1
    122     cdef int mzd_cmp(mzd_t *, mzd_t *)
    123 
    124     # transpose
    125     cdef mzd_t *mzd_transpose(mzd_t *, mzd_t *)
    126 
    127     # density with given resolution
    128     cdef double mzd_density(mzd_t *, int resolution)
    129 
    130     ########################
    131     # LAPACK Level Functions
    132     ########################
    133 
    134     # cubic Gaussian elimination
    135     cdef int mzd_echelonize_naive(mzd_t *, int full)
    136 
    137     # row echelon form using Gray codes
    138     cdef int mzd_echelonize_m4ri(mzd_t *m, int full, int k)
    139 
    140     # reduced row echelon form from upper triangular form
    141     cdef void mzd_top_echelonize_m4ri(mzd_t *m, int k)
    142 
    143     # matrix inversion using Gray codes
    144     cdef mzd_t *mzd_invert_m4ri(mzd_t *m, mzd_t *identity, int k)
    145 
    146     # asymptotically fast PLUQ factorization
    147     cdef long mzd_pluq(mzd_t *A, mzp_t *P, mzp_t *Q, int cutoff)
    148 
    149     # PLUQ factorization using Gray codes
    150     cdef long _mzd_pluq_mmpf(mzd_t *A, mzp_t *P, mzp_t *Q, int k)
    151 
    152     # cubic PLUQ factorization
    153     cdef long _mzd_pluq_naive(mzd_t *A, mzp_t *P, mzp_t *Q)
    154 
    155     # asymptotically fast LQUP factorization
    156     cdef long mzd_lqup(mzd_t *A, mzp_t *P, mzp_t *Q, int cutoff)
    157 
    158     # LQUP factorization using Gray codes
    159     cdef long _mzd_lqup_mmpf(mzd_t *A, mzp_t *P, mzp_t *Q, int k)
    160 
    161     # cubic LQUP factorization
    162     cdef long _mzd_lqup_naive(mzd_t *A, mzp_t *P, mzp_t *Q)
    163 
    164     # reduced row echelon form using PLUQ factorization
    165     cdef long mzd_echelonize_pluq(mzd_t *A, int full)
    166 
    167     # reduced row echelon form using PLUQ factorization
    168     cdef mzd_t *mzd_kernel_left_pluq(mzd_t *A, int cutoff)
    169 
    170     ##################################
    171     # Internal Functions for Debugging
    172     ##################################
    173 
    174     cdef void mzd_write_zeroed_bits(mzd_t *m, int x, int y, int n, word values)
    175      
    176     cdef void mzd_clear_bits(mzd_t *m, int x, int y, int n)
    177 
    178 
    179 cimport matrix_dense
     5from sage.libs.m4ri cimport *
    1806
    1817cdef class Matrix_mod2_dense(matrix_dense.Matrix_dense):
    1828    cdef mzd_t *_entries
  • sage/matrix/matrix_mod2_dense.pyx

    diff -r 7e790bf43247 -r b29fdccf32e1 sage/matrix/matrix_mod2_dense.pyx
    a b  
    33
    44AUTHOR: Martin Albrecht <malb@informatik.uni-bremen.de>
    55
    6 EXAMPLES:
     6EXAMPLES::
     7
    78    sage: a = matrix(GF(2),3,range(9),sparse=False); a
    89    [0 1 0]
    910    [1 0 1]
     
    8889
    8990##############################################################################
    9091#       Copyright (C) 2004,2005,2006 William Stein <wstein@gmail.com>
    91 #       Copyright (C) 2007,2008 Martin Albrecht <M.R.Albrecht@rhul.ac.uk>
     92#       Copyright (C) 2007,2008,2009 Martin Albrecht <M.R.Albrecht@rhul.ac.uk>
    9293#  Distributed under the terms of the GNU General Public License (GPL)
    9394#  The full text of the GPL is available at:
    9495#                  http://www.gnu.org/licenses/
     
    100101include '../ext/random.pxi'
    101102
    102103cimport matrix_dense
    103 from sage.structure.element cimport Matrix
     104from sage.structure.element cimport Matrix, Vector
    104105from sage.structure.element cimport ModuleElement, Element
    105106
    106107from sage.misc.functional import log
     
    108109from sage.misc.misc import verbose, get_verbose, cputime
    109110
    110111from sage.modules.free_module import VectorSpace
     112from sage.modules.vector_mod2_dense cimport Vector_mod2_dense
    111113
    112114cdef extern from "gd.h":
    113115    ctypedef struct gdImagePtr "gdImagePtr":
     
    158160    ########################################################################
    159161    def __cinit__(self, parent, entries, copy, coerce, alloc=True):
    160162        """
    161         Creates a new dense matrix over GF(2).
     163        Dense matrix over GF(2) constructor.
    162164
    163165        INPUT:
    164             parent -- MatrixSpace (always)
    165             entries -- ignored
    166             copy -- ignored
    167             coerce -- ignored
    168             alloc -- if True a zero matrix is allocated (default:True)
     166
     167        - ``parent`` - MatrixSpace.
     168        - ``entries`` - may be list or 0 or 1
     169        - ``copy`` - ignored, elements are always copied
     170        - ``coerce`` - ignored, elements are always coerced to ints % 2
     171
     172        EXAMPLES::
     173
     174            sage: type(random_matrix(GF(2),2,2))
     175            <type 'sage.matrix.matrix_mod2_dense.Matrix_mod2_dense'>
     176
     177            sage: Matrix(GF(2),3,3,1) # indirect doctest
     178            [1 0 0]
     179            [0 1 0]
     180            [0 0 1]
    169181        """
    170182        matrix_dense.Matrix_dense.__init__(self, parent)
    171183
     
    186198        Dense matrix over GF(2) constructor.
    187199
    188200        INPUT:
    189             parent -- MatrixSpace.
    190             entries -- may be list or 0 or 1
    191             copy -- ignored, elements are always copied
    192             coerce -- ignored, elements are always coerced to ints % 2
    193201
    194         EXAMPLES:
     202        - ``parent`` - MatrixSpace.
     203        - ``entries`` - may be list or 0 or 1
     204        - ``copy`` - ignored, elements are always copied
     205        - ``coerce`` - ignored, elements are always coerced to ints % 2
     206
     207        EXAMPLES::
     208
    195209            sage: type(random_matrix(GF(2),2,2))
    196210            <type 'sage.matrix.matrix_mod2_dense.Matrix_mod2_dense'>
    197211
     
    213227            sage: Matrix(GF(2),1,1, [1/3])
    214228            [1]
    215229           
    216         TESTS:
     230        TESTS::
     231
    217232            sage: Matrix(GF(2),0,0)
    218233            []
    219234            sage: Matrix(GF(2),2,0)
     
    262277            right -- a matrix
    263278            op -- comparison operation
    264279
    265         EXAMPLE:
     280        EXAMPLE::
     281
    266282            sage: A = random_matrix(GF(2),2,2)
    267283            sage: B = random_matrix(GF(2),3,3)
    268284            sage: A < B
     
    273289            sage: A < B
    274290            True
    275291
    276         TESTS:
     292        TESTS::
     293
    277294            sage: A = matrix(GF(2),2,0)
    278295            sage: B = matrix(GF(2),2,0)
    279296            sage: A < B
     
    286303        The has of a matrix is computed as $\oplus i*a_i$ where the $a_i$ are
    287304        the flattened entries in a matrix (by row, then by column).
    288305       
    289         EXAMPLE:
     306        EXAMPLE::
     307
    290308            sage: B = random_matrix(GF(2),3,3)
    291309            sage: B.set_immutable()
    292310            sage: {B:0} # indirect doctest
     
    327345
    328346        cdef unsigned long i, j, truerow
    329347        cdef unsigned long start, shift
    330         cdef word row_xor
    331         cdef word end_mask = ~(((<word>1)<<(RADIX - self._ncols%RADIX))-1)
    332         cdef word top_mask, bot_mask
    333         cdef word cur
    334         cdef word* row
     348        cdef m4ri_word row_xor
     349        cdef m4ri_word end_mask = ~(((<m4ri_word>1)<<(RADIX - self._ncols%RADIX))-1)
     350        cdef m4ri_word top_mask, bot_mask
     351        cdef m4ri_word cur
     352        cdef m4ri_word* row
    335353       
    336354        # running_xor is the xor of all words in the matrix, as if the rows
    337355        # in the matrix were written out consecutively, without regard to
    338356        # word boundaries.
    339         cdef word running_xor = 0
     357        cdef m4ri_word running_xor = 0
    340358        # running_parity is the number of extra words that must be xor'd.
    341359        cdef unsigned long running_parity = 0
    342360       
     
    348366            row = self._entries.rows[i]
    349367            start = (i*self._entries.ncols) >> 6
    350368            shift = (i*self._entries.ncols) & 0x3F
    351             bot_mask = (~(<word>0)) << shift
     369            bot_mask = (~(<m4ri_word>0)) << shift
    352370            top_mask = ~bot_mask
    353371
    354372            if self._entries.width > 1:
     
    418436
    419437    def str(self):
    420438        """
    421         EXAMPLE:
     439        EXAMPLE::
     440
    422441            sage: B = random_matrix(GF(2),3,3)       
    423442            sage: B # indirect doctest
    424443            [0 1 0]
     
    467486                s.insert(i, row_divider)
    468487       
    469488        return "\n".join(s)
    470    
     489
     490    def row(self, Py_ssize_t i, from_list=False):
     491        """
     492        Return the ``i``'th row of this matrix as a vector.
     493
     494        This row is a dense vector if and only if the matrix is a dense
     495        matrix.
     496
     497        INPUT:
     498
     499        - ``i`` - integer
     500
     501        - ``from_list`` - bool (default: ``False``); if ``True``,
     502          returns the ``i``'th element of ``self.rows()`` (see
     503          :func:`rows`), which may be faster, but requires building a
     504          list of all rows the first time it is called after an entry
     505          of the matrix is changed.
     506
     507        EXAMPLES::
     508
     509            sage: A = random_matrix(GF(2),10,10); A
     510            [0 1 0 1 1 0 0 0 1 1]
     511            [0 1 1 1 0 1 1 0 0 1]
     512            [0 0 0 1 0 1 0 0 1 0]
     513            [0 1 1 0 0 1 0 1 1 0]
     514            [0 0 0 1 1 1 1 0 1 1]
     515            [0 0 1 1 1 1 0 0 0 0]
     516            [1 1 1 1 0 1 0 1 1 0]
     517            [0 0 0 1 1 0 0 0 1 1]
     518            [1 0 0 0 1 1 1 0 1 1]
     519            [1 0 0 1 1 0 1 0 0 0]
     520            sage: A.row(0)
     521            (0, 1, 0, 1, 1, 0, 0, 0, 1, 1)
     522           
     523            sage: A.row(-1)
     524            (1, 0, 0, 1, 1, 0, 1, 0, 0, 0)
     525           
     526            sage: A.row(2,from_list=True)
     527            (0, 0, 0, 1, 0, 1, 0, 0, 1, 0)
     528        """
     529        if self._nrows == 0:
     530            raise IndexError("matrix has no rows")
     531        if i >= self._nrows or i < -self._nrows:
     532            raise IndexError("row index out of range")
     533        if i < 0:
     534            i = i + self._nrows
     535        if from_list:
     536            return self.rows(copy=False)[i]
     537        cdef Py_ssize_t j
     538        cdef Vector_mod2_dense z = PY_NEW(Vector_mod2_dense)
     539        z._init(self._ncols, VectorSpace(self.base_ring(),self._ncols))
     540        mzd_submatrix(z._entries, self._entries, i, 0, i+1, self._ncols)
     541        return z
    471542   
    472543    ########################################################################
    473544    # LEVEL 2 functionality
     
    488559        INPUT:
    489560            right -- matrix of dimension self.nrows() x self.ncols()
    490561
    491         EXAMPLES:
     562        EXAMPLES::
     563
    492564            sage: A = random_matrix(GF(2),10,10)
    493565            sage: A + A == Matrix(GF(2),10,10,0)
    494566            True
    495567
    496568            sage: A = random_matrix(GF(2),257,253)
    497             sage: A + A == Matrix(GF(2),257,253,0)
     569            sage: A + A == Matrix(GF(2),257,253,0) # indirect doctest
    498570            True
    499571
    500         TESTS:
     572        TESTS::
     573
    501574            sage: A = matrix(GF(2),2,0)
    502575            sage: A+A
    503576            []
     
    523596        INPUT:
    524597            right -- matrix of dimension self.nrows() x self.ncols()
    525598
    526         EXAMPLES:
     599        EXAMPLES::
     600
    527601            sage: A = random_matrix(GF(2),10,10)
    528             sage: A - A == Matrix(GF(2),10,10,0)
     602            sage: A - A == Matrix(GF(2),10,10,0)  # indirect doctest
    529603            True
    530604        """
    531605        return self._add_(right)
     606
     607    cdef Vector _matrix_times_vector_(self, Vector v):
     608        """
     609        EXAMPLES::
     610       
     611            sage: A = random_matrix(GF(2),10^4,10^4)
     612            sage: v0 = random_matrix(GF(2),10^4,1)
     613            sage: v1 = v0.column(0)
     614            sage: r0 = A*v0
     615            sage: r1 = A*v1
     616            sage: r0.column(0) == r1
     617            True
     618        """
     619        cdef mzd_t *tmp
     620        if not PY_TYPE_CHECK(v, Vector_mod2_dense):
     621            M = VectorSpace(self._base_ring, self._nrows)
     622            v = M(v)
     623        if self.ncols() != v.degree():
     624            raise ArithmeticError("number of columns of matrix must equal degree of vector")
     625
     626        VS = VectorSpace(self._base_ring, self._nrows)
     627        cdef Vector_mod2_dense c = PY_NEW(Vector_mod2_dense)
     628        c._init(self._nrows, VS)
     629        c._entries = mzd_init(1, self._nrows)
     630        tmp = mzd_init(self._nrows, 1)
     631        _mzd_mul_naive(tmp, self._entries, (<Vector_mod2_dense>v)._entries, 0)
     632        mzd_transpose(c._entries, tmp)
     633        mzd_free(tmp)
     634        return c
    532635       
    533636    cdef Matrix _matrix_times_matrix_(self, Matrix right):
    534637        """
    535638        Matrix multiplication.
    536639
    537640        ALGORITHM: Uses the 'Method of the Four Russians
    538         Multiplication', see self._multiply_m4rm.
     641        Multiplication', see :func:`_multiply_m4rm`.
    539642        """
    540643        if get_verbose() >= 2:
    541644            verbose('matrix multiply of %s x %s matrix by %s x %s matrix'%(
     
    559662                 suitable value is chosen by the function.
    560663                 ($0<= k <= 16$, default: 0)
    561664
    562         EXAMPLE:
     665        EXAMPLE::
     666
    563667              sage: A = Matrix(GF(2), 4, 3, [0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1] )
    564668              sage: B = Matrix(GF(2), 3, 4, [0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0] )
    565669              sage: A
     
    577681              [0 1 0 1]
    578682              [1 1 0 0]
    579683
    580         TESTS:
     684        TESTS::
     685
    581686            sage: A = random_matrix(GF(2),0,0)
    582687            sage: B = random_matrix(GF(2),0,0)
    583688            sage: A._multiply_m4rm(B, 0)
     
    639744        the other routines fall back to this implementation in that
    640745        case anyway.
    641746
    642         EXAMPLE:
     747        EXAMPLE::
     748
    643749              sage: A = Matrix(GF(2), 4, 3, [0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1] )
    644750              sage: B = Matrix(GF(2), 3, 4, [0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0] )
    645751              sage: A
     
    700806            cutoff -- matrix dimension where M4RM should be used
    701807                      instead of Strassen (default: let M4RI decide)
    702808
    703         EXAMPLE:
     809        EXAMPLE::
     810
    704811              sage: A = Matrix(GF(2), 4, 3, [0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1] )
    705812              sage: B = Matrix(GF(2), 3, 4, [0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0] )
    706813              sage: A
     
    768875
    769876    def __neg__(self):
    770877        """
    771         EXAMPLES:
     878        EXAMPLES::
     879
    772880            sage: A = random_matrix(GF(2),100,100)
    773881            sage: A - A == A - -A
    774882            True
     
    783891        If \code{self} is not invertible a \code{ZeroDivisionError} is
    784892        raised.
    785893
    786         EXAMPLE:
     894        EXAMPLE::
     895
    787896            sage: A = Matrix(GF(2),3,3, [0, 0, 1, 0, 1, 1, 1, 0, 1])
    788897            sage: MS = A.parent()
    789898            sage: A
     
    830939        r"""
    831940        Returns a copy of \code{self}.
    832941
    833         EXAMPLES:
     942        EXAMPLES::
     943
    834944             sage: MS = MatrixSpace(GF(2),3,3)
    835945             sage: A = MS(1)
    836946             sage: A.__copy__() == A
     
    865975        Returns list of the elements of \code{self} in row major
    866976        ordering.
    867977
    868         EXAMPLE:
     978        EXAMPLE::
     979
    869980            sage: A = Matrix(GF(2),2,2,[1,0,1,1])
    870981            sage: A
    871982            [1 0]
     
    9101021                  3/4 * log_2( min(nrows, ncols) ) as suggested in the M4RI paper.
    9111022            reduced -- return reduced row echelon form (default:True)
    9121023
    913         EXAMPLE:
     1024        EXAMPLE::
     1025
    9141026             sage: A = random_matrix(GF(2), 10, 10)
    9151027             sage: B = A.__copy__(); B.echelonize() # fastest
    9161028             sage: C = A.__copy__(); C.echelonize(k=2) # force k
     
    9181030             sage: B == C == E
    9191031             True
    9201032
    921         TESTS:
     1033        TESTS::
     1034
    9221035             sage: VF2 = VectorSpace(GF(2),2)
    9231036             sage: WF2 = VF2.submodule([VF2([1,1])])
    9241037             sage: WF2
     
    10071120        Returns the pivot columns of \code{self} if \code{self} is in
    10081121        row echelon form.
    10091122
    1010         EXAMPLE:
     1123        EXAMPLE::
     1124
    10111125            sage: A = matrix(GF(2),5,5,[0,1,0,1,0,0,1,0,1,1,0,1,0,1,0,0,0,0,1,0,0,0,1,0,1])
    10121126            sage: E = A.echelon_form()
    10131127            sage: E
     
    10401154        Randomize density proportion of the entries of this matrix,
    10411155        leaving the rest unchanged.
    10421156
    1043         EXAMPLES:
     1157        EXAMPLES::
     1158
    10441159            sage: A = matrix(GF(2), 5, 5, 0)
    10451160            sage: A.randomize(0.5); A
    10461161            [0 0 0 1 1]
     
    10841199        cdef int i, j, k
    10851200        cdef int nc
    10861201        cdef unsigned int low, high
    1087         cdef word mask = 1
     1202        cdef m4ri_word mask = 1
    10881203       
    10891204        if density == 1:
    1090             assert(sizeof(word) == 8)
     1205            assert(sizeof(m4ri_word) == 8)
    10911206            mask = ~((mask<<(RADIX - self._entries.ncols%RADIX))-1)
    10921207            for i from 0 <= i < self._nrows:
    10931208                for j from 0 <= j < self._entries.width:
     
    11091224
    11101225    cdef rescale_row_c(self, Py_ssize_t row, multiple, Py_ssize_t start_col):
    11111226        """
    1112         EXAMPLE:
     1227        EXAMPLE::
     1228
    11131229            sage: A = random_matrix(GF(2),3,3); A
    11141230            [0 1 0]
    11151231            [0 1 1]
     
    11251241    cdef add_multiple_of_row_c(self,  Py_ssize_t row_to, Py_ssize_t row_from, multiple,
    11261242                               Py_ssize_t start_col):
    11271243        """
    1128         EXAMPLE:
     1244        EXAMPLE::
     1245
    11291246            sage: A = random_matrix(GF(2),3,3); A
    11301247            [0 1 0]
    11311248            [0 1 1]
     
    11401257   
    11411258    cdef swap_rows_c(self, Py_ssize_t row1, Py_ssize_t row2):
    11421259        """
    1143         EXAMPLE:
     1260        EXAMPLE::
     1261
    11441262            sage: A = random_matrix(GF(2),3,3)
    11451263            sage: A
    11461264            [0 1 0]
     
    11551273
    11561274    cdef swap_columns_c(self, Py_ssize_t col1, Py_ssize_t col2):
    11571275        """
    1158         EXAMPLE:
     1276        EXAMPLE::
     1277
    11591278            sage: A = random_matrix(GF(2),3,3)
    11601279            sage: A
    11611280            [0 1 0]
     
    11901309        Returns a string of self in \Magma form. Does not return \Magma
    11911310        object but string.
    11921311
    1193         EXAMPLE:
     1312        EXAMPLE::
     1313
    11941314            sage: A = random_matrix(GF(2),3,3)
    11951315            sage: A._magma_init_(magma)                             # optional - magma
    11961316            'Matrix(GF(2),3,3,StringToIntegerSequence("0 1 0 0 1 1 0 0 0"))'
     
    11991319            sage: magma(A*B) == magma(A) * magma(B)                 # optional - magma
    12001320            True
    12011321
    1202         TESTS:
     1322        TESTS::
     1323
    12031324            sage: A = random_matrix(GF(2),0,3)
    12041325            sage: magma(A)                          # optional - magma
    12051326            Matrix with 0 rows and 3 columns
     
    12181339        """
    12191340        Return the determinant of this matrix over GF(2).
    12201341
    1221         EXAMPLES:
     1342        EXAMPLES::
     1343
    12221344            sage: matrix(GF(2),2,[1,1,0,1]).determinant()
    12231345            1
    12241346            sage: matrix(GF(2),2,[1,1,1,1]).determinant()
     
    12321354        """
    12331355        Returns transpose of self and leaves self untouched.
    12341356
    1235         EXAMPLE:
     1357        EXAMPLE::
     1358
    12361359            sage: A = Matrix(GF(2),3,5,[1, 0, 1, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 1, 0])
    12371360            sage: A
    12381361            [1 0 1 0 0]
     
    12711394        r"""
    12721395        Augments \code{self} with \code{right}.
    12731396
    1274         EXAMPLE:
     1397        EXAMPLE::
     1398
    12751399            sage: MS = MatrixSpace(GF(2),3,3)
    12761400            sage: A = MS([0, 1, 0, 1, 1, 0, 1, 1, 1]); A
    12771401            [0 1 0]
     
    13331457        r"""
    13341458        Stack \code{self} on top of \code{other}.
    13351459
    1336         EXAMPLE:
     1460        EXAMPLE::
     1461
    13371462            sage: A = matrix(GF(2),2,2,[1,0,0,1])
    13381463            sage: B = matrix(GF(2),2,2,[0,1,1,0])
    13391464            sage: A.stack(B)
     
    13471472            [1 0]
    13481473            [0 1]
    13491474
    1350         TESTS:
     1475        TESTS::
     1476
    13511477            sage: A = random_matrix(GF(2),0,3)
    13521478            sage: B = random_matrix(GF(2),3,3)
    13531479            sage: A.stack(B)
     
    13941520            nrows -- number of rows of submatrix
    13951521            ncols -- number of columns of submatrix
    13961522       
    1397         EXAMPLES:
     1523        EXAMPLES::
     1524
    13981525             sage: A = random_matrix(GF(2),200,200)
    13991526             sage: A[0:2,0:2] == A.submatrix(0,0,2,2)
    14001527             True
     
    14421569        r"""
    14431570        Serialize \code{self}.
    14441571
    1445         EXAMPLE:
     1572        EXAMPLE::
     1573
    14461574            sage: A = random_matrix(GF(2),10,10)
    14471575            sage: f,s = A.__reduce__()
    14481576            sage: f(*s) == A
     
    15141642        INPUT:
    15151643            approx -- return floating point approximation (default: False)
    15161644
    1517         EXAMPLE:
     1645        EXAMPLE::
     1646
    15181647            sage: A = random_matrix(GF(2),1000,1000)
    15191648            sage: d = A.density(); d
    15201649            62483/125000
     
    15441673
    15451674        - ``algorithm`` - either "lqup" or "m4ri"
    15461675
    1547         EXAMPLE:
     1676        EXAMPLE::
     1677
    15481678            sage: A = random_matrix(GF(2), 1000, 1000)
    15491679            sage: A.rank()
    15501680            999
     
    16791809                           
    16801810# gmp's ULONG_PARITY may use special
    16811811# assembly instructions, could be faster
    1682 cpdef inline unsigned long parity(word a):
     1812cpdef inline unsigned long parity(m4ri_word a):
    16831813    """
    16841814    Returns the parity of the number of bits in a.
    16851815   
    1686     EXAMPLES:
     1816    EXAMPLES::
     1817
    16871818        sage: from sage.matrix.matrix_mod2_dense import parity
    16881819        sage: parity(1)
    16891820        1L
     
    16921823        sage: parity(0x10000101011)
    16931824        1L
    16941825    """
    1695     if sizeof(word) == 8:
     1826    if sizeof(m4ri_word) == 8:
    16961827        a ^= a >> 32
    16971828    a ^= a >> 16
    16981829    a ^= a >> 8
    16991830    return parity_table[a & 0xFF]
    17001831
    1701 cdef inline unsigned long parity_mask(word a):
     1832cdef inline unsigned long parity_mask(m4ri_word a):
    17021833    return -parity(a)
    17031834
    17041835
     
    17121843        s -- a string
    17131844        size -- length of the string s
    17141845
    1715     EXAMPLE:
     1846    EXAMPLE::
     1847
    17161848        sage: A = random_matrix(GF(2),100,101)
    17171849        sage: _,(r,c,s,s2) = A.__reduce__()
    17181850        sage: from sage.matrix.matrix_mod2_dense import unpickle_matrix_mod2_dense_v1
     
    17601892    INPUT:
    17611893        filename -- a string
    17621894
    1763     EXAMPLE:
     1895    EXAMPLE::
     1896
    17641897        sage: from sage.matrix.matrix_mod2_dense import from_png, to_png
    17651898        sage: A = random_matrix(GF(2),10,10)
    17661899        sage: fn = tmp_filename()
     
    18021935        A -- a matrix over GF(2)
    18031936        filename -- a string for a file in a writable position
    18041937
    1805     EXAMPLE:
     1938    EXAMPLE::
     1939
    18061940        sage: from sage.matrix.matrix_mod2_dense import from_png, to_png
    18071941        sage: A = random_matrix(GF(2),10,10)
    18081942        sage: fn = tmp_filename()
     
    18431977        param -- either k for 'mmpf' is chosen or matrix multiplication
    18441978                 cutoff for 'standard' (default: 0)
    18451979
    1846     EXAMPLE:
     1980    EXAMPLE::
     1981
    18471982        sage: from sage.matrix.matrix_mod2_dense import pluq
    18481983        sage: A = random_matrix(GF(2),4,4); A
    18491984        [0 1 0 1]
     
    19012036        param -- either k for 'mmpf' is chosen or matrix multiplication
    19022037                 cutoff for 'standard' (default: 0)
    19032038
    1904     EXAMPLE:
     2039    EXAMPLE::
     2040
    19052041        sage: from sage.matrix.matrix_mod2_dense import lqup
    19062042        sage: A = random_matrix(GF(2),4,4); A
    19072043        [0 1 0 1]
  • sage/modules/free_module.py

    diff -r 7e790bf43247 -r b29fdccf32e1 sage/modules/free_module.py
    a b  
    60666066        sage: sage.modules.free_module.element_class(FF, is_sparse=True)
    60676067        <type 'sage.modules.free_module_element.FreeModuleElement_generic_sparse'>
    60686068        sage: sage.modules.free_module.element_class(FF, is_sparse=False)
     6069        <type 'sage.modules.vector_mod2_dense.Vector_mod2_dense'>
     6070        sage: sage.modules.free_module.element_class(GF(7), is_sparse=False)
    60696071        <type 'sage.modules.vector_modn_dense.Vector_modn_dense'>
    60706072        sage: sage.modules.free_module.element_class(P, is_sparse=True)
    60716073        <type 'sage.modules.free_module_element.FreeModuleElement_generic_sparse'>
     
    60836085        return Vector_rational_dense
    60846086    elif sage.rings.integer_mod_ring.is_IntegerModRing(R) and not is_sparse:
    60856087        from vector_modn_dense import Vector_modn_dense, MAX_MODULUS
     6088        from vector_mod2_dense import Vector_mod2_dense
     6089        if R.order() == 2:
     6090            return Vector_mod2_dense
    60866091        if R.order() < MAX_MODULUS:
    60876092            return Vector_modn_dense
    60886093        else:
  • new file sage/modules/vector_mod2_dense.pxd

    diff -r 7e790bf43247 -r b29fdccf32e1 sage/modules/vector_mod2_dense.pxd
    - +  
     1cimport free_module_element
     2import  free_module_element
     3
     4from sage.libs.m4ri cimport mzd_t
     5
     6cdef class Vector_mod2_dense(free_module_element.FreeModuleElement):
     7    cdef mzd_t *_entries
     8    cdef object _base_ring
     9
     10    cdef _new_c(self)
     11    cdef _init(self, Py_ssize_t degree, parent)
  • new file sage/modules/vector_mod2_dense.pyx

    diff -r 7e790bf43247 -r b29fdccf32e1 sage/modules/vector_mod2_dense.pyx
    - +  
     1"""
     2Vectors with elements in GF(2).
     3
     4AUTHOR:
     5
     6- Martin Albrecht (2009-12): initial implementation
     7
     8EXAMPLES::
     9
     10    sage: VS = GF(2)^3
     11    sage: e = VS.random_element(); e
     12    (1, 0, 0)
     13    sage: f = VS.random_element(); f
     14    (0, 1, 1)
     15    sage: e + f
     16    (1, 1, 1)
     17"""
     18
     19##############################################################################
     20#       Copyright (C) 2009 Martin Albrecht <M.R.Albrecht@rhul.ac.uk>
     21#  Distributed under the terms of the GNU General Public License (GPL)
     22#  The full text of the GPL is available at:
     23#                  http://www.gnu.org/licenses/
     24##############################################################################
     25
     26include '../ext/interrupt.pxi'
     27include '../ext/stdsage.pxi'
     28
     29from sage.rings.integer_mod cimport IntegerMod_int, IntegerMod_abstract
     30from sage.rings.integer cimport Integer
     31
     32from sage.structure.element cimport Element, ModuleElement, RingElement, Vector
     33
     34cimport free_module_element
     35from free_module_element import vector
     36
     37from sage.libs.m4ri cimport *
     38
     39cdef class Vector_mod2_dense(free_module_element.FreeModuleElement):
     40    cdef _new_c(self):
     41        """
     42        EXAMPLE::
     43
     44            sage: VS = VectorSpace(GF(2),3)
     45            sage: VS([0,0,1])
     46            (0, 0, 1)
     47            sage: type(_)
     48            <type 'sage.modules.vector_mod2_dense.Vector_mod2_dense'>
     49        """
     50        cdef Vector_mod2_dense y
     51        y = PY_NEW(Vector_mod2_dense)
     52        y._init(self._degree, self._parent)
     53        return y
     54
     55    cdef bint is_dense_c(self):
     56        """
     57        EXAMPLE::
     58
     59            sage: VS = VectorSpace(GF(2),3)
     60            sage: VS([0,0,1]).is_dense()
     61            True
     62        """
     63        return 1
     64   
     65    cdef bint is_sparse_c(self):
     66        """
     67        EXAMPLE::
     68
     69            sage: VS = VectorSpace(GF(2),3)
     70            sage: VS([0,0,1]).is_sparse()
     71            False
     72        """
     73        return 0
     74
     75    def __copy__(self):
     76        """
     77        EXAMPLE::
     78
     79            sage: VS = VectorSpace(GF(2),10^4)
     80            sage: v = VS.random_element()
     81            sage: w = copy(v)
     82            sage: w == v
     83            True
     84            sage: v[:10]
     85            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     86            sage: w[:10]
     87            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     88        """
     89        cdef Vector_mod2_dense y = self._new_c()
     90        mzd_copy(y._entries, self._entries)
     91        return y
     92
     93    cdef _init(self, Py_ssize_t degree, parent):
     94        """
     95        EXAMPLE::
     96
     97            sage: VS = VectorSpace(GF(2),3)
     98            sage: VS([0,0,1])
     99            (0, 0, 1)
     100            sage: type(_)
     101            <type 'sage.modules.vector_mod2_dense.Vector_mod2_dense'>
     102        """
     103        self._degree = degree
     104        self._parent = parent
     105        self._base_ring = parent.base_ring()
     106        self._entries = mzd_init(1, degree)
     107        if self._entries == NULL:
     108            raise MemoryError("Allocation of Vector_mod2_dense failed.")
     109       
     110    def __cinit__(self, parent=None, x=None, coerce=True, copy=True):
     111        """
     112        EXAMPLE::
     113
     114            sage: VS = VectorSpace(GF(2),3)
     115            sage: VS((0,0,1/3))
     116            (0, 0, 1)
     117            sage: type(_)
     118            <type 'sage.modules.vector_mod2_dense.Vector_mod2_dense'>
     119        """
     120        self._entries = NULL
     121        self._is_mutable = 1
     122        if not parent is None:
     123            self._init(parent.degree(), parent)
     124
     125    def __init__(self, parent, x, coerce=True, copy=True):
     126        """
     127        EXAMPLE::
     128
     129            sage: VS = VectorSpace(GF(2),3)
     130            sage: VS((0,0,1/3))
     131            (0, 0, 1)
     132            sage: type(_)
     133            <type 'sage.modules.vector_mod2_dense.Vector_mod2_dense'>
     134            sage: VS((0,0,int(3)))
     135            (0, 0, 1)
     136            sage: VS((0,0,3))
     137            (0, 0, 1)
     138            sage: VS((0,0,GF(2)(1)))
     139            (0, 0, 1)
     140        """
     141        cdef Py_ssize_t i
     142        cdef int xi
     143        if isinstance(x, (list, tuple)):
     144            if len(x) != self._degree:
     145                raise TypeError("x must be a list of the right length")
     146            for i from 0 <= i < self._degree:
     147                if PY_TYPE_CHECK(x[i],IntegerMod_int) or PY_TYPE_CHECK(x[i],int) or PY_TYPE_CHECK(x[i],Integer):
     148                    xi = x[i]
     149                    mzd_write_bit(self._entries, 0, i, xi%2)
     150                else:
     151                    mzd_write_bit(self._entries, 0, i, x[i]%2)
     152            return
     153        if x != 0:
     154            raise TypeError("can't initialize vector from nonzero non-list")
     155        else:
     156            for i from 0 <= i < self._degree:
     157                mzd_set_ui(self._entries, 0)
     158               
     159    def __dealloc__(self):
     160        """
     161        EXAMPLE::
     162
     163        sage: VS = VectorSpace(GF(2),10^3)
     164        sage: import gc
     165        sage: for i in range(10):
     166        ...      v = VS.random_element()
     167        ...      del v
     168        ...      _ = gc.collect()
     169        """
     170        if self._entries:
     171            mzd_free(self._entries)
     172
     173    cdef int _cmp_c_impl(left, Element right) except -2:
     174        """
     175        EXAMPLES::
     176            sage: v = vector(GF(2), [0,0,0,0])
     177            sage: v == 0
     178            True
     179            sage: v == 1
     180            False
     181            sage: v == v
     182            True
     183            sage: w = vector(GF(2), [1,0,0,0])
     184            sage: w < v
     185            False
     186            sage: w > v
     187            True
     188        """
     189        return mzd_cmp(left._entries, (<Vector_mod2_dense>right)._entries)
     190
     191    def __len__(self):
     192        """
     193        EXAMPLES::
     194
     195            sage: len(vector(GF(2),[0,0,1,1,1]))
     196            5
     197        """
     198        return self._degree
     199
     200    def __setitem__(self, Py_ssize_t i, x):
     201        """
     202        EXAMPLES::
     203
     204            sage: VS = VectorSpace(GF(2),4)
     205            sage: v = VS.random_element(); v
     206            (1, 0, 0, 0)
     207            sage: v[0] = 0; v
     208            (0, 0, 0, 0)
     209            sage: v[4] = 0
     210            Traceback (most recent call last):
     211            ...
     212            IndexError: Index '4' out of bound.
     213        """
     214        if not self._is_mutable:
     215            raise ValueError("vector is immutable; please change a copy instead (use copy())")
     216        cdef IntegerMod_int n
     217        n = self.base_ring()(x)
     218        if i < 0 or i >= self._degree:
     219            raise IndexError("Index '%s' out of bound."%(i))
     220        else:
     221            mzd_write_bit(self._entries, 0, i, n)
     222           
     223    def __getitem__(self, Py_ssize_t i):
     224        """
     225        Return the ith entry of self.
     226
     227        EXAMPLES::
     228
     229            sage: v = vector(GF(2), [1,2,3]); v
     230            (1, 0, 1)
     231            sage: v[0]
     232            1
     233            sage: v[2]
     234            1
     235            sage: v[-2]
     236            0
     237            sage: v[5]
     238            Traceback (most recent call last):
     239            ...
     240            IndexError: index '5' out of range
     241
     242            sage: v[-5]
     243            Traceback (most recent call last):
     244            ...
     245            IndexError: index '-2' out of range
     246        """
     247        if i < 0:
     248            i += self._degree
     249           
     250        if i < 0 or i >= self._degree:
     251            raise IndexError("index '%s' out of range"%(i,))
     252
     253        return self._base_ring(mzd_read_bit(self._entries, 0, i))
     254
     255    def __reduce__(self):
     256        """
     257        EXAMPLE::
     258
     259            sage: VS = VectorSpace(GF(2),10^4)
     260            sage: e = VS.random_element()
     261            sage: loads(dumps(e)) == e
     262            True
     263        """
     264        return unpickle_v0, (self._parent, self.list(), self._degree, self._is_mutable)
     265
     266    cpdef ModuleElement _add_(self, ModuleElement right):
     267        """
     268        EXAMPLE::
     269
     270            sage: VS = VectorSpace(GF(2),10)
     271            sage: e = VS([0,0,1,1,0,0,1,1,0,0])
     272            sage: f = VS([0,1,0,1,0,1,0,1,0,1])
     273            sage: e + f #indirect doctest
     274            (0, 1, 1, 0, 0, 1, 1, 0, 0, 1)
     275        """
     276        cdef Vector_mod2_dense z = self._new_c()
     277        mzd_add(z._entries, self._entries, (<Vector_mod2_dense>right)._entries)
     278        return z
     279
     280    cpdef ModuleElement _sub_(self, ModuleElement right):
     281        """
     282        EXAMPLE::
     283
     284            sage: VS = VectorSpace(GF(2),10)
     285            sage: e = VS([0,0,1,1,0,0,1,1,0,0])
     286            sage: f = VS([0,1,0,1,0,1,0,1,0,1])
     287            sage: e - f #indirect doctest
     288            (0, 1, 1, 0, 0, 1, 1, 0, 0, 1)
     289        """
     290        cdef Vector_mod2_dense z = self._new_c()
     291        mzd_add(z._entries, self._entries, (<Vector_mod2_dense>right)._entries)
     292        return z
     293
     294    cpdef Element _dot_product_(self, Vector right):
     295        """
     296        EXAMPLES::
     297           sage: VS = VectorSpace(GF(2),3)
     298           sage: v = VS([1,1,1]); w = VS([0,0,0])
     299           sage: v * w, w * v #indirect doctest
     300           (0, 0)
     301           sage: v = VS([1,1,1]); w = VS([0,1,0])
     302           sage: v * w, w * v
     303           (1, 1)
     304           sage: v = VS([1,1,1]); w = VS([0,1,1])
     305           sage: v * w, w * v
     306           (0, 0)
     307           sage: v = VS([1,1,1]); w = VS([1,1,1])
     308           sage: v * w, w * v
     309           (1, 1)
     310
     311           sage: VS = VectorSpace(GF(2),10^4)
     312           sage: v = VS(0); w = VS(0)
     313           sage: v[1337] = 1; w[1337] = 1
     314           sage: v * w, w * v
     315           (1, 1)
     316           sage: v[9881] = 1; w[9881] = 1
     317           sage: v * w, w * v
     318           (0, 0)
     319           sage: v[5172] = 1; w[6178] = 1
     320           sage: v * w, w * v
     321           (0, 0)
     322        """
     323        cdef Py_ssize_t i
     324        cdef IntegerMod_int n
     325        cdef Vector_mod2_dense r = right
     326        cdef m4ri_word tmp = 0
     327        n =  IntegerMod_int.__new__(IntegerMod_int)
     328        IntegerMod_abstract.__init__(n, self.base_ring())
     329        n.ivalue = 0
     330       
     331        for i from 0 <= i < self._entries.width:
     332            tmp ^= self._entries.rows[0][i] & r._entries.rows[0][i]
     333
     334        for i in range(64):
     335            n.ivalue ^= <int>(tmp & 1)
     336            tmp = tmp >> 1
     337       
     338        return n
     339
     340    cpdef Vector _pairwise_product_(self, Vector right):
     341        """
     342        EXAMPLE::
     343
     344            sage: VS = VectorSpace(GF(2),10)
     345            sage: e = VS.random_element(); e
     346            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     347            sage: f = VS.random_element(); f
     348            (1, 1, 0, 1, 1, 1, 0, 0, 0, 1)
     349            sage: e.pairwise_product(f) #indirect doctest
     350            (1, 0, 0, 0, 1, 1, 0, 0, 0, 1)
     351        """
     352        cdef Vector_mod2_dense z, r
     353        r = right
     354        z = self._new_c()
     355        cdef Py_ssize_t i
     356        for i from 0 <= i < self._entries.width:
     357            z._entries.rows[0][i] = (self._entries.rows[0][i] & r._entries.rows[0][i])
     358        return z
     359       
     360    cpdef ModuleElement _rmul_(self, RingElement left):
     361        """
     362        EXAMPLE::
     363
     364            sage: VS = VectorSpace(GF(2),10)
     365            sage: e = VS.random_element(); e
     366            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     367            sage: 0 * e
     368            (0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
     369            sage: 1 * e
     370            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     371            sage: 2 * e #indirect doctest
     372            (0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
     373        """
     374        cdef IntegerMod_int a
     375       
     376        if left:
     377            return self.__copy__()
     378        else:
     379            return self._new_c()
     380
     381
     382    cpdef ModuleElement _lmul_(self, RingElement right):
     383        """
     384        EXAMPLE::
     385
     386            sage: VS = VectorSpace(GF(2),10)
     387            sage: e = VS.random_element(); e
     388            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     389            sage: e * 0 #indirect doctest
     390            (0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
     391            sage: e * 1
     392            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     393            sage: e * 2
     394            (0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
     395        """
     396        return self._rmul_(right)
     397
     398    cpdef ModuleElement _neg_(self):
     399        """
     400        EXAMPLE::
     401
     402            sage: VS = VectorSpace(GF(2),10)
     403            sage: e = VS.random_element()
     404            sage: -e == e
     405            True
     406        """
     407        return self.__copy__()
     408
     409    def n(self, *args, **kwargs):
     410        """
     411        Returns a numerical approximation of ``self`` by calling the
     412        :meth:`n()` method on all of its entries.
     413
     414        EXAMPLES::
     415
     416            sage: v = vector(GF(2), [1,2,3])
     417            sage: v.n()
     418            (1.00000000000000, 0.000000000000000, 1.00000000000000)
     419            sage: _.parent()
     420            Vector space of dimension 3 over Real Field with 53 bits of precision
     421            sage: v.n(prec=75)
     422            (1.000000000000000000000, 0.0000000000000000000000, 1.000000000000000000000)
     423            sage: _.parent()
     424            Vector space of dimension 3 over Real Field with 75 bits of precision
     425        """
     426        return vector( [e.n(*args, **kwargs) for e in self] )
     427
     428    def list(self, copy=True):
     429        """
     430        Return a list of entries in ``self``.
     431
     432        INPUT:
     433
     434        - ``copy`` - always ``True
     435
     436        EXAMPLE::
     437
     438            sage: VS = VectorSpace(GF(2),10)
     439            sage: e = VS.random_element(); e
     440            (1, 0, 0, 0, 1, 1, 1, 0, 0, 1)
     441            sage: e.list()
     442            [1, 0, 0, 0, 1, 1, 1, 0, 0, 1]
     443        """
     444        cdef Py_ssize_t d = self._degree
     445        cdef Py_ssize_t i
     446        cdef list v = [0]*d
     447        K = self.base_ring()
     448        z = K.zero_element()
     449        o = K.one_element()
     450        cdef list switch = [z,o]
     451        for i in range(d):
     452            v[i] = switch[mzd_read_bit(self._entries, 0, i)]
     453        return v
     454
     455def unpickle_v0(parent, entries, degree, is_mutable):
     456    """
     457    EXAMPLE::
     458
     459        sage: from sage.modules.vector_mod2_dense import unpickle_v0
     460        sage: VS = VectorSpace(GF(2),10)
     461        sage: unpickle_v0(VS, [0,1,2,3,4,5,6,7,8,9], 10, 0)
     462        (0, 1, 0, 1, 0, 1, 0, 1, 0, 1)
     463    """
     464    # If you think you want to change this function, don't.
     465    cdef Vector_mod2_dense v
     466    v = PY_NEW(Vector_mod2_dense)
     467    v._init(degree, parent)
     468    cdef int xi
     469
     470    for i from 0 <= i < degree:
     471        if PY_TYPE_CHECK(entries[i],IntegerMod_int) or PY_TYPE_CHECK(entries[i],int) or PY_TYPE_CHECK(entries[i],Integer):
     472            xi = entries[i]
     473            mzd_write_bit(v._entries, 0, i, xi%2)
     474        else:
     475            mzd_write_bit(v._entries, 0, i, entries[i]%2)
     476    v._is_mutable = int(is_mutable)
     477    return v
     478
  • sage/rings/polynomial/polynomial_gf2x.pyx

    diff -r 7e790bf43247 -r b29fdccf32e1 sage/rings/polynomial/polynomial_gf2x.pyx
    a b  
    2222
    2323from sage.libs.all import pari
    2424
    25 from sage.matrix.matrix_mod2_dense cimport mzd_write_bit, mzd_read_bit, Matrix_mod2_dense, word
     25from sage.libs.m4ri cimport mzd_write_bit, mzd_read_bit
     26from sage.matrix.matrix_mod2_dense cimport Matrix_mod2_dense
    2627
    2728cdef class Polynomial_GF2X(Polynomial_template):
    2829    """