SIMD
From Wiki.conus.info
Single instruction, multiple data.
Как можно судить по названию, это обработка сразу пачки данных исполняя только одну инструкцию.
Это, как и FPU, также выглядит отдельным процессором внутри x86.
В x86 все началось с MMX. Появилось 8 64-битных регистров MM0-MM7.
Каждый регистр может содержать 2 32-битных значения, 4 16-битных или же 8 байт. Можно, например, складывая значения в двух MMX-регистрах, складывать одновременно 8 значений, каждое из которых 8-битное.
Простой пример, это некий графический редактор, который хранит открытое изображение как двумерный массив. Когда пользователь меняет яркость изображения, редактору нужно, например, прибавить некий коэффициент ко всем пикселям, или отнять. Если гипотетически представить, что изображение у нас бело-серо-черное и каждый пиксель занимает один байт, то с помощью MMX можно менять яркость сразу у восьми пикселей.
Когда MMX только появилось, эти регистры на самом деле распологались в FPU. Можно было использовать либо FPU или MMX в одно и то же время. Можно подумать что Intel решило немного сэкономить на транзисторах, но на самом деле причина проще - операционная система не знающая о дополнительных регистрах процессора не будет сохранять их во время переключения задач, а вот регистры FPU сохранять будет. Таким образом, процессор с MMX + старая операционная система + задача использующая MMX = все это может работать вместе.
SSE - это расширение регистров до 128 бит, теперь уже отдельно от FPU.
Имеются также планы расширения регистров до 256 бит (AVX).
Немного о практическом применении.
Конечно же, копирование блоков в памяти (memcpy), сравнение (memcmp), и подобное.
Еще пример: имеется алгоритм шифрования DES, который берет 64-битный блок, 56-битный ключ, шифрует блок с ключем и образуется 64-битный результат. Алгоритм DES можно легко представить в виде очень большой электронной цифровой схемы, с проводами, элементами И, ИЛИ, НЕ.
Идея bitslice DES - это обработка сразу группы блоков-ключей одновременно. Скажем, на x86 перменная типа unsigned int вмещает в себе 32 бита, так что там можно хранить промежуточные результаты сразу для 32-х блоков-ключей.
http://www.darkside.com.au/bitslice/
Я написал утилиту для перебора паролей/хешей Oracle RDBMS (которые основаны на алгоритме DES), только переделал алгоритм bitslice DES для SSE2 - и теперь возможно шифровать одновременно 128 блоков-ключей, там и исходники:
http://conus.info/utils/#ops_sse2
Contents |
1
Теперь чуть сложнее.
Насчет векторизации.
http://en.wikipedia.org/wiki/Vectorization_(computer_science)
Векторизация это когда у вас есть цикл, который берет на вход несколько массивов и выдает, например, один массив данных. Тело цикла берет некоторые элементы из входных массивов, что-то делает с ними и кладет в выходной. Важно что операция применяемая ко всем элементам одна и та же. Векторизация это обрабатывать одновременно несколько элементов.
Например:
for (i = 0; i < 1024; i++) { C[i] = A[i]*B[i]; }
Этот кусок кода берет элемент из A и B, перемножает и отправляет результат в C.
Если представить что каждый элемент массива - это 32-битный integer, то их можно загружать сразу по 4 из А в 128-битный XMM-регистр, из B в другой XMM-регистр и выполнив инстукцию PMULLD (Multiply Packed Signed Dword Integers and Store Low Result) и PMULHW (Multiply Packed Signed Integers and Store High Result), можно получить 4 64-битных результата умножения.
Таким образом, тело цикла можно будет исполнять не 1024 раза, а 1024/4 раза, что в 4 раза меньше, и, возможно, быстрее.
Некоторые компиляторы немного умеют делать автоматическую векторизацию, например Intel C++.
Я написал очень простую функцию:
int f (int sz, int *ar1, int *ar2, int *ar3) { for (int i=0; i<sz; i++) ar3[i]=ar1[i]+ar2[i]; return 0; };
Intel C++
Компилю при помощи Intel C++ 11.1.051 win32:
icl intel.cpp /QaxSSE2 /Faintel.asm /Ox
Имеем такое (в IDA):
; int __cdecl f(int, int *, int *, int *) public ?f@@YAHHPAH00@Z ?f@@YAHHPAH00@Z proc near var_10 = dword ptr -10h sz = dword ptr 4 ar1 = dword ptr 8 ar2 = dword ptr 0Ch ar3 = dword ptr 10h push edi push esi push ebx push esi mov edx, [esp+10h+sz] test edx, edx jle loc_15B mov eax, [esp+10h+ar3] cmp edx, 6 jle loc_143 cmp eax, [esp+10h+ar2] jbe short loc_36 mov esi, [esp+10h+ar2] sub esi, eax lea ecx, ds:0[edx*4] neg esi cmp ecx, esi jbe short loc_55 loc_36: ; CODE XREF: f(int,int *,int *,int *)+21�j cmp eax, [esp+10h+ar2] jnb loc_143 mov esi, [esp+10h+ar2] sub esi, eax lea ecx, ds:0[edx*4] cmp esi, ecx jb loc_143 loc_55: ; CODE XREF: f(int,int *,int *,int *)+34�j cmp eax, [esp+10h+ar1] jbe short loc_67 mov esi, [esp+10h+ar1] sub esi, eax neg esi cmp ecx, esi jbe short loc_7F loc_67: ; CODE XREF: f(int,int *,int *,int *)+59�j cmp eax, [esp+10h+ar1] jnb loc_143 mov esi, [esp+10h+ar1] sub esi, eax cmp esi, ecx jb loc_143 loc_7F: ; CODE XREF: f(int,int *,int *,int *)+65�j mov edi, eax ; edi = ar1 and edi, 0Fh ; is ar1 16-byte aligned? jz short loc_9A ; yes test edi, 3 jnz loc_162 neg edi add edi, 10h shr edi, 2 loc_9A: ; CODE XREF: f(int,int *,int *,int *)+84�j lea ecx, [edi+4] cmp edx, ecx jl loc_162 mov ecx, edx sub ecx, edi and ecx, 3 neg ecx add ecx, edx test edi, edi jbe short loc_D6 mov ebx, [esp+10h+ar2] mov [esp+10h+var_10], ecx mov ecx, [esp+10h+ar1] xor esi, esi loc_C1: ; CODE XREF: f(int,int *,int *,int *)+CD�j mov edx, [ecx+esi*4] add edx, [ebx+esi*4] mov [eax+esi*4], edx inc esi cmp esi, edi jb short loc_C1 mov ecx, [esp+10h+var_10] mov edx, [esp+10h+sz] loc_D6: ; CODE XREF: f(int,int *,int *,int *)+B2�j mov esi, [esp+10h+ar2] lea esi, [esi+edi*4] ; is ar2+i*4 16-byte aligned? test esi, 0Fh jz short loc_109 ; yes! mov ebx, [esp+10h+ar1] mov esi, [esp+10h+ar2] loc_ED: ; CODE XREF: f(int,int *,int *,int *)+105�j movdqu xmm1, xmmword ptr [ebx+edi*4] movdqu xmm0, xmmword ptr [esi+edi*4] ; ar2+i*4 is not 16-byte aligned, so load it to xmm0 paddd xmm1, xmm0 movdqa xmmword ptr [eax+edi*4], xmm1 add edi, 4 cmp edi, ecx jb short loc_ED jmp short loc_127 ; --------------------------------------------------------------------------- loc_109: ; CODE XREF: f(int,int *,int *,int *)+E3�j mov ebx, [esp+10h+ar1] mov esi, [esp+10h+ar2] loc_111: ; CODE XREF: f(int,int *,int *,int *)+125�j movdqu xmm0, xmmword ptr [ebx+edi*4] paddd xmm0, xmmword ptr [esi+edi*4] movdqa xmmword ptr [eax+edi*4], xmm0 add edi, 4 cmp edi, ecx jb short loc_111 loc_127: ; CODE XREF: f(int,int *,int *,int *)+107�j ; f(int,int *,int *,int *)+164�j cmp ecx, edx jnb short loc_15B mov esi, [esp+10h+ar1] mov edi, [esp+10h+ar2] loc_133: ; CODE XREF: f(int,int *,int *,int *)+13F�j mov ebx, [esi+ecx*4] add ebx, [edi+ecx*4] mov [eax+ecx*4], ebx inc ecx cmp ecx, edx jb short loc_133 jmp short loc_15B ; --------------------------------------------------------------------------- loc_143: ; CODE XREF: f(int,int *,int *,int *)+17�j ; f(int,int *,int *,int *)+3A�j ... mov esi, [esp+10h+ar1] mov edi, [esp+10h+ar2] xor ecx, ecx loc_14D: ; CODE XREF: f(int,int *,int *,int *)+159�j mov ebx, [esi+ecx*4] add ebx, [edi+ecx*4] mov [eax+ecx*4], ebx inc ecx cmp ecx, edx jb short loc_14D loc_15B: ; CODE XREF: f(int,int *,int *,int *)+A�j ; f(int,int *,int *,int *)+129�j ... xor eax, eax pop ecx pop ebx pop esi pop edi retn ; --------------------------------------------------------------------------- loc_162: ; CODE XREF: f(int,int *,int *,int *)+8C�j ; f(int,int *,int *,int *)+9F�j xor ecx, ecx jmp short loc_127 ?f@@YAHHPAH00@Z endp
Инструкции которые имеют отношение к SSE2 это:
MOVDQU (Move Unaligned Double Quadword) - она просто загружает 16 байт из памяти в XMM-регистр.
PADDD (Add Packed Integers) - складывает сразу 4 пары чисел и оставит в первом операнде результат. Кстати, если произойдет переполнение, то исключения не произойдет, запишутся просто младшие 32 бита результата. Эти инструкции не выставляют никаких флагов. Если один из операндов PADDD - адрес в памяти, то требуется чтобы адрес был выровнен по 16-байтной границе. Если он не выровнен, произойдет исключение.
MOVDQA (Move Aligned Double Quadword) - тоже что и MOVDQU, только подразумевает что адрес в памяти выровнен по 16-байтной границе. Если он не выровнен, произойдет исключение. MOVDQA работает быстрее чем MOVDQU, но требует вышеозначенного.
См. также: http://en.wikipedia.org/wiki/Data_structure_alignment
Итак, SSE2-инструкции исполнятся только в том случае если осталось просуммировать 4 пары плюс если указатель ar3 выровнен по 16-байтной границе.
Более того, если еще и ar2 выровнен по 16-байтной границе, то будет выполняться этот кусок:
movdqu xmm0, xmmword ptr [ebx+edi*4] ; ar1+i*4 paddd xmm0, xmmword ptr [esi+edi*4] ; ar2+i*4 movdqa xmmword ptr [eax+edi*4], xmm0 ; ar3+i*4
А иначе, значение из ar2 загрузится в xmm0 используя инструкцию MOVDQU, которая не требует выровненного указателя в памяти, зато может работать чуть медленнее:
movdqu xmm1, xmmword ptr [ebx+edi*4] ; ar1+i*4 movdqu xmm0, xmmword ptr [esi+edi*4] ; ar2+i*4 is not 16-byte aligned, so load it to xmm0 paddd xmm1, xmm0 movdqa xmmword ptr [eax+edi*4], xmm1 ; ar3+i*4
А во всех остальных случаях, будет исполняться код, который был бы как если бы не была включена поддержка SSE2.
Еще о том как Intel C++ умеет автоматически векторизировать циклы: http://www.intel.com/intelpress/sum_vmmx.htm (Excerpt: Effective Automatic Vectorization)
GCC
Но и GCC умеет кое-что векторизировать, если компилировать с опциями -O3 и включить поддержку SSE2: -msse2.
Вот что вышло (GCC 4.4.1):
; f(int, int *, int *, int *) public _Z1fiPiS_S_ _Z1fiPiS_S_ proc near var_18 = dword ptr -18h var_14 = dword ptr -14h var_10 = dword ptr -10h arg_0 = dword ptr 8 arg_4 = dword ptr 0Ch arg_8 = dword ptr 10h arg_C = dword ptr 14h push ebp mov ebp, esp push edi push esi push ebx sub esp, 0Ch mov ecx, [ebp+arg_0] mov esi, [ebp+arg_4] mov edi, [ebp+arg_8] mov ebx, [ebp+arg_C] test ecx, ecx jle short loc_80484D8 cmp ecx, 6 lea eax, [ebx+10h] ja short loc_80484E8 loc_80484C1: ; CODE XREF: f(int,int *,int *,int *)+4B�j ; f(int,int *,int *,int *)+61�j ... xor eax, eax nop lea esi, [esi+0] loc_80484C8: ; CODE XREF: f(int,int *,int *,int *)+36�j mov edx, [edi+eax*4] add edx, [esi+eax*4] mov [ebx+eax*4], edx add eax, 1 cmp eax, ecx jnz short loc_80484C8 loc_80484D8: ; CODE XREF: f(int,int *,int *,int *)+17�j ; f(int,int *,int *,int *)+A5�j add esp, 0Ch xor eax, eax pop ebx pop esi pop edi pop ebp retn ; --------------------------------------------------------------------------- align 8 loc_80484E8: ; CODE XREF: f(int,int *,int *,int *)+1F�j test bl, 0Fh jnz short loc_80484C1 lea edx, [esi+10h] cmp ebx, edx jbe loc_8048578 loc_80484F8: ; CODE XREF: f(int,int *,int *,int *)+E0�j lea edx, [edi+10h] cmp ebx, edx ja short loc_8048503 cmp edi, eax jbe short loc_80484C1 loc_8048503: ; CODE XREF: f(int,int *,int *,int *)+5D�j mov eax, ecx shr eax, 2 mov [ebp+var_14], eax shl eax, 2 test eax, eax mov [ebp+var_10], eax jz short loc_8048547 mov [ebp+var_18], ecx mov ecx, [ebp+var_14] xor eax, eax xor edx, edx nop loc_8048520: ; CODE XREF: f(int,int *,int *,int *)+9B�j movdqu xmm1, xmmword ptr [edi+eax] movdqu xmm0, xmmword ptr [esi+eax] add edx, 1 paddd xmm0, xmm1 movdqa xmmword ptr [ebx+eax], xmm0 add eax, 10h cmp edx, ecx jb short loc_8048520 mov ecx, [ebp+var_18] mov eax, [ebp+var_10] cmp ecx, eax jz short loc_80484D8 loc_8048547: ; CODE XREF: f(int,int *,int *,int *)+73�j lea edx, ds:0[eax*4] add esi, edx add edi, edx add ebx, edx lea esi, [esi+0] loc_8048558: ; CODE XREF: f(int,int *,int *,int *)+CC�j mov edx, [edi] add eax, 1 add edi, 4 add edx, [esi] add esi, 4 mov [ebx], edx add ebx, 4 cmp ecx, eax jg short loc_8048558 add esp, 0Ch xor eax, eax pop ebx pop esi pop edi pop ebp retn ; --------------------------------------------------------------------------- loc_8048578: ; CODE XREF: f(int,int *,int *,int *)+52�j cmp eax, esi jnb loc_80484C1 jmp loc_80484F8 _Z1fiPiS_S_ endp
Почти то же самое, хотя и не так дотошно как Intel C++.
Подробнее о GCC: http://gcc.gnu.org/projects/tree-ssa/vectorization.html
2
Прежде всего, следует заметить, что SIMD-инструкции можно вставлять в код при помощи специальных макросов.
http://msdn.microsoft.com/en-us/library/y0dh78ez(VS.80).aspx
В MSVC, часть находится в файле intrin.h.
Имеется возможность сделать strlen() при помощи SIMD-инструкций, работающий в 2-2.5 раза быстрее обычного. Эта функция будет загружать в XMM-регистр сразу 16 байт и проверять каждый на ноль.
size_t strlen_sse2(const char *str) { register size_t len = 0; const char *s=str; bool str_is_aligned=(((unsigned int)str)&0xFFFFFFF0) == (unsigned int)str; if (str_is_aligned==false) return strlen (str); __m128i xmm0 = _mm_setzero_si128(); __m128i xmm1; int mask = 0; for (;;) { xmm1 = _mm_load_si128((__m128i *)s); xmm1 = _mm_cmpeq_epi8(xmm1, xmm0); if ((mask = _mm_movemask_epi8(xmm1)) != 0) { unsigned long pos; _BitScanForward(&pos, mask); len += (size_t)pos; break; } s += sizeof(__m128i); len += sizeof(__m128i); }; return len; }
(пример базируется на исходнике отсюда).
Компилируем в MSVC 2010 с опцией /Ox:
_pos$75552 = -4 ; size = 4 _str$ = 8 ; size = 4 ?strlen_sse2@@YAIPBD@Z PROC ; strlen_sse2 push ebp mov ebp, esp and esp, -16 ; fffffff0H mov eax, DWORD PTR _str$[ebp] sub esp, 12 ; 0000000cH push esi mov esi, eax and esi, -16 ; fffffff0H xor edx, edx mov ecx, eax cmp esi, eax je SHORT $LN4@strlen_sse lea edx, DWORD PTR [eax+1] npad 3 $LL11@strlen_sse: mov cl, BYTE PTR [eax] inc eax test cl, cl jne SHORT $LL11@strlen_sse sub eax, edx pop esi mov esp, ebp pop ebp ret 0 $LN4@strlen_sse: movdqa xmm1, XMMWORD PTR [eax] pxor xmm0, xmm0 pcmpeqb xmm1, xmm0 pmovmskb eax, xmm1 test eax, eax jne SHORT $LN9@strlen_sse $LL3@strlen_sse: movdqa xmm1, XMMWORD PTR [ecx+16] add ecx, 16 ; 00000010H pcmpeqb xmm1, xmm0 add edx, 16 ; 00000010H pmovmskb eax, xmm1 test eax, eax je SHORT $LL3@strlen_sse $LN9@strlen_sse: bsf eax, eax mov ecx, eax mov DWORD PTR _pos$75552[esp+16], eax lea eax, DWORD PTR [ecx+edx] pop esi mov esp, ebp pop ebp ret 0 ?strlen_sse2@@YAIPBD@Z ENDP ; strlen_sse2
Итак, прежде всего, мы проверяем указатель str, выровнен ли он по 16-байтной границе. Если нет, то мы вызовем обычную strlen().
Далее мы загружаем по 16 байт в регистр XMM при помощи команды MOVDQA.
Наблюдательный читатель может спросить, почему в этом месте мы не можем использовать MOVDQU, которая может загружать откуда угодно не взирая на тот факт, выровнен ли указатель?
Да, можно было бы сделать вот как: если указатель выровнен, загружаем используя MOVDQA, иначе используем работающую чуть медленнее MOVDQU.
Однако здесь имеется не сразу заметная проблема, которая проявляется вот в чем:
В ОС линии Windows NT, и не только, память выделяется страницами по 4 кибибайта (4096). Каждый win32-процесс якобы имеет в наличии 4 гибибайта, но на самом деле, только некоторые части этого адресного пространства присоеденены к реальной физической памяти. Если процесс обратится к блоку памяти, которого не существует, сработает исключение. Так работает виртуальная память.
http://en.wikipedia.org/wiki/Page_(computer_memory)
Так вот, функция, читающая сразу по 16 байт, имеет возможность нечаянно вылезти за границу выделенного блока памяти. Предположим, ОС выделила программе 8192 (0x2000) байт по адресу 0x008c0000. Таким образом, блок занимает байты с адреса 0x008c0000 по 0x008c1fff включительно.
За этим блоком, то есть начиная с адреса 0x008c2000 нет вообще ничего, т.е., ОС не выделила для нас ничего. Обращение к памяти начиная с этого адреса вызовет исключение.
И предположим, что программа хранит некую строку из, скажем, пяти символов почти в самом конце блока, что не является преступлением:
0x008c1ff8 - 'h' 0x008c1ff9 - 'e' 0x008c1ffa - 'l' 0x008c1ffb - 'l' 0x008c1ffc - 'o' 0x008c1ffd - '\x00' 0x008c1ffe - здесь случайный мусор 0x008c1fff - здесь случайный мусор
Затем программа хочет вызвать strlen() передав ей указатель на строку 'hello' по адресу 0x008c1ff8. strlen() будет читать по одному байту до 0x008c1ffd, где ноль, и здесь он закончит.
Теперь, если мы напишем strlen() читающий сразу по 16 байт, с любого адреса, будь он выровнен по 16-байтной границе или нет, MOVDQU попытается загрузить 16 байт с адреса 0x008c1ff8 по 0x008c2008, и произойдет исключение. Это ситуация которой, конечно, хочется избежать.
Поэтому мы будем работать только с адресами выровненными по 16 байт, что в сочетании со знанием что размер страницы также как правило выровнен по 16 байт, даст некоторую гарантию что наша функция не будет пытаться читать из неожиданных мест.
Вернемся к нашей функции.
_mm_setzero_si128 это макрос, генерирующий pxor xmm0, xmm0 - это просто обнуляет регистр XMM0.
_mm_load_si128 - это макрос для MOVDQA, он просто загружает 16 байт по адресу из указателя в XMM1.
_mm_cmpeq_epi8 - это макрос для PCMPEQB, это инструкция которая побайтово сравнивает значения из двух XMM регистров. Иными словами, 16 раз сравниваются все байтовые части регистров. И если какой-то из байт равен другому, то в результирующем значении будет выставлено на месте этого байта 0xff, или 0, если байты не равны.
Например.
XMM1: 11223344556677880000000000000000 XMM0: 11ab3444007877881111111111111111
После исполнения pcmpeqb xmm1, xmm0, регистр XMM1 будет содержать:
XMM1: ff0000ff0000ffff0000000000000000
Эта инструкция в нашем случае, сравнивает каждый 16-байтный блок с блоком состоящим из 16-и нулей, выставленным в XMM0 при помощи pxor xmm0, xmm0.
Следующий макрос _mm_movemask_epi8 - это инструкция PMOVMSKB.
Она очень удобна как раз для использования в паре с PCMPEQB.
pmovmskb eax, xmm1
Эта инструкция выставит самый первый бит eax в еденицу, если старший бит первого байта в регистре XMM1 является еденицей. Иными словами, если первый байт в регистре XMM1 является 0xff, то первый бит в EAX будет также еденицей, иначе нулем. Если второй байт в регистре XMM1 является 0xff, то второй бит в EAX также будет еденицей. Иными словами, инструкция отвечает на вопрос, какие из байт в XMM1 являются 0xff, в результате приготовит 16 бит и запишет в EAX. Остальные биты в EAX обнулятся.
Кстати, не забывайте также вот о какой особенности нашего алгоритма.
На вход может прийти 16 байт вроде hello\x00garbage\x00ab
Это строка 'hello', после нее терминирующий ноль, затем немного мусора в памяти.
Если мы загрузим эти 16 байт в XMM1 сравним нулевым XMM0, то в итоге получим такое (я использую здесь порядок с MSB до LSB):
XMM1: 0000ff00000000000000ff0000000000
Это означает что инструкция сравнения обнаружила два нулевых байта, что и не удивительно.
PMOVMSKB в нашем случае подготовит EAX вот так (в двоичном представлении): 0010000000100000b.
Совершенно очевидно что далее наша функция должна учитывать только первый встретившийся ноль и игнорировать все остальное.
Следующая инструкция - BSF (Bit Scan Forward). Это инструкция находит самый младший бит во втором операнде и записывает его позицию в первый операнд.
EAX=0010000000100000b
bsf eax, eax
После исполнения инструкции, в EAX будет 5, что означает, что еденица найдена в пятой позиции (считая с нуля).
Для использования этой инструкции, в MSVC также имеется макрос _BitScanForward.
А дальше все просто. Если нулевой байт найден, его позиция прибавляется к тому что мы уже насчитали и возвращается результат.
Почти всё.
Кстати, следует также отметить, что компилятор MSVC закодировал два тела цикла сразу, для оптимизации.
Кстати, в SSE 4.2 (который появился в Intel Core i7) все эти манипуляции со строками могут быть еще проще: http://www.strchr.com/strcmp_and_strlen_using_sse_4.2
