genericmatrix.py

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

"""

     Copyright 2003 Mitsubishi Electric Research Laboratories All Rights
     Reserved.  Permission to use, copy and modify this software and its
     documentation without fee for educational, research and non-profit
     purposes, is hereby granted, provided that the above copyright
     notice and the following three paragraphs appear in all copies.      

     To request permission to incorporate this software into commercial
     products contact:  Vice President of Marketing and Business         
     Development;  Mitsubishi Electric Research Laboratories (MERL), 201
     Broadway, Cambridge, MA   02139 or <[e-mail]>.        
                                                                               
     IN NO EVENT SHALL MERL BE LIABLE TO ANY PARTY FOR DIRECT, INDIRECT, 
     SPECIAL, INCIDENTAL, OR CONSEQUENTIAL DAMAGES, INCLUDING LOST        
     PROFITS, ARISING OUT OF THE USE OF THIS SOFTWARE AND ITS 
     DOCUMENTATION, EVEN IF MERL HAS BEEN ADVISED OF THE POSSIBILITY OF
     SUCH DAMAGES.

     MERL SPECIFICALLY DISCLAIMS ANY WARRANTIES, INCLUDING, BUT NOT
     LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
     FOR A PARTICULAR PURPOSE.  THE SOFTWARE PROVIDED HEREUNDER IS ON AN
     ``AS IS'' BASIS, AND MERL HAS NO OBLIGATIONS TO PROVIDE MAINTENANCE,
     SUPPORT, UPDATES, ENHANCEMENTS OR MODIFICATIONS.
"""

# Copyright Emin Martinian 2002.  See below for license terms.
# Version Control Info: $Id: genericmatrix.py,v 1.2 2004/01/04 19:43:20 emin Exp $

"""
This package implements the GenericMatrix class to provide matrix
operations for any type that supports the multiply, add, subtract,
and divide operators.  For example, this package can be used to
do matrix calculations over finite fields using the ffield package
available at http://martinian.com.

The following docstrings provide detailed information on various topics:

  GenericMatrix.__doc__   Describes the methods of the GenericMatrix
                          class and how to use them.

  license_doc             Describes the license and lack of warranty
                          for this code.

  testing_doc             Describes some tests to make sure the code works.

"""

import operator

class GenericMatrix:

    """
    The GenericMatrix class implements a matrix with works with
    any generic type supporting addition, subtraction, multiplication,
    and division.  Matrix multiplication, addition, and subtraction
    are implemented as are methods for finding inverses,
    LU (actually LUP) decompositions, and determinants.  A complete
    list of user callable methods is:

    __init__
    __repr__
    __mul__
    __add__
    __sub__
    __setitem__
    __getitem__
    Size
    SetRow
    GetRow
    GetColumn
    Copy
    MakeSimilarMatrix
    SwapRows
    MulRow
    AddRow
    AddCol
    MulAddRow
    LeftMulColumnVec
    LowerGaussianElim
    Inverse
    Determinant
    LUP

    A quick and dirty example of how to use the GenericMatrix class
    for matricies of floats is provided below.
    
>>> import genericmatrix
>>> v = genericmatrix.GenericMatrix((3,3))
>>> v.SetRow(0,[0.0, -1.0, 1.0])
>>> v.SetRow(1,[1.0, 1.0, 1.0])
>>> v.SetRow(2,[1.0, 1.0, -1.0])
>>> v
<matrix
  0.0 -1.0  1.0
  1.0  1.0  1.0
  1.0  1.0 -1.0>
>>> vi = v.Inverse()
>>> vi
<matrix
  1.0  0.0  1.0
 -1.0  0.5 -0.5
 -0.0  0.5 -0.5>
>>> (vi * v) - v.MakeSimilarMatrix(v.Size(),'i')
<matrix
 0.0 0.0 0.0
 0.0 0.0 0.0
 0.0 0.0 0.0>

# See what happens when we try to invert a non-invertible matrix

>>> v[0,1] = 0.0
>>> v
<matrix
  0.0  0.0  1.0
  1.0  1.0  1.0
  1.0  1.0 -1.0>
>>> abs(v.Determinant())
0.0
>>> v.Inverse()
Traceback (most recent call last):
     ...
ValueError: matrix not invertible

# LUP decomposition will still work even if Inverse() won't.

>>> (l,u,p) = v.LUP()
>>> l
<matrix
 1.0 0.0 0.0
 0.0 1.0 0.0
 1.0 0.0 1.0>
>>> u
<matrix
  1.0  1.0  1.0
  0.0  0.0  1.0
  0.0  0.0 -2.0>
>>> p
<matrix
 0.0 1.0 0.0
 1.0 0.0 0.0
 0.0 0.0 1.0>
>>> p * v - l * u
<matrix
 0.0 0.0 0.0
 0.0 0.0 0.0
 0.0 0.0 0.0>

# Operate on some column vectors using v.
# The LeftMulColumnVec methods lets us do this without having
# to construct a new GenericMatrix to represent each column vector.
>>> v.LeftMulColumnVec([1.0,2.0,3.0])
[3.0, 6.0, 0.0]
>>> v.LeftMulColumnVec([1.0,-2.0,1.0])
[1.0, 0.0, -2.0]

# Most of the stuff above could be done with something like matlab.
# But, with this package you can do matrix ops for finite fields.
>>> XOR = lambda x,y: x^y
>>> AND = lambda x,y: x&y
>>> DIV = lambda x,y: x  
>>> m = GenericMatrix(size=(3,4),zeroElement=0,identityElement=1,add=XOR,mul=AND,sub=XOR,div=DIV)
>>> m.SetRow(0,[0,1,0,0])
>>> m.SetRow(1,[0,1,0,1])
>>> m.SetRow(2,[0,0,1,0])
>>> # You can't invert m since it isn't square, but you can still 
>>> # get the LUP decomposition or solve a system of equations.
>>> (l,u,p) = v.LUP()
>>> p*v-l*u
<matrix
 0.0 0.0 0.0
 0.0 0.0 0.0
 0.0 0.0 0.0>
>>> b = [1,0,1]
>>> x = m.Solve(b)
>>> b == m.LeftMulColumnVec(x)
1

    """

   
    def __init__(self, size=(2,2), zeroElement=0.0, identityElement=1.0,
                 add=operator.__add__, sub=operator.__sub__,
                 mul=operator.__mul__, div = operator.__div__,
                 eq = operator.__eq__, str=lambda x:`x`,
                 equalsZero = None,fillMode='z'):
        """
        Function:     __init__(size,zeroElement,identityElement,
                               add,sub,mul,div,eq,str,equalsZero fillMode)

        Description:  This is the constructor for the GenericMatrix
                      class.  All arguments are optional and default
                      to producing a 2-by-2 zero matrix for floats.
                      A detailed description of arguments follows:

             size: A tuple of the form (numRows, numColumns)
                   zeroElement: An object representing the additive
                   identity (i.e. 'zero') for the data
                   type of interest.
                   
             identityElement: An object representing the multiplicative
                              identity (i.e. 'one') for the data
                              type of interest.

             add,sub,mul,div: Functions implementing basic arithmetic
                              operations for the type of interest.

             eq: A function such that eq(x,y) == 1 if and only if x == y.

             str: A function used to produce a string representation of
                  the type of interest.

             equalsZero: A function used to decide if an element is
                         essentially zero.  For floats, you could use
                         lambda x: abs(x) < 1e-6.

             fillMode: This can either be 'e' in which case the contents
                       of the matrix are left empty, 'z', in which case
                       the matrix is filled with zeros, 'i' in which
                       case an identity matrix is created, or a two
                       argument function which is called with the row
                       and column of each index and produces the value
                       for that entry.  Default is 'z'.
        """
        if (None == equalsZero):
            equalsZero = lambda x: self.eq(self.zeroElement,x)

        self.equalsZero = equalsZero
        self.add = add
        self.sub = sub
        self.mul = mul
        self.div = div
        self.eq = eq
        self.str = str
        self.zeroElement = zeroElement
        self.identityElement = identityElement
        self.rows, self.cols = size
        self.data = []


        def q(x,y,z):
            if (x):
                return y
            else:
                return z

        if (fillMode == 'e'):
            return
        elif (fillMode == 'z'):
            fillMode = lambda x,y: self.zeroElement
        elif (fillMode == 'i'):
            fillMode = lambda x,y: q(self.eq(x,y),self.identityElement,
                                     self.zeroElement)

        for i in range(self.rows):
            self.data.append(map(fillMode,[i]*self.cols,range(self.cols)))

    def MakeSimilarMatrix(self,size,fillMode):
        """
        MakeSimilarMatrix(self,size,fillMode)

        Return a matrix of the given size filled according to fillMode
        with the same zeroElement, identityElement, add, sub, etc.
        as self.

        For example, self.MakeSimilarMatrix(self.Size(),'i') returns
        an identity matrix of the same shape as self.
        """
        return GenericMatrix(size=size,zeroElement=self.zeroElement,
                               identityElement=self.identityElement,
                               add=self.add,sub=self.sub,
                               mul=self.mul,div=self.div,eq=self.eq,
                               str=self.str,equalsZero=self.equalsZero,
                               fillMode=fillMode)
        

    def __repr__(self):
        m = 0
        # find the fattest element
        for r in self.data:
            for c in r:
                l = len(self.str(c))
                if l > m:
                    m = l
        f = '%%%ds' % (m+1)
        s = '<matrix'
        for r in self.data:
            s = s + '\n'
            for c in r:
                s = s + (f % self.str(c))
        s = s + '>'
        return s

    def __mul__(self,other):
        if (self.cols != other.rows):
            raise ValueError, "dimension mismatch"
        result = self.MakeSimilarMatrix((self.rows,other.cols),'z')
                               
        for i in range(self.rows):