Ни для кого не секрет, что самыми крупными простыми числами, известными на сегодня являются числа Мерсенна. На сегодня существует широкомасштабный открытый интернет-проект по поиску простых чисел Мерсенна Great Internet Mersenne Prime Search (GIMPS) , целью которого является поиск больших простых чисел Мерсенна путем распределенных вычислений на компьютерах любителей со всего земного шара. Впрочем, наряду с ним есть также и иные открытые интернет - проекты. Вот, только не всякий из них имеет практическую пользу. Например, проект Fermat, занимается объединенным поиском простых делителей чисел Ферма. Вот, только найденные в данном проекте факторы, вряд ли кому пригодятся, и как результат, будут всего лишь отражены на странице соответствующего сайта. Иное дело, если кому удастся найти делители мелких чисел Ферма не превышающих F(24). Тогда уже можно претендовать и на денежную премию. Все дело в том, что хотя сами делители Ферма никакой практической пользы не несут, но факт обнаружения любого неизвестного ранее делителя чисел от F(12) до F(24) означает, что появился еще один более совершенный алгоритм факторизации больших составных чисел и следует умаляюще пересматривать криптостойкость используемых алгоритмов шифрования. Но найти подобные факторы с помощью персоналки пока практически нереально. Числа же Мерсенна имеют большое прикладное значение, поскольку они применяются для построения генераторов псевдослучайных чисел с большими периодами, например Mersenne Twister. А посему и их поиск осуществляется не ради праздного любопытства и испытания пределов возможности современной вычислительной техники. И всякое алгоритмическое усовершенствование, способное упростить поиск простых чисел Мерсенна, имеет значительный интерес. В данной статье как раз и пойдет речь о таких усовершенствованиях, для которых, как собрана достаточная доказательная база в части теоретической, так и разработана алгоритмическая база и на ее основе и соответствующее ПО, с помощью которого теория была проверена на практике, дабы избежать возможных ошибок от коих никакой теоретик не застрахован.

КРИТЕРИЙ РЕШЕТОВА

Для проверки простоты чисел Мерсенна имеющих вид M(n) = 2^{n - 1}, где n - простое число, обычно применяют критерий Люка и теорему Люка-Лемера (1930 г.), благодаря которой получается довольно таки быстрый алгоритм:

S(0) = 4, S(m) = S(m - 1)² - 2 mod M(n);

Если S(n - 2) = 0 mod M(n), то число Мерсенна прошло проверку на простоту. Например для числа 2⁷ - 1 = 127.

S(0) = 4
S(1) = (4 * 4 - 2) mod 127 = 14
S(2) = (14 * 14 - 2) mod 127 = 67
S(3) = (67 * 67 - 2) mod 127 = 42
S(4) = (42 * 42 - 2) mod 127 = 111
S(5) = (111 * 111 - 2) mod 127 = 0

Известно, что тест основанный на теореме Люка - Лемера (доказательство приводить не стоит - оно многократно публиковалось, но лишь напомню, что принцип основан на том, что любое M(n) + 1 не имеет других делителей кроме 2 и ее целочисленных степеней) является самым быстрым и может быть выполнен даже на не очень производительном персональном компьютере, благодаря теореме Шенхаге-Штрассена, согласно которой умножение больших чисел посредством прямого и обратного быстрого преобразования Фурье можно выполнить со сложностью О( n * ln ₂ n) операций (а возведение в квадрат выполняется примерно вдвое быстрее). Наибольшую сложность представляет деление по модулю M(n).

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

Далее речь пойдет о критерии Решетова, с помощью которого можно доказать простоту числа Мерсенна опять же самыми примитивными вычислительными средствами, достаточными для тестирования по критерию Люка-Лемера.

Суть метода заключается в схожем итерационном процессе. И в результате алгоритм принимает форму:

S(0) = 3 или -3, S(m) = S(m - 1)² (mod M(n));

Если S(n - 1) = -3 (mod M(n)), то число Мерсенна прошло проверку на простоту по критерию Решетова. Проверим:

S(0) = -3
S(1) = -3 * -3 mod 127 = 9
S(2) = 9 * 9 mod 127 = 81
S(3) = 81 * 81 mod 127 = 84
S(4) = 84 * 84 mod 127 = 71
S(5) = 71 * 71 mod 127 = 88
S(6) = 88 * 88 mod 127 = 124 = 127 - 3 = -3



Как видно количество операций возведения в квадрат в данном методе на 1 больше по сравнению с методом Люка-Лемера, но отсутствуют вычитания. Следовательно, привлечение дополнительных вычислительных ресурсов не потребуется.

Теперь приведем доказательство простоты чисел Мерсенна, прошедших проверку по критерию Решетова. 
На самом деле, доказательство настолько тривиально, что просто удивительно, что его до сих пор никто не вывел.

Подсчитаем количество (период) между первой тройкой в алгоритме и последней. В результате получим (M(n) - 1) / 2. Это и есть количество всех квадратичных вычетов числа M(n), если оно является простым.

Теперь взглянем на усиленную теорему Ферма:

x^{(M(n) - 1) / 2} = ±1 (mod M(n)), если M(n) - простое

Поскольку всякое число Мерсенна есть ни что иное, как 4*а - 1, то 1 является квадратичным невычетом, и следовательно:

x^{(M(n) - 1) / 2} = -1 (mod M(n)) = x^{2^(n - 1) -1}  = -1 (mod 2ⁿ - 1)

Если x = 3, а 3 - квадратичный невычет для всех чисел Мерсенна M(n) (поскольку M(n) = 12*a - 5 и символ Якоби - Лежанра для J(3, 12*a + 5) = -1 => J(3, M(n)) = -1), то, подставив в вышеприведенное выражение, получим:

3^{(2^(n - 1) -1 )} = -1 (mod 2ⁿ - 1) => 3^{(2^(n - 1) )} = -3 (mod 2ⁿ - 1)

Правая часть в вышеприведенном выражении и есть критерий Решетова. А поскольку период от первой 3 до последней в этом критерии пробегает все квадратичные (не)вычеты, то значит данный критерий будет справедлив и для них (каждый из них может встречаться не более одного раза за весь период, что справедливо для простых чисел, коли длина их периодов не может превышать расчетную. Впрочем для составных чисел все равно при такой длине повторов быть не может, поскольку длина их периода почти в квадрат больше). Следовательно один единственный тест по критерию Решетова, в силу этих обстоятельств, если доказывает простоту M(n) для модуля 3, то одновременно также доказывает и простоту для всех модулей являющихся квадратичными (не)вычетами. Чтобы доказать простоту сразу для стольких модулей методом Миллера - Рабина, необходимо будет прогнать его по каждому модулю отдельно. С учетом того, что метод Миллера - Рабина гораздо более ресурсоемкий, получаем высочайшую эффективность по сравнению с дедовскими методами проверки простоты. Единственное, что следует сделать в случае сомнений, так это погонять тест Миллера - Рабина после критерия Решетова, только для тех модулей, что не являются квадратичными (не)вычетами для M(n).

ДЕЛЕНИЕ ПО МОДУЛЮ M(n) по упрощенному методу Решетова

Еще одно усовершенствование, которое было применено в алгоритме и позволило повысить скорость вычислений, это быстрое деление по модулю M(n). Деление всегда сложнее умножения и потому даже высокоскоростной тест Люка-Лемера больше ресурсов тратил именно на эту операцию, нежели на все остальное. Критерий Решетова построен по тому же самому принципу, а следовательно за счет умаления вычислительных ресурсов на деление, можно значительно усовершенствовать алгоритм.

Докажем, что деление по модулю M(n) для числа, битность которого в х раз больше, нежели битность самого M(n) можно выполнить за (х + 1) * ( сложение + (деление по модулю 2 - xor) + 2 * (cдвигов на n разрядов)). Для этой цели предлагается алгоритм:

1. Пока битность делимого - х больше битности M(n) перейти на пункт 2, в противном случае на п. 4
2. Получим частное - a и остаток b от деления х на 2ⁿ. Т.е. а = x / 2ⁿ; b = x % 2ⁿ
3. Возьмем в качестве нового x сумму частного - а и остатка - b полученных на предыдущем этапе (x = a + b) и перейти к п. 1
4. x является остатком деления по модулю M(n)

Пояснение к п. 1. В качестве проверки конечно же учитывается не битность, а значение делимого. Но поскольку, любое число превышающее M(n) будет, как минимум на один бит длиннее, то необходимой и достаточной является сравнение битности делимого и делителя.

Что делает данный алгоритм. Предположим, что у нас есть некое число M(n). Поскольку M(n) + 1 = 2ⁿ, то всякое вычитание из х числа 2^y*M(n), не превосходящего по значению самого х, соответствует операции х - (2ⁿ + 1) * 2^y. Вычитание любых чисел кратных M(n) из x не влияет на конечный результат, т.е. остаток:

x = x - a*M(n) = r mod M(n)

Собственно, алгоритм каждый раз берет сразу большое множество 2ⁿ и вычитает их одним махом на втором шаге алгоритма, а на третьем шаге добавляет 1 соответствующие младшим разрядам на ровно n бит меньше. Тем самым каждый шаг алгоритма уменьшает х минимум на n - 1 бит.

Поскольку частное и остаток от любого 2^m - есть ни что иное как операции побитовых сдвигов, то мы получаем высокоскоростной алгоритм деления по модулю M(n).

Поскольку после каждого возведения в квадрат для получения критерия Решетова, мы получаем число битность которого не превышает вдвое битности M(n), то следовательно алгоритм деления по модулю M(n) будет на каждом шаге циклически выполнен не более двух раз.

Вот собственно и все. Напоследок приведу пример программной реализации всего вышеописанного, дабы любой желающий мог убедиться на практике (Пример на языке Java, но в контексте понятном разработчикам на С):


/**
* Тестирование числа Мерсенна - M(n) по критерию Решетова
* @param n - степень
* @return простота числа M(n) по критерию Решетова
*/
private boolean primalityTestByYuryReshetov(int n) {
// y = 2^n
long y = (1 << n);
// Получаем M(n);
long m = y - 1;
// инициализация
long r = 3;
// Цикл на n -1 шагов
for (int i = 0; i < n - 1; i++) {
//r = r^2 mod M(n);
r = modByYuryReshetov(r * r, m, n);
}
// Число прошло тест, если остаток равен -3
return r == (m - 3);
// В противном случае оно является составным
}

/**
* Деление числа x по модулю M(n)
* @param x - делимое
* @param m - делитель
* @param n - степень
* @return остаток от деления
*/
private long modByYuryReshetov (long x, long m, int n) {
long r = x;
while (r > m) { /* Пока остаток больше M(n) */
long a = r >> n; /* a = r / 2^n */
long b = (a << n) ^ r; /* b = r % 2^n */
r = a + b;
}
return r;
}

Вышеприведенная программная реализация скорее для наглядности, нежели для сколь нибудь практичного применения. Ведь, если в Java числа примитивного типа long имеют 61 бит в положительной части, то предельно возможным для данной проверки будет число М(29).

Более боеспособная, где верхний предел ограничен лишь объемом используемого ОЗУ (плюс свопинг) версия на Pure Java уже выглядит так:

package sample;

/**
* <p>Reshetov testing </p>
* <p>Description: Reshetov primality test</p>
* <p>Copyright: Copyright (c) 2005</p>
* <p>Company: Reshetov & Co.</p>
* @autor Yury V. Reshetov
* @version 1.0
*/

import java.math.*;

public class PrimalityTest {

private BigInteger two = new BigInteger("2");
private BigInteger n3 = new BigInteger("3");

public PrimalityTest() {
BigInteger n = BigInteger.ONE;
while (true) {
n = n.add(two);
if (n.isProbablePrime(1) && primalityTestByYuryReshetov(n)) {
System.out.println("2^" + n + "-1 is prime");
}
}
}

/**
* Тестирование числа Мерсенна - M(n) по критерию Решетова
* @param n - степень
* @return простота числа M(n) по критерию Решетова
*/
private boolean primalityTestByYuryReshetov(BigInteger n) {
int nn = n.intValue();
BigInteger y = BigInteger.ONE.shiftLeft(nn);
// Получаем M(n);
BigInteger m = y.subtract(BigInteger.ONE);
// Инициализация
BigInteger r = n3;
// Цикл
for (int i = 0; i < nn - 1; i++) {
//r = r * r % M(n);
r = this.modByYuryReshetov(r.multiply(r), m, nn);
}
// Число прошло тест, если остаток равен -3
return r.equals(m.subtract(n3));
}

/**
* Деление числа x по модулю M(n)
* @param x - делимое
* @param m - делитель M(n)
* @param n - степень
* @return остаток от деления
*/
private BigInteger modByYuryReshetov(BigInteger x, BigInteger m, int n) {
BigInteger r = x;
// Пока битовая длина остатка больше битовой длины M(n)
while (r.bitLength() > m.bitLength()) {
// a = r / 2^n = r >> n
BigInteger a = r.shiftRight(n);
// b = r % 2^n = (r << n) xor r
BigInteger b = a.shiftLeft(n).xor(r);
// r = a + b
r = a.add(b);
}
return r;
}

public static void main(String[] args) {
new PrimalityTest();
}
}

Чтобы убедиться в том, насколько данная реализация эффективнее стандартной проверки на простоту, которая встроена в библиотеку Java, следует включить эту самую стандартную проверку? заменив строку:

if (n.isProbablePrime(1) && primalityTestByYuryReshetov(n)) {

на:

if (n.isProbablePrime(1) && two.pow(n.intValue()).subtract(BigInteger.ONE)) .isProbablePrime(1) {



Авторские права:

Автором и разработчиком алгоритмов тестирования простоты чисел Мерсена, а также упрощенного деления по модулю M(n) является (с) Решетов Юрий Вячеславович (xUSSR, г. Ташкент).

Поскольку вышеизложенные алгоритмы имеют коммерческую ценность, т.к. всякому, кто обнаружит простое число размером 10⁷ символов причитается денежный от Electronic Frontier Foundation (http://www.eff.org/ и http://www.mersenne.org/prize.htm) , то следует оговорить условия использования данных алгоритмов. 

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

Оценка вероятности простоты M(n)

Рассмотрим ряд чисел:

-3^x (mod M(n)), -3^{(x + 1)} (mod M(n)) ... -3^(y-1) (mod M(n)), -3^y(mod M(n))

где: x и y - натуральные числа. А M(n) - простое число Мерсенна - 2ⁿ - 1

Как было уже доказано, данный ряд должен повториться, и период повторяемости равен 2^{(n - 1)} - 1 = 2ⁿ / 2 - 1 = (M(n) - 1) / 2.
Значит следующее число ряда после 2^{(n - 1)} - 1, есть 2^{(n - 1)} и является степенью двойки. Следовательно каждое число в ряду - квадратичный (не)вычет по модулю M(n), а общее количество этих самых квадратичных (не)вычетов равно периоду повторения, т.е. = (M(n) - 1) / 2. Поскольку для всякого простого p количество квадратичных (не)вычетов равно (p - 1) / 2, то значит наш ряд за один период пробегает все квадратичные невычеты без исключения. А посему, если для алгоритма проверки простоты чисел Мерсенна посредством критерия Решетова можно взять любое число из вышеуказанного ряда -3^x(mod M(n)) при любом х, и если число М(n) является простым, то:

-3^x(mod M(n))^{(2^(n - 1))} = -3^x(mod M(n)) (mod M(n))

Следовательно, проверка по критерию Решетова для любого -3^{x(mod M(n))} (mod M(n)) при любом натуральном х, комлексно доказывает простоту и остальных -3^y(mod M(n)) (mod M(n)), для всякого y <> x, поскольку доказывает соответствие вышеуказанной периодичности.

Нам также известно, что общее количество модулей по которым одновременно за один проход теста производится комплексная проверка равно (M(n) - 1) / 2.

Для оценки вероятности того, что число M(n) пройдя тест по критерию Решетова, окажется простым, воспользуемся оценкой Леманна (Lehmann R. S.), cогласно которой если число a^{((p -1) / 2)} = ±1 (mod p), то вероятность простоты p равна 50%. Если провести t проверок по критерию Леманна, то вероятность простоты числа будет равна 1 - 1 / 2^t.

В этом случае, критерий Решетова за один проход теста доказывает простоту M(n) c вероятностью 

1 - 1 / 2^{((M(n) - 1) / 2)}. 

Например, для 2⁷ - 1 = 127, критерий Решетова доказывает простоту с вероятностью: 

1 - 1 / 2⁶³ = 0,99999999999999999989157978275145.

Чем больше число Мерсенна, тем выше вероятность.

ПРИМЕЧАНИЕ:

В действительности числа Мерсенна можно представить в виде 12*x - 5. По критерию Эйлера всякое число вида p = |5| (mod 12) для любого простого p имеет квадратичный невычет 3.

А доказательство того, что числа Мерсенна имеют вид 12*x - 5 тривиально. Возьмем, например одно из них, скажем 2³ - 1 = 7 = 12 * 1 - 5. Поскольку всякое 2ⁿ - 1 = 4 * ( 2^{(n -2)} - 1) + 3, для любого нечетного натурального n, то:

4 * (12 * х - 5) + 3 = 48 * х - 20 + 3 = 48 * х - 17 = 48 * х - 12 - 5 = 12 * (4 * х - 1) - 5

Поскольку у чисел Мерсенна M(n), n - нечетно, то все они вида 12*x - 5 и всякое простое из них в силу критерия Эйлера должно обязательно иметь квадратичный невычет 3.