Sage
Sage

Reputation: 53

Matrix multiplication using slicing. Python

I have the following code:

from numpy import *
a = random.rand(3,4)
b = random.rand(4,2)
c = linspace(0,0,6)
c.shape = (3,2)
for i in range(a.shape[0]): 
  for j in range(b.shape[1]): 
      for k in range(b.shape[0]): 
         c[i][j] += a[i][k] * b[k][j]
for r in c: 
  print "C = ", r

But I need to change the last (innermost) loop and I need to use slicing. As I understood I have to do something like this:

for i in range(a.shape[0]): 
 for j in range(b.shape[1]): 
     c[i][j] += a[i][0:l-1] * b[0:l-1][j]

But, unfortunately, it doesn't work. Could somebody help me and give a hint how to do it?

Upvotes: 0

Views: 2079

Answers (2)

gboffi
gboffi

Reputation: 25043

Let's start with an helper function that creates a ListOfLists (a lol) of r rows, each composed of c columns:

In [1]: def lol(r,c): return [[i*c+j for j in range(c)] for i in range(r)]

and create two list of lists

In [2]: a = lol(2,5) ; b = lol(5,4)

We want to verify that our code below, meant to do a matrix product using two lol's is working correctly, so we create two ndarrays from a and b and form their inner product

In [3]: from numpy import array

In [4]: aa = array(a) ; ab = array(b) ; aa.dot(ab)
Out[4]: 
array([[120, 130, 140, 150],
       [320, 355, 390, 425]])

Now, we can test our code for the inner, or matrix, product of two lols

In [5]: [[sum(x*y for x, y in zip(ar,bc)) for bc in zip(*b)] for ar in a]
Out[5]: [[120, 130, 140, 150], [320, 355, 390, 425]]

It's ok, isn't it? (I have to say that in a first iteration of the code I came out with the transpose of the result...).

Now that we have a bit of confidence, let's try with something more substantial

In [6]: a = lol(200,50) ; b = lol(50,400)

In [7]: aa = array(a) ; ab = array(b)

In [8]: %timeit c = aa.dot(ab)
100 loops, best of 3: 4.53 ms per loop

In [9]: %timeit c = [[sum(x*y for x, y in zip(ar,bc)) for bc in zip(*b)] for ar in a]
1 loops, best of 3: 469 ms per loop

As you can see, numpy is two order of magnitude faster than operating on lists, but it is more interesting, in the context of the OP question, to try our list code on the ndarrays:

In [10]: %timeit c = [[sum(x*y for x, y in zip(ar,bc)) for bc in zip(*ab)] for ar in aa] 
1 loops, best of 3: 1.32 s per loop

Oh, if you have numpy arrays it's way better to use the array methods rather than operating on the individual elements... but wait, we have a faster alternative to the inner zip:

In [11]: %timeit c = [[sum(x*y for x, y in zip(ar,bc)) for bc in ab.T] for ar in aa]
1 loops, best of 3: 1.34 s per loop

In [12]: 

no, even if we use the transpose property of a ndarray we have the same results.

Summary: never ever use individually accessed numpy array's elements to do heavy computational lifting...


Thank you to ipython and its %timeit magic that made this rambling easier and fun (for me).

Upvotes: 0

Marcus Müller
Marcus Müller

Reputation: 36402

What you're trying to do here is the dot product of the row vector from a, and the column vector from b:

c[i][j] += a[i][0:l-1] * b[0:l-1][j]

which would be

c[i][j] = np.dot(a[i], b[:][j]) 

which is the same as

sum([a_*b_ for a_,b_ in zip(a[i],b[:][j])])

or

sum(a[i]*b[:][j])

but faster.

However, if you're using np.dot, anyways:

c = np.dot(a,b)

is definitely faster.

Upvotes: 1

Related Questions