How does the performance of Python/Numpy array operations scale with increasing array dimensions?

Question

How do Python/Numpy arrays scale with increasing array dimensions?

This is based on some behaviour I noticed while benchmarking Python code for this question: How to express this complicated expression using numpy slices

The problem mostly involved indexing to populate an array. I found that the advantages of using (not-very-good) Cython and Numpy versions over a Python loop varied depending on the size of the arrays involved. Both Numpy and Cython experience an increasing performance advantage up to a point (somewhere broadly around $N=500$ for Cython and $N=2000$ for Numpy on my laptop), after which their advantages declined (the Cython function remained the fastest).

Is this hardware defined? In terms of working with large arrays, what are best practices that one should adhere to for code where performance is appreciated?

This question (Why isn't my Matrix-Vector Multiplication Scaling?) may be related, but I am interested in knowing more about how different ways of treating arrays in Python scale relative to each other.

Answer 1

I don't know how this benchmark was done, but probably with floating point numbers, that in Python default to doubles. The sizes correspond, respectively, to $4$ and $16 kB$. Those are reasonable values for L1 and L2 cache sizes of a (a bit old) AMD CPU. To be certain, I did my own benchmark:

def timeit(size):
     t0 = time.time()
     for _ in xrange(10):
         np.random.random(size)
     return time.time() - t0

sizes = np.logspace(1, 6, 40)
times = [timeit(s) for s in sizes]

There are a few things wrong with this benchmark, for a start, I am not disabling garbage collection and I am taking the sum, not the best time, but bear with me.

The time it takes is quite proportional to the size of the array, but there is a change in the slope around size $8000$. This is an array of $64 kB$, that is the size of the L1 cache of an i5 (my computer).

Should one worry about the cache size? As a general rule, I say no. Optimizing for it in Python means making the code much more complicated, for dubious performance gains. Don't forget that Python objects add several overheads that are difficult to track and predict. I can only think of two cases where this is a relevant factor:

• Basic operations on large arrays (like evaluate a polynomial), limited by memory bandwidth. Use Numexpr or (if the data is much bigger) Pytables. They are optimised to take the cache size into account amongst other optimisations.
• Performance critical code: if you want to squeeze every microsecond, you should not be using Python in the first place. Writing vectorized Cython and leaving the compiler do what it does best is probably the painless way to go.

In the comments, Evert mentioned CArray. Note that, even working, the development has stopped and it has been abandoned as a standalone project. The functionality will be included in Blaze instead, an ongoing project to make a "new generation Numpy".