Ticket #11574: m4ri_20110601.patch

File m4ri_20110601.patch, 10.0 KB (added by malb, 6 years ago)

rebased to 4.7.1.alpha4

  • module_list.py

    # HG changeset patch
    # User Martin Albrecht <martinralbrecht@googlemail.com>
    # Date 1306771873 -3600
    # Node ID 632260cc66e03df6c20be829b28bee185cae9c66
    # Parent  a48af224839e4eb2555175d15d703ecd81d0b4c1
    #11574: updated interface for M4RI version 20110601 and up.
    
    diff --git a/module_list.py b/module_list.py
    a b  
    886886    # TODO -- change to use BLAS at some point.
    887887    Extension('sage.matrix.matrix_integer_dense',
    888888              sources = ['sage/matrix/matrix_integer_dense.pyx'],
     889              extra_compile_args = ['-std=c99'],
    889890              # order matters for cygwin!!
    890891              libraries = ['iml', 'pari', 'm', 'gmp', BLAS, BLAS2]),
    891892
     
    896897    Extension('sage.matrix.matrix_mod2_dense',
    897898              sources = ['sage/matrix/matrix_mod2_dense.pyx'],
    898899              libraries = ['gmp','m4ri', 'gd', 'png12', 'z'],
     900              extra_compile_args = ['-std=c99'],
    899901              depends = [SAGE_INC + "png.h", SAGE_INC + "m4ri/m4ri.h"]),
    900902
    901903    Extension('sage.matrix.matrix_modn_dense',
     
    10971099    Extension('sage.modules.vector_mod2_dense',
    10981100              sources = ['sage/modules/vector_mod2_dense.pyx'],
    10991101              libraries = ['gmp','m4ri', 'png12', 'gd'],
     1102              extra_compile_args = ['-std=c99'],
    11001103              depends = [SAGE_INC + "png.h", SAGE_INC + "m4ri/m4ri.h"]),
    11011104   
    11021105    Extension('sage.modules.vector_rational_dense',
  • sage/libs/m4ri.pxd

    diff --git a/sage/libs/m4ri.pxd b/sage/libs/m4ri.pxd
    a b  
    11cdef extern from "m4ri/m4ri.h":
     2    ctypedef int rci_t
     3    ctypedef int wi_t
    24    ctypedef unsigned long long m4ri_word "word"
     5    ctypedef int BIT
    36
    47    ctypedef struct mzd_t:
    5         int nrows
    6         int ncols
    7         int width
     8        rci_t nrows
     9        rci_t ncols
     10        wi_t width
     11        int offset
    812        m4ri_word **rows
    913
    1014    ctypedef struct mzp_t:
    11         int *values
    12         int size
     15        rci_t *values
     16        rci_t size
    1317
    14     ctypedef int BIT
    15 
    16     cdef int RADIX
     18    cdef int m4ri_radix
    1719
    1820    ##############
    1921    # Maintainance
     
    2830    ##############
    2931 
    3032    # create empty matrix
    31     cdef mzd_t *mzd_init(int,int)
     33    cdef mzd_t *mzd_init(rci_t , rci_t)
    3234
    3335    # create the identity permutation
    34     cdef mzp_t *mzp_init(long)
     36    cdef mzp_t *mzp_init(rci_t)
    3537   
    3638    # free memory for the matrix
    3739    cdef void mzd_free(mzd_t *)
     
    4345    cdef void mzd_randomize(mzd_t *)
    4446
    4547    # identity matrix if i%2
    46     cdef void mzd_set_ui(mzd_t *, unsigned int i)
     48    cdef void mzd_set_ui(mzd_t *, unsigned int )
    4749
    4850    # [A],[B] -> [AB]
    49     cdef mzd_t *mzd_concat(mzd_t *C, mzd_t *A, mzd_t *B)
     51    cdef mzd_t *mzd_concat(mzd_t *, mzd_t *, mzd_t *)
    5052
    5153    # [A],[B] -> | A |
    5254    #            | B |
    53     cdef mzd_t *mzd_stack(mzd_t *C, mzd_t *A, mzd_t *B)
     55    cdef mzd_t *mzd_stack(mzd_t *, mzd_t *, mzd_t *)
    5456   
    5557    # 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)
     58    cdef mzd_t *mzd_submatrix(mzd_t *, mzd_t *, rci_t lowr, rci_t lowc, rci_t highr, rci_t highc)
    5759
    5860    # return a matrix window to A
    59     cdef mzd_t *mzd_init_window( mzd_t *A, int lowr, int lowc, int highr, int highc)
     61    cdef mzd_t *mzd_init_window(mzd_t *, rci_t lowr, rci_t lowc, rci_t highr, rci_t highc)
    6062
    6163    cdef void mzd_free_window(mzd_t *)
    6264
    6365    # deep copy
    6466    cdef mzd_t *mzd_copy(mzd_t *, mzd_t *)
    6567
     68    # printing
     69    cdef void mzd_print(mzd_t *)
     70
    6671    ##############
    6772    # Bit Level IO
    6873    ##############
    6974
    7075    # set BIT
    71     cdef void mzd_write_bit( mzd_t *m, int row, int col, BIT value)
     76    cdef void mzd_write_bit( mzd_t *m, rci_t row, rci_t col, BIT value)
    7277
    7378    # get BIT
    74     cdef BIT mzd_read_bit( mzd_t *m, int row, int col )
     79    cdef BIT mzd_read_bit( mzd_t *m, rci_t row, rci_t col )
    7580
    7681    # get BITs (n<=64)
    77     cdef m4ri_word mzd_read_bits( mzd_t *m, int row, int col, int n)
     82    cdef m4ri_word mzd_read_bits( mzd_t *m, rci_t row, rci_t col, int n)
    7883   
    7984    #####################
    8085    # Row/Column Based IO
    8186    #####################
    8287   
    83     cdef void mzd_row_swap(mzd_t *, int, int)
     88    cdef void mzd_row_swap(mzd_t *, rci_t, rci_t)
    8489
    85     cdef void mzd_col_swap(mzd_t *, int, int)
     90    cdef void mzd_col_swap(mzd_t *, rci_t, rci_t)
    8691
    87     cdef void mzd_row_clear_offset(mzd_t *m, int, int)
     92    cdef void mzd_row_clear_offset(mzd_t *m, rci_t, rci_t)
    8893
    89     cdef void mzd_row_add_offset(mzd_t *m, int, int, int)
     94    cdef void mzd_row_add_offset(mzd_t *m, rci_t, rci_t, rci_t)
    9095
    9196    ############
    9297    # Arithmetic
     
    171176    # reduced row echelon form using PLUQ factorization
    172177    cdef mzd_t *mzd_kernel_left_pluq(mzd_t *A, int cutoff)
    173178
     179    ########################
     180    # Bit operations
     181    ########################
     182
     183    cdef m4ri_word __M4RI_LEFT_BITMASK(int)
     184
     185    cdef m4ri_word m4ri_swap_bits(m4ri_word)
     186
    174187    ##################################
    175188    # Internal Functions for Debugging
    176189    ##################################
    177190
    178     cdef void mzd_write_zeroed_bits(mzd_t *m, int x, int y, int n, m4ri_word values)
    179      
    180191    cdef void mzd_clear_bits(mzd_t *m, int x, int y, int n)
    181192
  • sage/matrix/matrix_mod2_dense.pyx

    diff --git a/sage/matrix/matrix_mod2_dense.pyx b/sage/matrix/matrix_mod2_dense.pyx
    a b  
    350350        cdef unsigned long i, j, truerow
    351351        cdef unsigned long start, shift
    352352        cdef m4ri_word row_xor
    353         cdef m4ri_word end_mask = ~(((<m4ri_word>1)<<(RADIX - self._ncols%RADIX))-1)
     353        cdef m4ri_word end_mask = __M4RI_LEFT_BITMASK(self._ncols%m4ri_radix)
    354354        cdef m4ri_word top_mask, bot_mask
    355355        cdef m4ri_word cur
    356356        cdef m4ri_word* row
     
    370370            row = self._entries.rows[i]
    371371            start = (i*self._entries.ncols) >> 6
    372372            shift = (i*self._entries.ncols) & 0x3F
    373             bot_mask = (~(<m4ri_word>0)) << shift
     373            bot_mask = __M4RI_LEFT_BITMASK(m4ri_radix - shift)
    374374            top_mask = ~bot_mask
    375375
    376376            if self._entries.width > 1:
     
    379379
    380380                for j from 1 <= j < self._entries.width - 1:
    381381                    row_xor ^= row[j]
    382                     cur = ((row[j-1] << (63-shift)) << 1) ^ (row[j] >> shift)
     382                    cur = ((row[j-1] >> (63-shift)) >> 1) ^ (row[j] << shift)
    383383                    running_parity ^= (start+j) & parity_mask(cur)
    384384
    385385                running_parity ^= (start+j) & parity_mask(row[j-1] & top_mask)
     
    393393            running_parity ^= (start+j) & parity_mask(cur & bot_mask)
    394394            running_parity ^= (start+j+1) & parity_mask(cur & top_mask)
    395395
    396             start = (i*self._entries.ncols) & (RADIX-1)
    397             running_xor ^= (row_xor >> start) ^ ((row_xor << (RADIX-1-start)) << 1)
     396            running_xor ^= (row_xor << shift) ^ ((row_xor >> (63-shift)) >> 1)
    398397
    399398        cdef unsigned long bit_is_set
    400399        cdef unsigned long h
    401400
    402401        # Now we assemble the running_parity and running_xor to get the hash.
    403402        # Viewing the flattened matrix as a list of a_i, the hash is the xor
    404         # of the i for which a_i is non-zero. We split i into the lower RADIX
    405         # bits and the rest, so i = i1 << RADIX + i0. Now two matching i0
     403        # of the i for which a_i is non-zero. We split i into the lower m4ri_radix
     404        # bits and the rest, so i = i1 << m4ri_radix + i0. Now two matching i0
    406405        # would cancel, so we only need the parity of how many of each
    407406        # possible i0 occur. This is stored in the bits of running_xor.
    408407        # Similarly, running_parity is the xor of the i1 needed. It's called
     
    410409        # the number of i1 to add is equal to the number of set bits in that
    411410        # word (but because two cancel, we only need keep track of the
    412411        # parity.
    413         h = RADIX * running_parity
    414         for i from 0 <= i < RADIX:
    415             bit_is_set = (running_xor >> (RADIX-1-i)) & 1
    416             h ^= (RADIX-1) & ~(bit_is_set-1) & i
     412
     413        h = m4ri_radix * running_parity
     414        for i from 0 <= i < m4ri_radix:
     415            bit_is_set = (running_xor >> i) & 1
     416            h ^= (m4ri_radix-1) & ~(bit_is_set-1) & i
    417417
    418418        if h == -1:
    419419            h = -2
     
    665665        cdef Vector_mod2_dense c = PY_NEW(Vector_mod2_dense)
    666666        c._init(self._nrows, VS)
    667667        c._entries = mzd_init(1, self._nrows)
    668         tmp = mzd_init(self._nrows, 1)
    669         _mzd_mul_naive(tmp, self._entries, (<Vector_mod2_dense>v)._entries, 0)
    670         mzd_transpose(c._entries, tmp)
    671         mzd_free(tmp)
     668        if c._entries.nrows and c._entries.ncols:
     669            tmp = mzd_init(self._nrows, 1)
     670            _mzd_mul_naive(tmp, self._entries, (<Vector_mod2_dense>v)._entries, 0)
     671            mzd_transpose(c._entries, tmp)
     672            mzd_free(tmp)
    672673        return c
    673674       
    674675    cdef Matrix _matrix_times_matrix_(self, Matrix right):
     
    12661267        cdef int i, j, k
    12671268        cdef int nc
    12681269        cdef unsigned int low, high
    1269         cdef m4ri_word mask = 1
     1270        cdef m4ri_word mask = 0
    12701271       
    12711272        # Original code, before adding the ``nonzero`` option.
    12721273        if not nonzero:
    12731274            if density == 1:
    12741275                assert(sizeof(m4ri_word) == 8)
    1275                 mask = ~((mask<<(RADIX - self._entries.ncols%RADIX))-1)
     1276                mask = __M4RI_LEFT_BITMASK(self._entries.ncols % m4ri_radix)
    12761277                for i from 0 <= i < self._nrows:
    12771278                    for j from 0 <= j < self._entries.width:
    12781279                        # for portability we get 32-bit twice rather than 64-bit once
    12791280                        low = gmp_urandomb_ui(rstate.gmp_state, 32)
    12801281                        high = gmp_urandomb_ui(rstate.gmp_state, 32)
    1281                         self._entries.rows[i][j] = ((<unsigned long long>high)<<32) | (<unsigned long long>low)
     1282                        self._entries.rows[i][j] = m4ri_swap_bits( ((<unsigned long long>high)<<32) | (<unsigned long long>low) )
    12821283                    self._entries.rows[i][self._entries.width - 1] &= mask
    12831284            else:
    12841285                nc = self._ncols