Вычислите f(x)=2*(x^2)+5 с насыщением, используя набор инструкций MMX для 128 чисел размером 2 байта, загруженных из двоичного файла
У меня есть эта проблема, где мне нужно вычислить функцию f(x)=2*(x^2)+5 с набором инструкций MMX. У меня две проблемы. Это мой код сейчас:
section .data
print_fmt db '%d', 10, 0
my_loaded_data times 128 dw 0
fives times 4 dw 5
twos times 4 dw 2
my_loaded_data_file_name db 'test_numbers.bin', 0
mod_f db 'r', 0
section .text
extern printf, fopen, fread
global main
main:
PUSH rbp
mov rbp, rsp
mov rax, 0
mov rdi, my_loaded_data_file_name
mov rsi, mod_f
call fopen
cmp rax, 0
je .end
PUSH rax
PUSH rdi
PUSH rsi
PUSH rdx
PUSH rcx
mov rdi, my_loaded_data
mov rsi, 2
mov rdx, 128
mov rcx, rax
mov rax, 0
call fread
POP rcx
POP rdx
POP rsi
POP rdi
POP rax
mov rsi, my_loaded_data
mov rcx, 32
jmp .square_loop
.square_loop:
movq mm0, [rsi]
movq mm1, [rsi]
pmulhw mm0, mm1
movq [rsi], mm0
add rsi, 8
loop .square_loop
mov rcx, 32
mov rsi, my_loaded_data
movq mm1, [twos]
jmp .mult_with2_loop
.mult_with2_loop:
movq mm0, [rsi]
pmulhw mm0, mm1
movq [rsi], mm0
add rsi, 8
loop .mult_with2_loop
mov rcx, 32
mov rsi, my_loaded_data
movq mm1, [fives]
jmp .add_five_loop
.add_five_loop:
movq mm0, [rsi]
paddusw mm0, mm1
movq [rsi], mm0
add rsi, 8
loop .add_five_loop
jmp .write_it
.write_it:
mov r8, my_loaded_data
mov rcx, 128
.write_loop:
mov rax, 0
mov ax, [r8]
PUSH r8
PUSH rcx
PUSH rdi
PUSH rsi
PUSH rax
mov rdi, print_fmt
mov rsi, rax
mov rax, 0
call printf
POP rax
POP rsi
POP rdi
POP rcx
POP r8
add r8, 2
loop .write_loop
.end:
mov rax, 0
POP rbp
ret
Моя первая проблема - инструкция умножения. Какую инструкцию я использую для насыщения. Сначала я подумал, что будет такая инструкция pmulsw но, похоже, нет. pmulhw сохраняет верхние 16-битные 32-битные результаты. Я не могу найти инструкцию, которая дала бы 16-битный результат. Это единственный способ сохранить 32-битные результаты?
Вторая проблема связана с printf. Это продолжает давать ошибку сегментации, и я не знаю почему. Это из моего терминала:
Program received signal SIGSEGV, Segmentation fault.
__GI___tcgetattr (fd=1, termios_p=termios_p@entry=0x7ffffffed9a8) at ../sysdeps/unix/sysv/linux/tcgetattr.c:42
42 ../sysdeps/unix/sysv/linux/tcgetattr.c: No such file or directory.
Вот make-файл:
zad7.o: zad7.asm
nasm -f elf64 -g -F dwarf zad7.asm -o zad7.o
zad7: zad7.o
gcc -o zad7 zad7.o -no-pie -ggdb
А для вашего удобства есть небольшая программа на C, которая может генерировать двоичный файл для чтения:
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
void make_data_file(unsigned int number_of_data, unsigned int size_in_bytes, char* file_name)
{
FILE *write_ptr = fopen(file_name, "wb");
if(write_ptr==NULL)
{
printf("Error creating file '%s'.\n", file_name);
return;
}
double n_bits = size_in_bytes * 8;
unsigned int max_num = pow(2, n_bits);
unsigned short random_number;
for(int i=0; i< number_of_data; i++)
{
random_number = i;
fwrite(&random_number, size_in_bytes, 1, write_ptr);
}
fclose(write_ptr);
}
int main()
{
make_data_file(128, 2, "test_numbers.bin");
return 0;
}
2 ответа
Если вы заботитесь о производительности, избегайте loop инструкция по современным процессорам. Почему инструкция цикла медленная? Разве Intel не могла реализовать это эффективно?, Также используйте SSE2 вместо MMX; ваш размер массива кратен 16 и 8, и вы используете x86-64, который гарантированно будет иметь SSE2. MMX абсолютно бессмысленен, если вы не делаете 32-разрядную версию для Pentium III / Athlon-XP и более ранних версий.
(Весь код в моем ответе будет работать с 8-байтовыми регистрами MMX вместо 16-байтовых регистров XMM, потому что есть MMX-версии всех инструкций, которые я использовал. Согласно приложению B к руководству NASM, pmullw, pxor, pcmpgtw, а также paddusw были все доступны в оригинальном P5 Pentium MMX. Некоторые инструкции, которые в руководстве Intel указаны как "MMX" (например, pmulhuw а также pshufw) были добавлены только позже, как, например, с Pentium II или вместе с SSE в Pentium III, но это не относится к инструкциям, которые были здесь полезны.)
См. https://stackru.com/tags/x86/info для получения руководств по производительности / оптимизации, а также ссылок ABI / соглашений о вызовах, которые объясняют 16-байтовое выравнивание стека, необходимое для вызова функций.
mov rax, 0 / mov ax, [r8] это действительно глупо. использование movzx eax, word [r8] как нормальный человек. Вам также не нужно jmp на следующую строку источника, как jmp .square_loop / .square_loop: Выполнение всегда переходит к следующей строке само по себе, если нет инструкции перехода.
x86 SIMD не имеет умножающего насыщения, только насыщающее добавление со знаком / без знака и насыщающее уплотнение для более узких элементов. (MMX / SSE2 paddsw / paddusw). Поскольку вы печатаете с %d Может быть, вы хотите подписанную насыщенность? Но это только после распаковки в 32-битную версию, и ваша формула всегда будет давать положительный результат, поэтому вы можете использовать беззнаковое насыщение. Я вижу, это то, что ваш код использует paddusw,
Кроме того, использование 3 отдельных циклов, которые сохраняют / перезагружают ваши данные между каждым шагом формулы, действительно ужасно. Вы (почти) всегда хотите увеличить вычислительную интенсивность (объем ALU-работы, выполняемой с вашими данными на каждую полосу пропускания памяти / кэша). Вам также не нужна инструкция умножения, чтобы удвоить число: просто добавьте его к себе. padd* работает на большем количестве портов, чем pmul* и имеет лучшую задержку и (в данном случае, соответствующую) пропускную способность.
default rel
;;; Efficient but *doesn't* saturate the multiply
lea rcx, [rsi + length] ; rcx = end_pointer
movdqa xmm5, [fives]
.loop: ; do{
movdqu xmm0, [rsi] ; use movdqa if your data is aligned, faster on very old CPUs.
pmullw xmm0, xmm0 ; x*x ; does NOT saturate. will wrap.
paddusw xmm0, xmm0 ; *= 2 ; with unsigned saturation
paddusw xmm0, xmm5 ; += 5 ; xmm5 = _mm_set1_epi16(5) outside the loop
movdqu [rsi], xmm0
add rsi, 16
cmp rsi, rcx ; }while(p<endp);
jb .loop
...
section .rodata
align 16
fives: times 8 dw 5
Для насыщения вы можете использовать SSSE3 https://www.felixcloutier.com/x86/pmaddubsw, но для этого нужны только байтовые входы. Насыщает горизонтальную сумму пар i8 x u8 => i16 произведений.
В противном случае вам, вероятно, придется распаковать мечи и packssdw (подписано) или packusdw (без знака насыщения) вернуться к словам. Но умножение мечей происходит медленно с SSE4.1 pmulld (2 мопа на Haswell и позже). Однако на некоторых старых процессорах это всего 1 моп. И, конечно же, распаковка создает в два раза больше работы благодаря наличию более широких элементов.
В этом случае ваша формула монотонна с величиной входа, поэтому вы можете просто сравнить на входе и насытить вручную.
Если мы предположим, что ваши входные данные также не подписаны, нам не нужно вводить абсолютное значение перед сравнением. Но (до AVX512) у нас нет сравнения целых чисел без знака, только подписи больше-больше, поэтому большие входы без знака будут сравниваться как отрицательные.
2 * 0x00b5^2 + 5 = 0xfff7который умещается в 16 бит2 * 0x00b6^2 + 5 = 0x102cdкоторый не делает, и мы хотим, чтобы он был насыщенным0xffff
Точка отсечения переполнения является четным числом, поэтому мы могли бы решить проблему сравнения со знаком, сдвигая вправо. Это было бы беззнаковым делением на 2, что позволило бы безопасно использовать результат как неотрицательное целое число со знаком. 0xb6 >> 1 = 0x5b, Но pcmpgtw это сравнение для > не >=,
Если мы обратим операнды к pcmpgtw сравнивать для (x>>1) < (0xb6>>1) тогда нам придется movdqa константа, чтобы избежать его разрушения, но нам все еще нужно сдвинуть вправо x с movdqa+psrlw. И более эффективно иметь вектор, который 0xffff когда необходимо насыщение (иначе 0), потому что мы можем применить это непосредственно с POR или PADDUSW.
Таким образом, наш лучший выбор - просто сдвинуть диапазон как x, так и 0xb5 подписать и сделать (x-0x8000) > (0xb5 - 0x8000) с использованием pcmpgtw SIMD подписано сравнить.
Другие худшие варианты включают в себя:
- Проверка переполнения при умножении с
pmulhuwрассчитать верхнюю половину (и проверить, если она не равна нулю). Мы были бы в опасности узкого места в умножении пропускной способности и проверки на ненулевое значение сpcmpeqwявляется обратным условию, которое мы хотим. psubusw x, 0xb5и проверьте, что это == 0.pcmpeqwдаст нам перевернутую маску, но мы не можем использоватьpcmpgtwпроверить наusat16(x-0xb5) > 0потому что большие входные данные (с установленным старшим битом) останутся "отрицательными" после вычитания только небольшого числа, например0xb5,padduswи проверить на== 0xffff: только небольшие входные данные оставляют конечный результат ненасыщенным. Некоторые процессоры могут работатьpxorна большем количестве портов, чемpadd*и это не требует меньшего числа ненулевых векторных констант, так что это никак не лучше. Но это одинаково хорошо на Skylake.
default rel
;;; With a check on the input to saturate the output
lea rcx, [rsi + length] ; rcx = end_pointer
movdqa xmm4, [input_saturation_cutoff]
movdqa xmm5, [fives]
pcmpeqd xmm3, xmm3
psllw xmm3, 15 ; xmm3 = set1_epi16(0x8000) for range-shifting to signed
.loop:
movdqu xmm0, [rsi]
movdqa xmm1, xmm0
; if x>0xb5 (unsigned), saturate output to 0xffff
pxor xmm1, xmm3 ; x - 0x8000 = range shift to signed for compare
pcmpgtw xmm1, xmm4 ; xmm1 = (x > 0xb5) ? -1 : 0
pmullw xmm0, xmm0 ; x*x
paddusw xmm0, xmm0 ; *= 2 ; saturating add or not doesn't matter here
por xmm1, xmm5 ; 0xffff (saturation needed) else 5. Off the critical path to give OoO exec an easier time.
paddusw xmm0, xmm1 ; += 5 or += 0xffff with unsigned saturation.
movdqu [rsi], xmm0
add rsi, 16
cmp rsi, rcx
jb .loop
...
section .rodata
align 16
input_saturation_cutoff: times 8 dw (0x00b5 - 0x8000) ; range-shifted to signed for pcmpgtw
fives: times 8 dw 5
; 5 = 0xb6 >> 5 or 0xb6 >> 5 but we'd need to knock off the high bit from input_saturation_cutoff
; Or we could materialize constants like this:
; movdqa xmm4, [set1w_0xb5]
; pcmpeqd xmm3, xmm3
; psllw xmm3, 15 ; rangeshift_mask = set1(0x8000)
; movdqa xmm5, xmm4
; psrlw xmm5, 5 ; fives = set1(5)
; pxor xmm4, xmm3 ; input_sat_cutoff = set1(0x80b5)
;; probably not worth it since we already need to load 1 from memory.
Я проверил это, и paddusw делает 0x2 + 0xffff = 0xffff например.
Мы могли бы просто POR окончательный результат с 0 или 0xffff, чтобы оставить его неизменным или установить 0xffff, но изменив ввод до последнего paddusw создает больше параллелизма на уровне команд за одну итерацию. Таким образом, выполнение вне порядка не должно перекрывать столько независимых итераций, чтобы скрыть задержку тела цикла. (Если бы мы действительно планировали это для Atom или P5 Pentium-MMX по порядку, мы бы чередовали больше из двух цепочек dep.)
На самом деле, смещение вправо на 1 сработало бы: нам нужно только сравнение, чтобы поймать входные данные, настолько большие, что умножение оборачивается до небольшого результата. 0xb6 * 0xb6 не оборачивается, поэтому он сам по себе просто насыщается paddubsw,
Хорошо, если мы проверим (x>>1) > (0xb6>>1) с pcmpgtw (вместо >=) чтобы поймать входы как 0xffff (где pmullw с 0xffff дает нам 0x0001). Таким образом, мы могли бы сохранить 1 векторную константу, но в остальном это не лучше.
pxor + pcmpgtw по крайней мере так же дешево, как psrlw xmm, 1 + pcmpgtw, если только мы не настраиваемся на семейство Intel P6 (Core2/Nehalem) и работаем в стойках регистрации-чтения-порта. Но это маловероятно: xmm0, xmm1 и rsi всегда должны быть горячими (недавно записанными и, таким образом, прочитанными из ROB, а не из файла постоянного реестра). Мы читаем только 2 векторных константы в первой группе из 4 инструкций в цикле, а затем 1 позже.
Как я скажу ниже, на многих процессорах Intel, psrlw может работать только на том же порту, что и pmullw, на vec-int shift+ умножить модуль выполнения. Вероятно, здесь не узкое место в пропускной способности.
Но pcmp а также padd работать на ограниченных портах (на Intel до Skylake), пока pxor может работать на любом векторном порте АЛУ. Смесь чисто padd / pcmp / pmul / psrl` uops оставил бы один векторный порт ALU неиспользованным.
Альтернативная идея проверки насыщенности
(Я написал эту часть, забывая о *2 в формуле, когда искал самый высокий ввод, который не переполняется.)
Если бы формула была (0x00ff)^2 + 5 проверка насыщенности будет проще.
Мы могли бы просто проверить битовые позиции.
(0x00ff)^2 + 5 = 0xfe06который умещается в 16 бит(0x0100)^2 + 5 = 0x10005который не делает, и мы хотим, чтобы он был насыщенным0xffff
Таким образом, мы должны проверить, что все старшие 16 бит равны нулю, или что x&0xFF == x, или это x>>8 == 0,
Это требует меньше констант, но на самом деле это хуже, чем смещение диапазона для подписи с PXOR, потому что на некоторых процессорах исполнительные блоки с векторным смещением и умножением вектора находятся на одном и том же порту. (И поэтому psrlw а также pmullw конкурировать друг с другом за пропускную способность. Всего этого достаточно, чтобы не было узкого места в порту 0 на Nehalem / Sandybridge / Haswell, но это не повредит.)
lea rcx, [rsi + length] ; rcx = end_pointer
movq xmm5, [fives]
punpcklqdq xmm5, xmm5 ; or with SSE3, movddup xmm5, [fives] to broadcast-load
pxor xmm4, xmm4 ; xmm4 = 0
.loop:
movdqu xmm0, [rsi]
movdqa xmm1, xmm0
; if x>0xffU, i.e. if the high byte >0, saturate output to 0xffff
psrlw xmm1, 8 ; x>>8 (logical right shift)
pcmpgtw xmm1, xmm4 ; xmm1 = ((x>>8) > 0) ? -1 : 0
pmullw xmm0, xmm0 ; x*x ; does NOT saturate. will wrap.
paddusw xmm0, xmm0 ; *= 2 ; with unsigned saturation
por xmm1, xmm5 ; 0xffff (saturation needed) or 5 (normal). Off the critical path to give OoO exec an easier time.
paddusw xmm0, xmm1 ; += 5 or += 0xffff with unsigned saturation.
movdqu [rsi], xmm0
add rsi, 16
cmp rsi, rcx
jb .loop
С AVX512BW (Skylake-X) для сравнения без знака с использованием регистров маски
Наконец, мы можем сравнивать целое число без знака с AVX512F, а по размеру элемента word - с AVX512BW. Но результат находится в регистре маски вместо вектора, поэтому мы не можем просто vpor это с вектором set1(5) создать вход для насыщения добавить.
Вместо этого мы можем смешать между вектором 5 а также 0xffff по маске сравнения.
;; AVX512BW version
;; On a Skylake Xeon you probably only want to use YMM regs unless you spend a lot of time in this
;; to avoid reducing max turbo much.
;; Even with xmm or ymm regs (AVX512VL + BW), this demonstrates
;; that we gain even more efficiency than just widening the vectors
;; Just having non-destructive AVX encodings saves the `movdqa xmm1,xmm0` in the SSE2 version.
;; With YMM or XMM regs, most of these instructions can still use shorter VEX encoding (AVX), not the longer EVEX (AVX512)
;; (Use vmovdqa/u instead of vmovdqu64. The 64 is element-size, irrelevant when not masking.)
;;;;;;;;;;; UNTESTED ;;;;;;;;;;;;;;;;;
mov eax, 0xb5 ;; materialize vector constants from an immediate
vpbroadcastd zmm4, eax ; largest input that won't overflow/saturate
vpsrlw zmm5, zmm4, 5 ; fives = 0xb5 >> 5 = 5
;vpcmpeqd xmm3, xmm3 ; all-ones: result on saturation
vpternlogd zmm3,zmm3,zmm3, 0xff ; alternative for ZMM regs, where there's no compare-into-vector
.loop:
; alignment recommended for 512-bit vectors, but `u` lets it still work (slower) on unaligned.
vmovdqu64 zmm0, [rsi]
;; if x>0xb5 (unsigned), saturate output to 0xffff
vpcmpuw k1, zmm0, zmm4, 2 ; k1 = x <= 0xb5. 2 = LE predicate
; k1 set for elements that WON'T overflow
vpmullw zmm0, zmm0 ; x*x
vpaddusw zmm0, zmm0 ; *= 2 ; saturating add or not doesn't matter here
vmovdqa64 zmm1, zmm3 ; 0xffff
vpaddusw zmm1{k1}, zmm0, zmm5 ; 0xffff or 2*x^2 + 5 merge masking
vmovdqu64 [rsi], zmm1
add rsi, 64
cmp rsi, rcx
jb .loop
(NASM позволяет vpmullw a, b как ярлык для vpaddusw a, a, b когда вы не хотите использовать неразрушающее целевое 3-операндное кодирование, как это делается для imul eax, 123.)
Более ранняя идея применения насыщенности заключалась в использовании vpblendmw выбирать между векторами 5 а также 0xffff по маске.
vpcmpuw k1, xmm4, xmm0, 1 ; k1 = 0xb5<x = x>0xb5. 1 = LT predicate numerically because NASM doesn't seem to accept vpcmpltuw the way it accepts vpcmpeqb
; k1 = 1 for elements that WILL overflow.
multiply and add as usual
...
vpblendmw xmm1{k1}, xmm5, xmm3 ; select (k1==0) ? 5 : 0xffff
vpaddusw xmm0, xmm1 ; += 5 or += 0xffff with unsigned saturation.
Копирование регистра по-прежнему занимает входную операцию, но не внутреннюю операцию ALU. Так (особенно для 512-битных регистров, где порт 1 отключается для векторных мопов на SKX), это vpblendmb Идея хуже, чем копировать и слить маску.
Кроме того, МАКА считает, vpblendmw xmm1{k1}, xmm3, xmm5 имеет выходную зависимость от XMM1, хотя на самом деле это только запись. (Я проверил, поместив 8 из этого в цикле, с / без разрыва vpxor). Инструкции смешивания маски представляют собой особый случай: для неустановленных битов маски означает, что он копирует из src1 (или ноль для маскирования нуля), а для установленных битов маски он копирует из src2.
Но машинное кодирование использует маскирование слиянием, поэтому возможно, что HW будет рассматривать его как любую другую операцию ALU с маскированием слиянием. (Где выходной вектор является третьей входной зависимостью, vpaddw xmm1{k1}, xmm2, xmm3: если маска имеет какие-либо нули, результатом в XMM1 будет входное значение этого элемента.)
Вероятно, это не проблема: в соответствии с IACA, SKX может выполнять это за одну итерацию в течение 2,24 циклов (узкое место на внешнем интерфейсе), так что перенос цепи через XMM1 не является проблемой, если только 1 задержка цикла (При развертывании с целью уменьшения узких мест в верхней / внешней части цикла вы должны использовать отдельный вектор для мест назначения смешивания, чтобы разделить итерации, даже если вы не можете получить его где-либо рядом с 1 циклом за итерацию.)
(И версия с использованием копирования + слияния-маскирования в вектор 0xffff также работает с такой пропускной способностью, даже для векторов ZMM. Но IACA считает, что версия vpblendmb будет медленнее с ZMM, даже несмотря на то, что на переднем конце написано, что оба узких места...)
Итак, я нашел решение. Инструкция по использованию pmullw, инструкция pmullw mm0, mm1 рассчитает произведение 4 слов в регистрах по порядку и сохранит их в mm0, И для printf проблема, я просто нажал другой регистр, rdx раньше звонил и теперь работает. Я предполагаю, что это как-то связано с упомянутым смещением стека. Если кто-то может объяснить мне более подробно, как это работает, было бы здорово.