Проверка параметров умножения на константу в 64 бит

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

gcc реализует беззнаковый 64-битный divmod с константой 0x28f5c28f5c28f5c3. Я не рассматривал подробно генерирование констант для деления, но есть алгоритмы их генерирования, которые дают заведомо хорошие результаты (поэтому исчерпывающее тестирование не требуется).

Код на самом деле имеет несколько важных отличий: он использует константу иначе, чем константа OP.

См. комментарии для анализа того, что это на самом деле делает: сначала разделите на 4, чтобы он мог использовать константу, которая работает для деления на 25 только тогда, когда делимое достаточно мало. Это также позволяет избежать необходимости добавления позже.

#include <stdint.h>

// rem, quot ordering takes one extra instruction
struct divmod { uint64_t quotient, remainder; }
 div_by_100(uint64_t x) {
    struct divmod retval = { x%100, x/100 };
    return retval;
}

компилируется в (gcc 5.3 -O3 -mtune=haswell):

    movabs  rdx, 2951479051793528259
    mov     rax, rdi            ; Function arg starts in RDI (SysV ABI)
    shr     rax, 2
    mul     rdx
    shr     rdx, 2
    lea     rax, [rdx+rdx*4]    ; multiply by 5
    lea     rax, [rax+rax*4]    ; multiply by another 5
    sal     rax, 2              ; imul rax, rdx, 100 is better here (Intel SnB).
    sub     rdi, rax
    mov     rax, rdi
    ret
; return values in rdx:rax

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

Часть умножения на 100.

gcc использует указанную выше последовательность lea/lea/shl, как и в вашем вопросе. В вашем ответе используется последовательность mov imm/mul.

Каждый из ваших комментариев говорит, что версия, которую они выбрали, быстрее. Если да, то это из-за тонкого выравнивания инструкций или другого вторичного эффекта: в Intel SnB-семействе это одинаковое количество мопов (3) и такая же задержка критического пути (mov imm находится за пределами критического пути, а mul составляет 3 цикла).

clang использует то, что я считаю лучшим выбором (imul rax, rdx, 100). Я подумал об этом до того, как увидел, что это выбрал лязг, но это не имеет значения. Это 1 uop слитого домена (который может выполняться только на p0), все еще с задержкой 3c. Поэтому, если вы ограничены задержкой, используя эту процедуру для множественной точности, это, вероятно, не поможет, но это лучший выбор. (Если вы привязаны к задержке, встраивание вашего кода в цикл вместо передачи одного из параметров через память может сэкономить много циклов.)

imul работает, потому что вы только с использованием младших 64b результата. Для mul не существует 2-х или 3-х операндной формы, потому что младшая половина результата одинакова независимо от подписанной или беззнаковой интерпретации входных данных.

Кстати, clang с -march=native использует mulx для 64x64->128 вместо mul, но ничего от этого не выигрывает. Согласно таблицам Агнера Фога, задержка на один цикл меньше, чем mul.

AMD имеет задержку хуже, чем 3c для imul r,r,i (особенно версия 64b), возможно, поэтому gcc ее избегает. IDK, сколько труда сопровождающие gcc вложили в настройку затрат, чтобы такие настройки, как -mtune=haswell, работали хорошо, но много кода не скомпилировано ни с одним параметром -mtune (даже с тем, который подразумевается -march), так что я не удивлен, когда gcc делает выбор, который был оптимален для старых процессоров или для AMD.

clang по-прежнему использует imul r64, r64, imm с -mtune=bdver1 (Bulldozer), что экономит m-ops, но требует задержки на 1c больше, чем при использовании lea/lea/shl. (значение со шкалой>1 соответствует задержке 2c на Bulldozer).

Я нашел решение с помощью libdivide.h. Вот немного более сложная часть для Win64:

{$ELSE WIN64}
asm
        .NOFRAME

        MOV     RAX,[RCX]
        MOV     R8,RAX
        XOR     RDX,RDX
        MOV     R9,Div100Const       // New: $47AE147AE147AE15
        MUL     R9                   // Preliminary result Q in RDX

        // Additional part: add/shift

        ADD     RDX,R8               // Q := Q + X shr 1;
        RCR     RDX,1

        SHR     RDX,Div100PostShift  // Q := Q shr 6;
        MOV     [RCX],RDX            // X := Q;

        // Faster than LEA and SHL

        MOV     RAX,RDX
        MOV     R9D,100
        MUL     R9
        SUB     R8,RAX
        MOV     RAX,R8         // Remainder
end;
{$ENDIF WIN32}

Код в ответе @Rudy является результатом следующих шагов:

1. Запишите 1/100 в двоичной форме: 0.000000(10100011110101110000);
2. Подсчет ведущих нулей после запятой: S = 6;
3. 72 первых значащих бита:

1010 0011 1101 0111 0000 1010 0011 1101 0111 0000 1010 0011 1101 0111 0000 1010 0011 1101

4. Округлить до 65 бит; есть какая-то магия в том, как выполняется это округление; путем обратного проектирования константы из ответа Руди правильное округление:

1010 0011 1101 0111 0000 1010 0011 1101 0111 0000 1010 0011 1101 0111 0000 1010 1

5. Удалите начальный бит 1:

0100 0111 1010 1110 0001 0100 0111 1010 1110 0001 0100 0111 1010 1110 0001 0101

6. Напишите в шестнадцатеричном виде (получив обратно отомщенную константу):

A = 47 AE 14 7A E1 47 AE 15

7. X div 100 = (((uint128(X) * uint128(A)) shr 64) + X) shr 7 (7 = 1 + S)