Оптимизация умножения по модулю небольшой прайм

Question

Мне нужно сделать следующую операцию много раз:

1. Возьмите два целых числа a, b
2. Compute a * b mod p, где p = 1000000007 и a, b одного и того же порядка, что и p

Мое ощущение кишки наивное

result = a * b 
result %= p 

неэффективен. Могу ли я оптимизировать умножение по модулю p так же, как возведение в степень по модулю p оптимизирован с pow(a, b, p)?

Answer 1

Чтобы сделать этот расчет в сборке, но иметь его отозван из Python, я бы попробовать inline assembly из Python module written in C. Оба компилятора GCC и MSVC имеют встроенную сборку только с различным синтаксисом.

Обратите внимание, что наш модуль p = 1000000007 просто подходит для 30-бит. Результат желаемый (a*b)%p может быть вычислен в регистрах Intel 80x86 с учетом некоторых слабых ограничений на a,b не намного больше, чем p.

Ограничения по размеру a,b

(1) a,b являются 32-разрядные целые числа без знака

(2) a*b меньше p << 32, т.е.p раз 2^32

В частности, если a,b каждого меньше, чем 2*p, переполнение будет предотвращено. Учитывая (1), также достаточно, чтобы либо один из них был меньше p.

Инструкция Intel 80x86 MUL может умножить два 32-разрядных целых числа без знака и сохранить 64-разрядный результат в пачке регистра аккумуляторов EDX: EAX. Некоторые детали и причуды MUL обсуждаются в Разделе 10.2.1 этого полезного summary.

Инструкция DIV может затем разделить этот 64-разрядный результат на 32-битную константу (модуль p), сохраняя значение в EAX и остальную часть EDX. См. Раздел 10.2.2 последней ссылки. Результат, который мы хотим, это остаток.

Это такое разделение команд DIV, что влечет за собой риск переполнения, должны 64-битный продукт в числитель EDX: EAX дают фактор больше, чем 32 бита , будучи не удовлетворяют условию (2) выше.

Я работаю над фрагментом кода в сборке C/inline для «доказательства концепции». Однако максимальное преимущество в скорости будет зависеть от массивов дозаций данных a,b для обработки, амортизации служебных вызовов функций и т. Д. В Python (если это целевая платформа).

Answer 2

Хотя это тривиально просто, вы можете попробовать и сэкономить время на mod p этапе путем создания списка продуктов на основе 1000000007 (размер списка зависит от размера a и b). Испытайте по модулю по каждому из них (начиная с самого высокого). Конечно, это только помогает, если a & b >= sqrt(p) * 2.

Answer 3

Вы упомянули, что "a, b имеют такой же порядок величины, что и p." Часто в криптографии это означает, что a,b - большие числа рядом с p, но строго меньше p.

Если это так, то вы могли бы использовать простую идентичность

превратить ваш расчет в

result = ((a-p)*(b-p))%p 

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

Answer 4

Это не отвечает на вопрос напрямую, но я бы рекомендовал не делать этого в чистом Python, если вы ищете производительность. Некоторые опции:

• сделать небольшую библиотеку в C, что делает ваши вычисления и использовать в Python ctypes, чтобы поговорить с ним.
• numpy; вероятно, лучший вариант, если вы хотите не вмешиваться в компиляцию материала самостоятельно. Выполнение операций по одному не будет быстрее, чем собственные операторы Python, но если вы можете поместить несколько в массив numpy, вычисления на них будут намного быстрее, чем эквивалент в Python.
• Используйте cython, чтобы объявить переменные как целые числа; снова, так же, как numpy, вы выиграете от этого больше всего, если будете делать это партиями (потому что тогда вы также можете оптимизировать цикл).

Answer 5

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

Say неоптимизированный цикл был:

p = 1000000007 
b = 123456789 
a = 0 
while a < p: 
    result = (a * b) % p 
    dosomething(a, b, result) 
    a += 1 

вы могли бы оптимизировать из * и% от высокой петли частоты:

p = 1000000007 
b = 123456789 
a = 0 
result = (a * b) % p 
while a < p: 
    dosomething(a, b, result) 
    a += 1 
    result += b 
    if result >= p: 
     result -= p