[PYTHON MATRIX-SIG] Matrix multiplication
carlos fonseca
fonseca@acse.shef.ac.uk
Wed, 19 Jun 1996 19:38:40 +0100
Hello,
I may have missed any previous discussion on the implementation of
matrix multiplication in Numerical Python. Anyway, I understand that
the current implementation of .matrixMultiply() requires arrays to be
2-dimensional (has this changed?).
While trying to port some matlab code to Numerical Python and make it
more general (eliminate for loops, add functionality, etc), I realised
that it would be nice to be able to:
1) treat (e.g.) the last 2 dimensions of an array as matrices, i.e.,
having arrays of matrices. For example,
A[m-by-n-by-p] * B[p-by-q] = A[m-by-n-by-p] * B[1-p-by-q]
= C[m-by-n-by-q]
In this case, each of the m n-by-p matrices in A is multiplied by B.
Another example would be
A[m-by-n-by-p] * D[m-p-by-q] = E[m-by-n-by-q]
Now, each of the m n-by-p matrices in A is multiplied by the corresponding
p-by-q matrix in B
2) be able to do the sort of thing Yorick does where, for example,
A[n-by-m-by-p] * B[p-by-q-by-r] = C[n-by-m-by-q-by-r]
This looks more like (2-dimensional) matrix mutiplication.
A simple way of achieving both would be to define a .dot(B,n) method
such that A.dot(B,n) would compute the dot product of the rows of A by
the corresponding rows of B (or any other dimension n, for that
matter), reducing the number of dimensions by one like add.reduce does
for a single array. The arrays should be conformable as for addition
and multiplication.
Essentially, A.dot(B,n) = add.reduce(A*B,n), but memory for the
intermediate result A*B would not need to be allocated (a saving when
casting of dimensions occurs, as in the example below).
For two-dimensional arrays, A.matrixMultiply(B) could be written
something like A[...,NewAxis].dot(B,-2). The automatic casting of the
dimension introduced into A to match that of B would produce the
desired results. Note that add.reduce(A[...,NewAxis]*B) would imply
storing a 3-dimensional array. This would generalize well for any
number of dimensions, I believe.
For example, case 1) above could be written
C = A[...,NewAxis].dot(B,-2)
E = A[...,NewAxis].dot(D[...,NewAxis,:,:],-2)
and case 2) could be written
A[...,NewAxis,NewAxis].dot(B,-3)
How do others feel about this issue?
Carlos Fonseca
P.S.: Why doesn't A[RubberIndex,NewAxis,NewAxis] work?
>>> a
0 1 2 3 4
5 6 7 8 9
10 11 12 13 14
>>> a[RubberIndex,NewAxis]
Traceback (innermost last):
File "<stdin>", line 1, in ?
ArrayError: each subindex must be either a slice, an integer,
RubberIndex, or None
^^^^^^^^^^^ !??
=================
MATRIX-SIG - SIG on Matrix Math for Python
send messages to: matrix-sig@python.org
administrivia to: matrix-sig-request@python.org
=================