python 深度模仿 matlab 矩阵语法

来源:互联网 发布:学淘宝美工 编辑:程序博客网 时间:2024/05/25 01:35

直接上代码。MyMat (298行开始)类模仿matlab矩阵索引与赋值。从1开始计数,倒数从0开始。可以用单指标索引:1为第一个元素,2 为第一列第二行元素……以后还会加以改进。网友们也可以提出一些建议。(代码以pypi为准


If you like matlab as well as python, then this is your choice. Class MyMat imitates the reference and assignment of matlab matrices. The index is form 1, instead of 0, and the last index is 0, instead of -1. I will add more functions or classes in the future to improve the modules. You can give me any advice under the article.

You can download the codes from pypi. https://pypi.python.org/pypi/mymat


截止10月1日(接近一个月),至少已有九千人下载代码。



</pre><p></p><p><pre name="code" class="python">import numpy as npimport numpy.linalg as LA'''notation (abbreviation):num: complex numberind: index (int, slice or list, or tuple of index sometimes)sglind: single index (int, slice or list)vec: vectorcolvec: column vectorslc: sliceran: rangelst: listmat: matrixr(m): rowc(n): column'''# python for my matrix type# special matricesdef O(m, n=None):    # 0 matrix    if n==None:        n = m    return MyMat(np.zeros((m, n)))def One(m, n=None):    # 1 matrix    if n==None:        n = m    return MyMat(np.ones((m, n)))def I(m, n=None):    # identity matrix    if n==None:        n = m    return MyMat(np.eye(m, n))def E(i,j,m,n=None):    if n==None: n = m    A = O(m,n)    A[i,j] = 1    return Adef R(m, n=None, randfun='rand'):    # random matrix    if n==None:        n = m    return MyMat(getattr(np.random, randfun)(m, n))def H(m, n=None):    # Hilbert matrix    if n==None:        n = m    return MyMat([[1 / (k + l - 1) for l in range(1, n+1)] for k in range(1, m+1)])def TestMat(m, n=None):    if n==None:        n = m    return MyMat([[l + (k-1)*n for l in range(1, n+1)] for k in range(1, m+1)])# mapsdef mat_map(f,*As):    # matrices (np.matrix) in As have the same size    rn, cn = As[0].shape    # size (shape) of matrix    return np.mat([[f(*(A[k, l] for A in As)) for l in range(cn)] for k in range(rn)])    def mymat_map(f,*As):    # matrices (MyMat) in As have the same size    rn, cn = As[0].shape    # size (shape) of matrix    return MyMat([[f(*(A[k, l] for A in As)) for l in range(1, cn+1)] for k in range(1, rn+1)])# relations# will be used in testdef equal(A, *Ms):    if myMat(A).isempty():        for M in Ms:            if not MyMat(M).isempty():                return False        return True    for M in Ms:        M = np.mat(M)        if A.shape != M.shape or (A != M).any():            return False    return Truedef equal_tol(A, *Ms, tol=0.001):    for k, M in enumerate(Ms):        M = np.mat(M)        for N in Ms[k+1:]:            N = np.mat(N)            if A.shape != M.shape or LA.norm(A - M)>=tol:                return False    return True''' # Index class# it may be used in the futuredef indadd(ind, x):    if isinstance(ind, int):        return ind + x    elif isinstance(ind, slice):        return slice(ind.start + x, ind.stop + x, ind.step)    else:        return [a + x for a in ind]def indinv(ind):    if isinstance(ind, int):        return ind    elif isinstance(ind, slice):        return slice(ind.stop, ind.start, -ind.step)    else:        ind.reverse()        return inddef indneg(ind):    if isinstance(ind, int):        return -ind    elif isinstance(ind, slice):        return slice(-ind.stop, -ind.start, -ind.step)    else:        return [-a for a in ind]    class Index:    # value: int, slice, list    def __init__(self, value=0):        self.value=value    def __add__(self, other):        if isinstance(other, int):            return Index(indadd(self.value, other))    def __neg__(self):        return Index(indneg(self.value))        def __invert__(self):        return Index(indinv(self.value))    def __str__(self):        return ','.join(map(str, ind2iter(self.value)))'''# MyMat Classdef isreal(x):    return isinstance(x, (int, float))def isnum(x):    return np.isscalar(x)def isscalar(x):    return isnum(x) or isinstance(x, MyMat) and x.isscalar()def issgl(x):    # single index used in python    return isinstance(x, (int, slice, list))def isemp(x):    # empty index    return isinstance(x, list) and x==[] or isinstance(x, tuple) and (x[0]==[] or x[1]==[])def iscrd(x):    # coordinate type index    # isinstance(x, tuple) and len(x) == 2    return isinstance(x, tuple) and isinstance(x[0], list) and isinstance(x[1], list)COLON = slice(None)# transform of the indicesimport itertoolsdef times(ind):    '''>>> times(([1,2],[3,4,6]))([1, 1, 1, 2, 2, 2], [3, 4, 6, 3, 4, 6])'''    if iscrd(ind):        lst = itertools.product(ind[0],ind[1])        b=([], [])        for a in lst:            b[0].append(a[0])            b[1].append(a[1])        return b    else:        return inddef ind2iter(ind, N):    # slice to range, int and list to list    if isinstance(ind, list):        return ind    elif isinstance(ind, slice): # slc is a slice        a = ind.start if ind.start else 0        b = ind.stop if ind.stop else N # b <= N        return range(a, b, c) if ind.step else range(a,b)    else:        return [ind]def ind2list(ind, N):    # slice to range, int and list to list    if isinstance(ind, list):        return ind    elif isinstance(ind, slice): # slc is a slice        a = ind.start if ind.start else 0        b = ind.stop if ind.stop else N # b <= N        return list(range(a, b, c)) if ind.step else list(range(a,b))    else:        return [ind]    def ind2ind(ind):    '''index of matlab type to index of python typeThis is the main gimmick to define MyMatalso see slice2slice'''    if isinstance(ind, tuple):        return tuple(ind2ind(a) for a in ind)    if isinstance(ind, list):        return [a - 1 for a in ind]    elif isinstance(ind, slice):        return slice2slice(ind)       # see its definition    else:        return ind - 1def slice2slice(slc):    '''slice of matlab type to slice of python typefor example: M[1:3] := M[0:3] >>> slice2slice(slice(1, 6, 2))slice(0, 6, 2)>>> slice2slice(slice(-3, 0))slice(-4, None, None)'''    return slice(slc.start-1 if slc.start else None, None if slc.stop == 0 else slc.stop, slc.step)def int2pair(ind, r, c=1):    '''>>> int2pair(ind2ind([3,4,7,10]),5,5)([2, 3, 1, 4], [0, 0, 1, 1])'''    if isinstance(ind, int):        return np.mod(ind, r), np.floor_divide(ind, r)    else:        ind=ind2iter(ind, r*c)        return [np.mod(a, r) for a in ind], [np.floor_divide(a, r) for a in ind]def pair2int(ind, r, c=None):    '''inverse of int2pair>>> pair2int(int2pair([4,7,12,22],5,5),5,5)[4, 7, 12, 22]'''    if isinstance(ind[0], int) and isinstance(ind[1], int):        return ind[0] + r * ind[1]    else:        if c==None:            c=max(ind[1])        return [a + r * b for a, b in zip(ind2iter(ind[0], r), ind2iter(ind[1], c))]"""# int2pair(ind2ind) == ind2ind(int2pair2)def int2pair2(ind, r, c=1):    '''single index (int) => pair index (tuple)example:>>> int2pair(7, 3)(1, 3)    if isinstance(ind, int):        return (r, np.floor_divide(ind, r)) if np.mod(ind, r)==0 else (np.mod(ind, r), np.floor_divide(ind, r)+1)    else:  # extend to lists or slices        ind=ind2iter(ind, r * c)        p=([],[])        for a in ind:            if np.mod(a, r)==0:                p[0].append(r); p[1].append(np.floor_divide(a, r))            else:                p[0].append(np.mod(a, r)); p[1].append(np.floor_divide(a, r)+1)        return p"""def maxind(ind):    # if issgl(ind):    if isinstance(ind, int):        return ind    elif isinstance(ind, list):        return max(ind)    elif isinstance(ind, slice):        return ind.stop if ind.stop else 0    else:        return maxind(ind[0]), maxind(ind[1])# definition of MyMat Classclass MyMat(np.matrix):    # python for imitating matrices (2D) in matlab    # matrix on complex numbers    def __new__(cls, data, *args, **kw):        if isinstance(data, list):            data=np.array(data)        elif isinstance(data, str):            data = np.array(np.matrix(data))        elif hasattr(data, 'tolist'):            data = np.array(data.tolist())        return super(MyMat, cls).__new__(cls, data, *args, **kw)    def __mul__(self, other):        if isinstance(other, MyMat):            if other.isscalar():                return super(MyMat, self).__mul__(other.toscalar())            else:                return super(MyMat, self).__mul__(other)        else:            if self.isscalar():                return MyMat([self[1,1] * other])            else:                return super(MyMat, self).__mul__(other)                def __rmul__(self, other):        return self * other       def __imul__(self, other):        self[:,:] = self * other        return self    def __pow__(self, other):        # A**B := A.*B        return MyMat(np.multiply(self, other))    def __rpow__(self, other):        # B**A := B.*A        return MyMat(np.multiply(other, self))    def __ipow__(self, other):        self[:,:] = self ** other        return self    def __floordiv__(self, other):        # A//B := A./B        return MyMat(np.divide(self, other))    def __rfloordiv__(self, other):        # B//A := B./A        return MyMat(np.divide(other, self))    def __ifloordiv__(self, other):        self[:,:] = self // other        return self    def __truediv__(self, other):        # A/B := A/B (A*inv(B))        if isinstance(other, MyMat) and other.isscalar():            return MyMat(np.divide(self, other.toscalar()))        elif isnum(other):            return self * (1 / other)  # = 1/other * self               return self * other.I    def __rtruediv__(self, other):        # r/A := r/A (r * inv(A))        return other * self.I    def __itruediv__(self, other):        self[:,:] = self / other        return self        def __xor__(self, num):        # A^n := A^n, n is int        if isinstance(num, MyMat) and num.isscalar():            return super(MyMat, self).__pow__(num.toscalar())        return super(MyMat, self).__pow__(num)    def __rxor__(self, num):        pass    def __ixor__(self, num):        self[:,:] = self ^ other        return self    def __lshift__(self, other):        # A<<B := A.^B        return MyMat(np.power(self, other))    def __rlshift__(self, other):        # r<<A := r.^A        return MyMat(np.power(other, self))    def __ilshift__(self, other):        self[:,:] = self << other        return self    @property    def I(self):        return self.getI()    def __getitem__(self, ind):        r, c = self.shape        if iscrd(ind):            return self.get((ind[0], COLON)).get((COLON, ind[1]))        else:            return self.get(ind)    def __setitem__(self, ind, val):        r, c = self.shape        if issgl(ind):            self.set(ind, val)        else:            if iscrd(ind):                if isinstance(val, MyMat):                    if val.isscalar():                        self.set(times(ind), val.toscalar())                    else:                        self.set(times(ind), val.tovec(dim=2).tolist())                elif isnum(val):                    self.set(times(ind), val)                else: # isinstance(val, list)                    self.set(times(ind), val)            else:                self.set(ind, val)    def get(self, ind):        r, c = self.shape        if issgl(ind):            return super(MyMat, self).__getitem__(int2pair(ind2ind(ind), r, c))        return super(MyMat, self).__getitem__(ind2ind(ind))    def set(self, ind, val):        r, c = self.shape        if issgl(ind):   # single index            mx=maxind(ind)            if mx > r*c:                mx1, mx2=int2pair(mx-1, r, c)                self.just(max(r, mx1+1), max(c, mx2+1))            if isnum(val):                super(MyMat, self).__setitem__(int2pair(ind2ind(ind), r, c), val)            elif isinstance(val, list):                super(MyMat, self).__setitem__(int2pair(ind2ind(ind), r, c), val)            elif isinstance(val, MyMat):                if val.isscalar():                    super(MyMat, self).__setitem__(int2pair(ind2ind(ind), r, c), val[1, 1])                elif val.isvec():                    super(MyMat, self).__setitem__(int2pair(ind2ind(ind), r, c), val.tolist())                elif val.iscolvec():                    super(MyMat, self).__setitem__(int2pair(ind2ind(ind), r, c), val.T.tolist())            elif isinstance(val, np.matrix):                self.set(ind, MyMat(val))        else:      # double index            mx1, mx2 = maxind(ind)            if mx1 > r or mx2 > c:                self.just(max(r, mx1), max(c, mx2))            if isinstance(ind[0], slice) and isinstance(ind[0], slice):                super(MyMat, self).__setitem__(ind, val)            else:                super(MyMat, self).__setitem__(ind2ind(ind), val)    ''' def set_(self, ind, val):        # extension of setitem        r, c = self.shape        if issgl(ind):            mx=maxind(ind)            if mx > r*c:                mx1, mx2=int2pair(mx-1, r, c)                self.just(max(r, mx1+1), max(c, mx2+1))        else:            mx1, mx2 = maxind(ind)            if mx1 > r or mx2 > c:                self.just(max(r, mx1), max(c, mx2))        self.set(ind, val)        return self'''    def update(self, other):        self.resize(other.shape, refcheck=False)        super(MyMat, self).__setitem__((COLON, COLON), other)    def __call__(self, *ind):        if len(ind)==1:            return self[ind[0]]        return self[ind]    def __str__(self):        return '['+self.join(ch2=';\n ')+']'    def join(self, ch1=', ', ch2='; '):        r, c = self.shape        return ch2.join(ch1.join(str(super(MyMat, self).__getitem__((k, l))) for l in range(c)) for k in range(r))    def cat(self, *others, dim=1):        if self.isempty():            return np.concatenate(others, axis=dim-1)        elif not others:            return self        else:            return np.concatenate((self,)+others, axis=dim-1)    def __or__(self, other):        return np.concatenate((self, other), 1)    def __ior__(self, other):        self.update(self | other)        return self    def __and__(self, other):        return np.concatenate((self, other))    def __iand__(self, other):        self.update(self & other)        return self    def wequal(self, *others):        # weakly equal        if self.isempty():            for other in others:                if not MyMat(other).isempty():                    return False            return True        for other in others:            M = MyMat(other)            if self.shape != M.shape or (self != M).any():                return False        return True    def equal(self, *others):        # others are MyMat        if self.isempty():            for other in others:                if not other.isempty():                    return False            return True        for other in others:            if self.shape != other.shape or (self != M).any():                return False        return True        def toscalar(self):        return super(MyMat, self).__getitem__((0,0))        def tolist(self):        if self.isvec():            return super(MyMat, self).tolist()[0]        else:            return super(MyMat, self).tolist()    def toarray(self):        if self.isvec():            return self.A1        else:            return self.A    def tomatrix(self):        return np.matrix(self.tolist())    def tomat(self):        # == tomatrix        return np.mat(self.tolist())    def tovec(self, dim=1):        if dim==1:            return MyMat([[self[k, l]] for l in range(1, self.shape[1]+1) for k in range(1, self.shape[0]+1)])        elif dim==2:            lst=[]            for k in range(1, self.shape[0]+1):                lst.extend(self[k,:].tolist())             return MyMat(lst)    def just(self, m, n=None):        if n == None: n = m        r, c = self.shape        if m>0 and n>0:            if m<=r:                temp = super(MyMat, self).__getitem__((slice(m), slice(r)))            else:                temp = self | O(m - r, c)            if n<=c:                temp = super(MyMat, temp).__getitem__((slice(c), slice(n)))            else:                temp = temp & O(m, n - c)        elif m==0 or n==0:            temp = Empty        else:            temp = O(abs(m), abs(n))        self.update(temp)        # return self    def delete(self, ind1=[], ind2=[]):        r, c = self.shape        if not isemp(ind1):            ind1=ind2list(ind2ind(ind1), r)            if ind1 == list(range(r)):                return Empty            self=np.concatenate(tuple(super(MyMat, self).__getitem__((slice(a+1,b), COLON)) for a, b in zip([-1]+ind1, ind1+[r]) if b-a!=1))        if not isemp(ind2):            ind2=ind2list(ind2ind(ind2), c)            if ind2 == list(range(c)):                return Empty            return np.concatenate(tuple(super(MyMat, self).__getitem__((COLON, slice(a+1,b))) for a, b in zip([-1]+ind2, ind2+[c]) if b-a!=1), 1)        else:            return self    '''def delete(self, ind1, ind2):        ......        if not isemtpy(ind1):            self=super(PyMat, self).__getitem__((slice(ind1[0]), COLON)).cat(                *(tuple(super(PyMat, self).__getitem__((slice(a+1,b), COLON)) for a, b in zip(ind1[:-1], ind1[1:]) if b-a!=1)                +(super(PyMat, self).__getitem__((slice(ind1[-1]+1,None), COLON)),)))            if self.isempty():                return Empty        if not isemtpy(ind2):            return super(PyMat, self).__getitem__((COLON, slice(ind2[0]))).cat(                *(tuple(super(PyMat, self).__getitem__((COLON, slice(a+1,b))) for a, b in zip(ind2[:-1], ind2[1:]) if b-a!=1)                +(super(PyMat, self).__getitem__((COLON, slice(ind2[-1]+1,None))),)), dim=2)        else:            return self'''                def norm(self, p=2):        if self.isvec():            return LA.norm(self.toarray(), p)        elif self.iscolvec():            return LA.norm(self.T.toarray(), p)        else:            return LA.norm(self, p)    def repmat(self, m=1):        return np.tile(self, m-1)            def apply(self, fun):        return mymat_map(fun, self)    def isempty(self):        return self.shape[1]==0 or self.shape[0]==0    def isscalar(self):        return self.shape == (1,1)    def isvec(self):        return self.shape[0] == 1    def iscolvec(self):        return self.shape[1] == 1    @property    def size(self):        return self.shape    def poly(self, p):        # self is a square matrix        r, c = self.shape  #  r == c        if len(p)==1:            return p[0]*I(r)        return p[0]*I(r)+self.poly(p[1:])*self    def expm(self, N=10):        p=[1/np.prod(range(1, n+1)) for n in range(N)]        return self.poly(p)# matrix instanceMMat = myMat = MyMat # aliasEmpty = MyMat([])    # empty setim = MyMat([[0, 1], [-1, 0]])def c(a=0, b=0):    # linear represenation of complex number    return MyMat([[a, b], [-b, a]])def q(a=0, b=0, c=0, d=0):    # linear represenation of quaternions    return MyMat([[a, b, c, d], [-b, a, -d, c], [-c, d, a, -b], [-d, -c, b, a]])


0 0
原创粉丝点击