多核机上单精度与双精度阵列矩阵乘法的性能退化

Question

更新

不幸的是，由于我的疏忽，我有一个旧版本的MKL（11.1）与numpy相关联。更新版本的MKL（11.3.1）在C语言中以及从python调用时都具有相同的性能。在

使事情变得模糊不清的是，即使将编译好的共享库显式地与新的MKL链接起来，并通过LD*variables指向它们，然后在python中执行import numpy，也会使python调用旧的MKL库。只有在python lib文件夹中替换所有libmkl_*.so，我才能够匹配python和C调用的性能。在

背景/库信息。

矩阵乘法是通过sgemm（单精度）和dgemm（双精度）Intel的MKL库调用，通过纽比.dot功能。库函数的实际调用可以用例如oprof进行验证。在

这里使用2x18核的CPU E5-2699 v3，因此总共有36个物理内核。 KMP_AFFINITY=分散。在linux上运行。在

TL；DR

1）为什么纽比.dot，即使它调用相同的MKL库函数，但与C编译的代码相比最多慢两倍？在

2）为什么选择via纽比.dot随着核心数量的增加，性能会下降，而在C代码中（调用相同的库函数）却没有观察到相同的效果。在

问题

我观察到，做单精度/双精度的矩阵乘法纽比.dot，以及直接从编译的C共享库调用cblas sgemm/dgemm，与从纯C代码内部调用相同的MKL cblas sgemm/dgemm函数相比，性能明显较差。在

import numpy as np
import mkl
n = 10000
A = np.random.randn(n,n).astype('float32')
B = np.random.randn(n,n).astype('float32')
C = np.zeros((n,n)).astype('float32')

mkl.set_num_threads(3); %time np.dot(A, B, out=C)
11.5 seconds
mkl.set_num_threads(6); %time np.dot(A, B, out=C)
6 seconds
mkl.set_num_threads(12); %time np.dot(A, B, out=C)
3 seconds
mkl.set_num_threads(18); %time np.dot(A, B, out=C)
2.4 seconds
mkl.set_num_threads(24); %time np.dot(A, B, out=C)
3.6 seconds
mkl.set_num_threads(30); %time np.dot(A, B, out=C)
5 seconds
mkl.set_num_threads(36); %time np.dot(A, B, out=C)
5.5 seconds

执行与上述完全相同的操作，但使用双精度A、B和C，可以得到： 3芯20s，6芯10s，12芯5s，18芯4.3s，24芯3s，30芯2.8s，36芯2.8s

单精度浮点的速度增加似乎与缓存未命中有关。 对于28核运行，这里是perf的输出。 对于单精度：

^{pr2}$

双精度：

93,087,703 cache-misses # 5.164 % of all cache refs

C共享库，用

/opt/intel/bin/icc -o comp_sgemm_mkl.so -openmp -mkl sgem_lib.c -lm -lirc -O3 -fPIC -shared -std=c99 -vec-report1 -xhost -I/opt/intel/composer/mkl/include

#include <stdio.h>
#include <stdlib.h>
#include "mkl.h"

void comp_sgemm_mkl(int m, int n, int k, float *A, float *B, float *C);

void comp_sgemm_mkl(int m, int n, int k, float *A, float *B, float *C)
{
    int i, j;
    float alpha, beta;
    alpha = 1.0; beta = 0.0;

    cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,
                m, n, k, alpha, A, k, B, n, beta, C, n);
}

Python包装器函数，调用上面编译的库：

def comp_sgemm_mkl(A, B, out=None):
    lib = CDLL(omplib)
    lib.cblas_sgemm_mkl.argtypes = [c_int, c_int, c_int, 
                                 np.ctypeslib.ndpointer(dtype=np.float32, ndim=2), 
                                 np.ctypeslib.ndpointer(dtype=np.float32, ndim=2),
                                 np.ctypeslib.ndpointer(dtype=np.float32, ndim=2)]
    lib.comp_sgemm_mkl.restype = c_void_p
    m = A.shape[0]
    n = B.shape[0]
    k = B.shape[1]
    if np.isfortran(A):
        raise ValueError('Fortran array')
    if m != n:
        raise ValueError('Wrong matrix dimensions')
    if out is None:
        out = np.empty((m,k), np.float32)
    lib.comp_sgemm_mkl(m, n, k, A, B, out)

然而，从一个调用MKL的cblas_sgemm/cblas_ggemm的C编译二进制文件中显式调用，以及在C中通过malloc分配的数组，与python代码相比，性能提高了近2倍，即纽比.dot打电话。此外，还没有观察到随核心数量增加而导致性能下降的影响。单精度矩阵乘法的最佳性能为900ms，通过mkl_set_num_cores使用所有36个物理内核，并使用numactl--interleave=all运行C代码时，性能达到最佳。在

也许有什么新奇的工具或建议来进一步分析/检查/了解这种情况？任何阅读材料也很受欢迎。在

更新 遵循@Hristo Iliev的建议，运行numactl--interleave=all./ipython没有更改计时（在noise中），但改进了纯C二进制运行时。在

Answer 1

我怀疑这是由于不幸的线程调度。我能复制出一个和你相似的效果。Python的运行速度是2.2秒，而C版本的运行速度从1.4到2.2秒有很大的变化

申请： KMP_AFFINITY=scatter,granularity=thread 这样可以确保28个线程始终在同一个处理器线程上运行。在

将这两个运行时都减少到更稳定的状态：C为~1.24 s，python为~1.26 s。在

这是一个28核双插座Xeon E5-2680 v3系统。在

有趣的是，在一个非常相似的24核双套接字Haswell系统上，python和C的性能几乎完全相同，即使没有线程亲和力/固定。在

为什么python会影响调度？我假设有更多的运行时环境。底线是，没有固定的性能结果将是不确定的。在

另外，您需要考虑的是，“英特尔OpenMP运行时”会产生一个额外的管理线程，这会使调度程序感到困惑。对于固定有更多的选择，例如KMP_AFFINITY=compact-但由于某些原因，这在我的系统上完全是一团糟。您可以将,verbose添加到变量中，以查看运行时如何固定线程。在

likwid-pin是一个有用的替代方案，提供了更方便的控制。在

一般来说，单精度应至少与双精度一样快。双精度可能较慢，因为：

• 您需要更多的内存/缓存带宽来实现双精度。在
• 您可以构建单精度吞吐量更高的ALU，但这通常不适用于CPU，而是适用于GPU。在

我认为，一旦你摆脱了性能异常，这将反映在你的数字。在

当您扩大MKL/*gemm的线程数量时，请考虑

• 内存/共享缓存带宽可能成为瓶颈，限制可伸缩性
• Turbo模式将在提高利用率的同时有效地降低核心频率。这一点即使在正常频率下运行也适用：在Haswell EP处理器上，AVX指令将施加较低的“AVX基频”——但当使用较少的内核/可用的热余量时，处理器允许超过该频率，而且通常在短时间内甚至会更高。如果你想要完美的中性效果，你必须使用AVX的基频，这是1.9ghz。它被记录在here，并在one picture中解释。在

我不认为有一个真正简单的方法来衡量你的应用程序是如何受到不良调度的影响的。你可以用perf trace -e sched:sched_switch来展示这个，还有一个{a4}来可视化它，但是这会带来一个很高的学习曲线。再说一次-对于并行性能分析，您应该将线程固定。在