Ticket #6452: trac_6452-ring-codes.patch

File trac_6452-ring-codes.patch, 18.0 KB (added by wdj, 10 years ago)

ignore the previous patch; this applies to 4.1.2.rc2

  • module_list.py

    # HG changeset patch
    # User David Joyner <wdjoyner@gmail.com>
    # Date 1255895149 14400
    # Node ID 36ee0da8e15a83a72633998fe18933b4a6072949
    # Parent  3c9b44598ccfcb526239bb0954bbe77d105614a6
    new ring codes patch - wdj
    
    diff -r 3c9b44598ccf -r 36ee0da8e15a module_list.py
    a b  
    154154    Extension('sage.coding.binary_code',
    155155              sources = ['sage/coding/binary_code.pyx']),
    156156
     157    Extension('sage.coding.ring_code',
     158              sources = ['sage/coding/ring_code.pyx']),
     159
    157160    ################################
    158161    ##
    159162    ## sage.combinat
  • sage/coding/all.py

    diff -r 3c9b44598ccf -r 36ee0da8e15a sage/coding/all.py
    a b  
    5555
    5656from sd_codes import self_dual_codes_binary
    5757
     58from ring_code import RingCode
  • new file sage/coding/ring_code.pyx

    diff -r 3c9b44598ccf -r 36ee0da8e15a sage/coding/ring_code.pyx
    - +  
     1"""
     2Codes over rings of the form ZZ/mZZ.
     3
     4This module constructs codes over rings of the form ZZ/mZZ, that is, submodules
     5of FreeModule(IntegerModRing(m), n).
     6
     7    AUTHORS:
     8        -- Cesar Agustin Garcia-Vazquez
     9        -- Carlos A. Lopez-Andrade
     10        -- David Joyner
     11
     12TODO:
     13  - Rewrite nearest ngbr decoder completely using pure cython with no GF(q) classes
     14    Coerce answer at last step.
     15  - Generalize this code to allow for other finite rings.
     16
     17"""
     18
     19    #*****************************************************************************
     20    #       Copyright (C) 2008 William Stein <wstein@gmail.com>
     21    #                     2008 Cesar A. Garcia-Vazquez <cesarnda@gmail.com>
     22    #                     2008 Carlos A. Lopez-Andrade <calopez@cs.buap.mx>
     23    #
     24    #  Distributed under the terms of the GNU General Public License (GPL),
     25    #  version 2 or later (at your option).
     26    #
     27    # This code is distributed in the hope that it will be useful,
     28    # but WITHOUT ANY WARRANTY; without even the implied warranty of
     29    # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
     30    # General Public License for more details.
     31    #
     32    #  The full text of the GPL is available at:
     33    #
     34    #                  http://www.gnu.org/licenses/
     35    #*****************************************************************************
     36
     37import sage.modules.free_module as fm
     38import sage.modules.module as module
     39from sage.rings.integer_mod_ring import IntegerModRing
     40
     41#cython
     42cdef extern from *:
     43    void* malloc(int)
     44    int memset(void*, int, int)
     45    void free(int*)
     46    double ceil( double )
     47    double floor( double )
     48    double pow( double, double)
     49
     50from sage.modules.free_module_element import vector
     51from sage.modules.free_module import FreeModule
     52
     53cdef class RingCode: #(module.Module):
     54    r"""
     55    This class generates code over a ring ZZ/mZZ generated by a "generator matrix".
     56
     57    CALL FORMATS:
     58        1. CS = RingCode(matrixGen)
     59
     60    INPUT:
     61        matrixGen -- matrix to be used as a generator matrix
     62
     63    OUTPUT:
     64        list -- list to be used as a set.
     65
     66    EXAMPLES:
     67        sage: M = Matrix(IntegerModRing(25), [[1, 1, 6],[1, 6, 1],[6, 1, 1]])
     68        sage: CS = RingCode(M)
     69        sage: CS
     70        (3, 625, 1)-code over the Ring of integers modulo 25
     71        sage: M = Matrix(IntegerModRing(25), [[6, 6, 6],[6, 6, 6],[6, 6, 6]])
     72        sage: CS = RingCode(M)
     73        sage: CS
     74        (3, 25, 3)-code over the Ring of integers modulo 25
     75        sage: M = Matrix(IntegerModRing(9), [[5,4,8,1],[1,5,4,8],[8,1,5,4],[4,8,1,5]])
     76        sage: CS = RingCode(M)
     77        sage: CS
     78        (4, 81, 2)-code over the Ring of integers modulo 9
     79        sage: M = Matrix(IntegerModRing(25), [[1, 0, 4],[0, 4, 1],[4, 1, 0]])
     80        sage: CS = RingCode(M)
     81        sage: CS
     82        (3, 3125, 1)-code over the Ring of integers modulo 25
     83        sage: M = Matrix(IntegerModRing(9), [[5,1,2,1],[1,5,1,2],[2,1,5,1],[1,2,1,5]])
     84        sage: CS = RingCode(M)
     85        sage: CS
     86        (4, 81, 2)-code over the Ring of integers modulo 9
     87        sage: M = Matrix(IntegerModRing(9), [[5,0,0,0,4],[4,5,0,0,0],[0,4,5,0,0],[0,0,4,5,0],[0,0,0,4,5]])
     88        sage: CS = RingCode(M)
     89        sage: CS
     90        (5, 6561, 2)-code over the Ring of integers modulo 9
     91       
     92    NOTES:
     93        This function uses a d-heap to store the codewords, this way, it is as fast as Magma.
     94
     95    AUTHORS:
     96        -- Cesar A. Garcia-Vazquez
     97        -- Carlos A. Lopez-Andrade
     98    """
     99    cdef int codewords      #Cardinality of the set
     100    cdef int size
     101    cdef int sizeCols       #number of columns of the matrix
     102    cdef int modulus       
     103    cdef int minimum        #minimum Hamming distance
     104    cdef int count          #to iterate on self
     105    cdef long sizeArray
     106    cdef int** heap_matrix
     107    cdef _base_ring, ambient_module
     108    cdef gen_matrix, gen_mat, minimum_weight_codeword, minimum_distance, mincodeword
     109    cdef next, spanning_codewords, characteristic, _length
     110    cdef list codeSet
     111
     112    def __new__(self, matrixGen):
     113        cdef int size = <int>matrixGen.nrows()
     114        cdef int sizeCols = <int>matrixGen.ncols()
     115        #cdef int length = sizeCols
     116        cdef int lamb = 2*size
     117        self.size = size
     118        self.sizeCols = sizeCols
     119        self._length = sizeCols
     120        self._base_ring = matrixGen.base_ring()
     121        self.modulus = <int>matrixGen.base_ring().order()
     122        if sizeCols <= size + lamb:
     123            self._use_heap_matrix(matrixGen)
     124        self.gen_matrix = matrixGen
     125
     126    cdef void _use_heap_matrix(self, matrixGen):
     127        cdef int modulus = self.modulus
     128        cdef int size = self.size
     129        cdef int sizeCols = self.sizeCols
     130        cdef int* vectorGen =  <int*>malloc(sizeof(int) * size)
     131        cdef int* newVector =  <int*>malloc(sizeof(int) * sizeCols)
     132       
     133        #Number of columns of the matrix
     134        cdef long sizeArray = (<int>pow(modulus,sizeCols) - 1)/(modulus - 1) - 1   
     135        self.sizeArray = sizeArray
     136        self.count = 0
     137        #Limit used to determinate how many codewords are going to be generated.
     138        cdef int limit = <int>pow(modulus,size)
     139       
     140        cdef int localLimit     #Local limit to modify the values of the vector that multiplies the matrixGen.
     141        cdef int localIndex     #Local index to determinate which index is going to change.
     142        cdef int value      #Takes the first element of the vector to check if the second has to change.
     143        cdef int element        #Element that must be in the matrixGen.
     144        cdef int head       #Second element of the vector generated by multiplying the matrixGen and the vector.
     145        cdef int index1     #Index used in some fors.
     146        cdef int index2     #Another index used in some fors
     147        cdef int index3
     148        cdef int mul
     149             
     150        cdef long index
     151        cdef int** heap_matrix = <int**>malloc(sizeof(int*)*modulus)
     152       
     153        for index1 from 0 <= index1 < modulus:
     154            heap_matrix[index1]= <int*>malloc(sizeof(int)*sizeArray)
     155            memset(heap_matrix[index1],0,sizeof(int)*sizeArray)
     156       
     157        cdef int** matrixGenerator = <int**>malloc(sizeof(int*)*size)
     158       
     159       
     160        for index1 from 0 <= index1 < size:
     161            matrixGenerator[index1]= <int*>malloc(sizeof(int)*sizeCols)
     162            memset(matrixGenerator[index1],0,sizeof(int)*sizeCols)
     163       
     164        for index1 from 0 <= index1 < size:
     165            for index2 from 0 <= index2 < sizeCols:
     166                matrixGenerator[index1][index2] = matrixGen[index1][index2]
     167       
     168        for index1 from 0 <= index1 < modulus:
     169            for index2 from 0 <= index2 < sizeArray:
     170                heap_matrix[index1][index2] = -1
     171       
     172        for index1 from 0 <= index1 < size:
     173            vectorGen[index1] = 0
     174       
     175        vectorGen[0] = -1
     176       
     177        #main for for obtaining the whole set
     178        for index from 0 <= index < limit:         
     179            localLimit = <int>ceil( <double>index / <double>modulus )
     180            value = vectorGen[0] + 1
     181            vectorGen[0] = value % modulus
     182            localIndex = 1
     183       
     184            for index2 from 0 <= index2 < localLimit:
     185                if value == modulus:
     186                    value = vectorGen[localIndex] + 1
     187                    vectorGen[localIndex] = value % modulus
     188                    localIndex = localIndex + 1
     189                else:
     190                    break
     191       
     192            for index1 from 0 <= index1 < sizeCols:
     193                newVector[index1] = 0
     194            # multiplying vector and matrix
     195            for index1 from  0 <= index1 < sizeCols:
     196                for index2  from 0 <= index2 < size:
     197                    mul = matrixGenerator[index2][index1] * vectorGen[index2]
     198                    newVector[index1] = newVector[index1] + mul
     199                newVector[index1] = newVector[index1] % modulus
     200            head = newVector[0]
     201            element = -1    #To start filling
     202            index1 = 1
     203       
     204            #Filling While
     205            while index1 < sizeCols:
     206                element = (element + 1)* modulus + newVector[index1]
     207                heap_matrix[head][element] = newVector[index1]
     208                index1 = index1 + 1
     209        #End of main for
     210       
     211        self.heap_matrix = heap_matrix
     212        free(vectorGen)
     213        self._remove_zero()
     214        self._generate_codewords()
     215       
     216    cdef _remove_zero(self):
     217        cdef int sizeCols = self.sizeCols
     218        cdef int modulus = self.modulus
     219        cdef int* vectorGen =  <int*>malloc(sizeof(int) * sizeCols)
     220        cdef int index
     221        cdef int element = -1
     222        cdef int head = 0
     223        index = 1
     224        while index < sizeCols:
     225            element = (element + 1)* modulus
     226            self.heap_matrix[head][element] = -1
     227            index = index + 1
     228       
     229    cdef _generate_codewords(self):
     230        cdef int index1
     231        cdef int index2
     232        cdef int index3
     233        cdef int index
     234        cdef int head
     235        cdef int codewords = 0
     236        cdef int sizeCols = self.sizeCols
     237        cdef int modulus = self.modulus
     238        cdef int limit = sizeCols - 1
     239        cdef long sizeArray = self.sizeArray
     240        cdef long beginArray = sizeArray - <int>pow(modulus,(sizeCols - 1))
     241        cdef long endArray = <long>pow(modulus,(sizeCols - 2))
     242        cdef long index2Start = beginArray
     243        cdef int minimum = sizeCols #minimum distance
     244        cdef int localMin    #temp minium distance
     245        cdef int* vectorGen =  <int*>malloc(sizeof(int) * sizeCols)
     246        cdef list codeSet = []
     247        #cdef int** heap_matrix = self.heap_matrix
     248        for index1 from 0 <= index1 < modulus:
     249            for index2 from 0 <= index2 < endArray:
     250                for index from 0 <= index < modulus:
     251                    if self.heap_matrix[index1][beginArray + index] != -1:
     252                        addVector = ([])
     253                        #mincodeword = ([])
     254                        codewords = codewords + 1
     255                        vectorGen[limit] = index
     256                        localMin = 0
     257                        if index != 0:
     258                            localMin = localMin + 1
     259                        head = beginArray
     260                        #print "Before searching for fathers"
     261                        for index3 from limit > index3 >= 1:
     262                            head = <int> floor( (head - 1) / modulus)
     263                            vectorGen[index3] = self.heap_matrix[index1][ head ]
     264                            if vectorGen[index3] != 0:
     265                                localMin = localMin + 1
     266                            head = head - vectorGen[index3]
     267                        vectorGen[0] = index1
     268                        if index1 != 0:
     269                            localMin = localMin + 1
     270                        #print "Every father found"
     271                        for index3 from 0 <= index3 < sizeCols:
     272                            addVector.append(vectorGen[index3])
     273                        addVector = vector(IntegerModRing(modulus),addVector)
     274                        codeSet.append(addVector)
     275                        #print "Every vector in the set"
     276                        if localMin < minimum:              # Compares the min dist with a probable new min dist
     277                            minimum = localMin              # A new minimum distance has been found
     278                            #print "A new minimum distance has been found", minimum, addVector
     279                            mincodeword = addVector
     280                beginArray = beginArray + modulus
     281            beginArray = index2Start
     282        #free(vectorGen)
     283        addVector = ([])
     284        for index3 from 0 <= index3 < sizeCols:
     285            addVector.append(0)
     286        codeSet.append(vector(IntegerModRing(modulus), addVector))
     287        self.minimum = minimum
     288        self.mincodeword = mincodeword
     289        self.codewords = codewords + 1
     290        self.codeSet = codeSet
     291
     292    def base_ring(self):
     293        """
     294        Returns base ring.
     295
     296        EXAMPLES:
     297            sage: M = Matrix(IntegerModRing(12), [[1, 1, 6],[1, 6, 1],[6, 1, 1]])
     298            sage: C = RingCode(M)
     299            sage: C.base_ring()
     300            Ring of integers modulo 12
     301
     302        """
     303        return self._base_ring
     304
     305    def ambient_module(self):
     306        """
     307        Returns ambient module.
     308
     309        EXAMPLES:
     310            sage: M = Matrix(IntegerModRing(12), [[1, 1, 6],[1, 6, 1],[6, 1, 1]])
     311            sage: C = RingCode(M)
     312            sage: C.ambient_module()
     313            Ambient free module of rank 3 over Ring of integers modulo 12
     314        """
     315        return FreeModule(self.base_ring(),self._length)
     316
     317    def characteristic(self):
     318        """
     319        Returns characteristic of base ring.
     320
     321        EXAMPLES:
     322            sage: M = Matrix(IntegerModRing(12), [[1, 1, 6],[1, 6, 1],[6, 1, 1]])
     323            sage: C = RingCode(M)
     324            sage: C.characteristic()
     325            12
     326        """
     327        return (self.base_ring()).characteristic()
     328
     329    def __contains__(self,v):
     330        C = self.codeSet
     331        return C.__contains__(v)
     332
     333    def gen_mat(self):
     334        """
     335        Returns a generator matrix.
     336
     337        EXAMPLES:
     338            sage: M = Matrix(IntegerModRing(12), [[1, 1, 6, -1],[1, 6, 1, 2],[6, 1, 1, 0]])
     339            sage: C = RingCode(M)
     340            sage: C.gen_mat()
     341            [ 1  1  6 11]
     342            [ 1  6  1  2]
     343            [ 6  1  1  0]
     344
     345        """
     346        return self.gen_matrix
     347
     348    def __getitem__(self, item):
     349        """
     350        Iterator method.
     351
     352        EXAMPLES:
     353            sage: M = Matrix(IntegerModRing(6), [[2, 2, 2],[2, 2, 2],[2, 2, 2]])
     354            sage: CS = RingCode(M); CS
     355            (3, 3, 3)-code over the Ring of integers modulo 6
     356            sage: [c for c in CS]
     357            [(2, 2, 2), (4, 4, 4), (0, 0, 0)]
     358        """
     359        return self.codeSet[item]
     360   
     361    def __latex__(self):
     362        s = 'Linear Code over the %s \n(%d, %d, %d)'%(self.base_ring, self.sizeCols, self.codewords, self.minimum)
     363        return s
     364
     365    def length(self):
     366        """
     367        Returns the usual length of code.
     368
     369        EXAMPLES:
     370            sage: M = Matrix(IntegerModRing(12), [[1, 1, 6],[1, 6, 1],[6, 1, 1]])
     371            sage: C = RingCode(M)
     372            sage: C.length()
     373            3
     374        """
     375        return self._length
     376
     377    def list(self):
     378        """
     379        Returns list of all codewords.
     380
     381        EXAMPLES:
     382            sage: M = Matrix(IntegerModRing(2), [[1, 1, 0],[1, 0, 1],[0, 1, 1]])
     383            sage: C = RingCode(M)
     384            sage: C
     385            (3, 4, 2)-code over the Ring of integers modulo 2
     386            sage: C.list()
     387            [(0, 1, 1), (1, 0, 1), (1, 1, 0), (0, 0, 0)]
     388        """
     389        return self.codeSet
     390
     391    def minimum_distance(self):
     392        """
     393        Returns the usual minimum distance.
     394
     395        EXAMPLES:
     396            sage: M = Matrix(IntegerModRing(12), [[0, 1, 6, -1],[1, 6, 1, 2],[6, 1, 1, 0]])
     397            sage: C = RingCode(M)
     398            sage: C
     399            (4, 1728, 2)-code over the Ring of integers modulo 12
     400            sage: C.minimum_distance()
     401            2
     402
     403        """
     404        return self.minimum
     405
     406    def minimum_weight_codeword(self):
     407        """
     408        Returns a non-zero codeword of minimum weight.
     409
     410        EXAMPLES:
     411            sage: M = Matrix(IntegerModRing(12), [[0, 1, 6, -1],[1, 6, 1, 2],[6, 1, 1, 0]])
     412            sage: C = RingCode(M)
     413            sage: C
     414            (4, 1728, 2)-code over the Ring of integers modulo 12
     415            sage: c = C.minimum_weight_codeword(); c
     416            (0, 1, 0, 5)
     417            sage: c in C
     418            True
     419
     420        """
     421        return self.mincodeword
     422
     423    def next(self):
     424        """
     425        Returns an iterator fr the codewords.
     426
     427        EXAMPLES:
     428            sage: M = Matrix(IntegerModRing(12), [[0, 1, 6, -1],[1, 6, 1, 2],[6, 1, 1, 0]])
     429            sage: C = RingCode(M)
     430            sage: C[0]
     431            (0, 11, 1, 1)
     432            sage: C.next()
     433            (0, 11, 1, 1)
     434            sage: C[1]
     435            (0, 11, 2, 2)
     436            sage: C.next()
     437            (0, 11, 2, 2)
     438
     439        """
     440        cdef int c
     441        if self.count >= len(self.codeSet):
     442            self.count = 0
     443            raise StopIteration
     444        else:
     445            c = self.count
     446            self.count = c + 1
     447            return self.codeSet[c]
     448
     449    def __repr__(self):
     450        s = '(%d, %d, %d)-code over the %s'%(self.length(), self.codewords, self.minimum, self.base_ring())
     451        return s
     452       
     453    def spanning_codewords(self):
     454        """
     455        Returns a list of codewords which span the code.
     456
     457        EXAMPLES:
     458            sage: M = Matrix(IntegerModRing(12), [[0, 1, 6, -1],[1, 6, 1, 2],[6, 1, 1, 0]])
     459            sage: C = RingCode(M)
     460            sage: C.spanning_codewords()
     461            [(0, 1, 6, 11), (1, 6, 1, 2), (6, 1, 1, 0)]
     462
     463        """
     464        return self.gen_matrix.rows()
     465 No newline at end of file