formlp.py

著名的ldpc编解码的资料及源代码

        # parity checck sum to 0 modulo 2.

        As Martin Wainwright explained to me there are two ways to
        represent these constraints using 'configuration variables'
        as discussed in the paper

        J. Feldman, D. Karger and M. Wainwright
        Using linear programming to decode LDPC codes.
        Conference on Information Sciences and Systems, Baltimore, March 2003.

        While this form is conceptually simpler, the LP relaxation can
        be more compactly written by projecting out the configuration
        variables to get an LP as in example 1 of the paper

        J. Feldman, D. Karger and M. J. Wainwright,
        LP Decoding  (Invited Paper)
        Allerton Conference on Communication, Control, and Computing
        October 1--3, 2003; Urbana-Champaign, IL

        The rule for making the compact LP relaxation is as follows.
        Let x = x0, x1, ..., xC be the C variables involved in a parity
        check constraint.  let S be a subset of these C variables which
        are set to 1.  For example S = {0,2} means x0=x2=1 while the
        other xi are 0.  For each such subset with odd cardinality, the
        corresponding constraint in the LP is

           sum_{i in S} xi - sum_{i not in S} xi <= |S| - 1

        We obtain the required constraints for the LP by summing doing
        this for all possible subsets, S, with odd cardinality.
        
        For example, if x0,x1,x2 are connected to a check, then the
        required subsets of odd cardinality are
        S = {0,1,2}, S' = {0}, S'' = {1}, S'' = {2}
        yielding the constraints

        x0 + x1 + x2 <= 2
        x0 - x1 - x2 <= 0        
        -x0 + x1 - x2 <= 0
        -x0 - x1 + x2 <= 0                        
        
        """
        assert checkOn == 0 or checkOn == 1
        result = []
        varsForCheck = len(varIndexes)

        for configuration in range(1 << varsForCheck):
            if ( not ((CountOnes(varsForCheck,configuration)+checkOn) % 2)):
                continue
            constraintList = []                
            constraintList.append(
                's.t. c' + `checkNum` + 'c' + `checkOn` + '_' +
                `configuration` + ': ' + q(checkOn,'','-') + 
                'relaxedChecks[' + `checkNum` + '] ')
            for i in range(varsForCheck):
                if (configuration & (1 << i)):
                    constraintList.append('+')
                else:
                    constraintList.append('-')
                constraintList.append('relaxedVars[' + `varIndexes[i]` + ']')
            constraintList.append('<= ' +
                                  `(q(checkOn,0,-1) +
                                    CountOnes(varsForCheck,configuration))`)
            constraintList.append(';\n')
            result.append(string.join(constraintList))
        return string.join(result)

    def FormLPPreamble(self):

        self.output.append('set Vars := {0..' + `self.K-1` + '};\n')
        self.output.append('set Checks := {0..' + `self.N-1` + '};\n')
        self.output.append("""

var relaxedChecks{c in Checks};
var relaxedVars{v in Vars};

s.t. checksUB{c in Checks}: relaxedChecks[c] >= 0;
s.t. checksLB{c in Checks}: relaxedChecks[c] <= 1;
s.t. varsUB{v in Vars}: relaxedVars[v] <= 1;
s.t. varsLB{v in Vars}: relaxedVars[v] >= 0;

        """)
        

    def FormConstraintsConnectingChecksAndVars(self):

        for checkNum in range(len(self.codeMatrix)):
            self.output.append(self.MakeConstraintsConnectingCheckAndVars(
                checkNum,self.codeMatrix[checkNum],1))
            self.output.append(self.MakeConstraintsConnectingCheckAndVars(
                checkNum,self.codeMatrix[checkNum],0))

    def FormLPVarSection(self):
        """
        Form an linear program relaxation suitable for decoding the
        code codeMatrix.
        """

        self.FormConstraintsConnectingChecksAndVars()

        self.output.append(
            'param channel{c in Checks};\n\n' + 
            'minimize flipCost: sum{c in Checks} relaxedChecks[c]*channel[c];\n' +
            '/* a cost of +1 is incurred if we make check i a 1 when */\n' +
            '/* the channel evidence for it is 0 and a cost of -1 is */\n' +
            '/* incurred if we make check i a 0 when the channel */\n' +
            '/* evidence for it is 0.  Thus we try to match the channel */\n' +
            '/* evidence as well as possible. */\n')

    def FormErasureConstraints(self,channelData):
        """
        channelData: A description of the source values or channelData.
                     1's and 0's are represented as 1's and 0's while
                     erasures or don't cares are represented with 0.5.
        """
        for i in range(len(channelData)):
            if (0.5 != channelData[i]):
                self.output.append('s.t. e_' + `i` + '_constraint: ' +
                                   'relaxedChecks[' + `i` + ']=' +
                                   `channelData[i]` + ';\n')
        
class ECCLPFormer(LPFormer):
    """
    This class forms a linear programming relaxation to do error correction.
    The main methods intended to be called by the user are listed below,
    see the accompanying documentation for this class or the base class
    LPFormer for details.

    __init__
    FormLP
    FormErasureLP
    SolveLP
    IterSolveLP
    IterSolveLPDecimateAmbiguous
    GetLP
    """
    
    def FormLP(self,channelData,extraVarData=''):
        """
        channelData:   A description of the channel output values.
                       1's and 0's are represented as 1's and 0's while
                       erasures or don't cares are represented with 0.5.
        extraVarData:  A string representing extra stuff to put in the
                       var section of the linear program.

        This function forms the LP required to decode the channelData.
        After this function is called you can call one of the solve
        methods.

        Note that if you are doing erasure decoding you can either
        use this method or the FormErasureLP method but not both.
        """
        self.FormLPPreamble()
        self.FormLPVarSection()
        self.output.append(extraVarData)        
        self.FormLPDataSection()
        self.FormLPSourceData(channelData)        


    def FormErasureLP(self,channelData,extraVarData=''):
        """
        channelData:   A description of the recieved data from the channel
                       with 1's and 0's represented as 1's and 0's while
                       erasures or don't cares are represented with 0.5.
        extraVarData:  A string representing extra stuff to put in the
                       var section of the linear program.

        This function forms the LP required to decode the
        channelData using the erasure metric.  That is it attempts
        to exactly match all 0's and 1's and ignores 0.5's.  After
        this function is called you can call one of the solve methods.

        Note that if you are doing erasure decoding you can either
        use this method or the FormLP method but not both, but you
        *must* be doing erasure decoding to use this method.

        The main difference in this method is that it encodes the
        problem purely in the constraints of the linear program and
        does not have any objective function.
        """        
        self.FormLPPreamble()
        self.FormConstraintsConnectingChecksAndVars()
        self.output.append(extraVarData)
        self.FormErasureConstraints(channelData)        
        self.FormLPDataSection()

    def MakeConstraintsEnforcingParity(self,checkNum,varIndexes):
        """
        checkNum:   The index of the check to work on.
        varIndexes: List of indexes of the vars connected to check checkNum.

        This function adds constraints to the linear program requiring
        that the modulo-2 sum of the relaxedVars named by varIndexes
        equals 0.

        See the comment for MakeConstraintsEnforcingParity to see
        where the logic for this function comes from.
        """
        result = []
        varsForCheck = len(varIndexes)

        for configuration in range(1 << varsForCheck):
            if ( not (CountOnes(varsForCheck,configuration) % 2)):
                continue # only do stuff for subsets with odd number of 1's
            constraintList = []                
            constraintList.append(
                's.t. c' + `checkNum` + 'c' + '_' +
                `configuration` + ': ')
            for i in range(varsForCheck):
                if (configuration & (1 << i)):
                    constraintList.append('+')
                else:
                    constraintList.append('-')
                constraintList.append('relaxedVars[' + `varIndexes[i]` + ']')
            constraintList.append('<= ' +
                                  `CountOnes(varsForCheck,configuration)-1`)
            constraintList.append(';\n')
            result.append(string.join(constraintList))
        return string.join(result)

    def FormLPPreamble(self):

        self.output.append('set Vars := {0..' + `self.N-1` + '};\n')
        self.output.append("""

var relaxedVars{v in Vars};

s.t. varsUB{v in Vars}: relaxedVars[v] <= 1;
s.t. varsLB{v in Vars}: relaxedVars[v] >= 0;

        """)


    def FormLPVarSection(self):
        """
        Form an linear program relaxation suitable for decoding the
        code codeMatrix.
        """

        for checkNum in range(len(self.codeMatrix)):
            self.output.append(
                self.MakeConstraintsEnforcingParity(checkNum,
                                                    self.codeMatrix[checkNum]))
            

        self.output.append(
            'param channel{v in Vars};\n\n' + 
            'minimize flipCost: sum{v in Vars} relaxedVars[v]*channel[v];\n' +
            '/* a cost of +1 is incurred if we make var i a 1 when */\n' +
            '/* the channel evidence for it is 0 and a cost of -1 is */\n' +
            '/* incurred if we make var i a 0 when the channel */\n' +
            '/* evidence for it is 0.  Thus we try to match the channel */\n' +
            '/* evidence as well as possible. */\n')

    def FormErasureConstraints(self,channelData):
        """
        channelData: A description of the source values or channelData.
                     1's and 0's are represented as 1's and 0's while
                     erasures or don't cares are represented with 0.5.
        """
        for i in range(len(channelData)):
            if (0.5 != channelData[i]):
                self.output.append('s.t. e_' + `i` + '_constraint: ' +
                                   'relaxedVars[' + `i` + ']=' +
                                   `channelData[i]` + ';\n')

    

# The following code is used to make the doctest package
# check examples in docstrings.

def _test():
    import random
    random.setstate((1, (29245, 20096, 302), None)) # for consistent testing
    import doctest, FormLP
    return doctest.testmod(FormLP)

if __name__ == "__main__":
    _test()
    print 'Tests passed'