Numpy EIG approximately 4x slower than MATLAB EIG on Mac OS X 10.6. What am I doing wrong?

Question

I tried profiling the EIG function on MATLAB and NumPy to compare performance on my Macbook Pro (2 GHz, Quad core i7 running OS X 10.6). NumPy EIG appears to be quite slow compared to MATLAB.

Here's the code I profiled on NumPy:

s = '''\
x = numpy.random.random((2000,2000));
numpy.linalg.eig(x);
'''
t = timeit.Timer(stmt=s,setup="import numpy")
result = t.repeat(repeat=5,number=10)
result
Out[22]: 
[197.1737039089203,
 197.34872913360596,
 196.8160741329193,
 197.94081807136536,
 194.5740351676941]

That's around 19.5 secs/exec in NumPy.

Here's the same code in MATLAB:

clear all
tic;
for i = 1:50
    x = rand(2000,2000);
    eig(x);
end
toc;
Elapsed time is 267.976645 seconds.

That's around 5.36 secs/exec on MATLAB.

I guess something as simple as this should not depend much on JIT performance, so it probably boils down to BLAS and the routines that access the BLAS libraries. I know that MATLAB uses the Accelerate Framework on the Mac.

NumPy also appears to use the Accelerate Framework BLAS on my Macbook Pro; here's the output of numpy.show_config()

numpy.show_config()
lapack_opt_info:
    extra_link_args = ['-Wl,-framework', '-Wl,Accelerate']
    extra_compile_args = ['-msse3']
    define_macros = [('NO_ATLAS_INFO', 3)]
blas_opt_info:
    extra_link_args = ['-Wl,-framework', '-Wl,Accelerate']
    extra_compile_args = ['-msse3', '-I/System/Library/Frameworks/vecLib.framework/Headers']
    define_macros = [('NO_ATLAS_INFO', 3)]

I'm using Python 2.7.2 and NumPy 1.6 (both installed from MacPorts)

So here's my question to the NumPy folks: Why is NumPy slower in this instance? Have I left out some optimization while installing NumPy?

Answer 1

As far as I know, the MATLAB uses MKL libraries as BLAS, not the Accelerate Framework. My experience tells me, that Accelerate is significantly slower than MKL. To check it, you can try to get the academic version of the Enthought Python Distribution (EPD), where Numpy is compiled against MKL, and compare these timings. Moreover, by default, MATLAB uses all the threads (try running in single-threaded mode), but Numpy is not. In the EPD, it can be done with running

import mkl 
mkl.set_num_threads(4)