Использование SIMD для нахождения наибольшей разницы двух элементов

Я написал алгоритм для получения наибольшей разницы между двумя элементами в std::vector, где большее из двух значений должно иметь более высокий индекс, чем меньшее значение.

    unsigned short int min = input.front();
    unsigned short res = 0;

    for (size_t i = 1; i < input.size(); ++i)
    {
        if (input[i] <= min)
        {
            min = input[i];
            continue;
        }

        int dif = input[i] - min;
        res = dif > res ? dif : res;
    }

    return res != 0 ? res : -1;

Можно ли оптимизировать этот алгоритм с помощью SIMD? Я новичок в SIMD, и пока у меня не получилось с этим.


person Yum    schedule 24.09.2018    source источник


Ответы (1)


Вы не указали какую-либо конкретную архитектуру, поэтому я буду держать эту архитектуру в основном нейтральной с алгоритмом, описанным на английском языке. Но для этого требуется SIMD ISA, который может эффективно переходить к результатам сравнения SIMD для проверки обычно истинного условия, например, x86, но не ARM NEON.

Это не будет хорошо работать для NEON, потому что у него нет эквивалента маски перемещения, а SIMD -> integer вызывает зависание на многих микроархитектурах ARM.


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


Для каждого вектора v[0..7] входного массива (при условии, что 8 int16_t элементов на вектор (16 байтов), но это произвольно):

  • SIMD сравнивает vmin > v[0..7] и проверяет, верны ли какие-либо элементы. (например, x86 _mm_cmpgt_epi16 / if(_mm_movemask_epi8(cmp) != 0)) Если где-то есть новый min, у нас есть особый случай: старый min применяется к некоторым элементам, а новый min применяется к другим. И возможно, в векторе есть несколько обновлений new-min и кандидаты new-diff в любой из этих точек.

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

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

    В качестве альтернативы, SIMD-префикс-сумма (на самом деле префикс-мин) может создать vmin, где каждый элемент является минимумом до этого момента. (параллельная префиксная (кумулятивная) сумма с SSE). Вы можете всегда делать это, чтобы избежать каких-либо ветвлений, но если новые мини-кандидаты редки, то это дорого. Тем не менее, это может быть жизнеспособно на ARM NEON, где ветвление затруднено.

  • Если нового мин. нет, максимум diffmax[0..7] = max(diffmax[0..7], v[0..7]-vmin) SIMD упаковано. (Используйте вычитание с насыщением, чтобы не получить большую разницу без знака, если вы используете максимальное значение без знака для обработки всего диапазона.)

В конце цикла выполните горизонтальный максимум SIMD для вектора diffmax. Обратите внимание: поскольку нам не нужна позиция максимальной разницы, нам не нужно обновлять все элементы внутри цикла при обнаружении нового кандидата. Нам даже не нужно синхронизировать скалярный специальный случай diffmax и SIMD vdiffmax друг с другом, просто проверьте в конце, чтобы получить максимальное значение различий между скаляром и SIMD max.


SIMD min/max в основном такой же, как и горизонтальная сумма, за исключением того, что вы используете pack-max вместо packed-add. Для x86 см. самый быстрый способ сделать горизонтальная сумма векторов с плавающей запятой на x86.

Или на x86 с SSE4.1 для 16-битных целочисленных элементов phminposuw/_mm_minpos_epu16 можно использовать для минимального или максимального значения, со знаком или без знака, с соответствующими настройками ввода. max = -min(-diffmax). Вы можете рассматривать diffmax как неподписанный, поскольку известно, что он неотрицательный, но Горизонтальный минимум и максимум с использованием SSE показывает, как перевернуть знаковый бит на диапазон-сдвинуть знаковый на беззнаковый и обратно.


Вероятно, мы получаем неверное предсказание ветвления каждый раз, когда находим нового min кандидата, или же мы находим новых min кандидатов слишком часто, чтобы это было эффективно.

Если новые кандидаты min ожидаются часто, использование более коротких векторов может быть полезным. Или, обнаружив, что в текущем векторе есть новый-min, используйте более узкие векторы, чтобы использовать скаляр только для меньшего числа элементов. В x86 вы можете использовать bsf (битовое сканирование вперед), чтобы найти, у какого элемента был первый новый мин. Это дает вашему скалярному коду зависимость данных от векторной маски сравнения, но если ветвь к ней была неверно предсказана, маска сравнения будет готова. В противном случае, если предсказание ветвления может каким-то образом найти шаблон, в котором векторам нужен скалярный резерв, предсказание + спекулятивное выполнение нарушит эту зависимость данных.


Незавершенный/сломанный (мной) пример, адаптированный из удаленного ответа @harold полностью автономной версии, которая на лету строит вектор от минимума до этого элемента для x86 SSE2.

(@harold написал это с суффиксом-max вместо min, поэтому я думаю, почему он удалил его. Я частично преобразовал его из max в min.)

Версия встроенных функций без ответвлений для x86 может выглядеть примерно так. Но разветвленный, вероятно, лучше, если вы не ожидаете какого-то наклона или тренда, который делает частыми новые значения min.

// BROKEN, see FIXME comments.
// converted from @harold's suffix-max version

int broken_unfinished_maxDiffSSE(const std::vector<uint16_t> &input) {
    const uint16_t *ptr = input.data();

    // construct suffix-min
    // find max-diff at the same time
    __m128i min = _mm_set_epi32(-1);
    __m128i maxdiff = _mm_setzero_si128();

    size_t i = input.size();
    for (; i >= 8; i -= 8) {
        __m128i data = _mm_loadu_si128((const __m128i*)(ptr + i - 8));

   // FIXME: need to shift in 0xFFFF, not 0, for min.
   // or keep the old data, maybe with _mm_alignr_epi8
        __m128i d = data;
        // link with suffix
        d = _mm_min_epu16(d, _mm_slli_si128(max, 14));
        // do suffix-min within block.
        d = _mm_min_epu16(d, _mm_srli_si128(d, 2));
        d = _mm_min_epu16(d, _mm_shuffle_epi32(d, 0xFA));
        d = _mm_min_epu16(d, _mm_shuffle_epi32(d, 0xEE));
        max = d;

        // update max-diff
        __m128i diff = _mm_subs_epu16(data, min);  // with saturation to 0
        maxdiff = _mm_max_epu16(maxdiff, diff);
    }

    // horizontal max
    maxdiff = _mm_max_epu16(maxdiff, _mm_srli_si128(maxdiff, 2));
    maxdiff = _mm_max_epu16(maxdiff, _mm_shuffle_epi32(maxdiff, 0xFA));
    maxdiff = _mm_max_epu16(maxdiff, _mm_shuffle_epi32(maxdiff, 0xEE));
    int res = _mm_cvtsi128_si32(maxdiff) & 0xFFFF;

    unsigned scalarmin = _mm_extract_epi16(min, 7);  // last element of last vector
    for (; i != 0; i--) {
        scalarmin = std::min(scalarmin, ptr[i - 1]);
        res = std::max(res, ptr[i - 1] - scalarmin);
    }

    return res != 0 ? res : -1;
}

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

person Peter Cordes    schedule 25.09.2018
comment
Таким образом, я сделал эту ссылку с суффиксным шагом, который был прерван переключением на min, в любом случае был плохим подходом, действительно нужно было транслировать, а затем помещать его после в блоке, чтобы цикл переносил зависимость Быстрее - person harold; 25.09.2018
comment
@harold: о, хорошая мысль! Это даже не является узким местом при тасовке, а скорее в цепочке тасовки-›max-›... dep. Я хотел включить какой-то код в качестве примера того, как выглядят встроенные функции, чтобы OP мог больше узнать о Google. Так что я оставлю этот не очень хороший блок кода с предупредительными комментариями. Если вы улучшите / исправите его, не стесняйтесь редактировать мой ответ или лучше восстановить свой собственный. - person Peter Cordes; 25.09.2018