Подтвердить что ты не робот

Деградация производительности матричного умножения одиночных и двухточечных массивов на многоядерной машине

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.

4b9b3361

Ответ 1

Я подозреваю, что это связано с неудачным планированием потоков. Я смог воспроизвести эффект, подобный вашему. Python работал на ~ 2.2 с, а версия С показывала огромные изменения от 1.4-2.2 с.

Применение: KMP_AFFINITY=scatter,granularity=thread Это гарантирует, что 28 потоков всегда работают на одном и том же потоке процессора.

Сокращает обе время работы до более стабильной ~ 1,24 с для C и ~ 1,26 с для python.

Это на 28-ядерном двухъядерном процессоре Xeon E5-2680 v3.

Интересно, что на очень похожей 24-ядерной системе с двумя гнездами Haswell, как python, так и C выполняют почти идентичные даже без аффинности/привязки нитей.

Почему python влияет на планирование? Хорошо, я предполагаю, что вокруг есть среда выполнения. Итог, без фиксации результатов вашей работы, будет недетерминированным.

Также вам нужно учитывать, что среда Intel OpenMP запускает дополнительный поток управления, который может запутать планировщик. Есть больше вариантов для закрепления, например KMP_AFFINITY=compact - но по какой-то причине это полностью перепутано в моей системе. Вы можете добавить ,verbose к переменной, чтобы увидеть, как среда выполнения привязывает ваши потоки.

likwid-pin - полезная альтернатива, обеспечивающая более удобное управление.

В целом одиночная точность должна быть как минимум такой же быстрой, как двойная точность. Двойная точность может быть медленнее, потому что:

  • Для двойной точности вам потребуется больше пропускной способности памяти/кэша.
  • Вы можете создавать ALU с более высокой пропускной способностью для одиночной точности, но обычно это не относится к процессорам, а скорее к графическим процессорам.

Я бы подумал, что как только вы избавитесь от аномалии производительности, это будет отражено в ваших номерах.

При масштабировании количества потоков для MKL/* gemm рассмотрите

  • Полоса пропускания памяти/общего кэша может стать узким местом, ограничивая масштабируемость
  • Режим Turbo эффективно уменьшает частоту ядра при увеличении использования. Это применимо даже тогда, когда вы работаете на номинальной частоте: на процессорах Haswell-EP инструкции AVX будут устанавливать более низкую "базовую частоту AVX" - но процессору разрешено превышать это, когда используется меньшее количество ядер/имеется тепловой запас и вообще больше на короткое время. Если вы хотите получить абсолютно нейтральные результаты, вам придется использовать базовую частоту AVX, которая для вас составляет 1,9 ГГц. Здесь описано и объяснено в одно изображение.

Я не думаю, что есть очень простой способ оценить, как на ваше приложение влияет плохое планирование. Вы можете показать это с помощью perf trace -e sched:sched_switch и некоторое программное обеспечение, чтобы визуализировать это, но это будет связано с высокой кривой обучения. И снова - для параллельного анализа производительности вы должны иметь нитки, закрепленные в любом случае.