Часто от приложения требуется максимально возможная производительность, и нередко помочь ее достигнуть может переход на новую платформу. А дополнительный прирост производительности может быть получен не только за счет архитектурных изменений новой платформы, т.е. увеличения частоты процессора, пропускной способности памяти и прочих изменений, но и за счет использования новых наборов команд процессора. Этот подход и будет рассмотрен в данной статье на примере векторных инструкций, появлявшихся последовательно в процессорах Intel от SSE до нового набора инструкций AVX 3.1, опубликованного недавно в соответствующей спецификации на сайте компании Intel.
Автор данной статьи принимает участие в разработке библиотеки Intel IPP, поэтому в данной статье рассматривается прирост производительности на примере одной из функции из данной библиотеки. Хочется надеяться, что данная статья будет полезна следующим категориям читателей:
Помимо этих категорий читателей, разработчикам, заинтересованным в оптимизации своего приложения и желающим минимизировать свои затраты на низкоуровневую оптимизацию, данная статья может помочь познакомиться поближе с IPP библиотекой.
API данной функции выглядит следующим образом.
Функция читает исходное float изображение по указателю pSrc с шагом в srcStep между строками и размером width x height, переставляет каналы в каждом пикселе в соответствии с новым порядком dstOrder и записывает результирующее изображение по указателю pDst. Фукнция также может записать указанное пользователем фиксированное значение val в канал. Допускается дублирование каналов в выходном изображении. Схема работы функции представлена на рисунке 1.

Надеемся, что читатель уже знаком с понятием векторных регистров, если же нет, то можно ознакомиться и воспользоваться удобным навигатором интринсиков.
Итак, посмотрим, как можно было бы векторизовать, т.е. переписать код этой функции с использованием векторных регистров SSE. В соответствии с логикой функции, нам необходимо осуществить перестановку float элементов внутри векторного регистра, поэтому нам нужны соответствующие инструкции перестановки элементов. К таким инструкциям можно отнести следующие инструкции из набора SSE:
Данные инструкции очень полезны и для нашей задачи они тоже могут быть использованы, но возникает следующая проблема. Инструкция SHUFPS имеет непосредственный операнд imm8, формирующий перестановку элементов. Третий параметр данной инструкции это imm8 — значение, т.е. некая константа. При попытке же использовать в интринсике переменную компилятор выдаст ошибку. А инструкции UNPCKHPS/UNPCKLPS чередуют элементы в регистрах лишь единственным способом. Допустимых вариантов перестановки каналов, указываемых в dstOrder довольно много и поэтому для реализации всех возможных комбинаций пришлось бы разработать код со множеством веток и различным значением параметра mask на каждую возможную комбинацию каналов. Это довольно утомительная задача, к тому же на сегодняшний день изначальный SSE набор инструкций существенно расширен и код можно сделать значительно проще. Ближайшее к SSE поколение, которое содержит очень полезную для наших целей инструкцию, называется SSSE3. И речь конечно же об инструкции:
Эта инструкция может быть эффективно использована во множестве алгоритмов, где для вычисления одного выходного элемента, используются несколько соседних входных элементов. Данная инструкция переставляет байты приемника в абсолютно произвольном порядке, а также может обнулить байт при установленном старшем бите в регистре маске, см. Рис 2.

Несмотря на то, что инструкция переставляет байты, мы можем ее применить и для float данных.
Как видно, наш SSSE3 код использует довольно много SIMD инструкций, и в строке комментария поясняется, зачем они используются. Инструкции _mm_unpacklo_pd используются для того чтобы объединить в один XMM регистр, значения, прочитанные в разные регистры. Инструкции _mm_and_ps, _mm_or_ps используются для объединения по маске с немодифицированным каналом и val значением. Непосредственно для перестановки элементов мы и используем инструкцию _mm_shuffle_epi8. Измеряем, сколько времени уходит на выполнение и получаем P=0.4. Несмотря на то, что код выглядит более громоздко по сравнению с C кодом, прирост производительности составил примерно 2.5X. И это довольно неплохой результат.

Поэтому AVX вариант нашего кода может выглядеть примерно так:
Закономерно может возникнуть вопрос: почему же в нашем 256-битном AVX коде используется загрузка с помощью 128 битных инструкций? Дело в том, что в нашем коде была использована инструкция VPERMPS в строке
ymm0 = _mm256_permutevar_ps(ymm0, ymm8)
Как уже объяснялось выше, данная инструкция переставляет элементы внутри младшей и старшей части регистра (или по другому lane) независимо. Поэтому, несмотря на то, что 2 пикселя лежат в памяти рядом, мы вынуждены читать их таким образом, чтобы они попали в разные lane. Обратите внимание, как формируется регистр-маска ymm8, он содержит абсолютно одинаковые значения в старшей и младшей части для перестановки. Код написан таким образом, чтобы количество итераций в циклах сохранилось и не влияло на итоговую производительность. Измеряем производительность: P=0.31. Таким образом прирост AVX кода по отношению к C составил 3.2X, а у SSSE3 было 2.48X
Обратите внимание что маска перестановки в регистре ymm8 совершенно изменилась и является сквозной по всему регистру, т.е. индексы в ней относятся в целом ко всему 256 битному регистру.
Однако в этом коде есть проблема в том, как происходит чтение регистра ymm0. В нем загружаются лишние 2 float элемента. Так делать неправильно, поскольку может произойти чтение за границей изображения и исправить это нам поможет опять же появившаяся лишь в AVX2 инструкция чтения по маске
__m256i _mm256_maskload_epi32(int const *, __m256i);
Данная инструкция читает необходимые элементы из памяти и не читает ничего лишнего, что предотвращает возможную memory segmentation ошибку. Теперь наш правильный код будет выглядеть так.
Производительность этого кода составляет P=0.286. Ускорение таким образом изменилось от k=3.2X до k=3.5X. Помимо чтения по маске, мы можем использовать и сохранять немодифицированный канал без изменения, используя инструкцию
void _mm256_maskstore_ps(float *, __m256i, __m256);
Обратите внимание на то, что 128 битная предшественница данной инструкции
void _mm_maskmoveu_si128(__m128i, __m128i, char *)
сохраняет данные, используя non-temporal запись, а вот в AVX версиях _mm256_maskstore_ps и _mm_maskstore_ps от этого механизма было решено отказаться. О том, что такое non-temporal запись и как она применяется в IPP библиотеке будет рассказано в одной из следующих статей. Меряем производительность и получаем небольшой прирост, P=0.27, k=3.7X.
Как следует из названия длина векторных регистров увеличилась в два раза с 256 бит до 512. Lane концепция также сохранилась. Второе новшество заключается в том, что операции с масками приобрели массовый характер и затрагивают теперь большинство векторных операций. Для этого были введены специальные масочные регистры. Масочные регистры имеют свой тип данных __mmask8, __mmask16, __mmask32, __mmask64, но с ними можно оперировать как с целочисленными типами данных, поскольку они так и объявлены в заголовочном файле zmmintrin.h. Использование масок делает AVX512 обрабатывающий цикл функции довольно минималистичным.
Количество интринсиков в одной итерации цикла сократилось до трех, и каждый из них использует маску. Первый интринсик
z0 = _mm512_mask_loadu_ps(_mm512_setzero_ps(), 0xFFF, (float *)(pSrc+s))
использует маску 0xFFF для того чтобы прочитать из памяти 4 RGB пикселя или 12 float, в двоичном виде 0xFFF=1111111111111111b. Следующий интринсик
z0 = _mm512_mask_permutexvar_ps(zv, mskv, zi0, z0)
выполняет сразу 2 функции: переставляет каналы входного пикселя и объединяет их по маске со значением val. Маска mskv формируется перед началом цикла и биты маски используются следующим образом: значение бита равное ‘1’ означает «поместить результат операций в результирующий регистр», а ‘0’ – «cохранить значение того регистра, который указан до маски», в нашем случае zv. Третий интринсик
_mm512_mask_storeu_ps(pDst+d, mska, z0)
сохраняет полученные 4 пикселя в память, при этом он не пишет ничего в немодифицированный канал. Производительность кода составляет P=0.183, k=5.44X.
Читатель может сравнить код SSSE3 с AVX512 кодом, в котором используется всего 3 интринсика. Именно краткость avx512 кода и вдохновила автора на написание данной статьи.
Итоговая таблица роста производительности.

Ознакомившись с данной статьей, читатель может попробовать использовать SIMD код в своих приложениях, а читатель, уже имеющий такой опыт, может оценить по достоинству преимущества crosslane инструкций перестановки и инструкций с масками в новых процессорах Intel, поддерживающих набор инструкций AVX512.
Также хотелось бы отметить, что в данной статье приведена незначительная часть приемов оптимизации, используемых в библиотеке IPP и позволяющих достичь ее уникальной производительности.
Читатель может и сам попробовать использовать библиотеку IPP в своем приложении, освободив возможности для решения своих прикладных задач. Дополнительным эффектом использования IPP станет то, что после выхода на рынок процессоров с набором инструкции AVX512, библиотека IPP уже будет иметь оптимизированную на AVX512 версию наиболее востребованных функций.
Автор данной статьи принимает участие в разработке библиотеки Intel IPP, поэтому в данной статье рассматривается прирост производительности на примере одной из функции из данной библиотеки. Хочется надеяться, что данная статья будет полезна следующим категориям читателей:
- начинающим разработчикам, желающим познакомиться на конкретном примере, как с помощью SIMD инструкций можно повысить производительность кода;
- тем, кто уже имеет достаточный опыт в использовании SIMD инструкций и интересуется новым набором инструкций AVX3.1 для оптимизации своих приложений.
Помимо этих категорий читателей, разработчикам, заинтересованным в оптимизации своего приложения и желающим минимизировать свои затраты на низкоуровневую оптимизацию, данная статья может помочь познакомиться поближе с IPP библиотекой.
Функция ippiSwapChannels_32f_C3C4R
Библиотека IPP содержит оптимизированные 32-х и 64-х битные версии для процессоров, поддерживающих наборы инструкций SSSE3, SSE41, AVX, AVX2, AVX31. Поэтому и рассмотрим возможности оптимизации с использованием этих наборов данных. Вообще, IPP содержит много функций для различных областей обработки данных. Эти функции реализуют как довольно сложные для понимания математические алгоритмы, к примеру преобразование Фурье, так и алгоритмы, не представляющие особых сложностей для понимания. Именно такую довольно простую функцию ippiSwapChannels_32f_C3C4R и рассмотрим в данной статье, постараясь последовательно применить к ней наборы инструкций от SSSE3 до AVX31.API данной функции выглядит следующим образом.
int ippiSwapChannels_32f_C3C4R(
const float* pSrc, int srcStep, float* pDst, int dstStep,
int width, int height, const int dstOrder[4], float val )
Функция читает исходное float изображение по указателю pSrc с шагом в srcStep между строками и размером width x height, переставляет каналы в каждом пикселе в соответствии с новым порядком dstOrder и записывает результирующее изображение по указателю pDst. Фукнция также может записать указанное пользователем фиксированное значение val в канал. Допускается дублирование каналов в выходном изображении. Схема работы функции представлена на рисунке 1.

C код
Код функции на языке C довольно прост (см. ниже). В каждой итерации цикла формируется один пиксель выходного изображения при этом с каждым каналом выходного пикселя может произойти одно из 3-х возможных действий:- скопироваться канал пикселя исходного изображения, если индекс в dstOrder меньше 3;
- поместиться фиксированное значение val, если индекс в dstOrder равен 3;
- сохраниться старое немодифицированное значение, если индекс в dstOrder больше 3
Код
int s, d, h, ord0 = dstOrder[0], ord1 = dstOrder[1], ord2 = dstOrder[2], ord3 = dstOrder[3];
for( h = 0; h < height; h++ ) {
for( s = 0, d = 0; d < width*4; s += 3, d += 4 ) {
if( ord0 < 3 ) {
pDst[d] = pSrc[s + ord0];
} else if( ord0 == 3 ) {
pDst[d] = val;
}
if( ord1 < 3 ) {
pDst[d + 1] = pSrc[s + ord1];
} else if( ord1 == 3 ) {
pDst[d + 1] = val;
}
if( ord2 < 3 ) {
pDst[d + 2] = pSrc[s + ord2];
} else if( ord2 == 3 ) {
pDst[d + 2] = val;
}
if( ord3 < 3 ) {
pDst[d + 3] = pSrc[s + ord3];
} else if( ord3 == 3 ) {
pDst[d + 3] = val;
}
}
pSrc = (float*)((char*)pSrc + srcStep);
pDst = (float*)((char*)pDst + dstStep);
}
Первая версия SIMD кода с использованием SSSE3
А теперь, как и планировалось, попробуем оптимизировать данный код используя векторные SIMD инструкции вплоть до AVX3.1. Будем использовать Intel Composer, а производительность всех вариантов кода измерять на внутреннем симуляторе серверной версии SkyLake, поддерживающей набор AVX31 — SKX, поскольку нужно оценить именно преимущества новых инструкций без влияния архитектурных изменений. Для компилятора также указываем максимально возможный уровень оптимизации, одинаково для всех вариантов кода. Измерянную таким образом производительность C кода примем за некую условную величину P=1.Надеемся, что читатель уже знаком с понятием векторных регистров, если же нет, то можно ознакомиться и воспользоваться удобным навигатором интринсиков.
Итак, посмотрим, как можно было бы векторизовать, т.е. переписать код этой функции с использованием векторных регистров SSE. В соответствии с логикой функции, нам необходимо осуществить перестановку float элементов внутри векторного регистра, поэтому нам нужны соответствующие инструкции перестановки элементов. К таким инструкциям можно отнести следующие инструкции из набора SSE:
SHUFPS xmm1, xmm2/m128, imm8
UNPCKHPS xmm1, xmm2/m128
UNPCKLPS xmm1, xmm2/m128
__m128 _mm_shuffle_ps(__m128 a, __m128 b, unsigned int mask).
__m128 _mm_unpackhi_ps(__m128, __m128)
__m128 _mm_unpacklo_ps(__m128, __m128)
Данные инструкции очень полезны и для нашей задачи они тоже могут быть использованы, но возникает следующая проблема. Инструкция SHUFPS имеет непосредственный операнд imm8, формирующий перестановку элементов. Третий параметр данной инструкции это imm8 — значение, т.е. некая константа. При попытке же использовать в интринсике переменную компилятор выдаст ошибку. А инструкции UNPCKHPS/UNPCKLPS чередуют элементы в регистрах лишь единственным способом. Допустимых вариантов перестановки каналов, указываемых в dstOrder довольно много и поэтому для реализации всех возможных комбинаций пришлось бы разработать код со множеством веток и различным значением параметра mask на каждую возможную комбинацию каналов. Это довольно утомительная задача, к тому же на сегодняшний день изначальный SSE набор инструкций существенно расширен и код можно сделать значительно проще. Ближайшее к SSE поколение, которое содержит очень полезную для наших целей инструкцию, называется SSSE3. И речь конечно же об инструкции:
PSHUFB xmm1, xmm2/m128
__m128i _mm_shuffle_epi8 (__m128i a, __m128i b)
Эта инструкция может быть эффективно использована во множестве алгоритмов, где для вычисления одного выходного элемента, используются несколько соседних входных элементов. Данная инструкция переставляет байты приемника в абсолютно произвольном порядке, а также может обнулить байт при установленном старшем бите в регистре маске, см. Рис 2.

Несмотря на то, что инструкция переставляет байты, мы можем ее применить и для float данных.
Код
int s, d, h, ord0 = dstOrder[0], ord1 = dstOrder[1], ord2 = dstOrder[2], ord3 = dstOrder[3];
__m128 xmm0, xmm1, xmm2, xmm3, xmm4, xmm5,xmm6, xmm7;
__m128i xmm8 = _mm_set_epi8(4*ord3+3,4*ord3+2,4*ord3+1,4*ord3,
4*ord2+3,4*ord2+2,4*ord2+1,4*ord2,
4*ord1+3,4*ord1+2,4*ord1+1,4*ord1,
4*ord0+3,4*ord0+2,4*ord0+1,4*ord0);
__m128i xmm9 = _mm_set_epi32(0,0,0,0);//unmodified channel mask
__m128i xmm10 = _mm_set_epi32(0,0,0,0);//val
if(ord3 >= 3){
xmm8 = _mm_or_si128(xmm8, _mm_set_epi32(-1, 0, 0, 0));
if(ord3 > 3){ xmm9 = _mm_or_si128(xmm9, _mm_set_epi32(-1, 0, 0, 0));}
else { xmm10 = _mm_or_si128(xmm10, _mm_set_epi32(val, 0, 0, 0));};
}
if(ord2 >= 3){
xmm8 = _mm_or_si128(xmm8, _mm_set_epi32( 0,-1, 0, 0));
if(ord2 > 3){ xmm9 = _mm_or_si128(xmm9, _mm_set_epi32( 0, -1, 0, 0));}
else { xmm10 = _mm_or_si128(xmm10, _mm_set_epi32( 0, val,0, 0));}
}
if(ord1 >=3){
xmm8 = _mm_or_si128(xmm8, _mm_set_epi32( 0, 0,-1, 0));
if(ord1 > 3){xmm9 = _mm_or_si128(xmm9, _mm_set_epi32( 0, 0,-1, 0));}
else {xmm10 = _mm_or_si128(xmm10, _mm_set_epi32( 0, 0,val, 0));}
}
if(ord0 >= 3){
xmm8 = _mm_or_si128(xmm8, _mm_set_epi32( 0, 0, 0,-1));
if(ord0 > 3){xmm9 = _mm_or_si128(xmm9, _mm_set_epi32( 0, 0, 0,-1));}
else {xmm10 = _mm_or_si128(xmm10, _mm_set_epi32( 0, 0, 0,val));}
}
for( h = 0; h < height; h++ ) {
s = 0;
d = 0;
for(; d < width*4; s += 3*4, d += 4*4 ) {
xmm0 = _mm_load_sd(pSrc+s ); //читаем g0 r0
xmm1 = _mm_load_ss(pSrc+s+2 ); //читаем b0
xmm2 = _mm_load_sd(pSrc+s +3);
xmm3 = _mm_load_ss(pSrc+s+2+3);
xmm4 = _mm_load_sd(pSrc+s +6);
xmm5 = _mm_load_ss(pSrc+s+2+6);
xmm6 = _mm_load_sd(pSrc+s +9);
xmm7 = _mm_load_ss(pSrc+s+2+9);
xmm0 = _mm_unpacklo_pd(xmm0, xmm1); //получаем b0 g0 r0 в один регистр
xmm1 = _mm_unpacklo_pd(xmm2, xmm3);
xmm2 = _mm_unpacklo_pd(xmm4, xmm5);
xmm3 = _mm_unpacklo_pd(xmm6, xmm7);
xmm4 = _mm_loadu_ps(pDst+d+0); //читаем немодифицируемый канал
xmm5 = _mm_loadu_ps(pDst+d+4);
xmm6 = _mm_loadu_ps(pDst+d+8);
xmm7 = _mm_loadu_ps(pDst+d+12);
xmm0 = _mm_shuffle_epi8(xmm0, xmm8); //осуществляем перестановку в соответствии с маской
xmm1 = _mm_shuffle_epi8(xmm1, xmm8);
xmm2 = _mm_shuffle_epi8(xmm2, xmm8);
xmm3 = _mm_shuffle_epi8(xmm3, xmm8);
xmm4 = _mm_and_ps(xmm4,xmm9); //обнуляем ненужные значения прочитанные из dst
xmm5 = _mm_and_ps(xmm5,xmm9);
xmm6 = _mm_and_ps(xmm6,xmm9);
xmm7 = _mm_and_ps(xmm7,xmm9);
xmm0 = _mm_or_ps(xmm0,xmm4); //объединяем переставленные каналы с немодифицированным каналом
xmm1 = _mm_or_ps(xmm1,xmm5);
xmm2 = _mm_or_ps(xmm2,xmm6);
xmm3 = _mm_or_ps(xmm3,xmm7);
xmm0 = _mm_or_ps(xmm0,xmm10); //объединяем переставленные каналы с val значением
xmm1 = _mm_or_ps(xmm1,xmm10);
xmm2 = _mm_or_ps(xmm2,xmm10);
xmm3 = _mm_or_ps(xmm3,xmm10);
_mm_storeu_ps(pDst+d , xmm0); //записываем результат
_mm_storeu_ps(pDst+d+ 4, xmm1);
_mm_storeu_ps(pDst+d+ 8, xmm2);
_mm_storeu_ps(pDst+d+12, xmm3);
}
for(; d < width*4; s += 3, d += 4 ) {
__m128 xmm0, xmm1, xmm4;
xmm0 = _mm_load_sd(pSrc+s );
xmm1 = _mm_load_ss(pSrc+s+2);
xmm0 = _mm_unpacklo_pd(xmm0, xmm1);
xmm4 = _mm_loadu_ps(pDst+d);
xmm0 = _mm_shuffle_epi8(xmm0, xmm8);
xmm4 = _mm_and_ps(xmm4,xmm9); //apply unmodified mask
xmm0 = _mm_or_ps(xmm0,xmm4);
xmm0 = _mm_or_ps(xmm0,xmm10);
_mm_storeu_ps(pDst+d, xmm0);
}
pSrc = (float*)((char*)pSrc + srcStep);
pDst = (float*)((char*)pDst + dstStep);
}
Как видно, наш SSSE3 код использует довольно много SIMD инструкций, и в строке комментария поясняется, зачем они используются. Инструкции _mm_unpacklo_pd используются для того чтобы объединить в один XMM регистр, значения, прочитанные в разные регистры. Инструкции _mm_and_ps, _mm_or_ps используются для объединения по маске с немодифицированным каналом и val значением. Непосредственно для перестановки элементов мы и используем инструкцию _mm_shuffle_epi8. Измеряем, сколько времени уходит на выполнение и получаем P=0.4. Несмотря на то, что код выглядит более громоздко по сравнению с C кодом, прирост производительности составил примерно 2.5X. И это довольно неплохой результат.
Вторая версия SIMD кода с использованием AVX
На сегодняшний день процессоры Intel из топовой линейки уже давно поддерживают как минимум набор инструкций AVX. Векторные регистры, используемые в инструкциях AVX, имеют длину 256 бит и соответственно в 2 раза шире 128 битных SSE регистров. Однако AVX содержит два существенных для нас ограничения. Во-первых, AVX содержит 256-битные инструкции перестановки только для float и double типов данных и мы не можем использовать 256-битный аналог _mm_shuffle_epi8. Во-вторых, после выхода набора инструкций AVX была введена так называмая lane концепция для векторных инструкций. В рамках данной концепции один 256-битный AVX регистр представляет собой связку двух 128-битных регистров, именуемых словом lane. Большинство 256-битных расширений 128-битных инструкции с горизонтальным доступом работают только внутри каждой lane. Таким образом, зачастую AVX вариант кода, который уже имеет SSE версию выглядит так как это показано на рис. 3
Поэтому AVX вариант нашего кода может выглядеть примерно так:
Код
int s, d, h, ord0 = dstOrder[0], ord1 = dstOrder[1], ord2 = dstOrder[2], ord3 = dstOrder[3];
__m256 ymm0, ymm1, ymm2, ymm3, ymm4, ymm5,ymm6, ymm7;
__m128 xmm0, xmm1, xmm2,xmm3;
__m256 ymm10= _mm256_setzero_ps();//val
__m256i ymm8 = _mm256_set_epi32(
dstOrder[3],dstOrder[2],dstOrder[1],dstOrder[0],
dstOrder[3],dstOrder[2],dstOrder[1],dstOrder[0]);
__m256 ymm9 = _mm256_setzero_ps();
__m256 ymm11 = _mm256_setzero_ps();
__m256i YMM_N[4], YMM_V[4];
YMM_N[0]= _mm256_set_m128(_mm_set_epi32( 0, 0, 0,-1), _mm_set_epi32( 0, 0, 0,-1));
YMM_N[1]= _mm256_set_m128(_mm_set_epi32( 0, 0,-1, 0), _mm_set_epi32( 0, 0,-1, 0));
YMM_N[2]= _mm256_set_m128(_mm_set_epi32( 0,-1, 0, 0), _mm_set_epi32( 0,-1, 0, 0));
YMM_N[3]= _mm256_set_m128(_mm_set_epi32(-1, 0, 0, 0), _mm_set_epi32(-1, 0, 0, 0));
YMM_V[0]= _mm256_set_m128(_mm_set_epi32( 0 , 0, 0,val), _mm_set_epi32( 0 , 0, 0,val));
YMM_V[1]= _mm256_set_m128(_mm_set_epi32( 0 , 0,val, 0), _mm_set_epi32( 0 , 0,val, 0));
YMM_V[2]= _mm256_set_m128(_mm_set_epi32( 0 ,val, 0, 0), _mm_set_epi32( 0 ,val, 0, 0));
YMM_V[3]= _mm256_set_m128(_mm_set_epi32(val, 0, 0, 0), _mm_set_epi32(val, 0, 0, 0));
//генерируем маски
for(s=0;s<4;s++){
if(dstOrder[s] >= 3){
ymm9 = _mm256_or_ps(ymm9, YMM_N[s]);
if(dstOrder[s] > 3){ ymm11 = _mm256_or_si256(ymm11, YMM_N[s]);} //маска немодифицированного канала
else { ymm10 = _mm256_or_si256(ymm10, YMM_V[s]);}; //val значение
}
}
for( h = 0; h < height; h++ ) {
s = 0;
d = 0;
for(; d < (width&(~1))*4; s += 3*2, d += 4*2 ) {//обрабатываем по 2 пикселя за итерацию
xmm0 = _mm_loadu_ps (pSrc+s ); //читаем src
xmm2 = _mm_load_sd (pSrc+s +3);
xmm3 = _mm_load_ss (pSrc+s+2 +3);
xmm2 = _mm_unpacklo_pd (xmm2, xmm3 );
ymm0 = _mm256_set_m128 (xmm2, xmm0 ); //собираем их в один 256битный регистр
ymm4 = _mm256_loadu_ps (pDst+d ); //читаем выходное изображение
ymm0 = _mm256_permutevar_ps(ymm0, ymm8 ); //переставляем
ymm4 = _mm256_and_ps (ymm4, ymm11); //применяем маску немодифицированного канала
ymm0 = _mm256_andnot_ps (ymm9, ymm0 ); //обнуляем каналы где немодифицированное или val значение
ymm0 = _mm256_or_ps (ymm0, ymm4 ); //объединяем с немодифицированным каналом
ymm0 = _mm256_or_ps (ymm0, ymm10 ); //объединяем с val значением
_mm256_storeu_ps(pDst+d, ymm0);
}
//обрабатываем хвост по 1 пикселю
pSrc = (float *)((char*)pSrc + srcStep);
pDst = (float *)((char*)pDst + dstStep);
}
Закономерно может возникнуть вопрос: почему же в нашем 256-битном AVX коде используется загрузка с помощью 128 битных инструкций? Дело в том, что в нашем коде была использована инструкция VPERMPS в строке
ymm0 = _mm256_permutevar_ps(ymm0, ymm8)
Как уже объяснялось выше, данная инструкция переставляет элементы внутри младшей и старшей части регистра (или по другому lane) независимо. Поэтому, несмотря на то, что 2 пикселя лежат в памяти рядом, мы вынуждены читать их таким образом, чтобы они попали в разные lane. Обратите внимание, как формируется регистр-маска ymm8, он содержит абсолютно одинаковые значения в старшей и младшей части для перестановки. Код написан таким образом, чтобы количество итераций в циклах сохранилось и не влияло на итоговую производительность. Измеряем производительность: P=0.31. Таким образом прирост AVX кода по отношению к C составил 3.2X, а у SSSE3 было 2.48X
Третья версия SIMD кода с использованием AVX2
Развитием AVX стал набор AVX2. Длина векторных регистров в нем не изменилась, но появились новые инструкции, которые нам помогут еще ускорить код. В AVX2 появились acrosslane инструкции перестановки _mm256_permutevar8x32_ps. Эти инструкции позволяют переставить элементы внутри 256 битного регистра произвольным образом и наш код мог бы выглядеть уже так.Код
__m256i ymm8 = _mm256_set_epi32(
dstOrder[3]+3,dstOrder[2]+3,dstOrder[1]+3,dstOrder[0]+3,
dstOrder[3],dstOrder[2],dstOrder[1],dstOrder[0]);
for(; d < (width&(~1))*4; s += 3*2, d += 4*2 ) {
ymm0 = _mm256_loadu_ps (pSrc+s);
ymm4 = _mm256_loadu_ps (pDst+d );
ymm0 = _mm256_permutevar8x32_ps(ymm0, ymm8);
ymm4 = _mm256_and_ps (ymm4, ymm11); //apply unmodified mask
ymm0 = _mm256_andnot_ps (ymm9, ymm0 ); //apply unmodified and val mask with inverse
ymm0 = _mm256_or_ps (ymm0, ymm4 ); //join with unmodified
ymm0 = _mm256_or_ps (ymm0, ymm10 ); //join with val
_mm256_storeu_ps(pDst+d, ymm0);
}
Обратите внимание что маска перестановки в регистре ymm8 совершенно изменилась и является сквозной по всему регистру, т.е. индексы в ней относятся в целом ко всему 256 битному регистру.
Однако в этом коде есть проблема в том, как происходит чтение регистра ymm0. В нем загружаются лишние 2 float элемента. Так делать неправильно, поскольку может произойти чтение за границей изображения и исправить это нам поможет опять же появившаяся лишь в AVX2 инструкция чтения по маске
__m256i _mm256_maskload_epi32(int const *, __m256i);
Данная инструкция читает необходимые элементы из памяти и не читает ничего лишнего, что предотвращает возможную memory segmentation ошибку. Теперь наш правильный код будет выглядеть так.
Код
__m256i ld_mask = _mm256_set_epi32(0,0,-1,-1,-1,-1,-1,-1);
for(; d < (width&(~1))*4; s += 3*2, d += 4*2 ) {
ymm0 = _mm256_maskload_ps (pSrc+s, ld_mask);
ymm4 = _mm256_loadu_ps (pDst+d );
ymm0 = _mm256_permutevar8x32_ps(ymm0, ymm8);
ymm4 = _mm256_and_ps (ymm4, ymm11); //apply unmodified mask
ymm0 = _mm256_andnot_ps (ymm9, ymm0 ); //apply unmodified and val mask with inverse
ymm0 = _mm256_or_ps (ymm0, ymm4 ); //join with unmodified
ymm0 = _mm256_or_ps (ymm0, ymm10 ); //join with val
_mm256_storeu_ps(pDst+d, ymm0);
}
Производительность этого кода составляет P=0.286. Ускорение таким образом изменилось от k=3.2X до k=3.5X. Помимо чтения по маске, мы можем использовать и сохранять немодифицированный канал без изменения, используя инструкцию
void _mm256_maskstore_ps(float *, __m256i, __m256);
Код
__m256 ymmNA = _mm256_set_epi32( 0,-1,-1,-1, 0,-1,-1,-1);
for(; d < (width&(~1))*4; s += 3*2, d += 4*2 ) {
ymm0 = _mm256_maskload_ps (pSrc+s, ld_mask);
ymm0 = _mm256_permutevar8x32_ps(ymm0, ymm8);
ymm0 = _mm256_andnot_ps (ymm9, ymm0 ); //apply unmodified and val mask with inverse
ymm0 = _mm256_or_ps (ymm0, ymm10 ); //join with val
_mm256_maskstore_epi32(pDst+d, ymmNA, ymm0);
}
Обратите внимание на то, что 128 битная предшественница данной инструкции
void _mm_maskmoveu_si128(__m128i, __m128i, char *)
сохраняет данные, используя non-temporal запись, а вот в AVX версиях _mm256_maskstore_ps и _mm_maskstore_ps от этого механизма было решено отказаться. О том, что такое non-temporal запись и как она применяется в IPP библиотеке будет рассказано в одной из следующих статей. Меряем производительность и получаем небольшой прирост, P=0.27, k=3.7X.
Версия SIMD кода с использованием AVX512
Итак, переходим к самой интересной части нашей статьи – оптимизации кода на AVX512.Как следует из названия длина векторных регистров увеличилась в два раза с 256 бит до 512. Lane концепция также сохранилась. Второе новшество заключается в том, что операции с масками приобрели массовый характер и затрагивают теперь большинство векторных операций. Для этого были введены специальные масочные регистры. Масочные регистры имеют свой тип данных __mmask8, __mmask16, __mmask32, __mmask64, но с ними можно оперировать как с целочисленными типами данных, поскольку они так и объявлены в заголовочном файле zmmintrin.h. Использование масок делает AVX512 обрабатывающий цикл функции довольно минималистичным.
Код
int s, d, h, ord0 = dstOrder[0], ord1 = dstOrder[1], ord2 = dstOrder[2], ord3 = dstOrder[3];
__mmask16 mskv; //mask for joing with val
__mmask16 mska; //mask for unmodified channel
__m512 z0,z1,z2,z3,z4,z5;
__m512i zv = _mm512_set1_epi32(val);
__m512i zi0 = _mm512_set_epi32(ord3+9,ord2+9,ord1+9,ord0+9, ord3+6,ord2+6,ord1+6,ord0+6,
ord3+3,ord2+3,ord1+3,ord0+3, ord3+0,ord2+0,ord1+0,ord0+0);
mska = mskv = 0xFFFF;
for (int i = 0; i < 4; i++){
if(dstOrder[i] == 3) { mskv ^= (0x1111<<i); }
if(dstOrder[i] > 3) { mska ^= (0x1111<<i); }
}
for( h = 0; h < height; h++ ) {
s = 0; d = 0;
for(; d < ((width*4) & (~15)); s += 3*4, d += 4*4 ) {
z0 = _mm512_mask_loadu_ps(_mm512_setzero_ps(), 0xFFF,(float *)(pSrc+s));
z0 = _mm512_mask_permutexvar_ps(zv,mskv,zi0, z0);
_mm512_mask_storeu_ps(pDst+d, mska, z0);
}
pSrc = (float *)((char*)pSrc + srcStep);
pDst = (float *)((char *)pDst + dstStep);
}
Количество интринсиков в одной итерации цикла сократилось до трех, и каждый из них использует маску. Первый интринсик
z0 = _mm512_mask_loadu_ps(_mm512_setzero_ps(), 0xFFF, (float *)(pSrc+s))
использует маску 0xFFF для того чтобы прочитать из памяти 4 RGB пикселя или 12 float, в двоичном виде 0xFFF=1111111111111111b. Следующий интринсик
z0 = _mm512_mask_permutexvar_ps(zv, mskv, zi0, z0)
выполняет сразу 2 функции: переставляет каналы входного пикселя и объединяет их по маске со значением val. Маска mskv формируется перед началом цикла и биты маски используются следующим образом: значение бита равное ‘1’ означает «поместить результат операций в результирующий регистр», а ‘0’ – «cохранить значение того регистра, который указан до маски», в нашем случае zv. Третий интринсик
_mm512_mask_storeu_ps(pDst+d, mska, z0)
сохраняет полученные 4 пикселя в память, при этом он не пишет ничего в немодифицированный канал. Производительность кода составляет P=0.183, k=5.44X.
Читатель может сравнить код SSSE3 с AVX512 кодом, в котором используется всего 3 интринсика. Именно краткость avx512 кода и вдохновила автора на написание данной статьи.
Итоговая таблица роста производительности.

Заключение
Конечно, использование масок в данном алгоритме довольно очевидно, но они могут использованы и в других, более сложных алгоритмах, и об этом будет другая статья.Ознакомившись с данной статьей, читатель может попробовать использовать SIMD код в своих приложениях, а читатель, уже имеющий такой опыт, может оценить по достоинству преимущества crosslane инструкций перестановки и инструкций с масками в новых процессорах Intel, поддерживающих набор инструкций AVX512.
Также хотелось бы отметить, что в данной статье приведена незначительная часть приемов оптимизации, используемых в библиотеке IPP и позволяющих достичь ее уникальной производительности.
Читатель может и сам попробовать использовать библиотеку IPP в своем приложении, освободив возможности для решения своих прикладных задач. Дополнительным эффектом использования IPP станет то, что после выхода на рынок процессоров с набором инструкции AVX512, библиотека IPP уже будет иметь оптимизированную на AVX512 версию наиболее востребованных функций.