基准测试（Python vs. C++ 使用 BLAS）和（Numpy）

Question

我想写一个程序，广泛使用BLAS和LAPACK这两个线性代数库的功能。因为性能很重要，我做了一些基准测试，想知道我采用的方法是否合理。

我有三个“选手”，想通过简单的矩阵乘法来测试它们的性能。这三个选手是：

1. Numpy，只使用dot这个功能。
2. Python，通过共享对象调用BLAS的功能。
3. C++，同样通过共享对象调用BLAS的功能。

场景

我实现了不同维度i的矩阵乘法。i的范围是从5到500，每次增加5，矩阵m1和m2的设置如下：

m1 = numpy.random.rand(i,i).astype(numpy.float32)
m2 = numpy.random.rand(i,i).astype(numpy.float32)

1. Numpy

使用的代码如下：

tNumpy = timeit.Timer("numpy.dot(m1, m2)", "import numpy; from __main__ import m1, m2")
rNumpy.append((i, tNumpy.repeat(20, 1)))

2. Python，通过共享对象调用BLAS

使用的函数是

_blaslib = ctypes.cdll.LoadLibrary("libblas.so")
def Mul(m1, m2, i, r):

    no_trans = c_char("n")
    n = c_int(i)
    one = c_float(1.0)
    zero = c_float(0.0)

    _blaslib.sgemm_(byref(no_trans), byref(no_trans), byref(n), byref(n), byref(n), 
            byref(one), m1.ctypes.data_as(ctypes.c_void_p), byref(n), 
            m2.ctypes.data_as(ctypes.c_void_p), byref(n), byref(zero), 
            r.ctypes.data_as(ctypes.c_void_p), byref(n))

测试代码如下：

r = numpy.zeros((i,i), numpy.float32)
tBlas = timeit.Timer("Mul(m1, m2, i, r)", "import numpy; from __main__ import i, m1, m2, r, Mul")
rBlas.append((i, tBlas.repeat(20, 1)))

3. C++，通过共享对象调用BLAS

由于C++的代码自然会长一些，所以我把信息简化到最小。
我用以下方式加载函数：

void* handle = dlopen("libblas.so", RTLD_LAZY);
void* Func = dlsym(handle, "sgemm_");

我用gettimeofday来测量时间，方法如下：

gettimeofday(&start, NULL);
f(&no_trans, &no_trans, &dim, &dim, &dim, &one, A, &dim, B, &dim, &zero, Return, &dim);
gettimeofday(&end, NULL);
dTimes[j] = CalcTime(start, end);

其中j是一个循环，运行20次。我用以下方式计算经过的时间：

double CalcTime(timeval start, timeval end)
{
double factor = 1000000;
return (((double)end.tv_sec) * factor + ((double)end.tv_usec) - (((double)start.tv_sec) * factor + ((double)start.tv_usec))) / factor;
}

结果

结果在下面的图中展示：

问题

1. 你觉得我的方法公平吗？有没有什么不必要的开销可以避免？
2. 你是否预期C++和Python的结果会有如此大的差异？它们都是通过共享对象进行计算的。
3. 因为我更想用Python来写我的程序，有什么办法可以提高调用BLAS或LAPACK例程时的性能？

下载

完整的基准测试可以在这里下载。（J.F. Sebastian提供了这个链接^^）

Answer 1

这里有另一个基准测试（在Linux上，只需输入 make）：http://dl.dropbox.com/u/5453551/blas_call_benchmark.zip

http://dl.dropbox.com/u/5453551/blas_call_benchmark.png

我发现对于大矩阵，Numpy、Ctypes和Fortran这几种方法之间几乎没有什么区别。（Fortran代替C++——如果这很重要，你的基准测试可能有问题。）

你在C++中的 CalcTime 函数似乎有个符号错误。... + ((double)start.tv_usec)) 应该改成 ... - ((double)start.tv_usec))。 也许你的基准测试还有其他错误，比如在不同的BLAS库之间比较，或者在不同的BLAS设置（比如线程数量）之间比较，或者在实际时间和CPU时间之间比较？

编辑：在 CalcTime 函数中没有正确计算括号——没关系。

作为一个指导原则：如果你做基准测试，请务必把所有代码都发出来。评论基准测试，尤其是当结果令人惊讶时，如果没有完整的代码，通常是没有意义的。

---

要找出Numpy链接的是哪个BLAS，可以执行：

$ python
Python 2.7.2+ (default, Aug 16 2011, 07:24:41) 
[GCC 4.6.1] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> import numpy.core._dotblas
>>> numpy.core._dotblas.__file__
'/usr/lib/pymodules/python2.7/numpy/core/_dotblas.so'
>>> 
$ ldd /usr/lib/pymodules/python2.7/numpy/core/_dotblas.so
    linux-vdso.so.1 =>  (0x00007fff5ebff000)
    libblas.so.3gf => /usr/lib/libblas.so.3gf (0x00007fbe618b3000)
    libc.so.6 => /lib/x86_64-linux-gnu/libc.so.6 (0x00007fbe61514000)

更新：如果你无法导入 numpy.core._dotblas，说明你的Numpy正在使用它的内部备用BLAS版本，这个版本比较慢，不适合用于性能计算！下面@Woltan的回复表明，这就是他/她在Numpy和Ctypes+BLAS之间看到的差异的原因。

要解决这个问题，你需要安装ATLAS或MKL——请查看这些说明：http://scipy.org/Installing_SciPy/Linux 大多数Linux发行版都自带ATLAS，所以最好的选择是安装他们的 libatlas-dev 包（名称可能会有所不同）。

Answer 2

更新（2014年7月30日）：

我在我们的新高性能计算机（HPC）上重新进行了基准测试。硬件和软件的配置与之前的回答有所不同。

我把结果放在了一个谷歌表格中（也包含了原始回答的结果）。

硬件

我们的高性能计算机有两个不同的节点，一个是使用Intel Sandy Bridge处理器，另一个是使用更新的Ivy Bridge处理器：

Sandy（MKL，OpenBLAS，ATLAS）：

• CPU: 2个16核的Intel(R) Xeon(R) E2560 Sandy Bridge @ 2.00GHz（16个核心）
• 内存: 64 GB

Ivy（MKL，OpenBLAS，ATLAS）：

• CPU: 2个20核的Intel(R) Xeon(R) E2680 V2 Ivy Bridge @ 2.80GHz（20个核心，开启超线程后为40个核心）
• 内存: 256 GB

软件

两个节点的软件配置是一样的。使用OpenBLAS代替了GotoBLAS2，并且还有一个设置为8线程的多线程ATLAS BLAS。

• 操作系统: Suse
• Intel编译器: ictce-5.3.0
• Numpy: 1.8.0
• OpenBLAS: 0.2.6
• ATLAS:: 3.8.4

点积基准测试

基准测试的代码与下面的相同。不过在新机器上，我还进行了5000和8000大小矩阵的基准测试。
下面的表格包含了原始回答的基准测试结果（重命名：MKL -> Nehalem MKL，Netlib Blas -> Nehalem Netlib BLAS，等等）。

单线程性能：

多线程性能（8线程）：

线程数与矩阵大小（Ivy Bridge MKL）：

基准测试套件

单线程性能：

多线程（8线程）性能：

结论

新的基准测试结果与原始回答中的结果相似。OpenBLAS和MKL的性能大致相当，唯一的例外是特征值测试。特征值测试在OpenBLAS的单线程模式下表现尚可，但在多线程模式下性能较差。

“矩阵大小与线程数图表”也显示，虽然MKL和OpenBLAS在核心/线程数量上通常表现良好，但这取决于矩阵的大小。对于小矩阵，增加更多核心对性能的提升并不明显。

从Sandy Bridge到Ivy Bridge的性能提升约为30%，这可能是由于更高的时钟频率（+0.8 GHz）和/或更好的架构。

---

原始回答（2011年10月4日）：

之前我需要优化一些用Python编写的线性代数计算/算法，这些算法使用了numpy和BLAS，所以我对不同的numpy/BLAS配置进行了基准测试。

具体测试了：

• Numpy与ATLAS
• Numpy与GotoBlas2（1.13）
• Numpy与MKL（11.1/073）
• Numpy与加速框架（Mac OS X）

我进行了两个不同的基准测试：

1. 不同大小矩阵的简单点积
2. 可以在这里找到的基准测试套件。

以下是我的结果：

机器

Linux（MKL，ATLAS，无MKL，GotoBlas2）：

• 操作系统: Ubuntu Lucid 10.4 64位。
• CPU: 2个4核的Intel(R) Xeon(R) E5504 @ 2.00GHz（8个核心）
• 内存: 24 GB
• Intel编译器: 11.1/073
• Scipy: 0.8
• Numpy: 1.5

Mac Book Pro（加速框架）：

• 操作系统: Mac OS X Snow Leopard（10.6）
• CPU: 1个Intel Core 2 Duo 2.93 GHz（2个核心）
• 内存: 4 GB
• Scipy: 0.7
• Numpy: 1.3

Mac Server（加速框架）：

• 操作系统: Mac OS X Snow Leopard Server（10.6）
• CPU: 4个Intel(R) Xeon(R) E5520 @ 2.26 GHz（8个核心）
• 内存: 4 GB
• Scipy: 0.8
• Numpy: 1.5.1

点积基准测试

代码:

import numpy as np
a = np.random.random_sample((size,size))
b = np.random.random_sample((size,size))
%timeit np.dot(a,b)

结果:

    System        |  size = 1000  | size = 2000 | size = 3000 |
netlib BLAS       |  1350 ms      |   10900 ms  |  39200 ms   |    
ATLAS (1 CPU)     |   314 ms      |    2560 ms  |   8700 ms   |     
MKL (1 CPUs)      |   268 ms      |    2110 ms  |   7120 ms   |
MKL (2 CPUs)      |    -          |       -     |   3660 ms   |
MKL (8 CPUs)      |    39 ms      |     319 ms  |   1000 ms   |
GotoBlas2 (1 CPU) |   266 ms      |    2100 ms  |   7280 ms   |
GotoBlas2 (2 CPUs)|   139 ms      |    1009 ms  |   3690 ms   |
GotoBlas2 (8 CPUs)|    54 ms      |     389 ms  |   1250 ms   |
Mac OS X (1 CPU)  |   143 ms      |    1060 ms  |   3605 ms   |
Mac Server (1 CPU)|    92 ms      |     714 ms  |   2130 ms   |

基准测试套件

代码:
有关基准测试套件的更多信息，请参见这里。

结果:

    System        | eigenvalues   |    svd   |   det  |   inv   |   dot   |
netlib BLAS       |  1688 ms      | 13102 ms | 438 ms | 2155 ms | 3522 ms |
ATLAS (1 CPU)     |   1210 ms     |  5897 ms | 170 ms |  560 ms |  893 ms |
MKL (1 CPUs)      |   691 ms      |  4475 ms | 141 ms |  450 ms |  736 ms |
MKL (2 CPUs)      |   552 ms      |  2718 ms |  96 ms |  267 ms |  423 ms |
MKL (8 CPUs)      |   525 ms      |  1679 ms |  60 ms |  137 ms |  197 ms |  
GotoBlas2 (1 CPU) |  2124 ms      |  4636 ms | 147 ms |  456 ms |  743 ms |
GotoBlas2 (2 CPUs)|  1560 ms      |  3278 ms | 116 ms |  295 ms |  460 ms |
GotoBlas2 (8 CPUs)|   741 ms      |  2914 ms |  82 ms |  262 ms |  192 ms |
Mac OS X (1 CPU)  |   948 ms      |  4339 ms | 151 ms |  318 ms |  566 ms |
Mac Server (1 CPU)|  1033 ms      |  3645 ms |  99 ms |  232 ms |  342 ms |

安装

安装MKL需要安装完整的Intel编译器套件，这个过程相对简单。不过，由于一些错误/问题，配置和编译支持MKL的numpy有点麻烦。

GotoBlas2是一个小包，可以很容易地编译成共享库。不过，由于一个错误，你需要在构建后重新创建共享库才能与numpy一起使用。
此外，出于某种原因，针对多个目标平台的构建并没有成功。所以我必须为每个平台创建一个.so文件，以便获得优化后的libgoto2.so文件。

如果你从Ubuntu的仓库安装numpy，它会自动安装并配置numpy以使用ATLAS。从源代码安装ATLAS可能需要一些时间，并且需要一些额外的步骤（如Fortran等）。

如果你在Mac OS X机器上使用Fink或Mac Ports安装numpy，它会将numpy配置为使用ATLAS或Apple的加速框架。你可以通过运行ldd命令查看numpy.core._dotblas文件，或者调用numpy.show_config()来检查。

结论

MKL的性能最好，其次是GotoBlas2。
在特征值测试中，GotoBlas2的表现出乎意料地差，不知道为什么会这样。
Apple的加速框架在单线程模式下表现非常好（与其他BLAS实现相比）。

GotoBlas2和MKL在线程数量上扩展得很好。因此，如果你需要处理大矩阵，使用多线程会有很大帮助。

无论如何，不要使用默认的netlib blas实现，因为它对于任何严肃的计算工作来说都太慢了。

在我们的集群上，我还安装了AMD的ACML，性能与MKL和GotoBlas2相似，但我没有具体的数据。

我个人建议使用GotoBlas2，因为它更容易安装并且是免费的。

如果你想用C++/C编程，也可以看看Eigen3，它在某些情况下的性能应该优于MKL/GotoBlas2，而且使用起来也相对简单。

Answer 3

我运行了你的基准测试。在我的机器上，C++和numpy之间没有区别：

你觉得我的做法公平吗？还是说有一些不必要的开销可以避免？

看起来是公平的，因为结果没有差别。

你会期待C++和Python之间的结果会有这么大的差异吗？它们在计算时都使用了共享对象。

不会。

既然我更愿意用Python来写我的程序，那我可以做些什么来提高调用BLAS或LAPACK函数的性能呢？

确保numpy在你的系统上使用的是优化过的BLAS/LAPACK库版本。