Миллер-Рабин: лжецов не больше четверти

Если ты объяснишь, что именно тест ищет в цепочке квадратов, почему по простому модулю у единицы только два корня и откуда берётся вероятность ошибки 1/4ᵏ, - листай до «Теперь сам».

Как думаешь, что выдаёт число 561

Число 561 = 3 · 11 · 17 обманывает тест Ферма при любом взаимно простом основании: 2⁵⁶⁰ mod 561 = 1. Посмотрим, как эта единица получается. 560 = 35 · 16, поэтому 2⁵⁶⁰ - это 2³⁵, четыре раза возведённое в квадрат:

2³⁵ mod 561 = 263
263² mod 561 = 166
166² mod 561 = 67
67² mod 561 = 1
1² mod 561 = 1

Финал - единица, тест Ферма доволен. Прежде чем читать дальше, ответь себе: какая строка этой цепочки невозможна, если бы 561 было простым?

На что это похоже

Представь, что кассир проверяет не только итог чека, но и то, как он складывался. Итог сходится - но в середине чека стоит строка, какой в этом магазине быть не может: товар, которого нет в продаже. Итог подогнан, а дорога к нему выдаёт подделку.

Тест Ферма смотрит на итог. Миллер и Рабин предложили читать чек построчно.

Механизм: третий корень из единицы

Что такое квадратный корень из единицы по модулю n: число x, для которого x² ≡ 1. Два таких числа есть всегда - 1 и n − 1, то есть −1.

Пусть p простое и x² ≡ 1 (mod p). Тогда x² − 1 = (x − 1)(x + 1) делится на p. Простое число, которое делит произведение, обязано делить один из множителей: если p не делит x − 1, у x − 1 есть обратный по модулю p, и умножение на него оставляет x + 1 ≡ 0. Значит x ≡ 1 или x ≡ −1, и других корней по простому модулю нет.

У составного модуля множители могут поделить работу. В цепочке выше 67² ≡ 1 (mod 561):

(67 − 1) · (67 + 1) = 66 · 68 = 4488 = 8 · 561
66 делится на 3 и на 11,   68 делится на 17      ни один множитель не делится на 561 целиком

Число 67 - третий корень из единицы. Ответ на вопрос из начала: невозможна четвёртая строка. По простому модулю единица в цепочке квадратов может появиться только после −1 или стоять с самого начала.

Теперь сам тест. Нечётное n проверяется основанием a:

1. Записать n − 1 = d · 2ˢ, где d нечётное.         560 = 35 · 2⁴
2. Вычислить x = aᵈ mod n.
3. Если x = 1 или x = n − 1 - основание ничего не нашло.
4. Иначе возводить x в квадрат до s − 1 раза.
   Встретилось n − 1 - основание ничего не нашло.
5. Не встретилось - n составное, a - свидетель.

Свидетель здесь ловит один из двух случаев: либо в цепочке единица появилась не после −1 - третий корень, - либо единицы нет и в конце, а это провал уже по Ферма.

Правило

Тест Миллера-Рабина по основанию a объявляет нечётное n составным, если в цепочке из s чисел aᵈ, a²ᵈ, a⁴ᵈ, … по модулю n, где n − 1 = d · 2ˢ, нет ни начальной единицы, ни числа n − 1. Ответ «составное» всегда верен. Для любого нечётного составного n лжёт не больше четверти оснований (доказательство длиннее, чем для половины у Ферма, и здесь опущено), поэтому k случайных оснований дают ошибку не больше 1/4ᵏ. Простое число тест никогда не назовёт составным. Цена одного основания - одна степень по модулю.

Разбор: программа и три цепочки

def svidetel(a, n):
    d = n - 1
    s = 0
    while d % 2 == 0:
        d = d // 2
        s = s + 1
    x = pow(a, d, n)
    if x == 1 or x == n - 1:
        return False
    i = 1
    while i < s:
        x = pow(x, 2, n)
        if x == n - 1:
            return False
        i = i + 1
    return True

Функция возвращает True, если a доказывает составность. Квадрат считается через pow(x, 2, n), а не x * x % n: произведение двух чисел до n не должно выйти за точность.

Все числа в цепочках - остатки по модулю n:

n = 97, a = 5      96 = 3 · 2⁵        5³ = 28 → 8 → 64 → 22 → 96                 встретилось n − 1
n = 561, a = 2     560 = 35 · 2⁴      2³⁵ = 263 → 166 → 67 → 1                   n − 1 не было: свидетель
n = 91, a = 3      90 = 45 · 2¹       3⁴⁵ = 27, квадратов больше нет             свидетель

В первой строке простое число: цепочка дошла до 96 = −1, следующий квадрат дал бы единицу - всё законно. Вторая - число Кармайкла из начала урока. Третья особенно показательна: основание 3 для 91 было лжецом в тесте Ферма, 3⁹⁰ ≡ 1. Здесь оно свидетель: 27² = 729 = 8 · 91 + 1, число 27 - третий корень из единицы, а тест остановился на нём, не дойдя до финала.

Счёт лжецов для знакомых чисел:

n        оснований    лжецов Ферма    лжецов Миллера-Рабина    доля
91          90            36                  18               20 %
341        340           100                  50               15 %
561        560           320                  10                2 %
1729      1728          1296                 162                9 %

В счёт входят и основания 1 и n − 1, которые не свидетельствуют никогда. Больше всего лжецов у девятки: из восьми оснований лгут два, эти самые 1 и 8, - ровно четверть. Выше четверти доля не поднимается ни у одного нечётного составного числа.

Разбор: от вероятности к точному ответу

Двадцать случайных оснований - ошибка не больше 1/4²⁰, меньше одной триллионной. Сорок - 1/4⁴⁰, около 10⁻²⁴: сбой памяти от космической частицы за время работы теста вероятнее.

Для чисел ограниченного размера случайность не нужна вовсе: известно, какие наборы оснований не ошибаются ни разу.

основания                    первое составное, которое проходит все
2                                        2 047
2, 3                                 1 373 653
2, 3, 5                             25 326 001
2, 3, 5, 7                       3 215 031 751
первые 12 простых (до 37)        больше 3 · 10²³ - покрывает все 64-битные числа

def prostoe(n):
    if n < 2:
        return False
    for p in [2, 3, 5, 7]:
        if n % p == 0:
            return n == p
    for a in [2, 3, 5, 7]:
        if svidetel(a, n):
            return False
    return True

Первый цикл снимает чётные числа и сами основания. На тренажёре prostoe совпала с решетом на всех числах до ста тысяч. Цена против перебора делителей:

n                 перебор делителей    prostoe
9 973                    299              133
999 983                2 999              131
99 999 989            29 999              185

Перебор дорожает вдесятеро на каждые две цифры, тест - на десятки шагов, и не строго: его цена зависит от числа битов и от того, на каком шаге цепочка встретила n − 1. И граница: prostoe(3215031751) возвращает True, хотя 3 215 031 751 = 151 · 751 · 28 351. Это то самое первое число из таблицы; функция честна только ниже него, а выше ей нужны ещё основания.

Что ломается без этого

Первое: забытая проверка n − 1 до цикла. Если aᵈ уже равно −1, число нельзя объявлять составным. Без этой проверки тест называет составными простые.

Второе: лишний квадрат. Возводить в квадрат достаточно s − 1 раз. Ещё один квадрат дал бы aⁿ⁻¹, а оно не бывает равно −1 ни у одного нечётного числа: проверка впустую, ответ тот же, шаги потрачены.

Третье: чётные и маленькие числа. Для n = 2 и n = 3 диапазон оснований пуст, для чётных n − 1 нечётно и цепочка из одного элемента. Их отсекают до теста.

Четвёртое: известные основания против чужих чисел. Набор 2, 3, 5, 7 надёжен для чисел до трёх миллиардов. Для больших чисел фиксированный набор - приглашение: составное число, которое проходит все основания из известного списка, можно построить заранее. В 2018 году так нашли ошибки в нескольких распространённых библиотеках. Число, пришедшее извне, проверяют случайными основаниями.

Зачем это в работе

Так рождается каждый ключ RSA. Берётся случайное нечётное число из 1024 бит; оно делится на несколько сотен маленьких простых; уцелевшее проходит несколько десятков оснований Миллера-Рабина. Простым оказывается примерно один кандидат из 355 - это ln 2¹⁰²⁴ / 2. Вся генерация ключа занимает доли секунды.

В библиотеках тест спрятан за одной функцией: isProbablePrime в Java, mpz_probab_prime_p в GMP, isprime в SymPy. Часто он идёт в паре со вторым тестом другой природы - сочетание называют Бейли-ПСВ, и составного числа, которое проходит оба, не найдено ни одного.

Точный алгоритм без случайности тоже есть: в 2002 году Агравал, Каял и Саксена показали, что простота проверяется за полиномиальное время детерминированно. На практике он намного медленнее, и его значение в другом: вопрос «простое ли число» окончательно закрыт как лёгкий. Вопрос «на что оно делится» - нет, и про эту разницу следующий урок.

Теперь сам

1. n = 25, a = 7. Разложи 24 = d · 2ˢ и проведи тест: 7³ mod 25 = 18, 18² mod 25 = 24.
2. n = 15, a = 4. 14 = 7 · 2, 4⁷ mod 15 = 4. Что скажет тест? А что говорил тест Ферма?
3. Сколько случайных оснований нужно, чтобы вероятность ошибки стала меньше одной миллионной?
4. Тест по основаниям 2 и 3 назвал число 1 373 653 простым. Что можно утверждать?