ЭратосфенКак-то, ради интереса, я решил поискать в Сети материалы о нахождении всех простых чисел в диапазоне от 1 до N за линейное время. Нашлись некоторые статьи на английском, но описанные алгоритмы были несколько сложны для реализации. К сожалению, поиски на русском не увенчались успехом – были какие-то оптимизации, пытающиеся уменьшить объем памяти, константу в оценке времени, но ничего за O(N), так что я решил восполнить данный пробел.

Определение проблемы

Для начала, вспомним, что существует так называемое решето Эратосфена, позволяющее решать нашу задачу за O(N ln ln N):

bool notPrime[MAXN + 1];
int primes[MAXN];
 
int sieve(int n) {
	int primesCount = 0;
	int sqrt_n = (int)sqrt((double)n);
	for (int i = 2; i <= n; ++i) {
		if (!notPrime[i]) {
			primes[primesCount++] = i;
			if (i <= sqrt_n) {
				for (int j = i * i; j <= n; j += i) {
					notPrime[j] = true;
				}
			}
		}
	}
	return primesCount;
}

Почему решето работает не за линейное время? Потому что одно число в нем может вычеркиваться несколько раз (во внутреннем цикле) . Если добиться, чтобы каждое число вычеркивалось не более одного раза, то, очевидно, алгоритм будет линейным.

Решение

Чтобы сделать это, нам понадобится заменить массив notPrime[] на массив int leastPrime[MAXN + 1]. leastPrime[x] будет хранить индекс минимального простого числа в primes[], на которое делится x. Очевидно, что primes[leastPrime[x]] = x для простого x. Также, вместо того чтобы вычеркивать числа вида k * i, где i – простое, будем вычеркивать числа вида primes[k] * i, где i – это любое число, а primes[k] <= primes[leastPrime[i]]:

int leastPrime[MAXN + 1];
int primes[MAXN];
 
int advancedSieve(int n) {
	int primesCount = 0;
	for (int i = 2; i <= n; ++i) {
		if (!leastPrime[i]) {
			primes[primesCount++] = i;
			leastPrime[i] = primesCount;
		}
		for (int j = 0; j < leastPrime[i]; ++j) {
			int t = primes[j] * i;
			if (t <= n) leastPrime[t] = j + 1;
			else break;
		}
	}
	return primesCount;
}

Так как любое составное число однозначно представимо в виде primes[k] * i, где primes[k] <= primes[leastPrime[i]], то каждое из них в цикле во внутреннем цикле будет пройдено ровно один раз. Также ровно один раз в перед внутренним циклом будет обработано каждое простое число. Таким образом, каждое из чисел обрабатывается ровно один раз и, следовательно, получаем линейный алгоритм.

Заключение

На практике алгоритм работает примерно в 2 раза быстрее приведенного простого варианта решета и примерно столько же, сколько и оптимизированный вариант.

Преимущества алгоритма:

  • линейная оценка времени работы (против O(N ln ln N) у решета)
  • параллельно с простыми числами алгоритм вычисляет минимальный простой делитель каждого числа от 1 до N, что позволяет, если нужно, производить разложение чисел на простые множители за время, линейное от количества множителей.

Недостатки алгоритма:

  • большая константа из-за постоянного умножения во внутреннем цикле, избавится от которого пока не понятно как;
  • требования к памяти: алгоритм требует обязательного хранения списка простых чисел (в простом решете это необязательно), а также массива длины O(N), причем для каждого числа необходимо хранить количество битов, вмещающее sqrt(N) (против решета, в котором для одного числа достаточно хранить один бит).

Файлы:

P.S. Спасибо Андрею Лопатину, который рассказал об этом алгоритме на школьных сборах, где я и подслушал его.