Viewing file:  Matrix.py (6.05 KB) -rw-r--r--
Select action/file-type:
 (+) |  (+) |  (+) | Code (+) | Session (+) |  (+) | SDB (+) |  (+) |  (+) |  (+) |  (+) |  (+) |
import string
import types
import Numeric
import LinearAlgebra
from UserArray import UserArray, asarray
from Numeric import dot, identity, multiply
from MLab import squeeze
# make translation table
_table = [None]*256
for k in range(256):
_table[k] = chr(k)
_table = ''.join(_table)
_numchars = string.digits + ".-+jeEL"
_todelete = []
for k in _table:
if k not in _numchars:
_todelete.append(k)
_todelete = ''.join(_todelete)
def _eval(astr):
return eval(astr.translate(_table,_todelete))
def _convert_from_string(data):
data.find
rows = data.split(';')
newdata = []
count = 0
for row in rows:
trow = row.split(',')
newrow = []
for col in trow:
temp = col.split()
newrow.extend(map(_eval,temp))
if count == 0:
Ncols = len(newrow)
elif len(newrow) != Ncols:
raise ValueError, "Rows not the same size."
count += 1
newdata.append(newrow)
return newdata
_lkup = {'0':'000',
'1':'001',
'2':'010',
'3':'011',
'4':'100',
'5':'101',
'6':'110',
'7':'111'}
def _binary(num):
ostr = oct(num)
bin = ''
for ch in ostr[1:]:
bin += _lkup[ch]
ind = 0
while bin[ind] == '0':
ind += 1
return bin[ind:]
class Matrix(UserArray):
def __init__(self, data, typecode=None, copy=1, savespace=0):
if type(data) is types.StringType:
data = _convert_from_string(data)
UserArray.__init__(self,data,typecode, copy, savespace)
if len(self.array.shape) == 1:
self.array.shape = (1,self.array.shape[0])
self.shape = self.array.shape
def __getitem__(self, index):
out = self._rc(self.array[index])
try:
n = len(index)
if n > 1 and isinstance(index[0], types.SliceType) \
and isinstance(index[1], types.IntType):
sh = out.array.shape
out.array.shape = (sh[1],sh[0])
out.shape = out.array.shape
except TypeError: # Index not a sequence
pass
return out
def __mul__(self, other):
aother = asarray(other)
if len(aother.shape) == 0:
return self._rc(self.array*aother)
else:
return self._rc(dot(self.array, aother))
def __rmul__(self, other):
aother = asarray(other)
if len(aother.shape) == 0:
return self._rc(aother*self.array)
else:
return self._rc(dot(aother, self.array))
def __imul__(self,other):
aother = asarray(other)
self.array = dot(self.array, aother)
return self
def __pow__(self, other):
shape = self.array.shape
if len(shape)!=2 or shape[0]!=shape[1]:
raise TypeError, "matrix is not square"
if type(other) in (type(1), type(1L)):
if other==0:
return Matrix(identity(shape[0]))
if other<0:
result=Matrix(LinearAlgebra.inverse(self.array))
x=Matrix(result)
other=-other
else:
result=self
x=result
if other <= 3:
while(other>1):
result=result*x
other=other-1
return result
# binary decomposition to reduce the number of Matrix
# Multiplies for other > 3.
beta = _binary(other)
t = len(beta)
Z,q = x.copy(),0
while beta[t-q-1] == '0':
Z *= Z
q += 1
result = Z.copy()
for k in range(q+1,t):
Z *= Z
if beta[t-k-1] == '1':
result *= Z
return result
else:
raise TypeError, "exponent must be an integer"
def __rpow__(self, other):
raise TypeError, "x**y not implemented for matrices y"
def __invert__(self):
return Matrix(Numeric.conjugate(self.array))
def __setattr__(self,attr,value):
if attr=='shape':
if len(value) == 0:
value = (1,1)
if len(value) == 1:
value = (1,) + value
if len(value) != 2:
raise ValueError, "Can only reshape a Matrix as a 2-d array."
self.array.shape=value
self.__dict__[attr]=value
def __getattr__(self, attr):
if attr == 'A':
return squeeze(self.array)
elif attr == 'T':
return Matrix(Numeric.transpose(self.array))
elif attr == 'H':
if len(self.array.shape) == 1:
self.array.shape = (1,self.array.shape[0])
return Matrix(Numeric.conjugate(Numeric.transpose(self.array)))
elif attr == 'I':
return Matrix(LinearAlgebra.inverse(self.array))
elif attr == 'real':
return Matrix(self.array.real)
elif attr == 'imag':
return Matrix(self.array.imag)
elif attr == 'flat':
return Matrix(self.array.flat)
else:
raise AttributeError, attr + " not found."
def __setitem__(self, index, value):
value = asarray(value, self._typecode)
self.array[index] = squeeze(value)
def __setslice__(self, i, j, value): self.array[i:j] = asarray(squeeze(value),self._typecode)
def __float__(self):
return float(squeeze(self.array))
def m(self,b):
return Matrix(self.array * asarray(b))
def asarray(self,t=None):
if t is None:
return self.array
else:
return asarray(self.array,t)
if __name__ == '__main__':
from Numeric import *
m = Matrix( [[1,2,3],[11,12,13],[21,22,23]])
print m*m
print m.array*m.array
print transpose(m)
print m**-1
m = Matrix([[1,1],[1,0]])
print "Fibonacci numbers:",
for i in range(10):
mm=m**i
print mm[0][0],
print
|