kz-maxx

Решето Эратосфена

Быстрый алгоритм для поиска всех простых чисел до n

Что такое простое число?

Простое число — это число больше 1, которое делится только на 1 и на само себя.

Примеры простых чисел:

2 3 5 7 11 13 17 19

Например, 7 простое, потому что оно делится только на 1 и 7.

Что такое составное число?

Составное число — это число больше 1, у которого есть делители кроме 1 и самого себя.

Примеры составных чисел:

4 6 8 9 10 12 14 15

Например, 12 составное, потому что делится на:

1 2 3 4 6 12
Число 1 не считается ни простым, ни составным.

Задача

Дано число n. Нужно найти все простые числа от 1 до n.

Input

20

Output

2 3 5 7 11 13 17 19

Медленный способ

Можно для каждого числа отдельно проверять, простое оно или нет.

Проверка числа x:

bool is_prime(int x) {
    if (x < 2) {
        return false;
    }

    for (int d = 2; d * d <= x; d++) {
        if (x % d == 0) {
            return false;
        }
    }

    return true;
}

Если так проверять все числа от 1 до n, будет примерно:

O(n√n)
Для большого n это может быть слишком медленно.

Идея решета Эратосфена

Вместо того чтобы проверять каждое число отдельно, мы будем вычеркивать составные числа.

Идея:

1. Сначала считаем все числа простыми
2. Берем первое простое число p
3. Все кратные p помечаем как составные
4. Переходим к следующему непомеченному числу

Например, если p = 2, то вычеркиваем:

4 6 8 10 12 14 16 ...

Если p = 3, то вычеркиваем:

6 9 12 15 18 21 ...

Что будем хранить?

Создадим массив:

isPrime[i]

Значение:

  • true — число i пока считается простым
  • false — число i составное

Сначала все числа считаем простыми:

vector<bool> isPrime(n + 1, true);

Но 0 и 1 не простые:

isPrime[0] = false;
isPrime[1] = false;

Главный цикл

Перебираем число p. Если оно осталось простым, вычеркиваем его кратные.

for (int p = 2; p * p <= n; p++) {
    if (isPrime[p]) {
        for (int x = p * p; x <= n; x += p) {
            isPrime[x] = false;
        }
    }
}
Внутренний цикл начинается с p * p, потому что меньшие кратные уже были вычеркнуты раньше.

Почему начинаем с p * p?

Например, p = 5. Кратные 5:

10 15 20 25 30 35 ...

Но:

  • 10 уже вычеркнуло число 2
  • 15 уже вычеркнуло число 3
  • 20 уже вычеркнуло число 2

Первое число, которое точно нужно обработать для 5:

25 = 5 * 5
поэтому стартуем с p * p

Почему внешний цикл до sqrt(n)?

Если число x составное, то у него есть делитель не больше sqrt(x).

Поэтому, чтобы вычеркнуть все составные числа до n, достаточно рассматривать p до sqrt(n).

p * p <= n

После этого все оставшиеся непомеченные числа будут простыми.

Полная реализация решета

#include <bits/stdc++.h>
using namespace std;

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    int n;
    cin >> n;

    vector<bool> isPrime(n + 1, true);

    if (n >= 0) {
        isPrime[0] = false;
    }

    if (n >= 1) {
        isPrime[1] = false;
    }

    for (int p = 2; p * p <= n; p++) {
        if (isPrime[p]) {
            for (int x = p * p; x <= n; x += p) {
                isPrime[x] = false;
            }
        }
    }

    for (int i = 2; i <= n; i++) {
        if (isPrime[i]) {
            cout << i << ' ';
        }
    }

    return 0;
}

Input / Output

Input

30

Output

2 3 5 7 11 13 17 19 23 29

Как решето работает руками?

Пусть n = 20.

Сначала:

2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

Берем 2, вычеркиваем кратные:

4 6 8 10 12 14 16 18 20

Берем 3, вычеркиваем кратные:

9 15

Берем 5. Но 5 * 5 = 25, это уже больше 20.

Оставшиеся простые:

2 3 5 7 11 13 17 19

Код для подсчета количества простых чисел

Иногда нужно не вывести простые числа, а только посчитать их количество.

int cnt = 0;

for (int i = 2; i <= n; i++) {
    if (isPrime[i]) {
        cnt++;
    }
}

cout << cnt;

Полный код: количество простых до n

#include <bits/stdc++.h>
using namespace std;

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    int n;
    cin >> n;

    vector<bool> isPrime(n + 1, true);

    if (n >= 0) {
        isPrime[0] = false;
    }

    if (n >= 1) {
        isPrime[1] = false;
    }

    for (int p = 2; p * p <= n; p++) {
        if (isPrime[p]) {
            for (int x = p * p; x <= n; x += p) {
                isPrime[x] = false;
            }
        }
    }

    int cnt = 0;

    for (int i = 2; i <= n; i++) {
        if (isPrime[i]) {
            cnt++;
        }
    }

    cout << cnt;

    return 0;
}

Записать простые числа в vector

Часто удобно сохранить все простые числа в массив:

vector<int> primes;

for (int i = 2; i <= n; i++) {
    if (isPrime[i]) {
        primes.push_back(i);
    }
}

Потом можно использовать primes в других задачах.

Функция sieve

Лучше вынести решето в отдельную функцию.

vector<int> sieve(int n) {
    vector<bool> isPrime(n + 1, true);
    vector<int> primes;

    if (n >= 0) {
        isPrime[0] = false;
    }

    if (n >= 1) {
        isPrime[1] = false;
    }

    for (int p = 2; p * p <= n; p++) {
        if (isPrime[p]) {
            for (int x = p * p; x <= n; x += p) {
                isPrime[x] = false;
            }
        }
    }

    for (int i = 2; i <= n; i++) {
        if (isPrime[i]) {
            primes.push_back(i);
        }
    }

    return primes;
}

Безопасный вариант с long long

Если n большой, выражение p * p может переполнить int.

Поэтому безопаснее писать:

1LL * p * p <= n

И старт внутреннего цикла:

for (long long x = 1LL * p * p; x <= n; x += p)

Безопасная реализация

#include <bits/stdc++.h>
using namespace std;

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    int n;
    cin >> n;

    vector<bool> isPrime(n + 1, true);

    if (n >= 0) {
        isPrime[0] = false;
    }

    if (n >= 1) {
        isPrime[1] = false;
    }

    for (int p = 2; 1LL * p * p <= n; p++) {
        if (isPrime[p]) {
            for (long long x = 1LL * p * p; x <= n; x += p) {
                isPrime[x] = false;
            }
        }
    }

    for (int i = 2; i <= n; i++) {
        if (isPrime[i]) {
            cout << i << ' ';
        }
    }

    return 0;
}

Smallest Prime Factor

Иногда нужно не просто знать, простое число или нет, а быстро раскладывать числа на простые множители.

Тогда можно посчитать:

spf[x]

spf[x] — smallest prime factor, то есть минимальный простой делитель числа x.

spf[12] = 2
spf[15] = 3
spf[17] = 17

Код SPF

#include <bits/stdc++.h>
using namespace std;

const int N = 1000000;

int spf[N + 1];

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    for (int i = 0; i <= N; i++) {
        spf[i] = i;
    }

    for (int p = 2; p * p <= N; p++) {
        if (spf[p] == p) {
            for (int x = p * p; x <= N; x += p) {
                if (spf[x] == x) {
                    spf[x] = p;
                }
            }
        }
    }

    int x;
    cin >> x;

    while (x > 1) {
        cout << spf[x] << ' ';
        x /= spf[x];
    }

    return 0;
}

Пример SPF

Input

60

Output

2 2 3 5

Потому что:

60 = 2 * 2 * 3 * 5

Решето для проверки многих запросов

Часто есть много запросов: каждое число нужно проверить на простоту.

Тогда делаем решето один раз, а потом отвечаем за O(1).

if (isPrime[x]) {
    cout << "YES";
} else {
    cout << "NO";
}

Полный код с запросами

#include <bits/stdc++.h>
using namespace std;

const int N = 1000000;

bool isPrime[N + 1];

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    for (int i = 0; i <= N; i++) {
        isPrime[i] = true;
    }

    isPrime[0] = false;
    isPrime[1] = false;

    for (int p = 2; p * p <= N; p++) {
        if (isPrime[p]) {
            for (int x = p * p; x <= N; x += p) {
                isPrime[x] = false;
            }
        }
    }

    int q;
    cin >> q;

    while (q--) {
        int x;
        cin >> x;

        if (isPrime[x]) {
            cout << "YES\n";
        } else {
            cout << "NO\n";
        }
    }

    return 0;
}

Асимптотика

Метод Время Память
Проверять каждое число отдельно O(n√n) O(1)
Решето Эратосфена O(n log log n) O(n)
Ответ на запрос после решета O(1) O(n)
Решето очень быстро работает, когда нужно много простых чисел до n.

Когда использовать решето?

1. Нужно найти все простые числа до n
2. Нужно быстро отвечать на запросы isPrime[x]
3. Нужно разложение многих чисел на множители
4. Нужно посчитать количество простых чисел до n
5. Нужно использовать список primes в другой задаче

Типичные ошибки

  • Считать 1 простым числом
  • Забыть поставить isPrime[0] = false
  • Забыть поставить isPrime[1] = false
  • Начинать внутренний цикл с 2 * p, когда можно с p * p
  • Получить overflow в p * p
  • Использовать решето, когда нужен только один маленький isPrime
  • Сделать массив слишком маленьким

Короткий шаблон

vector<bool> isPrime(n + 1, true);

if (n >= 0) {
    isPrime[0] = false;
}

if (n >= 1) {
    isPrime[1] = false;
}

for (int p = 2; 1LL * p * p <= n; p++) {
    if (isPrime[p]) {
        for (long long x = 1LL * p * p; x <= n; x += p) {
            isPrime[x] = false;
        }
    }
}