Вычитание упакованных 8-битных целых чисел в 64-битное целое на 1 параллельно, SWAR без аппаратного SIMD

77

Если у меня есть 64-разрядное целое число, которое я интерпретирую как массив упакованных 8-разрядных целых чисел с 8 элементами. Мне нужно вычесть константу 1из каждого упакованного целого числа при обработке переполнения без влияния одного элемента на результат другого элемента.

У меня есть этот код на данный момент, и он работает, но мне нужно решение, которое выполняет вычитание каждого упакованного 8-битного целого числа параллельно и не осуществляет доступ к памяти. На x86 я мог бы использовать SIMD-инструкции, подобные psubbэтим, вычитая упакованные 8-битные целые числа параллельно, но платформа, для которой я кодирую, не поддерживает SIMD-инструкции. (RISC-V в этом случае).

Поэтому я пытаюсь сделать SWAR (SIMD внутри регистра), чтобы вручную отменить перенос переноса между байтами a uint64_t, делая что-то эквивалентное этому:

uint64_t sub(uint64_t arg) {
    uint8_t* packed = (uint8_t*) &arg;

    for (size_t i = 0; i < sizeof(uint64_t); ++i) {
        packed[i] -= 1;
    }

    return arg;
}

Я думаю, что вы могли бы сделать это с побитовыми операторами, но я не уверен. Я ищу решение, которое не использует инструкции SIMD. Я ищу решение на C или C ++, которое достаточно переносимо, или просто теорию, которая стоит за ним, чтобы я мог реализовать свое собственное решение.

Ж-белый
источник
5
Должны ли они быть 8-разрядными или могут быть 7-разрядными?
Тэдман
Они должны быть 8-битными извините :(
cam-white
12
Методы такого рода вещи называются Swar
Harold
1
Вы ожидаете, что байт содержит ноль, чтобы обернуть к 0xff?
Альнитак

Ответы:

75

Если у вас есть процессор с эффективными инструкциями SIMD, SSE / MMX paddb( _mm_add_epi8) также является жизнеспособным. В ответе Питера Кордеса также описывается векторный синтаксис GNU C (gcc / clang) и безопасность для строго псевдонимов UB. Я настоятельно рекомендую также рассмотреть этот ответ.

Делать это самостоятельно с uint64_tполностью переносимым, но все же требует осторожности, чтобы избежать проблем с выравниванием и строгим псевдонимом UB при доступе к uint8_tмассиву с помощью uint64_t*. Вы оставили эту часть вне вопроса, начав с ваших данных uint64_tуже, но для GNU C may_aliastypedef решает проблему (см. Ответ Питера или memcpy).

В противном случае вы можете выделить / объявить ваши данные как uint64_tи получить к ним доступ, uint8_t*когда вам нужны отдельные байты. unsigned char*разрешено создавать псевдонимы, чтобы обойти проблему для конкретного случая 8-битных элементов. (Если uint8_tсуществует вообще, вероятно, можно предположить, что это unsigned char.)


Обратите внимание, что это изменение от предыдущего неправильного алгоритма (см. Историю изменений).

Это возможно без зацикливания для произвольного вычитания и становится более эффективным для известной константы, как 1в каждом байте. Основная хитрость заключается в том, чтобы предотвратить перенос каждого байта путем установки старшего бита, а затем исправить результат вычитания.

Мы собираемся немного оптимизировать технику вычитания, приведенную здесь . Они определяют:

SWAR sub z = x - y
    z = ((x | H) - (y &~H)) ^ ((x ^~y) & H)

с Hопределенным как 0x8080808080808080U(то есть MSB каждого упакованного целого числа). Для декремента yесть 0x0101010101010101U.

Мы знаем, что yвсе его MSB очищены, поэтому мы можем пропустить один из шагов маски (т. Е. Такой y & ~Hже, как yв нашем случае). Расчет происходит следующим образом:

  1. Мы устанавливаем MSB каждого компонента равным x1, чтобы заем не мог распространяться за MSB до следующего компонента. Назовите это настроенным входом.
  2. Мы вычитаем 1 из каждого компонента, вычитая 0x01010101010101из скорректированного ввода. Это не вызывает межкомпонентные заимствования благодаря шагу 1. Назовите это скорректированным выходом.
  3. Теперь нам нужно исправить MSB результата. Мы скорректируем скорректированный вывод с инвертированными старшими битами исходного ввода, чтобы завершить исправление результата.

Операция может быть записана как:

#define U64MASK 0x0101010101010101U
#define MSBON 0x8080808080808080U
uint64_t decEach(uint64_t i){
      return ((i | MSBON) - U64MASK) ^ ((i ^ MSBON) & MSBON);
}

Предпочтительно, это указывается компилятором (используйте директивы компилятора, чтобы принудительно это сделать), или выражение записывается как часть другой функции.

Testcases:

in:  0000000000000000
out: ffffffffffffffff

in:  f200000015000013
out: f1ffffff14ffff12

in:  0000000000000100
out: ffffffffffff00ff

in:  808080807f7f7f7f
out: 7f7f7f7f7e7e7e7e

in:  0101010101010101
out: 0000000000000000

Детали исполнения

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

uint64t[rax] decEach(rcx):
    movabs  rcx, -9187201950435737472
    mov     rdx, rdi
    or      rdx, rcx
    movabs  rax, -72340172838076673
    add     rax, rdx
    and     rdi, rcx
    xor     rdi, rcx
    xor     rax, rdi
    ret

С некоторыми испытаниями IACA следующего фрагмента:

// Repeat the SWAR dec in a loop as a microbenchmark
uint64_t perftest(uint64_t dummyArg){
    uint64_t dummyCounter = 0;
    uint64_t i = 0x74656a6d27080100U; // another dummy value.
    while(i ^ dummyArg) {
        IACA_START
        uint64_t naive = i - U64MASK;
        i = naive + ((i ^ naive ^ U64MASK) & U64MASK);
        dummyCounter++;
    }
    IACA_END
    return dummyCounter;
}

мы можем показать, что на машине Skylake выполнение декремента, xor и сравнения + прыжок может быть выполнено всего за 5 циклов за итерацию:

Throughput Analysis Report
--------------------------
Block Throughput: 4.96 Cycles       Throughput Bottleneck: Backend
Loop Count:  26
Port Binding In Cycles Per Iteration:
--------------------------------------------------------------------------------------------------
|  Port  |   0   -  DV   |   1   |   2   -  D    |   3   -  D    |   4   |   5   |   6   |   7   |
--------------------------------------------------------------------------------------------------
| Cycles |  1.5     0.0  |  1.5  |  0.0     0.0  |  0.0     0.0  |  0.0  |  1.5  |  1.5  |  0.0  |
--------------------------------------------------------------------------------------------------

(Конечно, на x86-64 вы просто загрузите или movqв регистр XMM paddb, так что может быть интереснее посмотреть, как он компилируется для ISA, такого как RISC-V.)

нФ
источник
4
Мне нужен мой код для запуска на машинах RISC-V, у которых нет инструкций SIMD (пока), не говоря уже о поддержке MMX
cam-white
2
@ cam-white Понял - это, наверное, лучшее, что ты можешь сделать. Я прыгну на Годболт, чтобы проверить работоспособность сборки RISC. Редактировать: Нет поддержки RISC-V на Godbolt :(
нанофарад
7
На самом деле есть поддержка RISC-V для Godbolt, например, вот так (E: кажется, что компилятор слишком изобретателен при создании маски ...)
Гарольд
4
Дальнейшее чтение о том, как можно использовать трюк контроля четности (также называемый вектором выполнения) в различных ситуациях: emulators.com/docs/LazyOverflowDetect_Final.pdf
jpa
4
Я сделал еще одно редактирование; Собственные векторы GNU C фактически избегают проблем со строгим псевдонимом; вектор uint8_tдопускается для псевдонимов uint8_tданных. Вызывающие вашу функцию (которые должны получить uint8_tданные в uint64_t) должны беспокоиться о строгом псевдониме! Поэтому, вероятно, OP должен просто объявить / распределить массивы, uint64_tпотому что char*ему разрешено псевдоним в ISO C ++, но не наоборот.
Питер Кордес
16

Для RISC-V вы, вероятно, используете GCC / clang.

Интересный факт: GCC знает некоторые из этих хитростей SWAR-трюков (показанных в других ответах) и может использовать их для вас при компиляции кода с собственными векторами GNU C для целей без аппаратных инструкций SIMD. (Но clang для RISC-V просто наивно развернет его для скалярных операций, поэтому вам придется делать это самостоятельно, если вы хотите добиться хорошей производительности на всех компиляторах).

Одно из преимуществ нативного векторного синтаксиса заключается в том, что при нацеливании на компьютер с аппаратным SIMD он будет использовать его вместо автоматической векторизации вашего битхака или чего-то ужасного в этом роде.

Это облегчает написание vector -= scalarопераций; синтаксис Just Works, неявно транслирующий для вас скаляр.


Также обратите внимание, что uint64_t*загрузка из uint8_t array[]UB строго псевдонимов, так что будьте осторожны с этим. (См. Также Почему strlen glibc должен быть настолько сложным, чтобы быстро запускаться? Re: обеспечение безопасности строгих псевдонимов SWAR в чистом C). Возможно, вы захотите, чтобы что-то вроде этого объявляло, uint64_tчто вы можете приводить указатели для доступа к любым другим объектам, например, как это char*работает в ISO C / C ++.

используйте их, чтобы получить данные uint8_t в uint64_t для использования с другими ответами:

// GNU C: gcc/clang/ICC but not MSVC
typedef uint64_t  aliasing_u64 __attribute__((may_alias));  // still requires alignment
typedef uint64_t  aliasing_unaligned_u64 __attribute__((may_alias, aligned(1)));

Другой способ выполнить безопасные с точки зрения псевдонимов нагрузки заключается memcpyв использовании параметра uint64_t, который также устраняет alignof(uint64_tтребование выравнивания. Но на ISA без эффективных невыровненных нагрузок gcc / clang не встроен и не оптимизируется, memcpyкогда не может доказать, что указатель выровнен, что может иметь катастрофические последствия для производительности.

TL: DR: вам лучше всего объявить ваши данные какuint64_t array[...] или выделить их динамически как uint64_t, или, что желательно,alignas(16) uint64_t array[]; что обеспечивает выравнивание по крайней мере до 8 байтов, или 16, если вы укажете alignas.

Поскольку uint8_tэто почти наверняка unsigned char*, это безопасный доступ к байтам uint64_tvia uint8_t*(но не наоборот для массива uint8_t). Таким образом, для этого особого случая, когда тип элемента узкий unsigned char, вы можете обойти проблему строгого псевдонима, потому что charона особенная.


Пример синтаксиса собственного вектора GNU C:

GNU C родные векторы всегда может псевдоним с их базовым типом (например , int __attribute__((vector_size(16)))может безопасно псевдоним , intно не floatили uint8_tили что - нибудь еще.

#include <stdint.h>
#include <stddef.h>

// assumes array is 16-byte aligned
void dec_mem_gnu(uint8_t *array) {
    typedef uint8_t v16u8 __attribute__ ((vector_size (16), may_alias));
    v16u8 *vecs = (v16u8*) array;
    vecs[0] -= 1;
    vecs[1] -= 1;   // can be done in a loop.
}

Для RISC-V без какого-либо HW SIMD вы могли бы использовать vector_size(8)для выражения только степень детализации, которую вы можете эффективно использовать, и сделать в два раза больше меньших векторов.

Но vector_size(8)для x86 очень тупо компилируется как с GCC, так и clang: GCC использует битовые хаки SWAR в целочисленных регистрах GP, clang распаковывает в 2-байтовые элементы для заполнения 16-байтового регистра XMM, а затем перепаковывает. (MMX настолько устарел, что GCC / clang даже не потрудился использовать его, по крайней мере, для x86-64.)

Но с помощью vector_size (16)( Godbolt ) мы получаем ожидаемое movdqa/ paddb. (С вектором «все единицы», сгенерированным pcmpeqd same,same). При этом -march=skylakeмы по-прежнему получаем две отдельные операции XMM вместо одной YMM, поэтому, к сожалению, современные компиляторы также не «автоматически векторизуют» векторные операции в более широкие векторы: /

Для AArch64 это не так уж плохо в использовании vector_size(8)( Godbolt ); ARM / AArch64 может работать в 8- или 16-байтовых чанках с регистрами dили q.

Так что вы, вероятно, захотите vector_size(16)на самом деле скомпилировать, если вам нужна портативная производительность для x86, RISC-V, ARM / AArch64 и POWER . Однако некоторые другие ISA выполняют SIMD в 64-битных целочисленных регистрах, например, MIPS MSA, я думаю.

vector_size(8)облегчает просмотр asm (только один регистр данных): проводник компилятора Godbolt

# GCC8.2 -O3 for RISC-V for vector_size(8) and only one vector

dec_mem_gnu(unsigned char*):
        lui     a4,%hi(.LC1)           # generate address for static constants.
        ld      a5,0(a0)                 # a5 = load from function arg
        ld      a3,%lo(.LC1)(a4)       # a3 = 0x7F7F7F7F7F7F7F7F
        lui     a2,%hi(.LC0)
        ld      a2,%lo(.LC0)(a2)       # a2 = 0x8080808080808080
                             # above here can be hoisted out of loops
        not     a4,a5                  # nx = ~x
        and     a5,a5,a3               # x &= 0x7f... clear high bit
        and     a4,a4,a2               # nx = (~x) & 0x80... inverse high bit isolated
        add     a5,a5,a3               # x += 0x7f...   (128-1)
        xor     a5,a4,a5               # x ^= nx  restore high bit or something.

        sd      a5,0(a0)               # store the result
        ret

Я думаю, что это та же самая основная идея, что и другие нецикличные ответы; предотвращение переноса и исправление результата.

Это 5 инструкций ALU, хуже, чем главный ответ, я думаю. Но похоже, что задержка критического пути составляет всего 3 цикла, с двумя цепочками по 2 инструкции, каждая из которых ведет к XOR. Ответ @Reinstate Monica - ζ - компилируется в 4-тактную цепочку dep (для x86). Пропускная способность цикла с 5 циклами является узким местом, в том числе путем наивности subна критическом пути, а цикл создает узкое место с задержкой.

Тем не менее, это бесполезно с лязгом. Он даже не добавляет и не хранит в том же порядке, в котором загружен, поэтому он даже не выполняет хорошую программную конвейеризацию!

# RISC-V clang (trunk) -O3
dec_mem_gnu(unsigned char*):
        lb      a6, 7(a0)
        lb      a7, 6(a0)
        lb      t0, 5(a0)
...
        addi    t1, a5, -1
        addi    t2, a1, -1
        addi    t3, a2, -1
...
        sb      a2, 7(a0)
        sb      a1, 6(a0)
        sb      a5, 5(a0)
...
        ret
Питер Кордес
источник
13

Я хотел бы отметить, что код, который вы написали, на самом деле векторизируется, когда вы начинаете работать с более чем одним uint64_t.

https://godbolt.org/z/J9DRzd

robthebloke
источник
1
Не могли бы вы объяснить или дать ссылку на то, что там происходит? Это кажется довольно интересным.
n314159
2
Я пытался сделать это без SIMD инструкций, но тем не менее я нашел это интересным :)
cam-white
8
С другой стороны, этот код SIMD ужасен. Компилятор совершенно не понял, что здесь происходит. E: это пример "это было явно сделано компилятором, потому что ни один человек не был бы таким глупым"
Гарольд
1
@PeterCordes: я больше думал о __vector_loop(index, start, past, pad)конструкции, которую реализация может трактовать как for(index=start; index<past; index++)[то есть любая реализация может обрабатывать код, используя ее, просто путем определения макроса], но которая будет иметь более слабую семантику, чтобы пригласить компилятор обрабатывать вещи в любой размер фрагмента степени два pad, расширяющий начало вниз и конец вверх, если они еще не кратны размеру фрагмента. Побочные эффекты в каждом чанке не были бы последовательными, и если breakв цикле возникает a , другие представители ...
суперкат
1
@PeterCordes: Хотя restrictэто полезно (и было бы более полезно, если бы Стандарт признал концепцию «по крайней мере, потенциально основанную на», а затем определил «на основе» и «по крайней мере, потенциально основанную на» напрямую, без глупых и неработающих угловых случаев) Мое предложение также позволило бы компилятору выполнять больше циклов, чем запрошено, что значительно упростит векторизацию, но для этого в стандарте не предусмотрено никаких условий.
суперкат
11

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

uint64_t sub(uint64_t arg) {
    uint64_t x1 = arg | 0x80808080808080;
    uint64_t x2 = ~arg & 0x80808080808080;
    // or uint64_t x2 = arg ^ x1; to save one instruction if you don't have an andnot instruction
    return (x1 - 0x101010101010101) ^ x2;
}
Фальк Хюффнер
источник
Я думаю, что это работает для всех 256 возможных значений байта; Я поместил его на Godbolt (с помощью RISC-V clang) godbolt.org/z/DGL9aq, чтобы посмотреть результаты постоянного распространения для различных входных данных, таких как 0x0, 0x7f, 0x80 и 0xff (смещенных в середину числа). Выглядит хорошо. Я думаю, что главный ответ сводится к тому же, но объясняет это более сложным способом.
Питер Кордес
Компиляторы могли бы лучше создавать константы в регистрах здесь. clang тратит много инструкций на построение splat(0x01)и splat(0x80)вместо того, чтобы получать одну от другой за смену. Даже если написать так в исходном коде, godbolt.org/z/6y9v-u не помогает компилятору создавать лучший код; это просто делает постоянное распространение.
Питер Кордес
Интересно, почему он не просто загружает константу из памяти; это то, что делают компиляторы для Alpha (похожая архитектура).
Фальк Хюффнер
НКУ для RISC-V делает постоянные нагрузки из памяти. Похоже, что clang нуждается в некоторой настройке, если только не ожидаются ошибки в кеше данных и они дороги по сравнению с пропускной способностью команд. (Этот баланс, безусловно, мог измениться со времен Alpha, и предположительно разные реализации RISC-V различны. Компиляторы также могли бы добиться большего успеха, если бы поняли, что это повторяющийся шаблон, который они могли бы сместить / ИЛИ расширить после запуска с одного LUI / добавления. для 20 + 12 = 32 бита непосредственных данных. Битовые последовательности AArch64 могут даже использовать их в качестве непосредственных для AND / OR / XOR, интеллектуального декодирования и выбора плотности)
Питер Кордес
Добавлен ответ, показывающий исходный вектор GAR SWAR для RISC-V
Питер Кордес
7

Не уверен, что это то, что вам нужно, но он выполняет 8 вычитаний параллельно друг другу:

#include <cstdint>

constexpr uint64_t mask = 0x0101010101010101;

uint64_t sub(uint64_t arg) {
    uint64_t mask_cp = mask;
    for(auto i = 0; i < 8 && mask_cp; ++i) {
        uint64_t new_mask = (arg & mask_cp) ^ mask_cp;
        arg = arg ^ mask_cp;
        mask_cp = new_mask << 1;
    }
    return arg;
}

Объяснение: Битовая маска начинается с 1 в каждом из 8-битных чисел. Мы исправим это с помощью нашего аргумента. Если у нас была 1 в этом месте, мы вычли 1 и должны остановиться. Это делается путем установки соответствующего бита в 0 в new_mask. Если у нас был 0, мы устанавливаем его в 1 и должны выполнять перенос, поэтому бит остается на 1, и мы смещаем маску влево. Вы лучше сами убедитесь, что генерация новой маски работает, как задумано, я так думаю, но второе мнение не будет плохим.

PS: Я на самом деле не уверен, что проверка на mask_cpненулевое значение в цикле может замедлить работу программы. Без него код по-прежнему был бы корректным (поскольку маска 0 просто ничего не делает), и компилятору было бы намного проще выполнить развертывание цикла.

n314159
источник
forне будет работать параллельно, вас смущает for_each?
LTPCGO
3
@LTPCGO Нет, я не собираюсь распараллеливать это для цикла for, это фактически нарушит алгоритм. Но этот код работает с различными 8-битными целыми числами в 64-битном целом параллельно, то есть все 8 вычитаний выполняются одновременно, но им нужно до 8 шагов.
n314159
Я понимаю, что то, что я спрашивал, могло быть немного неразумно, но это было довольно близко к тому, что мне было нужно, спасибо :)
cam-white
4
int subtractone(int x) 
{
    int f = 1; 

    // Flip all the set bits until we find a 1 at position y
    while (!(x & f)) { 
        x = x^f; 
        f <<= 1; 
    } 

    return x^f; // return answer but remember to flip the 1 at y
} 

Вы можете сделать это с помощью побитовых операций, используя вышеизложенное, и вам просто нужно разделить ваше целое число на 8 битных частей, чтобы отправить 8 раз в эту функцию. Следующая часть была взята из раздела Как разбить 64-битное число на восемь 8-битных значений? со мной добавив в вышеупомянутой функции

uint64_t v= _64bitVariable;
uint8_t i=0,parts[8]={0};
do parts[i++] = subtractone(v&0xFF); while (v>>=8);

Это действительно C или C ++ независимо от того, как кто-то сталкивается с этим

LTPCGO
источник
5
Это не распараллеливает работу, хотя, это вопрос OP.
nickelpro
Да, @nickelpro прав, это будет делать каждое вычитание одно за другим, я хотел бы вычесть все 8-битные целые числа одновременно. Я ценю ответ Тхо спасибо братан
кулачковый-белый
2
@nickelpro, когда я начал отвечать, редактирование еще не было выполнено, в котором была указана параллельная часть вопроса, и поэтому я не заметил этого до тех пор, пока после отправки не уйду, если это будет полезно для других, так как, по крайней мере, отвечает на часть для выполнения побитовых операций, и его можно for_each(std::execution::par_unseq,...
заставить
2
Это плохо, я отправил вопрос, потом понял, что я не говорил, что это нужно для параллельной обработки, так что отредактировано
cam-white
2

Не буду пытаться придумать код, но для уменьшения на 1 вы можете уменьшить на группу 8 1, а затем проверить, чтобы убедиться, что младшие биты результатов «перевернулись». Любой LSB, который не переключался, указывает, что перенос произошел из соседних 8 битов. Должна быть возможность разработать последовательность операций AND / ORs / XOR, чтобы справиться с этим, без каких-либо ветвей.

Горячие лижет
источник
Это может сработать, но рассмотрим случай, когда перенос распространяется через одну группу из 8 битов в другую. Стратегия в хороших ответах (сначала установить MSB или что-то в этом роде), чтобы гарантировать, что перенос не распространяется, вероятно, по крайней мере настолько эффективен, насколько это возможно. Текущий целевой показатель (то есть хорошие непериодические ответы без ответвлений) - это 5 инструкций RISC-V asm ALU с параллелизмом на уровне инструкций, делающим критический путь всего 3 цикла и использующим две 64-битные константы.
Питер Кордес
0

Сосредоточьте работу на каждом байте в одиночку, затем верните его туда, где он был.

uint64_t sub(uint64_t arg) {
   uint64_t res = 0;

   for (int i = 0; i < 64; i+=8) 
     res += ((arg >> i) - 1 & 0xFFU) << i;

    return res;
   }
nonock
источник