2012-05-08 10 views
18

Sto provando a moltiplicare due matrici insieme usando puro python. Ingresso (X1 è un 3x3 e Xt è un 3x2):Moltiplicazione della matrice in python?

X1 = [[1.0016, 0.0, -16.0514], 
     [0.0, 10000.0, -40000.0], 
     [-16.0514, -40000.0, 160513.6437]] 
Xt = [(1.0, 1.0), 
     (0.0, 0.25), 
     (0.0, 0.0625)] 

dove Xt è la trasposta zip di un'altra matrice. Ora qui è il codice:

def matrixmult (A, B): 
    C = [[0 for row in range(len(A))] for col in range(len(B[0]))] 
    for i in range(len(A)): 
     for j in range(len(B[0])): 
      for k in range(len(B)): 
       C[i][j] += A[i][k]*B[k][j] 
    return C 

L'errore che Python mi dà è questo: indice di lista fuori campo: IndexError. Ora non sono sicuro che Xt sia riconosciuto come matrice ed è ancora un oggetto elenco, ma tecnicamente dovrebbe funzionare.

+5

perché non stai usando [numpy/scipy] (http://numpy.scipy.org/)? – ulmangt

+2

Se questo è compito, si prega di aggiungere il tag 'compiti a casa '. – agf

+3

@ulmangt: "usando puro python". Lui/lei vuole farlo senza moduli scaricabili, probabilmente per la sfida. – beary605

risposta

1

La forma della matrice C è errata; è la trasposizione di ciò che effettivamente vuoi che sia. (Ma sono d'accordo con ulmangt: la cosa giusta è quasi certamente di usare NumPy, davvero.)

2

il guasto si verifica qui:

C[i][j]+=A[i][k]*B[k][j] 

si blocca quando k = 2. Questo perché la tupla A[i] ha solo 2 valori e pertanto è possibile chiamarla solo su A [i] [1] prima che si verifichi un errore.

EDIT: Ascolta la risposta di Gerard, la tua C è sbagliata. Dovrebbe essere C=[[0 for row in range(len(A))] for col in range(len(A[0]))].

Solo un consiglio: si potrebbe sostituire il primo ciclo con una moltiplicazione, quindi sarebbe C=[[0]*len(A) for col in range(len(A[0]))]

+1

Vero se viene valutato matrixMult (Xt, X1) –

16

Questo è l'inizializzazione corretta. Hai scambiato fila con col!

C = [[0 for row in range(len(A))] for col in range(len(B[0]))] 

inizializzazione corretta sarebbe

C = [[0 for col in range(len(B[0]))] for row in range(len(A))] 

Anche io suggerirei di usare meglio convenzioni di denominazione. Ti aiuterà molto nel debugging. Per esempio:

def matrixmult (A, B): 
    rows_A = len(A) 
    cols_A = len(A[0]) 
    rows_B = len(B) 
    cols_B = len(B[0]) 

    if cols_A != rows_B: 
     print "Cannot multiply the two matrices. Incorrect dimensions." 
     return 

    # Create the result matrix 
    # Dimensions would be rows_A x cols_B 
    C = [[0 for row in range(cols_B)] for col in range(rows_A)] 
    print C 

    for i in range(rows_A): 
     for j in range(cols_B): 
      for k in range(cols_A): 
       C[i][j] += A[i][k] * B[k][j] 
    return C 

Si può fare molto di più, ma si ottiene l'idea ...

33

Se davvero non si vuole utilizzare numpy si può fare qualcosa di simile:

def matmult(a,b): 
    zip_b = zip(*b) 
    # uncomment next line if python 3 : 
    # zip_b = list(zip_b) 
    return [[sum(ele_a*ele_b for ele_a, ele_b in zip(row_a, col_b)) 
      for col_b in zip_b] for row_a in a] 

x = [[1,2,3],[4,5,6],[7,8,9],[10,11,12]] 
y = [[1,2],[1,2],[3,4]] 

import numpy as np # I want to check my solution with numpy 

mx = np.matrix(x) 
my = np.matrix(y)  

Risultato:

>>> matmult(x,y) 
[[12, 18], [27, 42], [42, 66], [57, 90]] 
>>> mx * my 
matrix([[12, 18], 
     [27, 42], 
     [42, 66], 
     [57, 90]]) 
+4

È possibile migliorare facilmente questo calcolo calcolando 'zip (* b)' una volta –

+0

@ gnibbler, questo è un buon punto, grazie. Ho modificato il codice per riflettere il tuo suggerimento. – Akavall

+0

Un bel consiglio per usare Numpy per controllare il tuo lavoro. –

4

Ecco un po 'di codice breve e semplice per le routine di matrix/Vector in puro Python t cappello Ho scritto molti anni fa:

'''Basic Table, Matrix and Vector functions for Python 2.2 
    Author: Raymond Hettinger  
''' 

Version = 'File MATFUNC.PY, Ver 183, Date 12-Dec-2002,14:33:42' 

import operator, math, random 
NPRE, NPOST = 0, 0     # Disables pre and post condition checks 

def iszero(z): return abs(z) < .000001 
def getreal(z): 
    try: 
     return z.real 
    except AttributeError: 
     return z 
def getimag(z): 
    try: 
     return z.imag 
    except AttributeError: 
     return 0 
def getconj(z): 
    try: 
     return z.conjugate() 
    except AttributeError: 
     return z 


separator = [ '', '\t', '\n', '\n----------\n', '\n===========\n' ] 

class Table(list): 
    dim = 1 
    concat = list.__add__  # A substitute for the overridden __add__ method 
    def __getslice__(self, i, j): 
     return self.__class__(list.__getslice__(self,i,j)) 
    def __init__(self, elems): 
     list.__init__(self, elems) 
     if len(elems) and hasattr(elems[0], 'dim'): self.dim = elems[0].dim + 1 
    def __str__(self): 
     return separator[self.dim].join(map(str, self)) 
    def map(self, op, rhs=None): 
     '''Apply a unary operator to every element in the matrix or a binary operator to corresponding 
     elements in two arrays. If the dimensions are different, broadcast the smaller dimension over 
     the larger (i.e. match a scalar to every element in a vector or a vector to a matrix).''' 
     if rhs is None:             # Unary case 
      return self.dim==1 and self.__class__(map(op, self)) or self.__class__([elem.map(op) for elem in self]) 
     elif not hasattr(rhs,'dim'):         # List/Scalar op 
      return self.__class__([op(e,rhs) for e in self]) 
     elif self.dim == rhs.dim:          # Same level Vec/Vec or Matrix/Matrix 
      assert NPRE or len(self) == len(rhs), 'Table operation requires len sizes to agree' 
      return self.__class__(map(op, self, rhs)) 
     elif self.dim < rhs.dim:          # Vec/Matrix 
      return self.__class__([op(self,e) for e in rhs] ) 
     return self.__class__([op(e,rhs) for e in self])   # Matrix/Vec 
    def __mul__(self, rhs): return self.map(operator.mul, rhs) 
    def __div__(self, rhs): return self.map(operator.div, rhs) 
    def __sub__(self, rhs): return self.map(operator.sub, rhs) 
    def __add__(self, rhs): return self.map(operator.add, rhs) 
    def __rmul__(self, lhs): return self*lhs 
    def __rdiv__(self, lhs): return self*(1.0/lhs) 
    def __rsub__(self, lhs): return -(self-lhs) 
    def __radd__(self, lhs): return self+lhs 
    def __abs__(self): return self.map(abs) 
    def __neg__(self): return self.map(operator.neg) 
    def conjugate(self): return self.map(getconj) 
    def real(self): return self.map(getreal ) 
    def imag(self): return self.map(getimag) 
    def flatten(self): 
     if self.dim == 1: return self 
     return reduce(lambda cum, e: e.flatten().concat(cum), self, []) 
    def prod(self): return reduce(operator.mul, self.flatten(), 1.0) 
    def sum(self): return reduce(operator.add, self.flatten(), 0.0) 
    def exists(self, predicate): 
     for elem in self.flatten(): 
      if predicate(elem): 
       return 1 
     return 0 
    def forall(self, predicate): 
     for elem in self.flatten(): 
      if not predicate(elem): 
       return 0 
     return 1 
    def __eq__(self, rhs): return (self - rhs).forall(iszero) 

class Vec(Table): 
    def dot(self, otherVec): return reduce(operator.add, map(operator.mul, self, otherVec), 0.0) 
    def norm(self): return math.sqrt(abs(self.dot(self.conjugate()))) 
    def normalize(self): return self/self.norm() 
    def outer(self, otherVec): return Mat([otherVec*x for x in self]) 
    def cross(self, otherVec): 
     'Compute a Vector or Cross Product with another vector' 
     assert len(self) == len(otherVec) == 3, 'Cross product only defined for 3-D vectors' 
     u, v = self, otherVec 
     return Vec([ u[1]*v[2]-u[2]*v[1], u[2]*v[0]-u[0]*v[2], u[0]*v[1]-u[1]*v[0] ]) 
    def house(self, index): 
     'Compute a Householder vector which zeroes all but the index element after a reflection' 
     v = Vec(Table([0]*index).concat(self[index:])).normalize() 
     t = v[index] 
     sigma = 1.0 - t**2 
     if sigma != 0.0: 
      t = v[index] = t<=0 and t-1.0 or -sigma/(t + 1.0) 
      v /= t 
     return v, 2.0 * t**2/(sigma + t**2) 
    def polyval(self, x): 
     'Vec([6,3,4]).polyval(5) evaluates to 6*x**2 + 3*x + 4 at x=5' 
     return reduce(lambda cum,c: cum*x+c, self, 0.0) 
    def ratval(self, x): 
     'Vec([10,20,30,40,50]).ratfit(5) evaluates to (10*x**2 + 20*x + 30)/(40*x**2 + 50*x + 1) at x=5.' 
     degree = len(self)/2 
     num, den = self[:degree+1], self[degree+1:] + [1] 
     return num.polyval(x)/den.polyval(x) 

class Matrix(Table): 
    __slots__ = ['size', 'rows', 'cols'] 
    def __init__(self, elems): 
     'Form a matrix from a list of lists or a list of Vecs' 
     Table.__init__(self, hasattr(elems[0], 'dot') and elems or map(Vec,map(tuple,elems))) 
     self.size = self.rows, self.cols = len(elems), len(elems[0]) 
    def tr(self): 
     'Tranpose elements so that Transposed[i][j] = Original[j][i]' 
     return Mat(zip(*self)) 
    def star(self): 
     'Return the Hermetian adjoint so that Star[i][j] = Original[j][i].conjugate()' 
     return self.tr().conjugate() 
    def diag(self): 
     'Return a vector composed of elements on the matrix diagonal' 
     return Vec([self[i][i] for i in range(min(self.size))]) 
    def trace(self): return self.diag().sum() 
    def mmul(self, other): 
     'Matrix multiply by another matrix or a column vector ' 
     if other.dim==2: return Mat(map(self.mmul, other.tr())).tr() 
     assert NPRE or self.cols == len(other) 
     return Vec(map(other.dot, self)) 
    def augment(self, otherMat): 
     'Make a new matrix with the two original matrices laid side by side' 
     assert self.rows == otherMat.rows, 'Size mismatch: %s * %s' % (`self.size`, `otherMat.size`) 
     return Mat(map(Table.concat, self, otherMat)) 
    def qr(self, ROnly=0): 
     'QR decomposition using Householder reflections: Q*R==self, Q.tr()*Q==I(n), R upper triangular' 
     R = self 
     m, n = R.size 
     for i in range(min(m,n)): 
      v, beta = R.tr()[i].house(i) 
      R -= v.outer(R.tr().mmul(v)*beta) 
     for i in range(1,min(n,m)): R[i][:i] = [0] * i 
     R = Mat(R[:n]) 
     if ROnly: return R 
     Q = R.tr().solve(self.tr()).tr()  # Rt Qt = At nn nm = nm 
     self.qr = lambda r=0, c=`self`: not r and c==`self` and (Q,R) or Matrix.qr(self,r) #Cache result 
     assert NPOST or m>=n and Q.size==(m,n) and isinstance(R,UpperTri) or m<n and Q.size==(m,m) and R.size==(m,n) 
     assert NPOST or Q.mmul(R)==self and Q.tr().mmul(Q)==eye(min(m,n)) 
     return Q, R 
    def _solve(self, b): 
     '''General matrices (incuding) are solved using the QR composition. 
     For inconsistent cases, returns the least squares solution''' 
     Q, R = self.qr() 
     return R.solve(Q.tr().mmul(b)) 
    def solve(self, b): 
     'Divide matrix into a column vector or matrix and iterate to improve the solution' 
     if b.dim==2: return Mat(map(self.solve, b.tr())).tr() 
     assert NPRE or self.rows == len(b), 'Matrix row count %d must match vector length %d' % (self.rows, len(b)) 
     x = self._solve(b) 
     diff = b - self.mmul(x) 
     maxdiff = diff.dot(diff) 
     for i in range(10): 
      xnew = x + self._solve(diff) 
      diffnew = b - self.mmul(xnew) 
      maxdiffnew = diffnew.dot(diffnew) 
      if maxdiffnew >= maxdiff: break 
      x, diff, maxdiff = xnew, diffnew, maxdiffnew 
      #print >> sys.stderr, i+1, maxdiff 
     assert NPOST or self.rows!=self.cols or self.mmul(x) == b 
     return x 
    def rank(self): return Vec([ not row.forall(iszero) for row in self.qr(ROnly=1) ]).sum() 

class Square(Matrix): 
    def lu(self): 
     'Factor a square matrix into lower and upper triangular form such that L.mmul(U)==A' 
     n = self.rows 
     L, U = eye(n), Mat(self[:]) 
     for i in range(n): 
      for j in range(i+1,U.rows): 
       assert U[i][i] != 0.0, 'LU requires non-zero elements on the diagonal' 
       L[j][i] = m = 1.0 * U[j][i]/U[i][i] 
       U[j] -= U[i] * m 
     assert NPOST or isinstance(L,LowerTri) and isinstance(U,UpperTri) and L*U==self 
     return L, U 
    def __pow__(self, exp): 
     'Raise a square matrix to an integer power (i.e. A**3 is the same as A.mmul(A.mmul(A))' 
     assert NPRE or exp==int(exp) and exp>0, 'Matrix powers only defined for positive integers not %s' % exp 
     if exp == 1: return self 
     if exp&1: return self.mmul(self ** (exp-1)) 
     sqrme = self ** (exp/2) 
     return sqrme.mmul(sqrme) 
    def det(self): return self.qr(ROnly=1).det() 
    def inverse(self): return self.solve(eye(self.rows)) 
    def hessenberg(self): 
     '''Householder reduction to Hessenberg Form (zeroes below the diagonal) 
     while keeping the same eigenvalues as self.''' 
     for i in range(self.cols-2): 
      v, beta = self.tr()[i].house(i+1) 
      self -= v.outer(self.tr().mmul(v)*beta) 
      self -= self.mmul(v).outer(v*beta) 
     return self 
    def eigs(self): 
     'Estimate principal eigenvalues using the QR with shifts method' 
     origTrace, origDet = self.trace(), self.det() 
     self = self.hessenberg() 
     eigvals = Vec([]) 
     for i in range(self.rows-1,0,-1): 
      while not self[i][:i].forall(iszero): 
       shift = eye(i+1) * self[i][i] 
       q, r = (self - shift).qr() 
       self = r.mmul(q) + shift 
      eigvals.append(self[i][i]) 
      self = Mat([self[r][:i] for r in range(i)]) 
     eigvals.append(self[0][0]) 
     assert NPOST or iszero((abs(origDet) - abs(eigvals.prod()))/1000.0) 
     assert NPOST or iszero(origTrace - eigvals.sum()) 
     return Vec(eigvals) 

class Triangular(Square): 
    def eigs(self): return self.diag() 
    def det(self): return self.diag().prod() 

class UpperTri(Triangular): 
    def _solve(self, b): 
     'Solve an upper triangular matrix using backward substitution' 
     x = Vec([]) 
     for i in range(self.rows-1, -1, -1): 
      assert NPRE or self[i][i], 'Backsub requires non-zero elements on the diagonal' 
      x.insert(0, (b[i] - x.dot(self[i][i+1:]))/self[i][i]) 
     return x 

class LowerTri(Triangular): 
    def _solve(self, b): 
     'Solve a lower triangular matrix using forward substitution' 
     x = Vec([]) 
     for i in range(self.rows): 
      assert NPRE or self[i][i], 'Forward sub requires non-zero elements on the diagonal' 
      x.append((b[i] - x.dot(self[i][:i]))/self[i][i]) 
     return x 

def Mat(elems): 
    'Factory function to create a new matrix.' 
    m, n = len(elems), len(elems[0]) 
    if m != n: return Matrix(elems) 
    if n <= 1: return Square(elems) 
    for i in range(1, len(elems)): 
     if not iszero(max(map(abs, elems[i][:i]))): 
      break 
    else: return UpperTri(elems) 
    for i in range(0, len(elems)-1): 
     if not iszero(max(map(abs, elems[i][i+1:]))): 
      return Square(elems) 
    return LowerTri(elems) 


def funToVec(tgtfun, low=-1, high=1, steps=40, EqualSpacing=0): 
    '''Compute x,y points from evaluating a target function over an interval (low to high) 
    at evenly spaces points or with Chebyshev abscissa spacing (default) ''' 
    if EqualSpacing: 
     h = (0.0+high-low)/steps 
     xvec = [low+h/2.0+h*i for i in range(steps)] 
    else: 
     scale, base = (0.0+high-low)/2.0, (0.0+high+low)/2.0 
     xvec = [base+scale*math.cos(((2*steps-1-2*i)*math.pi)/(2*steps)) for i in range(steps)] 
    yvec = map(tgtfun, xvec) 
    return Mat([xvec, yvec]) 

def funfit((xvec, yvec), basisfuns): 
    'Solves design matrix for approximating to basis functions' 
    return Mat([ map(form,xvec) for form in basisfuns ]).tr().solve(Vec(yvec)) 

def polyfit((xvec, yvec), degree=2): 
    'Solves Vandermonde design matrix for approximating polynomial coefficients' 
    return Mat([ [x**n for n in range(degree,-1,-1)] for x in xvec ]).solve(Vec(yvec)) 

def ratfit((xvec, yvec), degree=2): 
    'Solves design matrix for approximating rational polynomial coefficients (a*x**2 + b*x + c)/(d*x**2 + e*x + 1)' 
    return Mat([[x**n for n in range(degree,-1,-1)]+[-y*x**n for n in range(degree,0,-1)] for x,y in zip(xvec,yvec)]).solve(Vec(yvec)) 

def genmat(m, n, func): 
    if not n: n=m 
    return Mat([ [func(i,j) for i in range(n)] for j in range(m) ]) 

def zeroes(m=1, n=None): 
    'Zero matrix with side length m-by-m or m-by-n.' 
    return genmat(m,n, lambda i,j: 0) 

def eye(m=1, n=None): 
    'Identity matrix with side length m-by-m or m-by-n' 
    return genmat(m,n, lambda i,j: i==j) 

def hilb(m=1, n=None): 
    'Hilbert matrix with side length m-by-m or m-by-n. Elem[i][j]=1/(i+j+1)' 
    return genmat(m,n, lambda i,j: 1.0/(i+j+1.0)) 

def rand(m=1, n=None): 
    'Random matrix with side length m-by-m or m-by-n' 
    return genmat(m,n, lambda i,j: random.random()) 

if __name__ == '__main__': 
    import cmath 
    a = Table([1+2j,2,3,4]) 
    b = Table([5,6,7,8]) 
    C = Table([a,b]) 
    print 'a+b', a+b 
    print '2+a', 2+a 
    print 'a/5.0', a/5.0 
    print '2*a+3*b', 2*a+3*b 
    print 'a+C', a+C 
    print '3+C', 3+C 
    print 'C+b', C+b 
    print 'C.sum()', C.sum() 
    print 'C.map(math.cos)', C.map(cmath.cos) 
    print 'C.conjugate()', C.conjugate() 
    print 'C.real()', C.real() 

    print zeroes(3) 
    print eye(4) 
    print hilb(3,5) 

    C = Mat([[1,2,3], [4,5,1,], [7,8,9]]) 
    print C.mmul(C.tr()), '\n' 
    print C ** 5, '\n' 
    print C + C.tr(), '\n' 

    A = C.tr().augment(Mat([[10,11,13]]).tr()).tr() 
    q, r = A.qr() 
    assert q.mmul(r) == A 
    assert q.tr().mmul(q)==eye(3) 
    print 'q:\n', q, '\nr:\n', r, '\nQ.tr()&Q:\n', q.tr().mmul(q), '\nQ*R\n', q.mmul(r), '\n' 
    b = Vec([50, 100, 220, 321]) 
    x = A.solve(b) 
    print 'x: ', x 
    print 'b: ', b 
    print 'Ax: ', A.mmul(x) 

    inv = C.inverse() 
    print '\ninverse C:\n', inv, '\nC * inv(C):\n', C.mmul(inv) 
    assert C.mmul(inv) == eye(3) 

    points = (xvec,yvec) = funToVec(lambda x: math.sin(x)+2*math.cos(.7*x+.1), low=0, high=3, EqualSpacing=1) 
    basis = [lambda x: math.sin(x), lambda x: math.exp(x), lambda x: x**2] 
    print 'Func coeffs:', funfit(points, basis) 
    print 'Poly coeffs:', polyfit(points, degree=5) 
    points = (xvec,yvec) = funToVec(lambda x: math.sin(x)+2*math.cos(.7*x+.1), low=0, high=3) 
    print 'Rational coeffs:', ratfit(points) 

    print polyfit(([1,2,3,4], [1,4,9,16]), 2) 

    mtable = Vec([1,2,3]).outer(Vec([1,2])) 
    print mtable, mtable.size 

    A = Mat([ [2,0,3], [1,5,1], [18,0,6] ]) 
    print 'A:' 
    print A 
    print 'eigs:' 
    print A.eigs() 
    print 'Should be:', Vec([11.6158, 5.0000, -3.6158]) 
    print 'det(A)' 
    print A.det() 

    c = Mat([[1,2,30],[4,5,10],[10,80,9]])  # Failed example from Konrad Hinsen 
    print 'C:\n', c 
    print c.eigs() 
    print 'Should be:', Vec([-8.9554, 43.2497, -19.2943]) 

    A = Mat([ [1,2,3,4], [4,5,6,7], [2,1,5,0], [4,2,1,0] ]) # Kincaid and Cheney p.326 
    print 'A:\n', A 
    print A.eigs() 
    print 'Should be:', Vec([3.5736, 0.1765, 11.1055, -3.8556]) 

    A = rand(3) 
    q,r = A.qr() 
    s,t = A.qr() 
    print q is s    # Test caching 
    print r is t 
    A[1][1] = 1.1    # Invalidate the cache 
    u,v = A.qr() 
    print q is u    # Verify old result not used 
    print r is v 
    print u.mmul(v) == A  # Verify new result 

    print 'Test qr on 3x5 matrix' 
    a = rand(3,5) 
    q,r = a.qr() 
    print q.mmul(r) == a 
    print q.tr().mmul(q) == eye(3) 
+0

i collegamenti sono interrotti, i cuccioli piangono – Harvey

0
m=input("row") 
n=input("col") 
X=[] 
for i in range (m): 
    m1=[] 
    for j in range (n): 
     m1.append(input("num")) 
    X.append(m1) 
Y=[] 
for i in range (m): 
    n1=[] 
    for j in range (n): 
     n1.append(input("num")) 
    Y.append(n1) 


# result is 3x3 
result = [[0,0,0], 
     [0,0,0], 
     [0,0,0]] 

# iterate through rows of X 
for i in range(len(X)): 
    # iterate through columns of Y 
    for j in range(len(Y[0])): 
     # iterate through rows of Y 
     for k in range(len(Y)): 
      result[i][j] += X[i][k] * Y[k][j] 

for r in result: 
    print(r) 
1

moltiplicazione di matrici in puro Python.

def matmult(m1,m2): 
r=[] 
m=[] 
for i in range(len(m1)): 
    for j in range(len(m2[0])): 
     sums=0 
     for k in range(len(m2)): 
      sums=sums+(m1[i][k]*m2[k][j]) 
     r.append(sums) 
    m.append(r) 
    r=[] 
return m 
0

Tutte le risposte qui sotto sarebbero tornati la list.Your bisogno di convertirlo in matrice

def MATMUL(X, Y): 
    rows_A = len(X) 
    cols_A = len(X[0]) 
    rows_B = len(Y) 
    cols_B = len(Y[0]) 

    if cols_A != rows_B: 
     print "Matrices are not compatible to Multiply. Check condition C1==R2" 
     return 

    # Create the result matrix 
    # Dimensions would be rows_A x cols_B 
    C = [[0 for row in range(cols_B)] for col in range(rows_A)] 
    print C 

    for i in range(rows_A): 
     for j in range(cols_B): 
      for k in range(cols_A): 
       C[i][j] += A[i][k] * B[k][j] 

    C = numpy.matrix(C).reshape(len(A),len(B[0])) 

    return C 
1

in ritardo alla festa, ma quando ho dovuto fare un po 'di aritmetica matrice ho definito una nuova classe aiutare.All'interno di tale classe è possibile definire metodi come __add__ o, nel proprio caso d'uso, __matmul__, consentendo di eseguire a * b o a *= b anziché matrixMult(a,b).

Ho incluso un codice che implementa questa qui di seguito (esclusi il __init__ metodo di proibitivamente lungo, che crea essenzialmente una lista bidimensionale self.mat e una tupla self.order secondo quanto viene passato ad esso)

class matrix: 
    def __matmul__(self, multiplier): 
     if self.order[1] != multiplier.order[0]: 
      raise ValueError("The multiplier was non-conformable under multiplication.") 
     return [[sum(a*b for a,b in zip(srow,mcol)) for mcol in zip(*multiplier.mat)] for srow in self.mat] 

    def __imatmul__(self, multiplier): 
     self.mat = self @ multiplier 
     return self.mat 

    def __rmatmul__(self, multiplicand): 
     if multiplicand.order[1] != self.order[0]: 
      raise ValueError("The multiplier was non-conformable under multiplication.") 
     return [[sum(a*b for a,b in zip(mrow,scol)) for scol in zip(*self.mat)] for mrow in multiplicand.mat] 

Nota:

  • __rmatmul__ è necessario in modo che a @ b e b @ a sia il lavoro correttamente;
  • __imatmul__ è necessario per a @= b per funzionare correttamente;
  • Se una matrice non è conformabile rispetto alla moltiplicazione significa che non può essere moltiplicato, solitamente perché ha più o meno righe che ci sono colonne nella multiplicand

EDIT: Ho appena realizzato che __rmatmul__ non è richiesto in quanto viene utilizzato solo per valutare a @ b se a non ha metodo __matmul__. Dato che io sono solo moltiplicare le matrici da altre istanze di matrix, asarà hanno un metodo __matmul__, ma se mai modificarlo in modo che le operazioni di lavoro, per esempio, con le matrici 2D, ho bisogno di aggiungere nuovamente i __rmatmul__ in modo che list @ matrix opere così come matrix @ list.

+1

Per coloro che non lavorano sempre con le ultime versioni di Python: questo operatore di moltiplicazione di matrice è stato aggiunto in Python 3.5. – whydoubt

Problemi correlati