MAXimal

добавлено: 6 Sep 2011 1:03
редактировано: 23 Mar 2012 3:58

Решето Эратосфена с линейным временем работы

Дано число $n$ . Требуется найти все простые в отрезке $[2; n]$ .

Классический способ решения этой задачи — решето Эратосфена. Этот алгоритм очень прост, но работает за время $O (n \log \log n)$ .

Хотя в настоящий момент известно достаточно много алгоритмов, работающих за сублинейное время (т.е. за $o(n)$ ), описываемый ниже алгоритм интересен своей простотой — он практически не сложнее классического решета Эратосфена.

Кроме того, приводимый здесь алгоритм в качестве "побочного эффекта" фактически вычисляет факторизацию всех чисел в отрезке $[2; n]$ , что может быть полезно во многих практических применениях.

Недостатком приводимого алгоритма является то, что он использует больше памяти, чем классическое решето Эратосфена: требуется заводить массив из $n$ чисел, в то время как классическому решету Эратосфена достаточно лишь $n$ бит памяти (что получается в $32$ раза меньше).

Таким образом, описываемый алгоритм имеет смысл применять только до чисел порядка $10^7$ , не более.

Авторство алгоритма, по всей видимости, принадлежит Грайсу и Мисра (Gries, Misra, 1978 г. — см. список литературы в конце). (И, собственно говоря, называть данный алгоритм "решетом Эратосфена" некорректно: слишком отличаются эти два алгоритма.)

Описание алгоритма

Наша цель — посчитать для каждого числа $i$ от в отрезке $[2; n]$ его минимальный простой делитель $lp[i]$ .

Кроме того, нам потребуется хранить список всех найденных простых чисел — назовём его массивом $pr[]$ .

Изначально все величины $lp[i]$ заполним нулями, что означает, что мы пока предполагаем все числа простыми. В ходе работы алгоритма этот массив будет постепенно заполняться.

Будем теперь перебирать текущее число $i$ от $2$ до $n$ . У нас может быть два случая:

$lp[i] = 0$ — это означает, что число $i$ — простое, т.к. для него так и не обнаружилось других делителей.
Следовательно, надо присвоить $lp[i] = i$ и добавить $i$ в конец списка $pr[]$ .
$lp[i] \ne 0$ — это означает, что текущее число $i$ — составное, и его минимальным простым делителем является $lp[i]$ .

В обоих случаях дальше начинается процесс расстановки значений в массиве $lp[]$ : мы будем брать числа, кратные $i$ , и обновлять у них значение $lp[]$ . Однако наша цель — научиться делать это таким образом, чтобы в итоге у каждого числа значение $lp[]$ было бы установлено не более одного раза.

Утверждается, что для этого можно поступить таким образом. Рассмотрим числа вида:

$x_j = i \cdot p_j,$

где последовательность $p_j$ — это все простые, не превосходящие $lp[i]$ (как раз для этого нам понадобилось хранить список всех простых чисел).

У всех чисел такого вида проставим новое значение $lp[x_j]$ — очевидно, оно будет равно $p_j$ .

Почему такой алгоритм корректен, и почему он работает за линейное время — см. ниже, пока же приведём его реализацию.

Реализация

Решето выполняется до указанного в константе числа $N$ .

const int N = 10000000;
int lp[N+1];
vector<int> pr;
 
for (int i=2; i<=N; ++i) {
	if (lp[i] == 0) {
		lp[i] = i;
		pr.push_back (i);
	}
	for (int j=0; j<(int)pr.size() && pr[j]<=lp[i] && i*pr[j]<=N; ++j)
		lp[i * pr[j]] = pr[j];
}

Эту реализацию можно немного ускорить, избавившись от вектора $pr$ (заменив его на обычный массив со счётчиком), а также избавившись от дублирующегося умножения во вложенном цикле $for$ (для чего результат произведения надо просто запомнить в какой-либо переменной).

Доказательство корректности

Докажем корректность алгоритма, т.е. что он корректно расставляет все значения $lp[]$ , причём каждое из них будет установлено ровно один раз. Отсюда будет следовать, что алгоритм работает за линейное время — поскольку все остальные действия алгоритма, очевидно, работают за $O (n)$ .

Для этого заметим, что у любого числа $i$ единственно представление такого вида:

$i = lp[i] \cdot x,$

где $lp[i]$ — (как и раньше) минимальный простой делитель числа $i$ , а число $x$ не имеет делителей, меньших $lp[i]$ , т.е.:

$lp[i] \le lp[x].$

Теперь сравним это с тем, что делает наш алгоритм — он фактически для каждого $x$ перебирает все простые, на которые его можно домножить, т.е. простые до $lp[x]$ включительно, чтобы получить числа в указанном выше представлении.

Следовательно, алгоритм действительно пройдёт по каждому составному числу ровно один раз, поставив у него правильное значение $lp[]$ .

Это означает корректность алгоритма и то, что он работает за линейное время.

Время работы и требуемая память

Хотя асимптотика $O (n)$ лучше асимптотики $O (n \log \log n)$ классического решета Эратосфена, разница между ними невелика. На практике это означает лишь двукратную разницу в скорости, а оптимизированные варианты решета Эратосфена и вовсе не проигрывают приведённому здесь алгоритму.

Учитывая затраты памяти, которые требует этот алгоритм — массив чисел $lp[]$ длины $n$ и массив всех простых $pr[]$ длины примерно $n / \ln n$ — этот алгоритм кажется уступающим классическому решету по всем статьям.

Однако спасает его то, что массив $lp[]$ , вычисляемый этим алгоритмом, позволяет искать факторизацию любого числа в отрезке $[2; n]$ за время порядка размера этой факторизации. Более того, ценой ещё одного дополнительного массива можно сделать, чтобы в этой факторизации не требовались операции деления.

Знание факторизации всех чисел — очень полезная информация для некоторых задач, и этот алгоритм является одним из немногих, которые позволяют искать её за линейное время.

Литература

David Gries, Jayadev Misra. A Linear Sieve Algorithm for Finding Prime Numbers [1978]