Блог пользователя Forestmy17

Автор Forestmy17, история, 8 часов назад, перевод, По-русски

Возможно, вы всё ещё используете малую теорему Ферма, чтобы находить обратные элементы по модулю простого числа $$$p$$$.

Если p — простое число и x не делится на p, то по малой теореме Ферма:

$$$x^{p-1}\equiv 1\pmod p$$$

Следовательно,

$$$x^{p-2}\equiv x^{-1}\pmod p$$$

Поэтому стандартный способ найти обратный элемент — использовать бинарное возведение в степень:

int binpow(int a, int n) {
    int res = 1;
    while (n) {
        if (n&1) res = 1LL * res * a % mod;
        a = 1LL * a * a % mod;
        n >>= 1;
    }
    return res;
}
int inv(int x) {
    return binpow(x, mod - 2);
}

Этот способ отлично работает. Но для такой частой операции хочется иметь что-то короче.

Например, обратный элемент можно записать так:

int inv(int x) {
    return x == 1 ? 1 : mod - 1LL * (mod / x) * inv(mod % x) % mod;
}

Теперь докажем, почему эта формула работает.

Сначала запишем mod через деление с остатком:

$$$mod = (mod / x) * x + mod \% x$$$

Следовательно,

$$$mod \% x = mod - (mod / x) * x$$$

Так как мы работаем по модулю mod, значение mod сравнимо с 0, поэтому:

$$$mod \% x\equiv -(mod/x) * x \pmod{mod}$$$

Теперь умножим обе части на $$$inv(mod\%x)$$$ $$$\big ($$$ мы знаем, что $$$(mod \% x) * inv(mod \% x) \equiv 1 \pmod{mod}\big)$$$:

$$$1\equiv -(mod / x) * x * inv(mod \% x) \pmod{mod}$$$

Так как mod — простое число и x не делится на mod, то деление на x корректно по модулю mod.

Разделив обе части на x, получаем:

$$$x^{-1} \equiv -(mod / x) * inv(mod \% x) \pmod{mod}$$$

Но $$$x^{-1}$$$ — это обратный элемент к x, поэтому:

$$$inv(x) \equiv -(mod / x) * inv(mod \% x) \pmod{mod}$$$

Конечно, мы хотим получить неотрицательное значение, поэтому вместо того, чтобы возвращать отрицательное число напрямую, добавим mod:

$$$inv(x) \equiv mod - (mod / x) * inv(mod \% x) \pmod{mod}$$$

Так мы получили рекурсивную формулу:

int inv(int x) {
    return x == 1 ? 1 : mod - 1LL * (mod / x) * inv(mod % x) % mod;
}

Рекурсия останавливается при x = 1, потому что inv(1) = 1.

Сложность — $$$O(\log mod)$$$, как и в первом решении, но реализация получается гораздо короче, чем бинарное возведение в степень.

  • Проголосовать: нравится
  • +14
  • Проголосовать: не нравится

»
8 часов назад, скрыть # |
 
Проголосовать: нравится 0 Проголосовать: не нравится

Автокомментарий: текст был переведен пользователем Forestmy17 (оригинальная версия, переведенная версия, сравнить).

»
8 часов назад, скрыть # |
 
Проголосовать: нравится 0 Проголосовать: не нравится

Really Helpful

»
8 часов назад, скрыть # |
Rev. 2  
Проголосовать: нравится 0 Проголосовать: не нравится

Is this intuition your own or have you learned from anywhere else?

»
8 часов назад, скрыть # |
 
Проголосовать: нравится 0 Проголосовать: не нравится

Note that the time complexity is O(log mod) or more precisely O(log x) (x < mod) asx/2 > mod%x >= 1 so at every iteration step you are dividing the number by a factor more than 2. (T(x) <= T(x/2) + O(1))

»
7 часов назад, скрыть # |
 
Проголосовать: нравится +6 Проголосовать: не нравится

i prefer to do binary exponentiation with a while loop to avoid the recursion overhead

»
5 часов назад, скрыть # |
 
Проголосовать: нравится 0 Проголосовать: не нравится

очень интересно, спасибо!

кажется что он особенно полезен при необходимости предподсчета обратных значений до некоторого $$$n$$$, улучшая асимптотику до строго линейной

в частности, при наличии определенной пряморукости (у меня отсутствующей), позволяет уложиться по времени при решении сегодняшней f на python

»
5 часов назад, скрыть # |
 
Проголосовать: нравится 0 Проголосовать: не нравится
const int mod = 998244353;
int mod_div(int x, int y) {
  int m = mod, u = 1, v = 0;
  while (m) swap(u -= y / m * v, v), swap(y %= m, m);
  assert(y == 1);
  return 1LL * x * (u + mod) % mod;
}
// int quotient = mod_div(x, y); // returns x * y^-1
»
5 часов назад, скрыть # |
 
Проголосовать: нравится 0 Проголосовать: не нравится
// Finds x such that a * x ≡ 1 (mod mod)
// Requires gcd(a, mod) = 1
long long inv(long long a) {
    long long b = mod, u = 1, v = 0;
    while (b) {
        long long t = a / b;
        a -= t * b; swap(a, b);
        u -= t * v; swap(u, v);
    }
    return (u % mod + mod) % mod;
}