3

I am trying to efficiently index a 2D array in Python and have the problem that it is really slow.

This is what I tried (simplified example):

xSize = veryBigNumber
ySize = veryBigNumber
a = np.ones((xSize,ySize))
N = veryBigNumber  
const = 1

for t in range(N):

  for i in range(xSize):
    for j in range(ySize):
      a[i,j] *= f(i,j)*const   # f(i,j) is an arbitrary function of i and j.

Now I would like to substitute the nested loop by something more efficient. How do I do this?

Improve this question
2
  • You need to tell us much more about what you are trying to do - should be able to avoid declaring a huge array for a simple formula. If you really need "an array whose size is close to the maximal available RAM", that's pretty restricting and you should state that requirement. You might also look into disk-backed data structures. And take a look at pandas. Commented Nov 25, 2014 at 10:05
  • @smci: of coarse it is used later on! The above piece of code would not make any sense if not, would it? Commented Mar 16, 2015 at 15:30

3 Answers 3

Reset to default
5

Your 2D array could be produced using the following addition:

np.arange(200)[:,np.newaxis] + np.arange(200)

This type of vectorised operation is likely to be very fast:

>>> %timeit np.arange(200)[:,np.newaxis] + np.arange(200)
1000 loops, best of 3: 178 µs per loop

This method in not limited to addition. We can use the two arrays in the above operation as the arguments of any universal function (commonly abbreviated to ufunc).

For example:

>>> np.multiply(np.arange(5)[:,np.newaxis], np.arange(5))
array([[ 0,  0,  0,  0,  0],
       [ 0,  1,  2,  3,  4],
       [ 0,  2,  4,  6,  8],
       [ 0,  3,  6,  9, 12],
       [ 0,  4,  8, 12, 16]])

NumPy has built in ufuncs for all the basic arithmetic operations and some more interesting ones too. If you need a more exotic function, NumPy allows you to make your own ufunc.

Edit: To quickly explain the broadcasting happening in this method; you can think of it like this...

np.arange(5) produces 1D array which looks like this:

array([0, 1, 2, 3, 4])

The code np.arange(5)[:,np.newaxis] adds a second dimension (columns) to the range, producing this 2D array:

array([[0],
       [1],
       [2],
       [3],
       [4]])

To create the final 5x5 array using np.multiply (although we could use any ufunc or binary arithmetic operation), NumPy takes the 0 in the second array and mutliplies it with each elements it the first array making a row like this:

[ 0,  0,  0,  0,  0]

It then takes the second element in the second array, 1, and multiplies it with the first array, producing this row:

[ 0,  1,  2,  3,  4]

This continues until we have the final 5x5 matrix.

Improve this answer
Sign up to request clarification or add additional context in comments.

Comments

4

You could use the indices routine:

b=np.indices(a.shape)
a=b[0]+b[1]

Timings:

%%timeit
    ...: b=np.indices(a.shape)
    ...: c=b[0]+b[1]
1000 loops, best of 3: 370 µs per loop


%%timeit
for i in range(200):
  for j in range(200):
     a[i,j] = i + j

100 loops, best of 3: 10.4 ms per loop
Improve this answer

1 Comment

Ok, this solution works. However, it needs three times the memory neccessary for completing the task. What I will need is that the size of array "a" is close to the maximal available RAM. Then, I think this solution will not work. Sorry, cannot upvote your answer because I have not enough reputation :(
1

Since your output matrix a is the element-wise power of N of a matrix F with elements f_ij = f(i,j) * const your code can simplify to

F = np.empty((xSize, ySize))
for i in range(xSize):
    for j in range(ySize):
        F[i,j] = f(i,j) * const

a = F ** n

For even more speed you can exchange the creation of the F matrix with something more efficient, given that the function f(i,j) is vectorized:

xmap, ymap = numpy.meshgrid(range(xSize), range(ySize))
F = f(xmap, ymap) * const
Improve this answer

Comments