Ticket #11574: m4ri_20110601.patch

File m4ri_20110601.patch, 10.0 KB (added by malb, 3 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