Changeset 7761:04d1da14e7cc


Ignore:
Timestamp:
11/13/07 11:49:55 (6 years ago)
Author:
Robert L Miller <rlm@…>
Branch:
default
Message:

4th half... datastructures mostly done, algorithms remain

Location:
sage/coding
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • sage/coding/binary_code.pxd

    r7760 r7761  
    22include '../ext/cdefs.pxi' 
    33 
     4cdef int *hamming_weights() 
     5 
    46cdef class BinaryCode: 
    5     cdef unsigned int *basis 
    6 #    cdef unsigned int *columns 
    7     cdef unsigned int *words 
     7    cdef int *basis 
     8    cdef int *words 
    89    cdef int ncols 
    910    cdef int nrows 
    1011    cdef int radix 
    11     cdef unsigned int nwords 
    12     cdef int is_one(self, unsigned int, int) 
    13     cdef int is_automorphism(self, int *, unsigned int *) 
     12    cdef int nwords 
     13 
     14    cdef int is_one(self, int, int) 
     15    cdef int is_automorphism(self, int *, int *) 
    1416 
    1517cdef class OrbitPartition: 
    16     cdef unsigned int nwords 
     18    cdef int nwords 
    1719    cdef int ncols 
    18     cdef unsigned int *wd_parent 
    19     cdef unsigned int *wd_rank 
    20     cdef unsigned int *wd_min_cell_rep 
    21     cdef unsigned int *wd_size 
     20    cdef int *wd_parent 
     21    cdef int *wd_rank 
     22    cdef int *wd_min_cell_rep 
     23    cdef int *wd_size 
    2224    cdef int *col_parent 
    2325    cdef int *col_rank 
     
    2527    cdef int *col_size 
    2628 
    27     cdef unsigned int wd_find(self, unsigned int) 
    28     cdef void wd_union(self, unsigned int, unsigned int) 
     29    cdef int wd_find(self, int) 
     30    cdef void wd_union(self, int, int) 
    2931    cdef int col_find(self, int) 
    3032    cdef void col_union(self, int, int) 
    31     cdef int merge_perm(self, int *, unsigned int *) 
     33    cdef int merge_perm(self, int *, int *) 
    3234 
    3335cdef class PartitionStack: 
    34     cdef unsigned int *wd_ents 
     36    cdef int *wd_ents 
    3537    cdef int *wd_lvls 
    3638    cdef int *col_ents 
    3739    cdef int *col_lvls 
    38     cdef unsigned int nwords 
     40    cdef int *basis_locations 
     41    cdef int nwords 
     42    cdef int nrows 
    3943    cdef int ncols 
    4044    cdef int radix 
    41     cdef unsigned int *col_degs   # 
    42     cdef int *col_counts          # 
    43     cdef int *col_output          # 
    44     cdef int *wd_degs             # 
    45     cdef unsigned int *wd_counts  # These are just for scratch space... 
    46     cdef unsigned int *wd_output  # 
     45    cdef int flag 
     46    cdef int *col_degs   # 
     47    cdef int *col_counts # 
     48    cdef int *col_output # 
     49    cdef int *wd_degs    # 
     50    cdef int *wd_counts  # These are just for scratch space... 
     51    cdef int *wd_output  # 
    4752 
    4853    cdef int is_discrete(self, int) 
    4954    cdef int num_cells(self, int) 
    5055    cdef int sat_225(self, int) 
    51     cdef unsigned int min_cell_reps(self, int) 
    52     cdef unsigned int fixed_cols(self, unsigned int, int) 
    53     cdef unsigned int first_smallest_nontrivial(self, int) 
     56    cdef int min_cell_reps(self, int) 
     57    cdef int fixed_cols(self, int, int) 
     58    cdef int first_smallest_nontrivial(self, int) 
    5459    cdef void col_percolate(self, int, int) 
    55     cdef void wd_percolate(self, unsigned int, unsigned int) 
     60    cdef void wd_percolate(self, int, int) 
    5661    cdef int split_column(self, int, int) 
    57     cdef unsigned int col_degree(self, BinaryCode, int, unsigned int, int) 
    58     cdef int wd_degree(self, BinaryCode, unsigned int, int, int) 
     62    cdef int col_degree(self, BinaryCode, int, int, int) 
     63    cdef int wd_degree(self, BinaryCode, int, int, int, int *) 
    5964    cdef int sort_cols(self, int, int) 
    60     cdef unsigned int sort_wds(self, unsigned int, int) 
    61  
    62 ################################################################################ 
    63 ################################################################################ 
    64 ################################################################################ 
    65  
    66     cdef unsigned int refine(self, int, unsigned int *, int *, BinaryCode) 
    67     cdef void get_permutation(self, PartitionStack, int *, int *) 
     65    cdef int sort_wds(self, int, int) 
     66    cdef int refine(self, int, int *, int, BinaryCode, int *) 
     67    cdef void clear(self, int) 
    6868    cdef int cmp(self, PartitionStack, BinaryCode) 
     69    cdef void find_basis(self, int *) 
     70    cdef void get_permutation(self, PartitionStack, int *, int *, int *) 
    6971 
    7072cdef class BinaryCodeClassifier: 
     
    7981 
    8082 
    81  
    82  
    83  
    84  
    85  
    86  
    87  
    88  
    89  
  • sage/coding/binary_code.pyx

    r7760 r7761  
    3838from sage.structure.element import is_Matrix 
    3939from sage.misc.misc import cputime 
    40 from math import log, ceil 
     40from math import log, floor 
    4141from sage.rings.integer import Integer 
    4242 
     
    4545## essentially only for doctesting and debugging, have underscores. 
    4646 
     47cdef int *hamming_weights(): 
     48    cdef int *ham_wts 
     49    ham_wts = <int *> sage_malloc( 65536 * sizeof(int) ) 
     50    if not ham_wts: 
     51        sage_free(ham_wts) 
     52        raise MemoryError("Memory.") 
     53    ham_wts[0] = 0 
     54    ham_wts[1] = 1 
     55    ham_wts[2] = 1 
     56    ham_wts[3] = 2 
     57    for i from 4 <= i < 16: 
     58        ham_wts[i] = ham_wts[i & 3] + ham_wts[(i>>2) & 3] 
     59    for i from 16 <= i < 256: 
     60        ham_wts[i] = ham_wts[i & 15] + ham_wts[(i>>4) & 15] 
     61    for i from 256 <= i < 65536: 
     62        ham_wts[i] = ham_wts[i & 255] + ham_wts[(i>>8) & 255] 
     63    return ham_wts 
     64    # This may seem like overkill, but the worst case for storing the words 
     65    # themselves is 65536- in this case, we are increasing memory usage by a 
     66    # factor of 2. Also, this function will only be called once ever. 
     67 
    4768cdef class BinaryCode: 
    4869    """ 
    4970    Minimal, but optimized, binary code object. 
    5071 
     72    EXAMPLE: 
     73        sage: import sage.coding.binary_code 
     74        sage: from sage.coding.binary_code import * 
     75        sage: M = Matrix(GF(2), [[1,1,1,1]]) 
     76        sage: B = BinaryCode(M)     # create from matrix 
     77        sage: C = BinaryCode(B, 60) # create using glue 
     78        sage: D = BinaryCode(C, 240) 
     79        sage: E = BinaryCode(D, 85) 
     80        sage: B 
     81        Binary [4,1] linear code, generator matrix 
     82        [1111] 
     83        sage: C 
     84        Binary [6,2] linear code, generator matrix 
     85        [111100] 
     86        [001111] 
     87        sage: D 
     88        Binary [8,3] linear code, generator matrix 
     89        [11110000] 
     90        [00111100] 
     91        [00001111] 
     92        sage: E 
     93        Binary [8,4] linear code, generator matrix 
     94        [11110000] 
     95        [00111100] 
     96        [00001111] 
     97        [10101010] 
     98 
    5199    """ 
    52     def __new__(self, arg1): 
    53         cdef unsigned int parity, word, combination 
     100    def __new__(self, arg1, arg2=None): 
    54101        cdef int nrows, i, j 
    55         cdef unsigned int *self_words, *self_basis 
    56         self.radix = 8*sizeof(unsigned int) 
    57         self.ncols = arg1.ncols() 
    58         self.nrows = arg1.nrows() 
    59         self.nwords = 1 << self.nrows 
    60  
    61         if self.nrows > self.radix or self.ncols > self.radix: 
     102        cdef int nwords, other_nwords, parity, word, combination, glue_word 
     103        cdef BinaryCode other 
     104        cdef int *self_words, *self_basis, *other_basis 
     105 
     106        self.radix = sizeof(int) << 3 
     107 
     108        if is_Matrix(arg1): 
     109            self.ncols = arg1.ncols() 
     110            self.nrows = arg1.nrows() 
     111            nrows = self.nrows 
     112            self.nwords = 1 << nrows 
     113            nwords = self.nwords 
     114        elif isinstance(arg1, BinaryCode): 
     115            other = arg1 
     116            self.nrows = other.nrows + 1 
     117            glue_word = arg2 
     118            self.ncols = max( other.ncols , floor(log(glue_word,2))+1 ) 
     119            other_nwords = other.nwords 
     120            self.nwords = 2 * other_nwords 
     121            nrows = self.nrows 
     122            nwords = self.nwords 
     123        else: raise NotImplementedError("!") 
     124 
     125        if self.nrows >= self.radix or self.ncols > self.radix: 
    62126            raise NotImplementedError("Columns and rows are stored as ints. This code is too big.") 
    63127 
    64         self.words = <unsigned int *> sage_malloc( self.nwords * sizeof(unsigned int) ) 
    65         self.basis = <unsigned int *> sage_malloc( self.nrows * sizeof(unsigned int) ) 
     128        self.words = <int *> sage_malloc( nwords * sizeof(int) ) 
     129        self.basis = <int *> sage_malloc( nrows * sizeof(int) ) 
    66130        if not self.words or not self.basis: 
    67131            if self.words: sage_free(self.words) 
    68132            if self.basis: sage_free(self.basis) 
    69133            raise MemoryError("Memory.") 
    70  
    71         nrows = self.nrows 
    72134        self_words = self.words 
    73135        self_basis = self.basis 
    74136 
    75         rows = arg1.rows() 
    76         for i from 0 <= i < nrows: 
     137        if is_Matrix(arg1): 
     138            rows = arg1.rows() 
     139            for i from 0 <= i < nrows: 
     140                word = 0 
     141                for j in rows[i].nonzero_positions(): 
     142                    word += (1<<j) 
     143                self_basis[i] = word 
     144 
    77145            word = 0 
    78             for j in rows[i].nonzero_positions(): 
    79                 word += (1<<j) 
    80             self_basis[i] = word 
    81  
    82         word = 0 
    83         parity = 0 
    84         combination = 0 
    85         while True: 
    86             self_words[combination] = word 
    87             parity ^= 1 
    88             j = 0 
    89             if not parity: 
    90                 while not combination & (1 << j): j += 1 
    91                 j += 1 
    92             if j == nrows: break 
    93             else: 
    94                 combination ^= (1 << j) 
    95                 word ^= self_basis[j] 
    96         # NOTE: when implementing construction from parent 
    97         # matrices, remember to grab this data from the parent 
     146            parity = 0 
     147            combination = 0 
     148            while True: 
     149                self_words[combination] = word 
     150                parity ^= 1 
     151                j = 0 
     152                if not parity: 
     153                    while not combination & (1 << j): j += 1 
     154                    j += 1 
     155                if j == nrows: break 
     156                else: 
     157                    combination ^= (1 << j) 
     158                    word ^= self_basis[j] 
     159 
     160        else: # isinstance(arg1, BinaryCode) 
     161            other_basis = other.basis 
     162            for i from 0 <= i < nrows-1: 
     163                self_basis[i] = other_basis[i] 
     164            i = nrows - 1 
     165            self_basis[i] = glue_word 
     166 
     167            memcpy(self_words, other.words, other_nwords*(self.radix>>3)) 
     168 
     169            for combination from 0 <= combination < other_nwords: 
     170                self_words[combination+other_nwords] = self_words[combination] ^ glue_word 
    98171 
    99172    def __dealloc__(self): 
     
    101174        sage_free(self.basis) 
    102175 
    103     def delete(self): 
    104         sage_free(self.words) 
    105         sage_free(self.basis) 
     176    def print_data(self): 
     177        """ 
     178        Print all data for self. 
     179 
     180        EXAMPLES: 
     181            sage: import sage.coding.binary_code 
     182            sage: from sage.coding.binary_code import * 
     183            sage: M = Matrix(GF(2), [[1,1,1,1]]) 
     184            sage: B = BinaryCode(M) 
     185            sage: C = BinaryCode(B, 60) 
     186            sage: D = BinaryCode(C, 240) 
     187            sage: E = BinaryCode(D, 85) 
     188            sage: B.print_data() # random - actually "print P.print_data()" 
     189            ncols: 4 
     190            nrows: 1 
     191            nwords: 2 
     192            radix: 32 
     193            basis: 
     194            1111 
     195            words: 
     196            0000 
     197            1111 
     198            sage: C.print_data() # random - actually "print P.print_data()" 
     199            ncols: 6 
     200            nrows: 2 
     201            nwords: 4 
     202            radix: 32 
     203            basis: 
     204            111100 
     205            001111 
     206            words: 
     207            000000 
     208            111100 
     209            001111 
     210            110011 
     211            sage: D.print_data() # random - actually "print P.print_data()" 
     212            ncols: 8 
     213            nrows: 3 
     214            nwords: 8 
     215            radix: 32 
     216            basis: 
     217            11110000 
     218            00111100 
     219            00001111 
     220            words: 
     221            00000000 
     222            11110000 
     223            00111100 
     224            11001100 
     225            00001111 
     226            11111111 
     227            00110011 
     228            11000011 
     229            sage: E.print_data() # random - actually "print P.print_data()" 
     230            ncols: 8 
     231            nrows: 4 
     232            nwords: 16 
     233            radix: 32 
     234            basis: 
     235            11110000 
     236            00111100 
     237            00001111 
     238            10101010 
     239            words: 
     240            00000000 
     241            11110000 
     242            00111100 
     243            11001100 
     244            00001111 
     245            11111111 
     246            00110011 
     247            11000011 
     248            10101010 
     249            01011010 
     250            10010110 
     251            01100110 
     252            10100101 
     253            01010101 
     254            10011001 
     255            01101001 
     256        """ 
     257        from sage.graphs.graph_fast import binary 
     258        cdef int ui 
     259        cdef int i 
     260        s = '' 
     261        s += "ncols:" + str(self.ncols) 
     262        s += "\nnrows:" + str(self.nrows) 
     263        s += "\nnwords:" + str(self.nwords) 
     264        s += "\nradix:" + str(self.radix) 
     265        s += "\nbasis:\n" 
     266        for i from 0 <= i < self.nrows: 
     267            b = list(binary(self.basis[i]).zfill(self.ncols)) 
     268            b.reverse() 
     269            b.append('\n') 
     270            s += ''.join(b) 
     271        s += "\nwords:\n" 
     272        for ui from 0 <= ui < self.nwords: 
     273            b = list(binary(self.words[ui]).zfill(self.ncols)) 
     274            b.reverse() 
     275            b.append('\n') 
     276            s += ''.join(b) 
    106277 
    107278    def __repr__(self): 
     
    158329        return self.is_one(word, col) 
    159330 
    160     cdef int is_one(self, unsigned int word, int column): 
     331    cdef int is_one(self, int word, int column): 
    161332        if self.words[word] & (1 << column): 
    162333            return 1 
     
    193364        cdef int i 
    194365        cdef int *_col_gamma 
    195         cdef unsigned int *_basis_gamma 
    196         _basis_gamma = <unsigned int *> sage_malloc(self.nrows * sizeof(unsigned int)) 
     366        cdef int *_basis_gamma 
     367        _basis_gamma = <int *> sage_malloc(self.nrows * sizeof(int)) 
    197368        _col_gamma = <int *> sage_malloc(self.ncols * sizeof(int)) 
    198369        if not (_col_gamma and _basis_gamma): 
     
    209380        return result 
    210381 
    211     cdef int is_automorphism(self, int *col_gamma, unsigned int *basis_gamma): 
     382    cdef int is_automorphism(self, int *col_gamma, int *basis_gamma): 
    212383        cdef int i, j, self_nrows = self.nrows, self_ncols = self.ncols 
    213384        for i from 0 <= i < self_nrows: 
     
    229400    """ 
    230401    def __new__(self, nrows, ncols): 
    231         cdef unsigned int nwords, word 
     402        cdef int nwords, word 
    232403        cdef int col 
    233404        nwords = (1 << nrows) 
    234405        self.nwords = nwords 
    235406        self.ncols = ncols 
    236         self.wd_parent =       <unsigned int *> sage_malloc( nwords * sizeof(unsigned int) ) 
    237         self.wd_rank =         <unsigned int *> sage_malloc( nwords * sizeof(unsigned int) ) 
    238         self.wd_min_cell_rep = <unsigned int *> sage_malloc( nwords * sizeof(unsigned int) ) 
    239         self.wd_size =         <unsigned int *> sage_malloc( nwords * sizeof(unsigned int) ) 
     407        self.wd_parent =       <int *> sage_malloc( nwords * sizeof(int) ) 
     408        self.wd_rank =         <int *> sage_malloc( nwords * sizeof(int) ) 
     409        self.wd_min_cell_rep = <int *> sage_malloc( nwords * sizeof(int) ) 
     410        self.wd_size =         <int *> sage_malloc( nwords * sizeof(int) ) 
    240411        self.col_parent =       <int *> sage_malloc( ncols * sizeof(int) ) 
    241412        self.col_rank =         <int *> sage_malloc( ncols * sizeof(int) ) 
     
    273444        sage_free(self.col_size) 
    274445 
    275     def delete(self): 
    276         sage_free(self.wd_parent) 
    277         sage_free(self.wd_rank) 
    278         sage_free(self.wd_min_cell_rep) 
    279         sage_free(self.wd_size) 
    280         sage_free(self.col_parent) 
    281         sage_free(self.col_rank) 
    282         sage_free(self.col_min_cell_rep) 
    283         sage_free(self.col_size) 
    284  
    285446    def __repr__(self): 
    286447        """ 
     
    299460 
    300461        """ 
    301         cdef unsigned int i 
     462        cdef int i 
    302463        cdef int j 
    303464        s = 'OrbitPartition on %d words and %d columns. Data:\nWords:\n'%(self.nwords, self.ncols) 
     
    324485            0,1,2,3,4,5,6,7 
    325486            sage: O._wd_find(12) 
    326             12L 
     487            12 
    327488 
    328489        """ 
    329490        return self.wd_find(word) 
    330491 
    331     cdef unsigned int wd_find(self, unsigned int word): 
     492    cdef int wd_find(self, int word): 
    332493        if self.wd_parent[word] == word: 
    333494            return word 
     
    358519            0,1,2,3,4,5,6,7 
    359520            sage: O._wd_find(10) 
    360             1L 
     521            1 
    361522 
    362523        """ 
    363524        self.wd_union(x, y) 
    364525 
    365     cdef void wd_union(self, unsigned int x, unsigned int y): 
    366         cdef unsigned int x_root, y_root 
     526    cdef void wd_union(self, int x, int y): 
     527        cdef int x_root, y_root 
    367528        x_root = self.wd_find(x) 
    368529        y_root = self.wd_find(y) 
     
    479640 
    480641        """ 
    481         cdef unsigned int i 
     642        cdef int i 
    482643        cdef int *_col_gamma 
    483         cdef unsigned int *_wd_gamma 
    484         _wd_gamma = <unsigned int *> sage_malloc(self.nwords * sizeof(unsigned int)) 
     644        cdef int *_wd_gamma 
     645        _wd_gamma = <int *> sage_malloc(self.nwords * sizeof(int)) 
    485646        _col_gamma = <int *> sage_malloc(self.ncols * sizeof(int)) 
    486647        if not (_col_gamma and _wd_gamma): 
     
    497658        return result 
    498659 
    499     cdef int merge_perm(self, int *col_gamma, unsigned int *wd_gamma): 
    500         cdef unsigned int i, gamma_i_root 
     660    cdef int merge_perm(self, int *col_gamma, int *wd_gamma): 
     661        cdef int i, gamma_i_root 
    501662        cdef int j, gamma_j_root, return_value = 0 
    502         cdef unsigned int *self_wd_parent = self.wd_parent 
     663        cdef int *self_wd_parent = self.wd_parent 
    503664        cdef int *self_col_parent = self.col_parent 
    504665        for i from 0 <= i < self.nwords: 
     
    523684    """ 
    524685    def __new__(self, arg1, arg2=None): 
    525         cdef int k 
    526         self.nwords = int(arg1) 
    527         self.ncols = int(arg2) 
    528         self.radix      = 8*sizeof(unsigned int) 
     686        cdef int k, nwords, ncols 
     687        cdef PartitionStack other 
     688        cdef int *wd_ents, *wd_lvls, *col_ents, *col_lvls 
     689 
     690        try: 
     691            self.nrows = int(arg1) 
     692            self.nwords = 1 << self.nrows 
     693            self.ncols = int(arg2) 
     694        except: 
     695            other = arg1 
     696            self.nrows = other.nrows 
     697            self.nwords = other.nwords 
     698            self.ncols = other.ncols 
     699        self.radix = 8*sizeof(int) 
     700        self.flag = (1 << (self.radix-1)) 
    529701 
    530702        # data 
    531         self.wd_ents    = <unsigned int *> sage_malloc( self.nwords * sizeof(unsigned int) ) 
    532         self.wd_lvls    = <int *> sage_malloc( self.nwords * sizeof(int) ) 
    533         self.col_ents   = <int *> sage_malloc( self.ncols * sizeof(int) ) 
    534         self.col_lvls   = <int *> sage_malloc( self.ncols * sizeof(int) ) 
     703        self.wd_ents =    <int *> sage_malloc( self.nwords * sizeof(int) ) 
     704        self.wd_lvls =    <int *> sage_malloc( self.nwords * sizeof(int) ) 
     705        self.col_ents =   <int *> sage_malloc( self.ncols * sizeof(int) ) 
     706        self.col_lvls =   <int *> sage_malloc( self.ncols * sizeof(int) ) 
    535707 
    536708        # scratch space 
    537         self.col_degs   = <unsigned int *> sage_malloc( self.ncols * sizeof(unsigned int) ) 
     709        self.col_degs =   <int *> sage_malloc( self.ncols  * sizeof(int) ) 
    538710        self.col_counts = <int *> sage_malloc( self.nwords * sizeof(int) ) 
    539         self.col_output = <int *> sage_malloc( self.ncols * sizeof(int) ) 
    540         self.wd_degs    = <int *> sage_malloc( self.nwords * sizeof(int) ) 
    541         self.wd_counts  = <unsigned int *> sage_malloc( self.ncols * sizeof(unsigned int) ) 
    542         self.wd_output  = <unsigned int *> sage_malloc( self.nwords * sizeof(unsigned int) ) 
    543  
    544         if not (self.wd_ents and self.wd_lvls and self.col_ents and self.col_lvls \ 
     711        self.col_output = <int *> sage_malloc( self.ncols  * sizeof(int) ) 
     712        self.wd_degs =    <int *> sage_malloc( self.nwords * sizeof(int) ) 
     713        self.wd_counts =  <int *> sage_malloc( (self.ncols+1)  * sizeof(int) ) 
     714        self.wd_output =  <int *> sage_malloc( self.nwords * sizeof(int) ) 
     715 
     716        if not (self.wd_ents  and self.wd_lvls    and self.col_ents   and self.col_lvls \ 
    545717            and self.col_degs and self.col_counts and self.col_output \ 
    546             and self.wd_degs and self.wd_counts and self.wd_output): 
    547             if self.wd_ents: sage_free(self.wd_ents) 
    548             if self.wd_lvls: sage_free(self.wd_lvls) 
    549             if self.col_ents: sage_free(self.col_ents) 
    550             if self.col_lvls: sage_free(self.col_lvls) 
    551             if self.col_degs: sage_free(self.col_degs) 
     718            and self.wd_degs  and self.wd_counts and self.wd_output): 
     719            if self.wd_ents:    sage_free(self.wd_ents) 
     720            if self.wd_lvls:    sage_free(self.wd_lvls) 
     721            if self.col_ents:   sage_free(self.col_ents) 
     722            if self.col_lvls:   sage_free(self.col_lvls) 
     723            if self.col_degs:   sage_free(self.col_degs) 
    552724            if self.col_counts: sage_free(self.col_counts) 
    553725            if self.col_output: sage_free(self.col_output) 
    554             if self.wd_degs: sage_free(self.wd_degs) 
    555             if self.wd_counts: sage_free(self.wd_counts) 
    556             if self.wd_output: sage_free(self.wd_output) 
     726            if self.wd_degs:    sage_free(self.wd_degs) 
     727            if self.wd_counts:  sage_free(self.wd_counts) 
     728            if self.wd_output:  sage_free(self.wd_output) 
    557729            raise MemoryError("Memory.") 
    558730 
    559         for k from 0 <= k < self.nwords-1: 
    560             self.wd_ents[k] = k 
    561             self.wd_lvls[k] = self.nwords 
    562         for k from 0 <= k < self.ncols-1: 
    563             self.col_ents[k] = k 
    564             self.col_lvls[k] = self.ncols 
    565         self.wd_ents[self.nwords-1] = self.nwords-1 
    566         self.wd_lvls[self.nwords-1] = -1 
    567         self.col_ents[self.ncols-1] = self.ncols-1 
    568         self.col_lvls[self.ncols-1] = -1 
     731        if other: 
     732            memcpy(self.wd_ents,  other.wd_ents, self.nwords*(self.radix>>3)) 
     733            memcpy(self.wd_lvls,  other.wd_lvls, self.nwords*(self.radix>>3)) 
     734            memcpy(self.col_ents, other.col_ents, self.ncols*(self.radix>>3)) 
     735            memcpy(self.col_lvls, other.col_lvls, self.ncols*(self.radix>>3)) 
     736        else: 
     737            wd_ents = self.wd_ents 
     738            wd_lvls = self.wd_lvls 
     739            col_ents = self.col_ents 
     740            col_lvls = self.col_lvls 
     741            nwords = self.nwords 
     742            ncols = self.ncols 
     743            for k from 0 <= k < nwords-1: 
     744                wd_ents[k] = k 
     745                wd_lvls[k] = nwords 
     746            for k from 0 <= k < ncols-1: 
     747                col_ents[k] = k 
     748                col_lvls[k] = ncols 
     749            wd_ents[nwords-1] = nwords-1 
     750            wd_lvls[nwords-1] = -1 
     751            col_ents[ncols-1] = ncols-1 
     752            col_lvls[ncols-1] = -1 
    569753 
    570754    def __dealloc__(self): 
     755        if self.basis_locations: sage_free(self.basis_locations) 
    571756        sage_free(self.wd_ents) 
    572757        sage_free(self.wd_lvls) 
     
    580765        sage_free(self.wd_output) 
    581766 
    582     def delete(self): 
    583         sage_free(self.wd_ents) 
    584         sage_free(self.wd_lvls) 
    585         sage_free(self.col_ents) 
    586         sage_free(self.col_lvls) 
    587         sage_free(self.col_degs) 
    588         sage_free(self.col_counts) 
    589         sage_free(self.col_output) 
    590         sage_free(self.wd_degs) 
    591         sage_free(self.wd_counts) 
    592         sage_free(self.wd_output) 
     767    def print_data(self): 
     768        """ 
     769        Prints all data for self. 
     770 
     771        EXAMPLE: 
     772            sage: import sage.coding.binary_code 
     773            sage: from sage.coding.binary_code import * 
     774            sage: P = PartitionStack(2, 6) 
     775            sage: P.print_data() # random - actually "print P.print_data()" 
     776            nwords: 4 
     777            nrows: 2 
     778            ncols: 6 
     779            radix: 32 
     780            wd_ents: 
     781            0 
     782            1 
     783            2 
     784            3 
     785            wd_lvls: 
     786            4 
     787            4 
     788            4 
     789            -1 
     790            col_ents: 
     791            0 
     792            1 
     793            2 
     794            3 
     795            4 
     796            5 
     797            col_lvls: 
     798            6 
     799            6 
     800            6 
     801            6 
     802            6 
     803            -1 
     804            col_degs: 
     805            -1209339024 
     806            145606688 
     807            135493408 
     808            3 
     809            -1210787264 
     810            -1210787232 
     811            col_counts: 
     812            -1209339024 
     813            145666744 
     814            40129536 
     815            21248 
     816            col_output: 
     817            -1209339024 
     818            145654064 
     819            0 
     820            0 
     821            0 
     822            0 
     823            wd_degs: 
     824            -1209339024 
     825            145166112 
     826            16 
     827            3 
     828            wd_counts: 
     829            -1209339024 
     830            146261160 
     831            0 
     832            0 
     833            0 
     834            0 
     835            wd_output: 
     836            -1209339024 
     837            146424680 
     838            135508928 
     839            3 
     840 
     841        """ 
     842        cdef int i, j 
     843        s = '' 
     844        s += "nwords:" + str(self.nwords) + '\n' 
     845        s += "nrows:" + str(self.nrows) + '\n' 
     846        s += "ncols:" + str(self.ncols) + '\n' 
     847        s += "radix:" + str(self.radix) + '\n' 
     848        s += "wd_ents:" + '\n' 
     849        for i from 0 <= i < self.nwords: 
     850            s += str(self.wd_ents[i]) + '\n' 
     851        s += "wd_lvls:" + '\n' 
     852        for i from 0 <= i < self.nwords: 
     853            s += str(self.wd_lvls[i]) + '\n' 
     854        s += "col_ents:" + '\n' 
     855        for i from 0 <= i < self.ncols: 
     856            s += str(self.col_ents[i]) + '\n' 
     857        s += "col_lvls:" + '\n' 
     858        for i from 0 <= i < self.ncols: 
     859            s += str(self.col_lvls[i]) + '\n' 
     860        s += "col_degs:" + '\n' 
     861        for i from 0 <= i < self.ncols: 
     862            s += str(self.col_degs[i]) + '\n' 
     863        s += "col_counts:" + '\n' 
     864        for i from 0 <= i < self.nwords: 
     865            s += str(self.col_counts[i]) + '\n' 
     866        s += "col_output:" + '\n' 
     867        for i from 0 <= i < self.ncols: 
     868            s += str(self.col_output[i]) + '\n' 
     869        s += "wd_degs:" + '\n' 
     870        for i from 0 <= i < self.nwords: 
     871            s += str(self.wd_degs[i]) + '\n' 
     872        s += "wd_counts:" + '\n' 
     873        for i from 0 <= i < self.ncols: 
     874            s += str(self.wd_counts[i]) + '\n' 
     875        s += "wd_output:" + '\n' 
     876        for i from 0 <= i < self.nwords: 
     877            s += str(self.wd_output[i]) + '\n' 
     878        if self.basis_locations: 
     879            s += "basis_locations:" + '\n' 
     880            j = 1 
     881            while (1 << j) < self.nwords: 
     882                j += 1 
     883            for i from 0 <= i < j: 
     884                s += str(self.basis_locations[i]) + '\n' 
     885        return s 
    593886 
    594887    def __repr__(self): 
     
    599892            sage: import sage.coding.binary_code 
    600893            sage: from sage.coding.binary_code import * 
    601             sage: P = PartitionStack(4, 6) 
     894            sage: P = PartitionStack(2, 6) 
    602895            sage: P 
    603896            ({0,1,2,3}) 
     
    638931        """ 
    639932        Returns whether the partition at level k is discrete. 
    640          
    641         EXAMPLE: 
    642             sage: import sage.coding.binary_code 
    643             sage: from sage.coding.binary_code import * 
    644             sage: P = PartitionStack(4, 6) 
     933 
     934        EXAMPLE: 
     935            sage: import sage.coding.binary_code 
     936            sage: from sage.coding.binary_code import * 
     937            sage: P = PartitionStack(2, 6) 
    645938            sage: [P._split_column(i,i+1) for i in range(5)] 
    646939            [0, 1, 2, 3, 4] 
     
    679972            sage: import sage.coding.binary_code 
    680973            sage: from sage.coding.binary_code import * 
    681             sage: P = PartitionStack(4, 6) 
     974            sage: P = PartitionStack(2, 6) 
    682975            sage: [P._split_column(i,i+1) for i in range(5)] 
    683976            [0, 1, 2, 3, 4] 
     
    7201013            sage: import sage.coding.binary_code 
    7211014            sage: from sage.coding.binary_code import * 
    722             sage: P = PartitionStack(4, 6) 
     1015            sage: P = PartitionStack(2, 6) 
    7231016            sage: [P._split_column(i,i+1) for i in range(5)] 
    7241017            [0, 1, 2, 3, 4] 
     
    7781071            sage: import sage.coding.binary_code 
    7791072            sage: from sage.coding.binary_code import * 
    780             sage: P = PartitionStack(4, 6) 
     1073            sage: P = PartitionStack(2, 6) 
    7811074            sage: [P._split_column(i,i+1) for i in range(5)] 
    7821075            [0, 1, 2, 3, 4] 
     
    7991092        return self.min_cell_reps(k) 
    8001093 
    801     cdef unsigned int min_cell_reps(self, int k): 
     1094    cdef int min_cell_reps(self, int k): 
    8021095        cdef int i 
    803         cdef unsigned int reps = 1 
     1096        cdef int reps = 1 
    8041097        cdef int *self_col_lvls = self.col_lvls 
    8051098        for i from 0 < i < self.ncols: 
     
    8161109            sage: import sage.coding.binary_code 
    8171110            sage: from sage.coding.binary_code import * 
    818             sage: P = PartitionStack(4, 6) 
     1111            sage: P = PartitionStack(2, 6) 
    8191112            sage: [P._split_column(i,i+1) for i in range(5)] 
    8201113            [0, 1, 2, 3, 4] 
     
    8371130        return self.fixed_cols(mcrs, k) 
    8381131 
    839     cdef unsigned int fixed_cols(self, unsigned int mcrs, int k): 
     1132    cdef int fixed_cols(self, int mcrs, int k): 
    8401133        cdef int i 
    841         cdef unsigned int fixed = 0 
     1134        cdef int fixed = 0 
    8421135        cdef int *self_col_lvls = self.col_lvls 
    8431136        for i from 0 <= i < self.ncols: 
     
    8481141    def _first_smallest_nontrivial(self, k): 
    8491142        """ 
    850         Returns an integer representing the first, smallest nontrivial cell. 
    851  
    852         EXAMPLE: 
    853             sage: import sage.coding.binary_code 
    854             sage: from sage.coding.binary_code import * 
    855             sage: P = PartitionStack(4, 6) 
     1143        Returns an integer representing the first, smallest nontrivial cell of columns. 
     1144 
     1145        EXAMPLE: 
     1146            sage: import sage.coding.binary_code 
     1147            sage: from sage.coding.binary_code import * 
     1148            sage: P = PartitionStack(2, 6) 
    8561149            sage: [P._split_column(i,i+1) for i in range(5)] 
    8571150            [0, 1, 2, 3, 4] 
     
    8741167        return self.first_smallest_nontrivial(k) 
    8751168 
    876     cdef unsigned int first_smallest_nontrivial(self, int k): 
    877         cdef unsigned int cell 
     1169    cdef int first_smallest_nontrivial(self, int k): 
     1170        cdef int cell 
    8781171        cdef int i = 0, j = 0, location = 0, ncols = self.ncols 
    8791172        cdef int *self_col_lvls = self.col_lvls 
     
    9041197            sage: import sage.coding.binary_code 
    9051198            sage: from sage.coding.binary_code import * 
    906             sage: P = PartitionStack(4, 6) 
     1199            sage: P = PartitionStack(2, 6) 
    9071200            sage: P 
    9081201            ({0,1,2,3}) 
     
    9301223 
    9311224        """ 
    932         cdef unsigned int i 
     1225        cdef int i 
    9331226        for i from 0 <= i < len(col_ents): 
    9341227            self.col_ents[i] = col_ents[i] 
     
    9471240            sage: import sage.coding.binary_code 
    9481241            sage: from sage.coding.binary_code import * 
    949             sage: P = PartitionStack(4, 6) 
     1242            sage: P = PartitionStack(2, 6) 
    9501243            sage: P._dangerous_dont_use_set_ents_lvls(range(5,-1,-1), [1,2,2,3,3,-1], range(3,-1,-1), [1,1,2,-1]) 
    9511244            sage: P 
     
    9891282        """ 
    9901283        Do one round of bubble sort on ents. 
    991          
    992         EXAMPLE: 
    993             sage: import sage.coding.binary_code 
    994             sage: from sage.coding.binary_code import * 
    995             sage: P = PartitionStack(4, 6) 
     1284 
     1285        EXAMPLE: 
     1286            sage: import sage.coding.binary_code 
     1287            sage: from sage.coding.binary_code import * 
     1288            sage: P = PartitionStack(2, 6) 
    9961289            sage: P._dangerous_dont_use_set_ents_lvls(range(5,-1,-1), [1,2,2,3,3,-1], range(3,-1,-1), [1,1,2,-1]) 
    9971290            sage: P 
     
    10231316        self.wd_percolate(start, end) 
    10241317 
    1025     cdef void wd_percolate(self, unsigned int start, unsigned int end): 
    1026         cdef unsigned int i, temp 
    1027         cdef unsigned int *self_wd_ents = self.wd_ents 
     1318    cdef void wd_percolate(self, int start, int end): 
     1319        cdef int i, temp 
     1320        cdef int *self_wd_ents = self.wd_ents 
    10281321        for i from end >= i > start: 
    10291322            if self_wd_ents[i] < self_wd_ents[i-1]: 
     
    10391332            sage: import sage.coding.binary_code 
    10401333            sage: from sage.coding.binary_code import * 
    1041             sage: P = PartitionStack(4, 6) 
     1334            sage: P = PartitionStack(2, 6) 
    10421335            sage: [P._split_column(i,i+1) for i in range(5)] 
    10431336            [0, 1, 2, 3, 4] 
     
    10531346            ({0},{1},{2},{3},{4,5}) 
    10541347            ({0},{1},{2},{3},{4},{5}) 
    1055             sage: P = PartitionStack(4, 6) 
     1348            sage: P = PartitionStack(2, 6) 
    10561349            sage: P._split_column(0,1) 
    10571350            0 
     
    10981391            sage: import sage.coding.binary_code 
    10991392            sage: from sage.coding.binary_code import * 
    1100             sage: P = PartitionStack(4, 6) 
     1393            sage: P = PartitionStack(2, 6) 
    11011394            sage: M = Matrix(GF(2), [[1,1,1,1,0,0],[0,0,1,1,1,1]]) 
    11021395            sage: B = BinaryCode(M) 
     
    11191412            ({0},{1},{2},{3},{4},{5}) 
    11201413            sage: P._col_degree(B, 2, 0, 2) 
    1121             2L 
     1414            2 
    11221415 
    11231416        """ 
    11241417        return self.col_degree(C, col, wd_ptr, k) 
    11251418 
    1126     cdef unsigned int col_degree(self, BinaryCode CG, int col, unsigned int wd_ptr, int k): 
    1127         cdef unsigned int i = 0 
    1128         cdef int *self_wd_lvls = self.wd_lvls 
     1419    cdef int col_degree(self, BinaryCode CG, int col, int wd_ptr, int k): 
     1420        cdef int i = 0 
     1421        cdef int *self_wd_lvls = self.wd_lvls, *self_wd_ents = self.wd_ents 
    11291422        col = self.col_ents[col] 
    11301423        while True: 
    1131             if CG.is_one(wd_ptr, col): i += 1 
     1424            if CG.is_one(self_wd_ents[wd_ptr], col): i += 1 
    11321425            if self_wd_lvls[wd_ptr] > k: wd_ptr += 1 
    11331426            else: break 
     
    11421435            sage: import sage.coding.binary_code 
    11431436            sage: from sage.coding.binary_code import * 
    1144             sage: P = PartitionStack(4, 6) 
     1437            sage: P = PartitionStack(2, 6) 
    11451438            sage: M = Matrix(GF(2), [[1,1,1,1,0,0],[0,0,1,1,1,1]]) 
    11461439            sage: B = BinaryCode(M) 
     
    11661459 
    11671460        """ 
    1168         return self.wd_degree(C, wd, col_ptr, k) 
    1169  
    1170     cdef int wd_degree(self, BinaryCode CG, unsigned int wd, int col_ptr, int k): 
    1171         cdef int i = 0 
    1172         wd = self.wd_ents[wd] 
     1461        cdef int *ham_wts = hamming_weights() 
     1462        result = self.wd_degree(C, wd, col_ptr, k, ham_wts) 
     1463        sage_free(ham_wts) 
     1464        return result 
     1465 
     1466    cdef int wd_degree(self, BinaryCode CG, int wd, int col_ptr, int k, int *ham_wts): 
     1467 
     1468        cdef int mask = (1 << col_ptr) 
    11731469        cdef int *self_col_lvls = self.col_lvls 
    1174         while True: 
    1175             if CG.is_one(wd, col_ptr): i += 1 
    1176             if self_col_lvls[col_ptr] > k: col_ptr += 1 
    1177             else: break 
    1178         return i 
     1470        while self_col_lvls[col_ptr] > k: 
     1471            col_ptr += 1 
     1472            mask += (1 << col_ptr) 
     1473        mask &= CG.words[self.wd_ents[wd]] 
     1474        return ham_wts[mask & 65535] + ham_wts[(mask >> 16) & 65535] 
    11791475 
    11801476    def _sort_cols(self, start, degrees, k): 
     
    11901486            sage: import sage.coding.binary_code 
    11911487            sage: from sage.coding.binary_code import * 
    1192             sage: P = PartitionStack(4, 6) 
     1488            sage: P = PartitionStack(2, 6) 
    11931489            sage: [P._split_column(i,i+1) for i in range(5)] 
    11941490            [0, 1, 2, 3, 4] 
     
    12151511    cdef int sort_cols(self, int start, int k): 
    12161512        cdef int i, j, max, max_location, self_ncols = self.ncols 
     1513        cdef int self_nwords = self.nwords, ii 
    12171514        cdef int *self_col_counts = self.col_counts 
    12181515        cdef int *self_col_lvls = self.col_lvls 
    1219         cdef unsigned int *self_col_degs = self.col_degs 
     1516        cdef int *self_col_degs = self.col_degs 
    12201517        cdef int *self_col_ents = self.col_ents 
    12211518        cdef int *self_col_output = self.col_output 
    1222         for i from 0 <= i < self_ncols: 
    1223             self_col_counts[i] = 0 
     1519        for ii from 0 <= ii < self_nwords: 
     1520            self_col_counts[ii] = 0 
    12241521        i = 0 
    12251522        while self_col_lvls[i+start] > k: 
     
    12311528        max = self_col_counts[0] 
    12321529        max_location = 0 
    1233         for j from 0 < j < self_ncols: 
    1234             if self_col_counts[j] > max: 
    1235                 max = self_col_counts[j] 
    1236                 max_location = j 
    1237             self_col_counts[j] += self_col_counts[j-1] 
     1530        for ii from 0 < ii < self_nwords: 
     1531            if self_col_counts[ii] > max: 
     1532                max = self_col_counts[ii] 
     1533                max_location = ii 
     1534            self_col_counts[ii] += self_col_counts[ii-1] 
    12381535 
    12391536        for j from i >= j >= 0: 
     
    12461543            self_col_ents[start+j] = self_col_output[j] 
    12471544 
    1248         j = 1 
    1249         while j < self_ncols and self_col_counts[j] <= i: 
    1250             if self_col_counts[j] > 0: 
    1251                 self_col_lvls[start + self_col_counts[j] - 1] = k 
    1252             self.col_percolate(start + self_col_counts[j-1], start + self_col_counts[j] - 1) 
    1253             j += 1 
     1545        ii = 1 
     1546        while ii < self_nwords and self_col_counts[ii] <= i: 
     1547            if self_col_counts[ii] > 0: 
     1548                self_col_lvls[start + self_col_counts[ii] - 1] = k 
     1549            self.col_percolate(start + self_col_counts[ii-1], start + self_col_counts[ii] - 1) 
     1550            ii += 1 
    12541551 
    12551552        return max_location 
     
    12671564            sage: import sage.coding.binary_code 
    12681565            sage: from sage.coding.binary_code import * 
    1269             sage: P = PartitionStack(8, 6) 
     1566            sage: P = PartitionStack(3, 6) 
    12701567            sage: P._sort_wds(0, [0,0,3,3,3,3,2,2], 1) 
    1271             4L 
     1568            4 
    12721569            sage: P 
    12731570            ({0,1,6,7,2,3,4,5}) 
     
    12921589        return self.sort_wds(start, k) 
    12931590 
    1294     cdef unsigned int sort_wds(self, unsigned int start, int k): 
    1295         cdef unsigned int i, j, max, max_location, self_nwords = self.nwords 
    1296         cdef unsigned int *self_wd_counts = self.wd_counts 
     1591    cdef int sort_wds(self, int start, int k): 
     1592        cdef int i, j, max, max_location, self_nwords = self.nwords 
     1593        cdef int ii, self_ncols = self.ncols 
     1594        cdef int *self_wd_counts = self.wd_counts 
    12971595        cdef int *self_wd_lvls = self.wd_lvls 
    12981596        cdef int *self_wd_degs = self.wd_degs 
    1299         cdef unsigned int *self_wd_ents = self.wd_ents 
    1300         cdef unsigned int *self_wd_output = self.wd_output 
    1301  
    1302         for i from 0 <= i < self_nwords: 
    1303             self_wd_counts[i] = 0 
     1597        cdef int *self_wd_ents = self.wd_ents 
     1598        cdef int *self_wd_output = self.wd_output 
     1599 
     1600        for ii from 0 <= ii < self_ncols+1: 
     1601            self_wd_counts[ii] = 0 
    13041602        i = 0 
    13051603        while self_wd_lvls[i+start] > k: 
     
    13111609        max = self_wd_counts[0] 
    13121610        max_location = 0 
    1313         for j from 0 < j < self_nwords: 
    1314             if self_wd_counts[j] > max: 
    1315                 max = self_wd_counts[j] 
    1316                 max_location = j 
    1317             self_wd_counts[j] += self_wd_counts[j-1] 
     1611        for ii from 0 < ii < self_ncols+1: 
     1612            if self_wd_counts[ii] > max: 
     1613                max = self_wd_counts[ii] 
     1614                max_location = ii 
     1615            self_wd_counts[ii] += self_wd_counts[ii-1] 
    13181616 
    13191617        for j from i >= j >= 0: 
    1320             if j > i: break # cython bug with unsigned ints... 
     1618            if j > i: break # cython bug with ints... 
    13211619            self_wd_counts[self_wd_degs[j]] -= 1 
    13221620            self_wd_output[self_wd_counts[self_wd_degs[j]]] = self_wd_ents[start+j] 
     
    13271625            self_wd_ents[start+j] = self_wd_output[j] 
    13281626 
    1329         j = 1 
    1330         while j < self_nwords and self_wd_counts[j] <= i: 
    1331             if self_wd_counts[j] > 0: 
    1332                 self_wd_lvls[start + self_wd_counts[j] - 1] = k 
    1333             self.wd_percolate(start + self_wd_counts[j-1], start + self_wd_counts[j] - 1) 
    1334             j += 1 
     1627        ii = 1 
     1628        while ii < self_ncols+1 and self_wd_counts[ii] <= i: 
     1629            if self_wd_counts[ii] > 0: 
     1630                self_wd_lvls[start + self_wd_counts[ii] - 1] = k 
     1631            self.wd_percolate(start + self_wd_counts[ii-1], start + self_wd_counts[ii] - 1) 
     1632            ii += 1 
    13351633 
    13361634        return max_location 
     1635 
     1636    def _refine(self, k, alpha, CG): 
     1637        """ 
     1638        EXAMPLE: 
     1639 
     1640            sage: import sage.coding.binary_code 
     1641            sage: from sage.coding.binary_code import * 
     1642            sage: M = Matrix(GF(2), [[1,1,1,1,0,0,0,0],[0,0,1,1,1,1,0,0],[0,0,0,0,1,1,1,1],[1,0,1,0,1,0,1,0]]) 
     1643            sage: B = BinaryCode(M) 
     1644            sage: P = PartitionStack(4, 8) 
     1645            sage: P._refine(1, [[0,0],[1,0]], B) 
     1646            304 
     1647            sage: P._split_column(0, 2) 
     1648            0 
     1649            sage: P._refine(2, [[0,0]], B) 
     1650            346 
     1651            sage: P._split_column(1, 3) 
     1652            1 
     1653            sage: P._refine(3, [[0,1]], B) 
     1654            558 
     1655            sage: P._split_column(2, 4) 
     1656            2 
     1657            sage: P._refine(4, [[0,2]], B) 
     1658            1713 
     1659            sage: P._split_column(3, 5) 
     1660            3 
     1661            sage: P._refine(5, [[0,3]], B) 
     1662            641 
     1663            sage: P._split_column(4, 6) 
     1664            4 
     1665            sage: P._refine(6, [[0,4]], B) 
     1666            1609 
     1667            sage: P._is_discrete(5) 
     1668            0 
     1669            sage: P._is_discrete(6) 
     1670            1 
     1671            sage: P 
     1672            ({0,4,6,2,13,9,11,15,10,14,12,8,7,3,1,5}) 
     1673            ({0},{4,6,2,13,9,11,15,10,14,12,8,7,3,1},{5}) 
     1674            ({0},{4,6,2,13,9,11,15},{10,14,12,8,7,3,1},{5}) 
     1675            ({0},{4,6,2},{13,9,11,15},{10,14,12,8},{7,3,1},{5}) 
     1676            ({0},{4},{6,2},{13,9},{11,15},{10,14},{12,8},{7,3},{1},{5}) 
     1677            ({0},{4},{6,2},{13,9},{11,15},{10,14},{12,8},{7,3},{1},{5}) 
     1678            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1679            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1680            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1681            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1682            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1683            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1684            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1685            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1686            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1687            ({0},{4},{6},{2},{13},{9},{11},{15},{10},{14},{12},{8},{7},{3},{1},{5}) 
     1688            ({0,1,2,3,4,7,6,5}) 
     1689            ({0,1,2,3,4,7,6,5}) 
     1690            ({0},{1,2,3,4,7,6,5}) 
     1691            ({0},{1},{2,3,4,7,6,5}) 
     1692            ({0},{1},{2},{3,4,7,6,5}) 
     1693            ({0},{1},{2},{3},{4,7,6,5}) 
     1694            ({0},{1},{2},{3},{4},{7},{6},{5}) 
     1695            ({0},{1},{2},{3},{4},{7},{6},{5}) 
     1696 
     1697        """ 
     1698        cdef int i, alpha_length = len(alpha) 
     1699        cdef int *_alpha = <int *> sage_malloc( (self.nwords + self.ncols) * sizeof(int) ) 
     1700        cdef int *ham_wts = hamming_weights() 
     1701        if not _alpha: 
     1702            sage_free(_alpha) 
     1703            raise MemoryError("Memory.") 
     1704        for i from 0 <= i < alpha_length: 
     1705            if alpha[i][0]: 
     1706                _alpha[i] = alpha[i][1] ^ self.flag 
     1707            else: 
     1708                _alpha[i] = alpha[i][1] 
     1709        result = self.refine(k, _alpha, alpha_length, CG, ham_wts) 
     1710        sage_free(_alpha) 
     1711        sage_free(ham_wts) 
     1712        return result 
     1713 
     1714    cdef int refine(self, int k, int *alpha, int alpha_length, BinaryCode CG, int *ham_wts): 
     1715        cdef int q, r, s, t, flag = self.flag, self_ncols = self.ncols 
     1716        cdef int t_w, self_nwords = self.nwords, invariant = 0, i, j, m = 0 
     1717        cdef int *self_wd_degs = self.wd_degs, *self_wd_lvls = self.wd_lvls, *self_col_lvls = self.col_lvls 
     1718        cdef int *self_col_degs = self.col_degs 
     1719        while not self.is_discrete(k) and m < alpha_length: 
     1720#            print "m:", m 
     1721#            print "alpha:", ','.join(['w'+str(alpha[i]^flag) if alpha[i]&flag else 'c'+str(alpha[i]) for i from 0 <= i < alpha_length]) 
     1722            invariant += 1 
     1723            j = 0 
     1724            if alpha[m] & flag: 
     1725                while j < self_ncols: 
     1726                    i = j; s = 0 
     1727                    invariant += 8 
     1728                    while True: 
     1729                        self_col_degs[i-j] = self.col_degree(CG, i, alpha[m]^flag, k) 
     1730                        if s == 0 and self_col_degs[i-j] != self_col_degs[0]: s = 1 
     1731                        i += 1 
     1732                        if self_col_lvls[i-1] <= k: break 
     1733                    if s: 
     1734                        invariant += 8 
     1735                        t = self.sort_cols(j, k) 
     1736                        invariant += t + self_col_degs[i-j-1] 
     1737                        q = m 
     1738                        while q < alpha_length: 
     1739                            if alpha[q] == j: 
     1740                                alpha[q] = t 
     1741                                break 
     1742                            q += 1 
     1743                        r = j 
     1744                        while True: 
     1745                            if r == 0 or self.col_lvls[r-1] == k: 
     1746                                if r != t: 
     1747                                    alpha[alpha_length] = r 
     1748                                    alpha_length += 1 
     1749                            r += 1 
     1750                            if r >= i: break 
     1751                        while self_col_lvls[j] > k: 
     1752                            j += 1 
     1753                        j += 1 
     1754                        invariant += (i-j) 
     1755                    else: j = i 
     1756            else: 
     1757                while j < self.nwords: 
     1758                    i = j; s = 0 
     1759                    invariant += 64 
     1760                    while True: 
     1761                        self_wd_degs[i-j] = self.wd_degree(CG, i, alpha[m], k, ham_wts) 
     1762                        if s == 0 and self_wd_degs[i-j] != self_wd_degs[0]: s = 1 
     1763                        i += 1 
     1764                        if self_wd_lvls[i-1] <= k: break 
     1765                    if s: 
     1766                        invariant += 64 
     1767                        t_w = self.sort_wds(j, k) 
     1768                        invariant += t_w + self_wd_degs[i-j-1] 
     1769                        q = m 
     1770                        j ^= flag 
     1771                        while q < alpha_length: 
     1772                            if alpha[q] == j: 
     1773                                alpha[q] = t_w ^ flag 
     1774                                break 
     1775                            q += 1 
     1776                        j ^= flag 
     1777                        r = j 
     1778                        while True: 
     1779                            if r == 0 or self.wd_lvls[r-1] == k: 
     1780                                if r != t_w: 
     1781                                    alpha[alpha_length] = r^flag 
     1782                                    alpha_length += 1 
     1783                            r += 1 
     1784                            if r >= i: break 
     1785                        while self_wd_lvls[j] > k: 
     1786                            j += 1 
     1787                        j += 1 
     1788                        invariant += (i-j) 
     1789                    else: j = i 
     1790            m += 1 
     1791        return invariant 
     1792 
     1793    def _clear(self, k): 
     1794        """ 
     1795        EXAMPLE: 
     1796            sage: import sage.coding.binary_code 
     1797            sage: from sage.coding.binary_code import * 
     1798            sage: P = PartitionStack(2, 6) 
     1799            sage: [P._split_column(i,i+1) for i in range(5)] 
     1800            [0, 1, 2, 3, 4] 
     1801            sage: P 
     1802            ({0,1,2,3}) 
     1803            ({0,1,2,3}) 
     1804            ({0,1,2,3}) 
     1805            ({0,1,2,3}) 
     1806            ({0,1,2,3,4,5}) 
     1807            ({0},{1,2,3,4,5}) 
     1808            ({0},{1},{2,3,4,5}) 
     1809            ({0},{1},{2},{3,4,5}) 
     1810            ({0},{1},{2},{3},{4,5}) 
     1811            ({0},{1},{2},{3},{4},{5}) 
     1812            sage: P._clear(2) 
     1813            sage: P 
     1814            ({0,1,2,3}) 
     1815            ({0,1,2,3}) 
     1816            ({0,1,2,3}) 
     1817            ({0,1,2,3}) 
     1818            ({0,1,2,3,4,5}) 
     1819            ({0},{1,2,3,4,5}) 
     1820            ({0},{1,2,3,4,5}) 
     1821            ({0},{1},{2,3,4,5}) 
     1822            ({0},{1},{2},{3,4,5}) 
     1823            ({0},{1},{2},{3},{4,5}) 
     1824 
     1825        """ 
     1826        self.clear(k) 
     1827 
     1828    cdef void clear(self, int k): 
     1829        cdef int i = 0, j = 0 
     1830        cdef int *wd_lvls = self.wd_lvls, *col_lvls = self.col_lvls 
     1831        while wd_lvls[i] != -1: 
     1832            if wd_lvls[i] >= k: 
     1833                wd_lvls[i] += 1 
     1834            if wd_lvls[i] < k: 
     1835                self.wd_percolate(j, i) 
     1836                j = i + 1 
     1837            i+=1 
     1838        i = 0 
     1839        j = 0 
     1840        while col_lvls[i] != -1: 
     1841            if col_lvls[i] >= k: 
     1842                col_lvls[i] += 1 
     1843            if col_lvls[i] < k: 
     1844                self.col_percolate(j, i) 
     1845                j = i + 1 
     1846            i+=1 
     1847 
     1848    def _cmp(self, other, C): 
     1849        """ 
     1850        EXAMPLE: 
     1851            sage: import sage.coding.binary_code 
     1852            sage: from sage.coding.binary_code import * 
     1853            sage: M = Matrix(GF(2), [[1,1,1,1,0,0,0,0],[0,0,1,1,1,1,0,0],[0,0,0,0,1,1,1,1],[1,0,1,0,1,0,1,0]]) 
     1854            sage: B = BinaryCode(M) 
     1855            sage: P = PartitionStack(4, 8) 
     1856            sage: P._refine(0, [[0,0],[1,0]], B) 
     1857            304 
     1858            sage: P._split_column(0, 1) 
     1859            0 
     1860            sage: P._refine(1, [[0,0]], B) 
     1861            346 
     1862            sage: P._split_column(1, 2) 
     1863            1 
     1864            sage: P._refine(2, [[0,1]], B) 
     1865            558 
     1866            sage: P._split_column(2, 3) 
     1867            2 
     1868            sage: P._refine(3, [[0,2]], B) 
     1869            1713 
     1870            sage: P._split_column(4, 4) 
     1871            4 
     1872            sage: P._refine(4, [[0,4]], B) 
     1873            1609 
     1874            sage: P._is_discrete(4) 
     1875            1 
     1876            sage: Q = PartitionStack(P) 
     1877            sage: Q._clear(4) 
     1878            sage: Q._split_column(5, 4) 
     1879            4 
     1880            sage: Q._refine(4, [[0,4]], B) 
     1881            1609 
     1882            sage: Q._is_discrete(4) 
     1883            1 
     1884            sage: Q._cmp(P, B) 
     1885            1 
     1886 
     1887        """ 
     1888        return self.cmp(other, C) 
     1889 
     1890    cdef int cmp(self, PartitionStack other, BinaryCode CG): 
     1891        cdef int *cur_span = self.col_counts # grab spare scratch space, size self.nwords. 
     1892        cdef int *self_wd_ents = self.wd_ents 
     1893        cdef int i, j, k, l, m, word, span = 1, ncols = self.ncols, nwords = self.nwords 
     1894        cur_span[0] = 0 
     1895        for i from 0 <= i < CG.nwords/2 + 1: 
     1896            word = self_wd_ents[i] 
     1897            k = 0 
     1898            while k < span and word != cur_span[k]: 
     1899                k += 1 
     1900            if k == span: 
     1901                if (span << 1) != nwords: 
     1902                    for j from 0 <= j < span: 
     1903                        cur_span[span+j] = cur_span[j] ^ word 
     1904                span << 1 
     1905                for j from 0 <= j < ncols: 
     1906                    l = CG.is_one(self.wd_ents[i], self.col_ents[j]) 
     1907                    m = CG.is_one(other.wd_ents[i], other.col_ents[j]) 
     1908                    if l != m: 
     1909                        return l - m 
     1910                if span == nwords: break 
     1911        return 0 
     1912 
     1913    def print_basis(self): 
     1914        """ 
     1915        EXAMPLE: 
     1916            sage: import sage.coding.binary_code 
     1917            sage: from sage.coding.binary_code import * 
     1918            sage: P = PartitionStack(4, 8) 
     1919            sage: P._dangerous_dont_use_set_ents_lvls(range(8), range(7)+[-1], [4,7,12,11,1,9,3,0,2,5,6,8,10,13,14,15], [0]*16) 
     1920            sage: P 
     1921            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1922            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1923            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1924            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1925            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1926            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1927            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1928            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1929            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1930            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1931            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1932            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1933            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1934            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1935            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1936            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1937            ({0},{1,2,3,4,5,6,7}) 
     1938            ({0},{1},{2,3,4,5,6,7}) 
     1939            ({0},{1},{2},{3,4,5,6,7}) 
     1940            ({0},{1},{2},{3},{4,5,6,7}) 
     1941            ({0},{1},{2},{3},{4},{5,6,7}) 
     1942            ({0},{1},{2},{3},{4},{5},{6,7}) 
     1943            ({0},{1},{2},{3},{4},{5},{6},{7}) 
     1944            ({0},{1},{2},{3},{4},{5},{6},{7}) 
     1945            sage: P._find_basis() 
     1946            sage: P.print_basis() 
     1947            basis_locations: 
     1948            4 
     1949            8 
     1950            0 
     1951            11 
     1952 
     1953        """ 
     1954        cdef int i, j 
     1955        if self.basis_locations: 
     1956            print "basis_locations:" 
     1957            j = 1 
     1958            while (1 << j) < self.nwords: 
     1959                j += 1 
     1960            for i from 0 <= i < j: 
     1961                print self.basis_locations[i] 
     1962 
     1963    def _find_basis(self): 
     1964        """ 
     1965        EXAMPLE: 
     1966            sage: import sage.coding.binary_code 
     1967            sage: from sage.coding.binary_code import * 
     1968            sage: P = PartitionStack(4, 8) 
     1969            sage: P._dangerous_dont_use_set_ents_lvls(range(8), range(7)+[-1], [4,7,12,11,1,9,3,0,2,5,6,8,10,13,14,15], [0]*16) 
     1970            sage: P 
     1971            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1972            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1973            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1974            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1975            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1976            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1977            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1978            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1979            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1980            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1981            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1982            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1983            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1984            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1985            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1986            ({4},{7},{12},{11},{1},{9},{3},{0},{2},{5},{6},{8},{10},{13},{14},{15}) 
     1987            ({0},{1,2,3,4,5,6,7}) 
     1988            ({0},{1},{2,3,4,5,6,7}) 
     1989            ({0},{1},{2},{3,4,5,6,7}) 
     1990            ({0},{1},{2},{3},{4,5,6,7}) 
     1991            ({0},{1},{2},{3},{4},{5,6,7}) 
     1992            ({0},{1},{2},{3},{4},{5},{6,7}) 
     1993            ({0},{1},{2},{3},{4},{5},{6},{7}) 
     1994            ({0},{1},{2},{3},{4},{5},{6},{7}) 
     1995            sage: P._find_basis() 
     1996            sage: P.print_basis() 
     1997            basis_locations: 
     1998            4 
     1999            8 
     2000            0 
     2001            11 
     2002 
     2003        """ 
     2004        cdef int i 
     2005        cdef int *ham_wts = hamming_weights() 
     2006        self.find_basis(ham_wts) 
     2007        sage_free(ham_wts) 
     2008 
     2009    cdef void find_basis(self, int *ham_wts): 
     2010        cdef int i = 0, j, k, nwords = self.nwords, weight, basis_elts = 0, nrows = self.nrows 
     2011        cdef int *self_wd_ents = self.wd_ents 
     2012        if not self.basis_locations: 
     2013            self.basis_locations = <int *> sage_malloc( nrows * sizeof(int) ) 
     2014        if not self.basis_locations: 
     2015            raise MemoryError("Memory.") 
     2016        while i < nwords: 
     2017            j = self_wd_ents[i] 
     2018            weight = ham_wts[j & 65535] + ham_wts[(j>>16) & 65535] 
     2019            if weight == 1: 
     2020                basis_elts += 1 
     2021                k = 0 
     2022                while not (1<<k) & j: 
     2023                    k += 1 
     2024                self.basis_locations[k] = i 
     2025                if basis_elts == nrows: break 
     2026            i += 1 
     2027 
     2028    def _get_permutation(self, other): 
     2029        """ 
     2030        EXAMPLE: 
     2031            sage: import sage.coding.binary_code 
     2032            sage: from sage.coding.binary_code import * 
     2033            sage: M = Matrix(GF(2), [[1,1,1,1,0,0,0,0],[0,0,1,1,1,1,0,0],[0,0,0,0,1,1,1,1],[1,0,1,0,1,0,1,0]]) 
     2034            sage: B = BinaryCode(M) 
     2035            sage: P = PartitionStack(4, 8) 
     2036            sage: P._refine(0, [[0,0],[1,0]], B) 
     2037            304 
     2038            sage: P._split_column(0, 1) 
     2039            0 
     2040            sage: P._refine(1, [[0,0]], B) 
     2041            346 
     2042            sage: P._split_column(1, 2) 
     2043            1 
     2044            sage: P._refine(2, [[0,1]], B) 
     2045            558 
     2046            sage: P._split_column(2, 3) 
     2047            2 
     2048            sage: P._refine(3, [[0,2]], B) 
     2049            1713 
     2050            sage: P._split_column(4, 4) 
     2051            4 
     2052            sage: P._refine(4, [[0,4]], B) 
     2053            1609 
     2054            sage: P._is_discrete(4) 
     2055            1 
     2056            sage: Q = PartitionStack(P) 
     2057            sage: Q._clear(4) 
     2058            sage: Q._split_column(5, 4) 
     2059            4 
     2060            sage: Q._refine(4, [[0,4]], B) 
     2061            1609 
     2062            sage: Q._is_discrete(4) 
     2063            1 
     2064            sage: P._get_permutation(Q) 
     2065            ([1, 2, 4, 8], [0, 1, 2, 3, 5, 4, 6, 7]) 
     2066 
     2067        """ 
     2068        cdef int i 
     2069        cdef int *ham_wts = hamming_weights() 
     2070        cdef int *basis_g = <int *> sage_malloc( self.nrows * sizeof(int) ) 
     2071        cdef int *col_g = <int *> sage_malloc( self.ncols * sizeof(int) ) 
     2072        if not (basis_g and col_g): 
     2073            if basis_g: sage_free(basis_g) 
     2074            if col_g: sage_free(col_g) 
     2075            raise MemoryError("Memory.") 
     2076        self.get_permutation(other, basis_g, col_g, ham_wts) 
     2077        sage_free(ham_wts) 
     2078        basis_l = [basis_g[i] for i from 0 <= i < self.nrows] 
     2079        col_l = [col_g[i] for i from 0 <= i < self.ncols] 
     2080        sage_free(basis_g) 
     2081        sage_free(col_g) 
     2082        return basis_l, col_l 
     2083 
     2084    cdef void get_permutation(self, PartitionStack other, int *basis_gamma, int *col_gamma, int *ham_wts): 
     2085        cdef int i 
     2086        cdef int *bas_loc, *self_wd_ents = self.wd_ents, *self_col_ents = self.col_ents, *other_col_ents = other.col_ents 
     2087        if not other.basis_locations: 
     2088            other.find_basis(ham_wts) 
     2089        bas_loc = other.basis_locations 
     2090        # basis_gamma[i] := image of the ith row as linear comb of rows 
     2091        for i from 0 <= i < self.nrows: 
     2092            basis_gamma[i] = self_wd_ents[bas_loc[i]] 
     2093        for i from 0 <= i < self.ncols: 
     2094            col_gamma[other_col_ents[i]] = self_col_ents[i] 
    13372095 
    13382096################################################################################ 
     
    13402098################################################################################ 
    13412099 
    1342     def _refine(self, k, col_alpha, wd_alpha, Code): 
    1343         """ 
    1344         Refines the partition at level k, using the list of cells alpha, and Code. 
    1345  
    1346         EXAMPLE: 
    1347  
    1348         """ 
    1349         cdef unsigned int i 
    1350         cdef int j 
    1351         cdef unsigned int *_col_a = <unsigned int *> sage_malloc(4*self.ncols*sizeof(unsigned int)) 
    1352         cdef int *_wd_a = <int *> sage_malloc(4*self.nwords*sizeof(int)) 
    1353         if not _col_a or not _wd_a: 
    1354             if _col_a: sage_free(_col_a) 
    1355             if _wd_a: sage_free(_wd_a) 
    1356             raise MemoryError("Memory.") 
    1357     # TODO       for i from  TODO 
    1358         result = self.refine(k, _col_a, _wd_a, Code) 
    1359         sage_free(_col_a) 
    1360         sage_free(_wd_a) 
    1361         return result 
    1362  
    1363     cdef unsigned int refine(self, int k, unsigned int *col_alpha, int *wd_alpha, BinaryCode CG): 
    1364         cdef int m = 0, j 
    1365         cdef int i, q, r, s, t 
    1366         cdef unsigned int t_w 
    1367         cdef unsigned int invariant 
    1368         cdef unsigned int *col_degrees = col_alpha + self.ncols 
    1369         cdef int *wd_degrees = wd_alpha + self.nwords 
    1370         while not self.is_discrete(k) and col_alpha[m] != -1: 
    1371             invariant += 1 
    1372             j = 0 
    1373             while j < self.nwords: 
    1374                 invariant += 8 
    1375                 i = j; s = 0 
    1376                 while True: 
    1377                     wd_degrees[i-j] = self.wd_degree(CG, i, col_alpha[m], k) 
    1378                     if s == 0 and wd_degrees[i-j] != wd_degrees[0]: s = 1 
    1379                     i += 1 
    1380                     if self.wd_lvls[i-1] <= k: break 
    1381                 if s: 
    1382                     invariant += 8 
    1383 #FIX                 t_w = self.sort_wds(j, wd_degrees, k) 
    1384                     invariant += t_w + wd_degrees[i-j-1] 
    1385                     q = m 
    1386                     while col_alpha[q] != -1: 
    1387                         if col_alpha[q] == j: col_alpha[q] = t_w 
    1388                         q += 1 
    1389                     r = j 
    1390                     while True: 
    1391                         if r == 0 or self.wd_lvls[r-1] == k: 
    1392                             if r != t_w: 
    1393                                 col_alpha[q] = r 
    1394                                 q += 1 
    1395                         r += 1 
    1396                         if r >= i: break 
    1397                     col_alpha[q] = -1 
    1398                     while self.wd_lvls[j] > k: 
    1399                         j += 1 
    1400                     j += 1 
    1401                     invariant += (i - j) 
    1402                 else: j = i 
    1403             m += 1 
    1404         m = 0 
    1405         while not self.is_discrete(k) and wd_alpha[m] != -1: 
    1406             invariant += 1 
    1407             j = 0 
    1408             while j < self.ncols: 
    1409                 invariant += 8 
    1410                 i = j; s = 0 
    1411                 while True: 
    1412                     col_degrees[i-j] = self.col_degree(CG, i, wd_alpha[m], k) 
    1413                     if s == 0 and col_degrees[i-j] != col_degrees[0]: s = 1 
    1414                     i += 1 
    1415                     if self.col_lvls[i-1] <= k: break 
    1416                 if s: 
    1417                     invariant += 8 
    1418 #FIX                 t = self.sort_cols(j, col_degrees, k) 
    1419                     invariant += t + col_degrees[i-j-1] 
    1420                     q = m 
    1421                     while wd_alpha[q] != -1: 
    1422                         if wd_alpha[q] == j: wd_alpha[q] = t 
    1423                         q += 1 
    1424                     r = j 
    1425                     while True: 
    1426                         if r == 0 or self.col_lvls[r-1] == k: 
    1427                             if r != t: 
    1428                                 wd_alpha[q] = r 
    1429                                 q += 1 
    1430                         r += 1 
    1431                         if r >= i: break 
    1432                     wd_alpha[q] = -1 
    1433                     while self.col_lvls[j] > k: 
    1434                         j += 1 
    1435                     j += 1 
    1436                     invariant += (i - j) 
    1437                 else: j = i 
    1438             m += 1 
    1439         return invariant 
    1440  
    1441     cdef void get_permutation(self, PartitionStack zeta, int *wd_gamma, int *col_gamma): 
    1442         cdef int i 
    1443         for i from 0 <= i < self.ncols: 
    1444             col_gamma[zeta.col_ents[i]] = self.col_ents[i] 
    1445         for i from 0 <= i < self.nwords: 
    1446             wd_gamma[zeta.wd_ents[i]] = self.wd_ents[i] 
    1447  
    1448     cdef int cmp(self, PartitionStack other, BinaryCode CG): 
    1449         # if CG(self) > CG(other): return 1 
    1450         # if CG(self) < CG(other): return -1 
    1451         # else: return 0 
    1452         cdef int i, j, k, l 
    1453         for i from 0 <= i < CG.nwords: 
    1454             for j from 0 <= j < CG.ncols: 
    1455                 k = CG.is_one(self.wd_ents[i], self.col_ents[j]) 
    1456                 l = CG.is_one(other.wd_ents[i], other.col_ents[j]) 
    1457                 if k != l: 
    1458                     return k - l 
    1459         return 0 
    1460  
    14612100cdef class BinaryCodeClassifier: 
    14622101 
    14632102    def __new__(self): 
    1464         cdef int *self_ham_wts 
    1465         self.ham_wts = <int *> sage_malloc( 65536 * sizeof(int) ) 
    1466         if not self.ham_wts: 
    1467             sage_free(self.ham_wts) 
    1468         self_ham_wts = self.ham_wts 
    1469         self_ham_wts[0]  = 0; self_ham_wts[1]  = 1; self_ham_wts[2]  = 1; self_ham_wts[3]  = 2 
    1470         self_ham_wts[4]  = 1; self_ham_wts[5]  = 2; self_ham_wts[6]  = 2; self_ham_wts[7]  = 3 
    1471         self_ham_wts[8]  = 1; self_ham_wts[9]  = 2; self_ham_wts[10] = 2; self_ham_wts[11] = 3 
    1472         self_ham_wts[12] = 2; self_ham_wts[13] = 3; self_ham_wts[14] = 3; self_ham_wts[15] = 4 
    1473         for i from 16 <= i < 256: 
    1474             self_ham_wts[i] = self_ham_wts[i & 15] + self_ham_wts[(i>>4) & 15] 
    1475         for i from 256 <= i < 65536: 
    1476             self_ham_wts[i] = self_ham_wts[i & 255] + self_ham_wts[(i>>8) & 255] 
    1477         # This may seem like overkill, but the worst case for storing the words 
    1478         # themselves is 65536- in this case, we are increasing memory usage by a 
    1479         # factor of 2. 
     2103        self.ham_wts = hamming_weights() 
    14802104 
    14812105    def __dealloc__(self): 
     
    15002124 
    15012125 
    1502 def classify(BinaryCode C, lab=True, verbosity=0): 
    1503     """ 
    1504     """ 
    1505     cdef int i, j # local variables 
    1506     cdef OrbitPartition Theta # keeps track of which vertices have been 
    1507                               # discovered to be equivalent 
    1508     cdef int index = 0, size = 1 
    1509     cdef int L = 100 
    1510     cdef int **Phi, **Omega 
    1511     cdef int l = -1 
    1512     cdef PartitionStack nu, zeta, rho 
    1513     cdef int k_rho 
    1514     cdef int k = 0 
    1515     cdef int h = -1 
    1516     cdef int hb 
    1517     cdef int hh = 1 
    1518     cdef int ht 
    1519     cdef mpz_t *Lambda_mpz, *zf_mpz, *zb_mpz 
    1520     cdef int hzf 
    1521     cdef int hzb = -1 
    1522     cdef unsigned int *basis_gamma 
    1523     cdef int *col_gamma 
    1524     cdef int *alpha 
    1525     cdef int *v 
    1526     cdef int *e 
    1527     cdef int state 
    1528     cdef int tvc, tvh, nwords = C.nwords, ncols = C.ncols, n = nwords + ncols, nrows = C.nrows 
    1529  
    1530     # trivial case 
    1531     if ncols == 0: 
    1532         return [], {} 
    1533     elif nwords == 0 and ncols == 1: 
    1534         return [], {0:0} 
    1535     elif nwords == 0: 
    1536         output1 = [] 
    1537         dd = {} 
    1538         for i from 0 <= i < ncols-1: 
    1539             dd[i] = i 
    1540             perm = range(ncols) 
    1541             perm[i] = i+1 
    1542             perm[i+1] = i 
    1543             output1.append(perm) 
    1544         dd[ncols-1] = ncols-1 
    1545         return output1, dd 
    1546  
    1547     # allocate int pointers 
    1548     Phi = <int **> sage_malloc(L * sizeof(int *)) 
    1549     Omega = <int **> sage_malloc(L * sizeof(int *)) 
    1550  
    1551     # allocate GMP int pointers 
    1552     Lambda_mpz = <mpz_t *> sage_malloc((ncols+2)*sizeof(mpz_t)) 
    1553     zf_mpz     = <mpz_t *> sage_malloc((ncols+2)*sizeof(mpz_t)) 
    1554     zb_mpz     = <mpz_t *> sage_malloc((ncols+2)*sizeof(mpz_t)) 
    1555  
    1556     # check for memory errors 
    1557     if not (Phi and Omega and Lambda_mpz and zf_mpz and zb_mpz): 
    1558         if Lambda_mpz: sage_free(Lambda_mpz) 
    1559         if zf_mpz: sage_free(zf_mpz) 
    1560         if zb_mpz: sage_free(zb_mpz) 
    1561         if Phi: sage_free(Phi) 
    1562         if Omega: sage_free(Omega) 
    1563         raise MemoryError("Error allocating memory.") 
    1564  
    1565     # allocate int arrays 
    1566     basis_gamma = <unsigned int *> sage_malloc(nrows*sizeof(unsigned int)) 
    1567     col_gamma = <int *> sage_malloc(ncols*sizeof(int)) 
    1568     Phi[0] = <int *> sage_malloc(L*ncols*sizeof(int)) 
    1569     Omega[0] = <int *> sage_malloc(L*ncols*sizeof(int)) 
    1570     alpha = <int *> sage_malloc(4*ncols*sizeof(int)) 
    1571     v = <int *> sage_malloc(ncols*sizeof(int)) 
    1572     e = <int *> sage_malloc(ncols*sizeof(int)) 
    1573  
    1574     # check for memory errors 
    1575     if not (basis_gamma and col_gamma and Phi[0] and Omega[0] and alpha and v and e): 
    1576         if basis_gamma: sage_free(basis_gamma) 
    1577         if col_gamma: sage_free(col_gamma) 
    1578         if Phi[0]: sage_free(Phi[0]) 
    1579         if Omega[0]: sage_free(Omega[0]) 
    1580         if alpha: sage_free(alpha) 
    1581         if v: sage_free(v) 
    1582         if e: sage_free(e) 
    1583         sage_free(Lambda_mpz) 
    1584         sage_free(zf_mpz) 
    1585         sage_free(zb_mpz) 
    1586         sage_free(Phi) 
    1587         sage_free(Omega) 
    1588         raise MemoryError("Error allocating memory.") 
    1589  
    1590     # setup double index arrays 
     2126#def classify(BinaryCode C, lab=True, verbosity=0): 
     2127#    """ 
     2128#    """ 
     2129#    cdef int i, j # local variables 
     2130#    cdef OrbitPartition Theta # keeps track of which vertices have been 
     2131#                              # discovered to be equivalent 
     2132#    cdef int index = 0, size = 1 
     2133#    cdef int L = 100 
     2134#    cdef int **Phi, **Omega 
     2135#    cdef int l = -1 
     2136#    cdef PartitionStack nu, zeta, rho 
     2137#    cdef int k_rho 
     2138#    cdef int k = 0 
     2139#    cdef int h = -1 
     2140#    cdef int hb 
     2141#    cdef int hh = 1 
     2142#    cdef int ht 
     2143#    cdef mpz_t *Lambda_mpz, *zf_mpz, *zb_mpz 
     2144#    cdef int hzf 
     2145#    cdef int hzb = -1 
     2146#    cdef int *basis_gamma 
     2147#    cdef int *col_gamma 
     2148#    cdef int *alpha 
     2149#    cdef int *v 
     2150#    cdef int *e 
     2151#    cdef int state 
     2152#    cdef int tvc, tvh, nwords = C.nwords, ncols = C.ncols, n = nwords + ncols, nrows = C.nrows 
     2153 
     2154#    # trivial case 
     2155#    if ncols == 0: 
     2156#        return [], {} 
     2157#    elif nwords == 0 and ncols == 1: 
     2158#        return [], {0:0} 
     2159#    elif nwords == 0: 
     2160#        output1 = [] 
     2161#        dd = {} 
     2162#        for i from 0 <= i < ncols-1: 
     2163#            dd[i] = i 
     2164#            perm = range(ncols) 
     2165#            perm[i] = i+1 
     2166#            perm[i+1] = i 
     2167#            output1.append(perm) 
     2168#        dd[ncols-1] = ncols-1 
     2169#        return output1, dd 
     2170 
     2171#    # allocate int pointers 
     2172#    Phi = <int **> sage_malloc(L * sizeof(int *)) 
     2173#    Omega = <int **> sage_malloc(L * sizeof(int *)) 
     2174 
     2175#    # allocate GMP int pointers 
     2176#    Lambda_mpz = <mpz_t *> sage_malloc((ncols+2)*sizeof(mpz_t)) 
     2177#    zf_mpz     = <mpz_t *> sage_malloc((ncols+2)*sizeof(mpz_t)) 
     2178#    zb_mpz     = <mpz_t *> sage_malloc((ncols+2)*sizeof(mpz_t)) 
     2179 
     2180#    # check for memory errors 
     2181#    if not (Phi and Omega and Lambda_mpz and zf_mpz and zb_mpz): 
     2182#        if Lambda_mpz: sage_free(Lambda_mpz) 
     2183#        if zf_mpz: sage_free(zf_mpz) 
     2184#        if zb_mpz: sage_free(zb_mpz) 
     2185#        if Phi: sage_free(Phi) 
     2186#        if Omega: sage_free(Omega) 
     2187#        raise MemoryError("Error allocating memory.") 
     2188 
     2189#    # allocate int arrays 
     2190#    basis_gamma = <int *> sage_malloc(nrows*sizeof(int)) 
     2191#    col_gamma = <int *> sage_malloc(ncols*sizeof(int)) 
     2192#    Phi[0] = <int *> sage_malloc(L*ncols*sizeof(int)) 
     2193#    Omega[0] = <int *> sage_malloc(L*ncols*sizeof(int)) 
     2194#    alpha = <int *> sage_malloc(4*ncols*sizeof(int)) 
     2195#    v = <int *> sage_malloc(ncols*sizeof(int)) 
     2196#    e = <int *> sage_malloc(ncols*sizeof(int)) 
     2197 
     2198#    # check for memory errors 
     2199#    if not (basis_gamma and col_gamma and Phi[0] and Omega[0] and alpha and v and e): 
     2200#        if basis_gamma: sage_free(basis_gamma) 
     2201#        if col_gamma: sage_free(col_gamma) 
     2202#        if Phi[0]: sage_free(Phi[0]) 
     2203#        if Omega[0]: sage_free(Omega[0]) 
     2204#        if alpha: sage_free(alpha) 
     2205#        if v: sage_free(v) 
     2206#        if e: sage_free(e) 
     2207#        sage_free(Lambda_mpz) 
     2208#        sage_free(zf_mpz) 
     2209#        sage_free(zb_mpz) 
     2210#        sage_free(Phi) 
     2211#        sage_free(Omega) 
     2212#        raise MemoryError("Error allocating memory.") 
     2213 
     2214#    # setup double index arrays 
     2215 
     2216 
     2217 
     2218 
     2219 
     2220 
     2221 
     2222 
     2223 
Note: See TracChangeset for help on using the changeset viewer.