Решето Аткина

Решето Аткина — алгоритм нахождения всех простых чисел дозаданного целого числа N. Алгоритм был создан и . Заявленная авторамиасимптотическая скорость работы алгоритма соответствует скорости лучшихранее известных алгоритмов просеивания, но в сравнении с ними решетоАткина требует меньше памяти.

Описание


 Основная идея алгоритма состоит в использовании неприводимыхквадратичных форм (представление чисел в виде ax²+by²). Предыдущиеалгоритмы в основном представляли собой различные модификации решетаЭратосфена, где использовалось представление чисел в виде редуцированныхформ (как правило — в виде произведения xy).
 В упрощённом виде алгоритм может быть представлен следующим образом:

  • Все числа, равные (по модулю 60) 0, 2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, 24, 26, 28, 30, 32, 34, 36, 38, 40, 42, 44, 46, 48, 50, 52, 54, 56 или 58, делятся на два и заведомо не простые. Все числа, равные (по модулю 60) 3, 9, 15, 21, 27, 33, 39, 45, 51 или 57, делятся на три и тоже не являются простыми. Все числа, равные (по модулю 60) 5, 25, 35 или 55, делятся на пять и также не простые. Все эти остатки (по модулю 60) игнорируются.
    • Все числа, равные (по модулю 60) 1, 13, 17, 29, 37, 41, 49 или 53, имеют остаток от деления на 4, равный 1. Эти числа являются простыми тогда и только тогда, когда количество решений уравнения 4x² + y² = n нечётно и само число не кратно никакому квадрату простого числа (:en:square-free integer).
    • Числа, равные (по модулю 60) 7, 19, 31, или 43, имеют остаток от деления на 6, равный 1. Эти числа являются простыми тогда и только тогда, когда количество решений уравнения 3x² + y² = n нечётно и само число не кратно никакому квадрату простого.
    • Числа, равные (по модулю 60) 11, 23, 47, или 59, имеют остаток от деления на 12, равный 11. Эти числа являются простыми тогда и только тогда, когда количество решений уравнения 3x² − y² = n (для x \textgreater y) нечётно и само число n не кратно никакому квадрату простого.
    • Отдельный шаг алгоритма вычеркивает числа, кратные квадратам простых чисел. Так как ни одно из рассматриваемых чисел не делится на 2, 3, или 5, то, соответственно, они не делятся и на их квадраты. Поэтому проверка, что число не кратно квадрату простого числа, не включает 2², 3², и 5².

Сегментация


 Для уменьшения требований к памяти «просеивание» производится порциями(сегментами, блоками), размер которых составляет примерно N.

Предпросеивание


 Для ускорения работы алгоритм игнорирует все числа, которые кратныодному из нескольких первых простых чисел (2, 3, 5, 7, \ldots). Этоделается путём использования стандартных структур данных и алгоритмов ихобработки, предложенных ранее Полом Притчардом . Они известны подназванием . Количество первых простых чисел выбирается в зависимости отзаданного числа N. Теоретически предлагается брать первые простыепримерно до logN. Это позволяет улучшить асимптотическуюоценку скорости алгоритма на множительO(1loglogN). При этом требуется дополнительнаяпамять, которая с ростом N ограничена как explogN.Увеличение требований к памяти оценивается какO(No(1)).
 Версия, представленная на сайте одного из авторов, оптимизирована дляпоиска всех простых чисел до миллиарда(log109log230=305,5),в ней исключаются из вычислений числа, кратные двум, трём, пяти и семи(2 × 3 × 5 × 7 = 210).

Оценкасложности


 По оценке авторов алгоритм имеет асимптотическую сложностьO(NloglogN) и требуетO(N12+o(1)) бит памяти. Ранее былиизвестны алгоритмы столь же асимптотически быстрые, но требующиесущественно больше памяти. Теоретически в данном алгоритме сочетаетсямаксимальная скорость работы при наименьших требованиях к памяти.Реализация алгоритма, выполненная одним из авторов, показываетдостаточно высокую практическую скорость.
 В алгоритме используется два вида оптимизации, которые существенноповышают его эффективность (по сравнению с упрощённой версией).
 Ниже представлена реализация упрощённой версии на языке программированияC, иллюстрирующая основную идею алгоритма — использование квадратичныхформ:
 \beginShaded\beginHighlighting[]\DataTypeTokint \NormalToklimit = \DecValTok1000\NormalTok;\DataTypeTokint \NormalToksqr_lim;\DataTypeTokbool \NormalTokis_prime[\DecValTok1001\NormalTok];\DataTypeTokint \NormalTokx2, y2;\DataTypeTokint \NormalToki, j;\DataTypeTokint \NormalTokn;
 \CommentTok// Инициализация решета\NormalToksqr_lim = (\DataTypeTokint\NormalTok)sqrt((\DataTypeToklong \DataTypeTokdouble\NormalTok)limit);\KeywordTokfor \NormalTok(i = \DecValTok0\NormalTok; i <= limit; i++) is_prime[i] = \KeywordTokfalse\NormalTok;\NormalTokis_prime[\DecValTok2\NormalTok] = \KeywordToktrue\NormalTok;\NormalTokis_prime[\DecValTok3\NormalTok] = \KeywordToktrue\NormalTok;
 \CommentTok// Предположительно простые — это целые с нечётным числом\CommentTok// представлений в данных квадратных формах.\CommentTok// x2 и y2 — это квадраты i и j (оптимизация).\NormalTokx2 = \DecValTok0\NormalTok;\KeywordTokfor \NormalTok(i = \DecValTok1\NormalTok; i <= sqr_lim; i++) \{ \NormalTokx2 += \DecValTok2 \NormalTok* i - \DecValTok1\NormalTok; \NormalToky2 = \DecValTok0\NormalTok; \KeywordTokfor \NormalTok(j = \DecValTok1\NormalTok; j <= sqr_lim; j++) \{ \NormalToky2 += \DecValTok2 \NormalTok* j - \DecValTok1\NormalTok; \NormalTokn = \DecValTok4 \NormalTok* x2 + y2; \KeywordTokif \NormalTok((n <= limit) && (n \NormalTokis_prime[n] = !is_prime[n]; \CommentTok// n = 3 * x2 + y2; \NormalTokn -= x2; \CommentTok// Оптимизация \KeywordTokif \NormalTok((n <= limit) && (n \NormalTokis_prime[n] = !is_prime[n]; \CommentTok// n = 3 * x2 - y2; \NormalTokn -= \DecValTok2 \NormalTok* y2; \CommentTok// Оптимизация \KeywordTokif \NormalTok((i > j) && (n <= limit) && (n \NormalTokis_prime[n] = !is_prime[n]; \NormalTok\}\NormalTok\}
 \CommentTok// Отсеиваем кратные квадратам простых чисел в интервале [5, sqrt(limit)].\CommentTok// (основной этап не может их отсеять)\KeywordTokfor \NormalTok(i = \DecValTok5\NormalTok; i <= sqr_lim; i++) \{ \KeywordTokif \NormalTok(is_prime[i]) \{ \NormalTokn = i * i; \KeywordTokfor \NormalTok(j = n; j <= limit; j += n) \{ \NormalTokis_prime[j] = \KeywordTokfalse\NormalTok; \NormalTok\} \NormalTok\}\NormalTok\}
 \CommentTok// Вывод списка простых чисел в консоль.\NormalTokprintf(\StringTok"2, 3, 5"\NormalTok);\KeywordTokfor \NormalTok(i = \DecValTok6\NormalTok; i <= limit; i++) \{ \CommentTok// добавлена проверка делимости на 3 и 5. В оригинальной версии алгоритма потребности в ней нет. \KeywordTokif \NormalTok((is_prime[i]) && (i \NormalTokprintf(\StringTok", \NormalTok\}\NormalTok\}\endHighlighting\endShaded

Версия алгоритма на языкеPascal


 \beginShaded\beginHighlighting[]\KeywordTokprogram \NormalTokatkin;\KeywordTokvar \NormalTokis_prime:\KeywordTokarray\NormalTok[\DecValTok1\NormalTok..\DecValTok10001\NormalTok] \KeywordTokof \DataTypeTokboolean\NormalTok; jj:\DataTypeTokint64\NormalTok;\KeywordTokprocedure \NormalTokdosieve(limit:\DataTypeTokint64\NormalTok);\KeywordTokvar \NormalToki,k,x,y:\DataTypeTokint64\NormalTok; n:\DataTypeTokint64\NormalTok;\KeywordTokbegin \KeywordTokfor \NormalToki:=\DecValTok5 \KeywordTokto \NormalToklimit \KeywordTokdo \NormalTokis_prime[i]:=\KeywordTokfalse\NormalTok; \KeywordTokfor \NormalTokx:=\DecValTok1 \KeywordTokto \NormalToktrunc(sqrt(limit)) \KeywordTokdo \KeywordTokfor \NormalToky:=\DecValTok1 \KeywordTokto \NormalToktrunc(sqrt(limit)) \KeywordTokdo \KeywordTokbegin \NormalTokn:=\DecValTok4\NormalTok*sqr(x)+sqr(y); \KeywordTokif \NormalTok(n<=limit) \KeywordTokand \NormalTok((n \KeywordTokmod \DecValTok12 \NormalTok= \DecValTok1\NormalTok) \KeywordTokor \NormalTok(n \KeywordTokmod \DecValTok12 \NormalTok= \DecValTok5\NormalTok)) \KeywordTokthen \NormalTokis_prime[n]:=\KeywordToknot \NormalTokis_prime[n]; \NormalTokn:=n-sqr(x); \KeywordTokif \NormalTok(n<=limit) \KeywordTokand \NormalTok(n \KeywordTokmod \DecValTok12 \NormalTok= \DecValTok7\NormalTok) \KeywordTokthen \NormalTokis_prime[n]:=\KeywordToknot \NormalTokis_prime[n]; \NormalTokn:=n\DecValTok-2\NormalTok*sqr(y); \KeywordTokif \NormalTok(x>y) \KeywordTokand \NormalTok(n<=limit) \KeywordTokand \NormalTok(n \KeywordTokmod \DecValTok12 \NormalTok= \DecValTok11\NormalTok) \KeywordTokthen \NormalTokis_prime[n]:=\KeywordToknot \NormalTokis_prime[n]; \KeywordTokend\NormalTok; \KeywordTokfor \NormalToki:=\DecValTok5 \KeywordTokto \NormalToktrunc(sqrt(limit)) \KeywordTokdo \KeywordTokbegin \KeywordTokif \NormalTokis_prime[i] \KeywordTokthen \KeywordTokbegin \NormalTokk:=sqr(i); \NormalTokn:=k; \KeywordTokwhile \NormalTokn<=limit \KeywordTokdo \KeywordTokbegin \NormalTokis_prime[n]:=\KeywordTokfalse\NormalTok; \NormalTokn:=n+k; \KeywordTokend\NormalTok; \KeywordTokend\NormalTok; \KeywordTokend\NormalTok; \NormalTokis_prime[\DecValTok2\NormalTok]:=\KeywordToktrue\NormalTok; \NormalTokis_prime[\DecValTok3\NormalTok]:=\KeywordToktrue\NormalTok;\KeywordTokend\NormalTok;\KeywordTokbegin \NormalTokdosieve(\DecValTok10000\NormalTok); \KeywordTokfor \NormalTokjj:=\DecValTok1 \KeywordTokto \DecValTok10000 \KeywordTokdo \KeywordTokif \NormalTokis_prime[jj] \KeywordTokthen \NormalTokwriteln(jj);\KeywordTokend\NormalTok.\endHighlighting\endShaded