Быстро Евклидово подразделение в C

Редактировать: Нет умножения или ветвления.

int euc(int i, int n)
{
    int r;

    r = i % n;
    r += n & (-(r < 0));

    return r;
}

Вот сгенерированный код. Согласно инструментальному профилировщику MSVC ++ (мое тестирование) и тестированию OP, они работают почти одинаково.

; Original post
00401000  cdq              
00401001  idiv        eax,ecx 
00401003  mov         eax,edx 
00401005  test        eax,eax 
00401007  jge         euc+0Bh (40100Bh) 
00401009  add         eax,ecx 
0040100B  ret              

; Mine
00401020  cdq              
00401021  idiv        eax,ecx 
00401023  xor         eax,eax 
00401025  test        edx,edx 
00401027  setl        al   
0040102A  neg         eax  
0040102C  and         eax,ecx 
0040102E  add         eax,edx 
00401030  ret              

Без ветки, но с небольшими битами:

int euc2(int i, int n)
{
    int r;
    r = i % n;
    r += (((unsigned int)r) >> 31) * n;
    return r;
}

Без умножения:

int euc2(int i, int n)
{
    int r;
    r = i % n;
    r += (r >> 31) & n;
    return r;
}

Это дает:

; _i$ = eax
; _n$ = ecx

cdq
idiv   ecx
mov eax, edx
sar eax, 31
and eax, ecx
add eax, edx

Целочисленное умножение намного быстрее деления. Для большого количества вызовов с известным N вы можете заменить деление на N умножением на псевдообратное значение N.

Я проиллюстрирую это на примере. Возьмем N = 29. Затем вычислите один раз псевдообратное значение 2 ^ 16 / N: K = 2259 (усечено с 2259,86 ...). Я предполагаю, что I положительно и I * K подходит для 32 бита.

Quo = (I*K)>>16;   // replaces the division, Quo <= I/N
Mod = I - Quo*N;   // Mod >= I%N
while (Mod >= N) Mod -= N;  // compensate for the approximation

В моем примере возьмем I = 753, мы получим Quo = 25 и Mod = 28. (компенсация не требуется)

EDIT.

В вашем примере трехмерной свертки большинство вызовов i% n будет с i в 0..n-1, поэтому в большинстве случаев первая строка вроде

if (i>=0 && i<n) return i;

позволит обойти дорогостоящий и здесь бесполезный idiv.

Кроме того, если у вас достаточно RAM, просто выровняйте все измерения до степени 2 и используйте битовые манипуляции (сдвиг и) вместо делений.

EDIT 2.

] Я действительно пробовал это на 10 ^ 9 звонках. i% n: 2,93 с, мой код: 1.38с. Просто имейте в виду, что это подразумевает ограничение на I (I * K должно соответствовать 32 битам).

Еще одна мысль: если ваши значения равны x + dx, с x в 0..n-1 и dx маленьким, тогда следующее будет охватывать все случаи:

if (i<0) return i+n; else if (i>=n) return i-n;
return i;

Я думаю, что 280Z28 и Кристофер справились с ассемблерным гольфом лучше, чем я, и это касается произвольного доступа.

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

В этом случае сначала проверьте подходящие границы. вы можете выполнить внутренний цикл в том, что я назову «тире». Убедитесь, что следующие k приращений не t приводит к переполнению в наименьшем измерении в любом массиве, а затем «разбивает» k шагов с использованием нового внутреннего цикла, который просто увеличивает «физический» индекс на 1 каждый раз вместо того, чтобы выполнять еще один idiv. Вы или компилятор можете развернуть этот цикл, использовать устройство Даффа и т. Д.

Если ядро ​​маленькое, и особенно если оно имеет фиксированный размер, то это (или несколько его с подходящим развертыванием, чтобы время от времени вычитать вместо добавления) , вероятно, является значением, которое следует использовать для длины «тире». Вероятно, лучше всего подойдет постоянная длина тире времени компиляции, поскольку тогда вы (или компилятор) можете полностью развернуть цикл тире и исключить условие продолжения. Пока это не делает код слишком большим, чтобы быть быстрым, он по существу заменяет всю операцию положительного модуля целочисленным приращением.

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

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

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

int euc(int i, int n)
{
    return (i % n) + (((i % n) < 0) * n);
}

Мне очень нравится выражение:

r = ((i%n)+n)%n; 

Разборка очень короткая:

r = ((i% n) + n)% n;

004135AC  mov         eax,dword ptr [i] 
004135AF  cdq              
004135B0  idiv        eax,dword ptr [n] 
004135B3  add         edx,dword ptr [n] 
004135B6  mov         eax,edx 
004135B8  cdq              
004135B9  idiv        eax,dword ptr [n] 
004135BC  mov         dword ptr [r],edx 

В ней нет скачков ( 2 idivs, что может быть дорогостоящим), и он может быть полностью встроен, избегая накладных расходов на вызов функции.

Что вы думаете?

Я рассчитал время для всех предложений в gcc -O3, используя TSC (за исключением предложения для постоянного N), и все они заняли одинаковое количество времени (в пределах 1%).

Моя мысль была такой. что либо ((i% n) + n)% n (без ветвления), либо (i + (n << 16))% n (очевидно, не работает для больших n или крайне отрицательных i) будет быстрее, но все они в то же время.

Если у вас достаточно низкий диапазон, создайте таблицу поиска - два тусклых массива. Также вы можете сделать функцию встроенной и убедиться в этом, посмотрев созданный код.

Если вы также можете гарантировать, что i никогда не меньше -n, ​​вы можете просто поставить необязательное сложение перед модулем. Таким образом, вам не нужна ветвь, и модуль modulo вырезает то, что вы добавили, если вам это не нужно.

int euc(int i, int n)
{
    return (i + n) % n;
}

Если i меньше -n, ​​вы все равно можете использовать этот метод. В таком случае вы, вероятно, точно знаете, в каком диапазоне будут находиться ваши значения. Поэтому вместо добавления n к i вы можете добавить x * n к i, где x - любое целое число, которое дает вам достаточный диапазон. Для увеличения скорости (на процессорах, которые не поддерживают однократное умножение) вы можете использовать левый сдвиг вместо умножения.

В своем ответе Эрику Бейнвиллу вы говорите, что в большинстве случаев 0 <= i и что у вас

if (i>=0 && i<n) return i;

int euc(int i, int n)
{
    if (n <= i)            return i % n;
    else if (i < 0)        return ((i + 1) % n) + n - 1;
    else /* 0 <= i < n */  return i;  // fastest possible response for common case
}