Comments 72
Можно поинтересоваться судьбой софтового рендера?
Сначала писал статью. Потом понял, что одна не получится. Потом понял, что сразу всё писать сложно. Эта статья как раз вышла из неё. Решил посмотреть как пойдёт и будет ли интерес.
Если кратко, то: Quake2 1600x1200, MT (multythreaded), SSE3, i7-3770, 30 fps.
Переделываю с нуля под AVX. Цикл статьей планируется по SSE варианту.
Вы настолько верите IACA? Тогда вас может ждать много удивительных открытий при погружении в реальный мир.
2. Про реальный мир вы зря: я давно уже использую SIMD и тестирую и так, и так. Тут решил обойтись.
3. Не могу протестировать всё. У меня максимум AVX (i7-3770).
4. Это будет синтетика. Я уже отмечал, что тысячами матрицы я не умножаю.
5. Однако вы правы: пункт 4 имеет место быть, но IACA не может учитывать реальную работу с данными, памятью. Поэтому стоит написать. Интересно будет не только вам.
В матрице 16 элементов, в сумме это 64 умножения и 48 сложений.
Small matrix multiplication kernels (SMMs) are generated for Intel SSE, Intel AVX, Intel AVX2, IMCI (KNCni) for Intel Xeon Phi coprocessors (KNC), and Intel AVX‑512 as found in the Intel Xeon Phi processor family (KNL, KNM) and Intel Xeon processors (SKX). Highly optimized code for small convolutions is targeting Intel AVX2 and Intel AVX‑512, whereas other targets can automatically leverage specialized SMMs to perform convolutions.
Интересно, удалось ли обогнать оптимизатор GCC с флагами "-m64 -O3 -march=native -funroll-loops -ffast-math"? Перемножение двух целочисленных uint32_t матриц 4096x4096 у меня заняло ~6.8 сек на одном ядре Intel Xeon E3-1220 v5, задача здесь. Алгоритм — слегка оптимизированный O(n^3), отданный на растерзание компилятору.
Статью стоит рассматривать как мини введение в векторные инструкции, а не как конкретную реализацию умножения матриц либо еще какой конкретной задачи.
Потому что мне кажется, что тот же ICC очень неплохо занимается векторизацией. Было бы интересно посмотреть на результаты бенчмарка.
А пока разве только разные сайты использовать, где можно онлайн скомпилить код любым компилятором, и посмотреть результат. Но кому интересно, уже наверняка это сделали. Код то я выложил.
ICC компилятор имеет уже давно бесплатную(возобновляемую каждые 90 дней) даже коммерческую лицензию. Просто заходите и качаете. К тому же Intel System Studio работает и на Windows, так что можете провести измерения msvc vs icc.
Как уже отметили, иногда сложно предсказать, какие способы будут работать быстрее, глядя просто на то, что там нагенерировал компилятор — а вот по бенчмаркам уже проще смотреть.
AVX+FMA вариант (финальный)… В данном случае нам требуется вычислить a*b + c*d + e*f + g*h.
Выполняется одной командой _mm_dp_ps из SSE4 — она гораздо быстрее
У вас покрасивее, чем у Intel-а код выглядит https://software.intel.com/en-us/articles/benefits-of-intel-avx-for-small-matrices.
А почему вы не используете __restrict__
и __builtin_assume_aligned
/__assume_aligned
? Может быть тогда и инлайниться лучше будет.
__forceinline __m128 __vectorcall operator + (__m128 const) const noexcept;
Однако здесь всё итак работает хорошо. Всегда проверяю сгенерированный ассемблерный код. Насчёт *assume_aligned я посмотрю что это. Скорее всего то, чего мне сильно и давно не хватает. Но если оно linux/intel_compiler/… specific, то мне, к сожалению, не подойдёт.
Когда SSE4.1 только рекламировали, я уже надеялся, что такая классная инструкция dpps всё упростит. И ждал. Но увы…
1. Она дорогая: 2 такта (и судя по всему не параллелится как mulps), то есть даже без учёта транспонирования, их надо 16 по 2 такта, получается 32!!! (IACA согласна)
1.1 dpps — 4 умножения, 3 сложения
1.2 mulps — 4 умножения, параллелится, это 8 умножений за такт, 16 за два такта, или, грубо, 8 умножений + 3 сложения за 2.5 такта
2. Одну матрицу надо предварительно трансформировать: dpps(строка M, столбец N)
3.1 дополнительная инструкция movss для извлечения результата и работа с памятью на уровне 1-го float
ИЛИ
3.2 дополнительный код для компоновки результата из 4 float в одном xmm
Не стоит думать об инструкции, как вещи в себе. Инструкции нужны данные, они должны быть в нужном представлении как на входе, так и на выходе.
По совокупности факторов её даже Intel не использует в своих примерах матричного умножения.
… Очень разумный вопрос. Статья по факту представляет чистый положительный опыт. Если собрать весь отрицательный, то она будет минимум в 3 раза больше…А может будет это вторая статья?
По моему опыту, dpps имеет приличный выигрыш, но случаи бывают всякие. И подготовка данных, полностью согласен, стоит времени ЦПУ
Процессоры Intel снижают частоту ниже базовой при большом числе AVX инструкций.
https://en.wikichip.org/wiki/intel/frequency_behavior
Например, разницы при сложении двух очень длинных векторов с использованием SSE и AVX при максимальном распараллеливании не будет, т.к. всё упрётся в пропускную способность памяти.
Зачем быть заложником одной компании?
www.agner.org/optimize/blog/read.php?i=838
AMD has four 128-bit units for floating point and vector operations. Two of these can do addition and two can do multiplication. Intel has two 256-bit units, both of which can do addition as well as multiplication. This means that floating point code with scalars or vectors of up to 128 bits will execute on the AMD processor at a maximum rate of four instructions per clock (two additions and two multiplications), while the Intel processor can do only two. With 256-bit vectors, AMD and Intel can both do two instructions per clock.
Т.е. тоже самое.
Но FMA быстрее на Интел.
Intel beats AMD on 256-bit fused multiply-and-add instructions, where AMD can do one while Intel can do two per clock. Intel is also better than AMD on 256-bit memory writes, where Intel has one 256-bit write port while the AMD processor has one 128-bit write port.
Очень интересная статья! Кстати, кроме как тем, что нужен свой формат хранения, матрички из GLM не такие же шустрые? Кажется, там была какая-то поддержка SIMD.
Пожалуй, лучше чем здесь я не отвечу.
По поводу разрядности. Я так понимаю что float в данном случае 32 бита?
А если точность нужна выше? 64 бита на Double. Этот же код сможет работать с double — данные влезут в регистры? Или надо будет модифицировать?
И возможен ли вариант на 80 бит (extended)? Или SIMD не поддерживает 80-битный float сопроцессора x86?
1 регистр AVX — это 8 float или 4 double
Код просто будет аналогичным.
80-битные числа SSE/AVX не поддерживаются
А вот насчет 80 бит жаль :-( Странно что Intel не дала возможности использовать его в своих же расширениях. Пусть бы и не так быстро как float32.
Странно что Intel не дала возможности использовать его в своих же расширениях.Нет в этом ничего странного. Операции с 80битами сильно дороже, чем с 64мя. А используются редко. Тратить драгоценные транзисторы на то, что не используется никто не будет.
Собственно изначально MMX работал с целыми числами, SSE — с 32-битными float'ами, SSE2 — с 64битными double'ами. И каждый раз собиралась статистика: кто будет использовать эти команды и зачем. Очевидно, что когда «выпекали» SSE3/4/etc толп людей, желающих воспользоваться 90битными флоатами у маркетиногово отдела Intel — не стояло…
Код просто будет аналогичным.
Если бы! AVX — это не в два раза больший SSE модуль, это два SSE модуля, работающих в паре. Если для сложения и умножения это ничего не меняет, то например _mm256_shuffle_ps ведёт себя совсем иначе.
это не в два раза больший SSE модуль
В 2 раза больший.
это два SSE модуля, работающих в паре.
Амд-проблемы.
Вы вообще не поняли о чем я, видимо ни разу не использовали AVX. Я не про какие-то реализации каких-то процессоров, я про весь набор команд и их архитектуру. Посмотрите как работает например _mm256_shuffle_ps — в dst[127:0]
могут попасть только биты из a[127:0]
и b[127:0]
, но не из a[255:128]
и b[255:128]
. То есть нижние 128 бит работают только с нижними 128 битами, старшие со старшими — это больше похоже на два модуля SSE, чем на один модуль SSE, но в два раза больше. И это касается большинства команд. То же самое можно сказать про AVX512 — это скорее 4 модуля SSE, работающих параллельно, чем один очень длинный модуль SSE.
Посмотрите как работает например _mm256_shuffle_ps — в dst[127:0] могут попасть только биты из a[127:0] и b[127:0], но не из a[255:128] и b[255:128]. То есть нижние 128 бит работают только с нижними 128 битами, старшие со старшими — это больше похоже на два модуля SSE, чем на один модуль SSE, но в два раза больше. И это касается большинства команд. То же самое можно сказать про AVX512 — это скорее 4 модуля SSE, работающих параллельно, чем один очень длинный модуль SSE.
Похоже вы ничего не знаете об микроархитектуре, странно как-то вы используете avx и ничего о ней не знаете.
Во-первых, что значит «например» — если это единственный пример(вернее там много похожих инструкций, но суть одна). Зачем вы делаете вид, будто бы у вас много примеров?
Во-вторых ваша аналитика ничего не имеет общего с реальность. И если бы вы действительно использовали тот самых avx, то вы бы знали в чём именно заключается проблема, а заключается проблема в 8битном селекторе.
Опять же, ваши рассуждения имели бы смысл, если данные операции работали в рамках 128 бит. Т.е. если бы из А можно было перемещать в рамках первой части dst, а из Б в рамках второй. Это было бы так, если бы они попросту склеили два 128битных шафла, но это не так. Они почему-то шафлят aa bb aa bb. Неужели там mmx?
И да, никаких «модулей» не существует. Это обывательская терминология, ничего общего с реальностью не имеющая.
Нет. Видимо, вы путаете AVX и AVX2. В AVX все команды работают с регистрами как с двумя 128-битными независимыми половинками (за исключением пары команд, которые обращаются с 128-битными частями как единым целым).
И если бы вы действительно использовали тот самых avx, то вы бы знали в чём именно заключается проблема, а заключается проблема в 8битном селекторе.
Когда 8-битного селектора не хватает, то в качестве аргумента команды может использоваться другой регистр, например, в команде AVX2 _mm256_permutevar8x32_ps.
Они почему-то шафлят aa bb aa bb
Потому что lower(result) = op(lower(aa), lower(bb)), higher(result) = op(higher(aa), higher(bb)).
SELECT4(src, control)
{
CASE(control[1:0])
0: tmp[31:0] := src[31:0]
1: tmp[31:0] := src[63:32]
2: tmp[31:0] := src[95:64]
3: tmp[31:0] := src[127:96]
ESAC RETURN tmp[31:0]
}
dst[31:0] := SELECT4(a[127:0], imm8[1:0])
dst[63:32] := SELECT4(a[127:0], imm8[3:2])
dst[95:64] := SELECT4(b[127:0], imm8[5:4])
dst[127:96] := SELECT4(b[127:0], imm8[7:6])
dst[159:128] := SELECT4(a[255:128], imm8[1:0])
dst[191:160] := SELECT4(a[255:128], imm8[3:2])
dst[223:192] := SELECT4(b[255:128], imm8[5:4])
dst[255:224] := SELECT4(b[255:128], imm8[7:6])
dst[MAX:256] := 0
Во-первых, что значит «например» — если это единственный пример(вернее там много похожих инструкций, но суть одна). Зачем вы делаете вид, будто бы у вас много примеров?
Ооох, да пожалуйста:
__m256 _mm256_unpacklo_ps (__m256 a, __m256 b)
INTERLEAVE_DWORDS(src1[127:0], src2[127:0]){
dst[31:0] := src1[31:0]
dst[63:32] := src2[31:0]
dst[95:64] := src1[63:32]
dst[127:96] := src2[63:32]
RETURN dst[127:0]
}
dst[127:0] := INTERLEAVE_DWORDS(a[127:0], b[127:0])
dst[255:128] := INTERLEAVE_DWORDS(a[255:128], b[255:128])
dst[MAX:256] := 0
__m256i _mm256_packs_epi32 (__m256i a, __m256i b)
dst[15:0] := Saturate_Int32_To_Int16 (a[31:0])
dst[31:16] := Saturate_Int32_To_Int16 (a[63:32])
dst[47:32] := Saturate_Int32_To_Int16 (a[95:64])
dst[63:48] := Saturate_Int32_To_Int16 (a[127:96])
dst[79:64] := Saturate_Int32_To_Int16 (b[31:0])
dst[95:80] := Saturate_Int32_To_Int16 (b[63:32])
dst[111:96] := Saturate_Int32_To_Int16 (b[95:64])
dst[127:112] := Saturate_Int32_To_Int16 (b[127:96])
dst[143:128] := Saturate_Int32_To_Int16 (a[159:128])
dst[159:144] := Saturate_Int32_To_Int16 (a[191:160])
dst[175:160] := Saturate_Int32_To_Int16 (a[223:192])
dst[191:176] := Saturate_Int32_To_Int16 (a[255:224])
dst[207:192] := Saturate_Int32_To_Int16 (b[159:128])
dst[223:208] := Saturate_Int32_To_Int16 (b[191:160])
dst[239:224] := Saturate_Int32_To_Int16 (b[223:192])
dst[255:240] := Saturate_Int32_To_Int16 (b[255:224])
dst[MAX:256] := 0
__m256 _mm256_dp_ps (__m256 a, __m256 b, const int imm8)
DP(a[127:0], b[127:0], imm8[7:0]) {
FOR j := 0 to 3
i := j*32
IF imm8[(4+j)%8]
temp[i+31:i] := a[i+31:i] * b[i+31:i]
ELSE
temp[i+31:i] := 0
FI
ENDFOR
sum[31:0] := (temp[127:96] + temp[95:64]) + (temp[63:32] + temp[31:0])
FOR j := 0 to 3
i := j*32
IF imm8[j%8]
tmpdst[i+31:i] := sum[31:0]
ELSE
tmpdst[i+31:i] := 0
FI
ENDFOR
RETURN tmpdst[127:0]
}
dst[127:0] := DP(a[127:0], b[127:0], imm8[7:0])
dst[255:128] := DP(a[255:128], b[255:128], imm8[7:0])
dst[MAX:256] := 0
И в все остальные команды арифметических и логических действий соблюдают это правило: результат работы с нижними 128 битами записывается в нижние 128 битов. А теперь ваша очередь показать все остальные инструкции спокойно перекидываются битами из нижних и верхних частей? Такие команды конечно есть, но работают они со всеми 128 битами одновременно, то есть являются инструкциями копирования, а никак не SIMD-операциями.
SELECT4(src, control){ CASE(control[1:0]) 0: tmp[63:0] := src[63:0] 1: tmp[63:0] := src[127:64] 2: tmp[63:0] := src[191:128] 3: tmp[63:0] := src[255:192] ESAC RETURN tmp[63:0] } dst[63:0] := SELECT4(a[255:0], imm8[1:0]) dst[127:64] := SELECT4(a[255:0], imm8[3:2]) dst[191:128] := SELECT4(a[255:0], imm8[5:4]) dst[255:192] := SELECT4(a[255:0], imm8[7:6]) dst[MAX:256] := 0
80 только FPU. SSE: 32*4=64*2=128bit, AVX: 64*4=32*8=256bit.
По сути точность несколько меньше, плата за производительность.
Здесь вскользь упомянут вариант, как делать для разных платформ. Пишем разные реализации, и потом после анализа процессора через cpuid подставляем нужный указателю на функцию. И потом работаем с данным функционалом через указатель. Данный подход полезен для больших, ёмких операций. Если использовать его для одной инструкции — точно будет только хуже.
Как пример: функция для умножения матрицы на вектор — станет медленнее, функция для умножения матрицы на массив векторов — имеет смысл.
Гибкий вариант есть здесь optimizing_cpp.pdf. Искать по: CriticalFunction_Dispatch.
Разные руководства по оптимизации рекомендуют лишний раз не дёргать кэш для массовых операций сохранения.
Это явно не ваш случай. Эти инструкции полезны, когда вам нужно обработать огромный массив данных (по крайней мере больше кеша) и вы уверены, что к нему не будет обращений до конца вашей работы. В вашем же случае вы оперируете 64 байтами, к которым с большой вероятностью в ближайшее время будет дальнейшая работа.
По опыту использования в обработке изображений стриминг давал незначительный прирост на больших изображениях и проигрыш в несколько раз на изображениях поменьше.
Не понял к чему вы это. Финальный код содержит стриминговые инструкции, которые в данном случае сильно вредны. Ваше нежелание или невозможность тестировать на реальных кейсах подводит вас.
У вас данные выровнены по линиям кеша? 4x4x4 — это как раз одна линия. Возможно, если данные пересекают линии, результат будет отличаться.
Вдруг будут идеи.
Ускоряем умножение матриц float 4x4 с помощью SIMD