UPDATE
К сожалению, из-за моего недосмотра у меня была более старая версия MKL (11.1), связанная с numpy. Более новая версия MKL (11.3.1) дает такую же производительность в C и при вызове из python.
То, что было скрыто, заключалось даже в том, что связь с компилируемыми разделяемыми библиотеками явно с новым MKL и указанием на них переменных LD_ *, а затем в python, выполняющих импорт numpy, каким-то образом превращала python в старые библиотеки MKL. Только заменяя в python lib папку все libmkl _ *. Так что с более новым MKL я смог сопоставить производительность в вызовах python и C.
Информация о фоновом режиме/библиотеке.
Матричное умножение выполнялось через вызовы библиотеки sgemm (single-precision) и dgemm (с двойной точностью) Intel MKL, через функцию numpy.dot. Фактический вызов функций библиотеки можно проверить, например. oprof.
Используя здесь 2x18 core CPU E5-2699 v3, следовательно, в общей сложности 36 физических ядер. KMP_AFFINITY = разброс. Работает на Linux.
TL; DR
1) Почему numpy.dot, хотя он вызывает одни и те же функции библиотеки MKL, в два раза медленнее в сравнении с C скомпилированным кодом?
2) Почему с помощью numpy.dot вы получаете снижение производительности с увеличением количества ядер, тогда как тот же эффект не наблюдается в коде C (вызывающий одни и те же функции библиотеки).
Проблема
Я заметил, что выполнение матричного умножения одиночных/двойных прецизионных поплавков в numpy.dot, а также вызов cblas_sgemm/dgemm непосредственно из скомпилированной C общей библиотеки дают заметно худшую производительность по сравнению с вызовом те же функции MKL cblas_sgemm/dgemm изнутри чистого кода C.
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 ядра: 20 с, 6 ядер: 10 с, 12 ядер: 5 с, 18 ядер: 4,3 с, 24 ядра: 3 с, 30 ядер: 2,8 с, 36 ядер: 2,8 с.
Поглощение скорости для плавающих точек с одинарной точностью, похоже, связано с промахами кеша. Для 28-го запуска, здесь выведено перформанс. Для одиночной точности:
perf stat -e task-clock,cycles,instructions,cache-references,cache-misses ./ptestf.py
631,301,854 cache-misses # 31.478 % of all cache refs
И двойная точность:
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)
Однако явные вызовы с C-компилируемого двоичного вызова MKL cblas_sgemm/cblas_dgemm с массивами, выделенными через malloc в C, дают почти 2x лучшую производительность по сравнению с кодом python, т.е. вызов numpy.dot. Кроме того, не наблюдается эффекта снижения производительности при увеличении количества ядер. Наибольшая производительность составила 900 мс для однократного матричного умножения и была достигнута при использовании всех 36 физических ядер через mkl_set_num_cores и запуска кода C с numactl --interleave = all.
Возможно, какие-нибудь причудливые инструменты или советы для профилирования/проверки/понимания этой ситуации? Любой материал для чтения также очень ценится.
UPDATE Следуя совету @Hristo Iliev, запуск numactl --interleave = all./ipython не менял таймингов (в пределах шума), но улучшает чистые времена выполнения C.