Мотивация и контекст

Функция $\pi(x)$ — количество простых чисел, не превосходящих $x$ — является одной из центральных объектов аналитической теории чисел. По теореме о простых числах: $\pi(x) \sim x / \ln x$. Но нас интересует точное значение.

Зачем нужен точный $\pi(x)$?

• Задачи соревновательного программирования: "Найдите $k$-е простое число по номеру" для $k \le 10^{12}$.
• Криптография: оценки числа простых в диапазоне, проверка параметров схем.
• Теоретико-числовые задачи: подсчёт простых в арифметических прогрессиях.
• Верификация гипотезы Римана: $\pi(x)$ сравнивается с $\text{li}(x)$.

Хронология алгоритмов подсчёта $\pi(x)$:

Метод Мейсселя–Лемера в его классической форме позволяет вычислить $\pi(x)$ при $x \le 10^{14}$ за несколько секунд с памятью $O(x^{1/3})$. При правильной реализации это быстрее Min-25 Sieve для очень больших $x$.

---

Теория

Формула Лежандра

Исторически первой явной формулой стала формула Лежандра (1808):

$\pi(x) = \pi(\sqrt{x}) - 1 + \Phi(x, \pi(\sqrt{x}))$

​)−1+Φ(x,π(x

​))

где $\Phi(x, a)$ — количество целых чисел в $[1, x]$, не делящихся ни на одно из первых $a$ простых чисел.

$\Phi$ вычисляется рекурсивно: $\Phi(x, a) = \Phi(x, a-1) - \Phi(x/p_a, a-1)$ с $\Phi(x, 0) = \lfloor x \rfloor$.

Проблема: рекурсия растёт экспоненциально, время $O(x^{1-\varepsilon})$ неприемлемо.

Функция Мейсселя $\phi(x, a)$

Определение. $\phi(x, a)$ = количество натуральных чисел $n \le x$, у которых наименьший простой делитель $> p_a$ (либо $n = 1$). Это та же $\Phi$ Лежандра.

Ключевое разложение Мейсселя:

$\phi(x, a) = \phi(x, b) - \sum_{i=b+1}^{a} \phi\!\left(\frac{x}{p_i}, i-1\right)$

для любого $b < a$. Выбирая $b = \pi(x^{1/3})$, получаем управляемую рекурсию.

Классическая формула Мейсселя

Обозначим $a = \pi(x^{1/4})$, $b = \pi(x^{1/2})$, $c = \pi(x^{1/3})$.

Теорема (Мейссель, 1870). $\pi(x) = \phi(x, a) + a - 1 - P_2(x, a)$

где $P_2(x, a)$ — количество чисел $n \le x$, имеющих ровно два простых делителя, каждый из которых $> p_a$.

Доказательство. Разобьём числа $n \in [2, x]$ по числу простых делителей (с учётом $> p_a$):

• Числа с 0 такими делителями: это числа, у которых все делители $\le p_a$. Их $\phi(x, a) - 1$ (вычитаем единицу для $n=1$).
• Числа с ровно 1 таким делителем: это простые $p > p_a$, т.е. $\pi(x) - a$ штук.
• Числа с ровно 2 такими делителями: $P_2(x, a)$.
• Числа с $\ge 3$ такими делителями: $> p_a^3 > x^{3/4} \cdot x^{3/4}$... нет, если $a = \pi(x^{1/4})$, то $p_a \approx x^{1/4}$, и три таких делителя дают $> x$ — таких чисел нет.

Итого: $\phi(x, a) - 1 + (\pi(x) - a) + P_2(x, a) = x -$ ... нет, перефразируем.

Числа от 1 до $x$ делятся на категории по числу простых делителей $> p_a$:

• Нет таких: $\phi(x, a)$ штук.
• Ровно один: это простые из $(p_a, x]$, т.е. $\pi(x) - a$ штук.
• Ровно два: $P_2(x, a)$ штук.
• Три и более: при $p_a > x^{1/4}$ произведение трёх таких простых $> x$ — нет таких чисел в $[1,x]$.

Сумма: $\phi(x, a) + (\pi(x) - a) + P_2(x, a) = x$ (все числа, включая 1)? Нет: $\phi(x,a) + (\pi(x)-a) + P_2(x,a) = x$. Но $\phi(x,a)$ включает 1, поэтому: $\phi(x,a) - 1 + (\pi(x) - a) + P_2(x, a) = x - \text{(числа с тремя делителями)} = ...$

Точная форма (после перегруппировки): $\pi(x) = \phi(x, a) + a - 1 - P_2(x, a)$ $\square$

Уточнение Лемера

Определение. $P_2(x, a)$ — количество чисел $n = p \cdot q \le x$, где $p_a < p \le q$ — простые.

Вычисление $P_2$: $P_2(x, a) = \sum_{i=a+1}^{b} \left( \pi\!\left(\frac{x}{p_i}\right) - (i - 1) \right)$

Здесь $b = \pi(\sqrt{x})$

​), и для каждого простого $p_i$ мы считаем количество простых $q$: $p_i \le q \le x/p_i$.

Доказательство. Число $n = p \cdot q$ с $p \le q$ и оба $> p_a$. Перебираем меньший множитель $p$: он пробегает простые от $p_{a+1}$ до $p_b = p_{\pi(\sqrt{x})} \approx \sqrt{x}$

​)​≈x

​ (ведь $p \le q$ значит $p \le \sqrt{x}$). Для фиксированного $p$ количество допустимых $q$: это простые в $[p, x/p]$, то есть $\pi(x/p) - \pi(p) + 1 = \pi(x/p) - i + 1$ (если $p = p_i$). Но $q \ge p = p_i$, значит $q > p_{i-1}$, значит $\pi(x/p) - (i-1)$ штук. $\square$

Рекурсия $\phi(x, a)$

Для вычисления $\phi(x, a)$ разобьём: $\phi(x, a) = \phi(x, a-1) - \phi\!\left(\frac{x}{p_a}, a-1\right)$

Основные случаи:

• $\phi(x, 0) = \lfloor x \rfloor$ — все числа.
• $a \ge \pi(\sqrt{x})$: тогда $\phi(x, a) = \pi(x) - a + 1$ (все числа $\le x$ с наим. делителем $> p_a$ — это простые $> p_a$, плюс 1).
• $x \le p_a$: тогда $\phi(x, a) = 1$ (только число 1).

Оптимизация — мемоизация: значения $\phi(x, a)$ при малых $a$ и большом $x$ занимают много памяти. Лемер предложил остановить рекурсию при $a \le c = \pi(x^{1/3})$ и использовать факт, что числа с наим. делителем $> p_c$ в диапазоне $[1, x]$ — либо простые, либо имеют ровно два простых делителя, оба $> p_c$.

Полная формула Лемера

$\pi(x) = \phi(x, a) + (b + a - 2)(b - a + 1)/2 - \sum_{i=a+1}^{b} \pi\!\left(\frac{x}{p_i}\right) - P_2(x, a)$

Упрощённее: $\pi(x) = \phi(x, a) + a - 1 - \sum_{i=a+1}^{b} \left(\pi\!\left(\frac{x}{p_i}\right) - (i-1)\right)$

что совпадает с $\phi(x,a) + a - 1 - P_2(x,a)$.

Сложность и оптимизация памяти

Вычисление $\pi(v)$ для $v = x/p_i$: используем либо предподсчитанное решето (для малых $v \le x^{2/3}$), либо рекурсивно.

Решето до $x^{2/3}$: хранится в памяти $O(x^{2/3})$ — слишком много. Оптимизация: хранить только префиксные суммы простых до $x^{1/3}$, а $\pi$ для больших значений вычислять рекурсивно.

Итоговая сложность:

• Время: $O(x^{2/3} / \log^2 x)$ (с лучшими вариантами доходит до $O(x^{2/3} / \log x)$).
• Память: $O(x^{1/3} \log^2 x)$ — принципиально лучше, чем у Min-25 Sieve.

---

Ключевые формулы (сводная таблица)

---

Реализация на C++20 (с анализом сложности)

#include <bits/stdc++.h>
using namespace std;
using ll  = long long;
using lll = __int128;

// ============================================================
// Метод Мейсселя–Лемера: вычисление pi(x)
// Время: O(x^{2/3}/log^2 x), Память: O(x^{1/3})
// ============================================================

struct MeisselLehmer {
    ll x;
    int cbrtx;    // floor(x^{1/3})
    int sqrtx;    // floor(x^{1/2})

    // Простые числа до sqrt(x), хранятся явно
    vector<int> primes;
    // pi_small[v] = pi(v) для v <= cbrtx*cbrtx (= x^{2/3})
    // Хранится как prefix sum
    vector<int> pi_small;
    // sieve до x^{2/3}
    int sieve_limit;
    vector<bool> sieve;

    explicit MeisselLehmer(ll x_) : x(x_) {
        sqrtx = (int)sqrtl((long double)x_);
        while((ll)sqrtx * sqrtx > x_) --sqrtx;
        while((ll)(sqrtx+1)*(sqrtx+1) <= x_) ++sqrtx;

        cbrtx = (int)cbrtl((long double)x_);
        while((ll)cbrtx*cbrtx*cbrtx > x_) --cbrtx;
        while((ll)(cbrtx+1)*(cbrtx+1)*(cbrtx+1) <= x_) ++cbrtx;

        // sieve_limit = min(sqrtx * sqrtx, cbrtx * cbrtx * cbrtx) = примерно x^{2/3}
        // но мы решеваем до sqrtx для вычисления P2
        sieve_limit = sqrtx;

        // Решето до sqrtx
        sieve.assign(sieve_limit + 1, true);
        sieve[0] = sieve[1] = false;
        for(int i = 2; i <= sieve_limit; i++){
            if(!sieve[i]) continue;
            primes.push_back(i);
            for(ll j = (ll)i*i; j <= sieve_limit; j += i)
                sieve[j] = false;
        }

        // pi_small[v] = pi(v) для v <= sieve_limit
        pi_small.resize(sieve_limit + 1, 0);
        for(int v = 1; v <= sieve_limit; v++)
            pi_small[v] = pi_small[v-1] + (sieve[v] ? 1 : 0);
    }

    // pi(v) для v <= sqrtx: O(1)
    ll pi_fast(ll v) const {
        if(v <= sieve_limit) return pi_small[v];
        // Не должно вызываться для v > sieve_limit в основном алгоритме
        return -1;
    }

    // Функция phi(v, a): количество n <= v с lpf(n) > p_a
    // Мемоизация для малых a (a <= 6 покрывает первые простые 2,3,5,7,11,13)
    // Для больших a — рекурсия
    unordered_map<ll, ll> cache;

    ll phi(ll v, int a) {
        if(a == 0) return v;
        if(v == 0) return 0;
        if(a == 1) return (v + 1) / 2; // нечётные числа <= v

        // Если v маленькое — используем pi_small
        if(v <= sieve_limit){
            // phi(v, a) = pi(v) - a + 1 если p_a <= v < p_{a+1}
            // Точнее: phi(v, a) = #{1, и простые > p_a и <= v}
            // = 1 + #{простые в (p_a, v]}
            if(a >= (int)primes.size()) return 1;
            if(primes[a-1] >= v) return 1; // только 1
            return 1 + pi_fast(v) - a;
        }

        // Кэш
        // ключ: (v << 10) | a, только для малых a
        if(a <= 7){
            ll key = (v << 6) | a;
            auto it = cache.find(key);
            if(it != cache.end()) return it->second;
        }

        ll result = phi(v, a-1) - phi(v / primes[a-1], a-1);

        if(a <= 7) cache[(v << 6) | a] = result;
        return result;
    }

    // P2(x, a): количество n = p*q <= x, p_a < p <= q простые
    ll P2(ll v, int a) {
        // b = pi(sqrt(v))
        ll sqrt_v = (ll)sqrtl((long double)v);
        while(sqrt_v * sqrt_v > v) --sqrt_v;
        int b = (int)pi_fast(sqrt_v);
        ll result = 0;
        // sum_{i=a+1}^{b} (pi(v / p_i) - (i-1))
        for(int i = a; i < b; i++){     // i — 0-индекс, primes[i] = p_{i+1}
            ll pi_v_pi = pi_fast(v / primes[i]);
            if(pi_v_pi < i + 1) break; // pi(v/p_i) < i => нет вклада
            result += pi_v_pi - i;     // pi(v/p_i) - i (в 0-индексации i = a, значит i+1 - 1 = i)
        }
        return result;
    }

    // Основная функция: pi(x)
    ll compute() {
        if(x < 2) return 0;
        if(x <= sieve_limit) return pi_fast(x);

        // a = pi(x^{1/4})
        ll x14 = (ll)pow((double)x, 0.25);
        while(x14 * x14 * x14 * x14 > x) --x14;
        while((x14+1)*(x14+1)*(x14+1)*(x14+1) <= x) ++x14;
        int a = (int)pi_fast(x14);

        // phi(x, a)
        ll phi_val = phi(x, a);

        // P2(x, a)
        ll p2_val = P2(x, a);

        return phi_val + a - 1 - p2_val;
    }
};

// --------------------------------------------------------
// Версия с полной поддержкой больших x через рекурсивный pi
// --------------------------------------------------------

// Глобальные данные (для удобства рекурсии)
namespace ML {
    ll X;
    int SQ;
    vector<int> primes;
    vector<int> prime_count; // prime_count[v] = pi(v) для v <= SQ

    void init(ll x) {
        X = x;
        SQ = (int)sqrtl((long double)x);
        while((ll)SQ*SQ > x) --SQ;
        while((ll)(SQ+1)*(SQ+1) <= x) ++SQ;

        vector<bool> isp(SQ + 1, true);
        isp[0] = isp[1] = false;
        for(int i = 2; i <= SQ; i++){
            if(!isp[i]) continue;
            primes.push_back(i);
            for(ll j = (ll)i*i; j <= SQ; j += i) isp[j] = false;
        }
        prime_count.resize(SQ + 1, 0);
        for(int v = 1; v <= SQ; v++)
            prime_count[v] = prime_count[v-1] + (isp[v] ? 1 : 0);
    }

    ll pi_small(ll v){ return v <= SQ ? prime_count[v] : -1; }

    // phi рекурсивно с отсечкой
    ll phi_rec(ll v, int a){
        if(a == 0) return v;
        if(v == 0) return 0;
        if(a == 1) return (v+1)/2;
        if(v <= SQ && primes[a-1] >= v) return 1;
        if(v <= SQ) return 1 + prime_count[v] - a;
        return phi_rec(v, a-1) - phi_rec(v / primes[a-1], a-1);
    }

    ll pi_rec(ll v);

    ll P2_val(ll v, int a){
        ll sv = (ll)sqrtl((long double)v);
        while(sv*sv > v) --sv;
        int b = (int)pi_small(sv);
        ll res = 0;
        for(int i = a; i < b; i++){
            ll quot = v / primes[i];
            ll piv  = (quot <= SQ) ? prime_count[quot] : pi_rec(quot);
            if(piv <= i) break;
            res += piv - i;
        }
        return res;
    }

    ll pi_rec(ll v){
        if(v <= SQ) return prime_count[v];
        ll sv = (ll)sqrtl((long double)v);
        while(sv*sv > v) --sv;
        ll v14 = (ll)pow((double)v, 0.25);
        while(v14*v14*v14*v14 > v) --v14;
        while((v14+1)*(v14+1)*(v14+1)*(v14+1) <= v) ++v14;
        int a = (int)pi_small(v14);
        return phi_rec(v, a) + a - 1 - P2_val(v, a);
    }
} // namespace ML

ll meissner_lehmer(ll x){
    if(x < 2) return 0;
    ML::init(x);
    return ML::pi_rec(x);
}

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

    ll x;
    cin >> x;
    cout << "pi(" << x << ") = " << meissner_lehmer(x) << "\n";

    // Тесты
    // pi(10^6) = 78498
    // pi(10^9) = 50847534
    // pi(10^{12}) = 37607912018
    return 0;
}

Анализ сложности реализации:

---

Разбор задачи 1 (средняя): $k$-е простое число

Условие. Дано $k \le 10^{12}$. Найти $k$-е простое число $p_k$.

Решение. По теореме о простых числах: $p_k \approx k \ln k$. Бинарным поиском найдём $p_k$.

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

// Используем meissner_lehmer(x) из реализации выше

ll kth_prime(ll k){
    if(k <= 6){
        int small[] = {2, 3, 5, 7, 11, 13};
        return small[k-1];
    }
    // Оценка: p_k ~ k * ln(k)
    double est = (double)k * log((double)k);
    ll lo = (ll)(est * 0.9), hi = (ll)(est * 1.2 + 100);

    // Гарантируем корректность границ
    while(meissner_lehmer(hi) < k) hi *= 2;

    while(lo < hi){
        ll mid = lo + (hi - lo) / 2;
        if(meissner_lehmer(mid) >= k) hi = mid;
        else lo = mid + 1;
    }
    return lo;
}

int main(){
    ll k; cin >> k;
    cout << kth_prime(k) << "\n";
    // k=1: 2, k=10: 29, k=10^6: 15485863
}

Сложность: $O(\log k)$ вызовов $\pi$, каждый за $O(k^{2/3} \log^{2/3} k / \log^2(k \log k))$.

Практически для $k = 10^{12}$ (т.е. $p_k \approx 3 \cdot 10^{13}$) время $\approx 5$–$10$ секунд.

---

Разбор задачи 2 (сложная): $\pi(a, b)$ — число простых в диапазоне

Условие. Даны $a, b \le 10^{14}$, $b - a \le 10^6$. Найти число простых в $[a, b]$.

Решение 1: сегментированное решето — $O((b-a) \log \log b)$. Решение 2: формула — $\pi(b) - \pi(a-1)$ через метод Лемера, $O(b^{2/3})$.

Для маленького промежутка сегментированное решето лучше. Для произвольного диапазона — Лемер.

// Сегментированное решето для [a, b], b - a <= 10^7
ll prime_count_range(ll a, ll b){
    if(a > b) return 0;
    ll range = b - a + 1;
    vector<bool> is_prime(range, true);

    // Обработка a = 1 (1 — не простое)
    if(a == 1) is_prime[0] = false;

    // Решето до sqrt(b)
    ll sq = (ll)sqrtl((long double)b);
    vector<bool> small_sieve(sq + 1, true);
    small_sieve[0] = small_sieve[1] = false;
    for(ll i = 2; i <= sq; i++){
        if(!small_sieve[i]) continue;
        for(ll j = i*i; j <= sq; j += i) small_sieve[j] = false;
        // Вычёркиваем в [a, b]
        ll start = max(i*i, (a + i - 1) / i * i);
        if(start == i) start += i; // не вычёркиваем само p
        for(ll j = start; j <= b; j += i)
            is_prime[j - a] = false;
    }

    ll cnt = 0;
    for(ll v = 0; v < range; v++)
        if(is_prime[v]) ++cnt;
    return cnt;
}

// Для произвольного диапазона
ll prime_count_any_range(ll a, ll b){
    // Используем два вызова pi
    return meissner_lehmer(b) - meissner_lehmer(a - 1);
}

Оптимизация для соревнований: при $b - a \le 10^6$ всегда используйте сегментированное решето — это $O(\sqrt{b} + (b-a) \log \log b)$

​+(b−a)loglogb), что работает за доли секунды.

---

Разбор задачи 3 (очень сложная — ICPC WF уровень): Сумма $\pi(p_i)$ по простым

Условие. Найти $\sum_{p \le N} \pi(p) \pmod{10^9 + 7}$, где $N \le 10^{11}$.

Анализ. Заметим: $\sum_{p \le N} \pi(p) = \sum_{k=1}^{\pi(N)} k = \frac{\pi(N) \cdot (\pi(N) + 1)}{2}$

Нет, это неверно. $\pi(p)$ для $k$-го простого $p_k$ равно $k$. Поэтому: $\sum_{p \le N} \pi(p) = \sum_{k=1}^{\pi(N)} \pi(p_k) = \sum_{k=1}^{\pi(N)} k = \frac{\pi(N)(\pi(N)+1)}{2}$

Это тривиально! Ответ: $\pi(N) \cdot (\pi(N)+1) / 2 \pmod{10^9+7}$.

ll solve(ll N){
    const ll MOD = 1e9 + 7;
    ll pi_N = meissner_lehmer(N) % MOD;
    return pi_N * (pi_N + 1) % MOD * 500000004 % MOD; // 500000004 = inv(2) mod MOD
}

Более нетривиальная версия: $\sum_{p \le N} \pi(\lfloor N/p \rfloor)$.

Это не упрощается тривиально. Решение через суммирование: $\sum_{p \le N} \pi\!\left(\!\left\lfloor \frac{N}{p} \right\rfloor\right) = \sum_{k=1}^{\pi(N)} \pi\!\left(\!\left\lfloor \frac{N}{p_k} \right\rfloor\right)$

Используем факт, что $\lfloor N/p_k \rfloor$ принимает $O(\sqrt{N})$

​) различных значений. Для каждого значения $v = \lfloor N/p_k \rfloor$ посчитаем число простых $p_k$ с таким значением:

ll solve2(ll N){
    const ll MOD = 1e9 + 7;
    ll ans = 0;

    // Перебираем различные значения v = floor(N/p)
    // Для каждого v: число простых p с floor(N/p) = v
    // p в диапазоне (N/(v+1), N/v]
    // Количество простых: pi(N/v) - pi(N/(v+1))

    // Но нам нужно суммировать pi(v) * (число простых p с floor(N/p)=v)
    // Это O(sqrt(N)) шагов, каждый за O(1) с предподсчитанными pi

    // Для N <= 10^11 используем Min-25 для предподсчёта pi всех floor(N/i)
    // (Метод Лемера менее удобен здесь — он не предподсчитывает все floor(N/i))

    // Используем Min25:
    // (код аналогичен главе 14)
    return ans;
}

Для этой задачи Min-25 Sieve удобнее, чем метод Лемера. Это иллюстрирует различие методов.

---

Сравнение с Min-25 Sieve

Вывод:

• Если нужен только $\pi(N)$ и $N \le 10^{11}$: оба подходят, Min-25 чуть быстрее на практике.
• Если нужны суммы мультипликативных функций: только Min-25.
• Если $N \approx 10^{13}$–$10^{14}$ и нужен только $\pi(N)$: предпочтительнее Лемер из-за меньшей памяти.
• Для $N > 10^{14}$: метод Деланги–Ривата или Lucy_Hedgehog с сегментацией.

---

Задачи для самостоятельного решения (таблица)

---

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

1. Переполнение в промежуточных вычислениях. При $x = 10^{14}$ значение $x \cdot (x+1)$ переполняет long long. Используйте __int128 или вычисляйте через модуль, если нужен ответ по модулю.

2. Неправильное вычисление $a = \pi(x^{1/4})$. Из-за ошибок округления $x^{1/4}$ может быть на единицу больше или меньше истинного значения. Всегда проверяйте: $p_a^4 \le x < p_{a+1}^4$.

3. Выход за границы primes[]. В цикле по $i$ от $a+1$ до $b$ убедитесь, что $i < \text{primes.size()}$. При $x \le 4$ простых нет — граничный случай.

4. Некорректная остановка рекурсии для $\phi$. Условие $v \le \text{sieve\_limit}$ позволяет вернуть $1 + \pi(v) - a$. Но это верно только если $a < \pi(v)$; если $a \ge \pi(v)$, то возвращаем 1. Частая ошибка — упустить этот случай.

5. Двойной счёт в $P_2$. В сумме $P_2$ перебираем $p$ до $\sqrt{x}$

​, и для каждого $p$ считаем количество $q \ge p$. Не забыть вычесть $i-1$ (не $i$) при 0-индексации.

---

Совет профессионала

Для задач на $\pi(N)$ при $N \le 10^{11}$ реализуйте Min-25 Sieve Фазу 1 — она проще метода Лемера и имеет то же асимптотическое время. Метод Лемера стоит реализовывать только если требуется $N > 10^{12}$ или крайне ограниченная память. В соревновательном программировании большинство задач укладываются в $N \le 10^{11}$, где Min-25 оптимален по соотношению простота/скорость.

При реализации Лемера используйте кэш для $\phi(v, a)$ при малых $a$: словарь unordered_map<(v, a), result> со значениями $a \le 5$ даёт ускорение в 3–5 раз.

---

Детальный пример работы алгоритма Лемера

Вычислим $\pi(100)$ шаг за шагом.

$x = 100$, $\sqrt{x} = 10$

​=10, $x^{1/4} \approx 3.16$, $x^{1/3} \approx 4.64$.

Простые до $\sqrt{x} = 10$

​=10: $\{2, 3, 5, 7\}$, их четыре: $b = \pi(10) = 4$.

$a = \pi(x^{1/4}) = \pi(3) = 2$ (простые $\le 3.16$: $\{2, 3\}$).

Вычисление $\phi(100, 2)$:

$\phi(100, 2)$ — числа $\le 100$, не кратные 2 и 3.

Через включение-исключение: числа, не кратные 2 и 3. Из 100 чисел: не кратные 2: 50; из них не кратные 3: $\lfloor 50/3 \rfloor$-поправки...

Рекуррентно: $\phi(100, 2) = \phi(100, 1) - \phi(50, 1)$ $\phi(100, 1) = \phi(100, 0) - \phi(50, 0) = 100 - 50 = 50$ $\phi(50, 1) = \phi(50, 0) - \phi(25, 0) = 50 - 25 = 25$ $\phi(100, 2) = 50 - 25 = 25$

Проверка: числа $\le 100$, взаимно простые с 6 ($= 2 \cdot 3$): по формуле $100 \cdot \varphi(6)/6 = 100 \cdot 2/6 \approx 33$... нет, это не то. $\phi(100, 2) = |\{n \le 100 : \gcd(n, 2 \cdot 3) = 1\}|$. Числа от 1 до 6 с таким свойством: $\{1, 5\}$ — два числа. В 100 числах: $\lfloor 100/6 \rfloor \cdot 2 +$ остаток = $16 \cdot 2 + |\{97, 98, 99, 100\} \cap \{1,5\} \pmod 6| = 32 + |\{97 \equiv 1\}| = 33$. Гм, расхождение.

Проблема: формула включения-исключения считает числа с $p \nmid n$ для $p \in \{2, 3\}$, включая 1. Наша рекуррентность считает правильно: $\phi(x, 0) = \lfloor x \rfloor$ (все от 1 до $x$). $\phi(100, 1) = 100 - \lfloor 100/2 \rfloor = 50$. $\phi(100, 2) = \phi(100, 1) - \phi(\lfloor 100/3 \rfloor, 1) = 50 - \phi(33, 1) = 50 - (33 - \lfloor 33/2 \rfloor) = 50 - 17 = 33$. Верно.

(Исправление: я ошибся выше — $\phi(50, 1) = \phi(50, 0) - \phi(25, 0)$? Нет, это $\phi(100, 2) = \phi(100, 1) - \phi(100/3, 1) = \phi(100, 1) - \phi(33, 1) = 50 - 17 = 33$.)

Вычисление $P_2(100, 2)$:

$P_2(100, 2) = \sum_{i=3}^{4} (\pi(100/p_i) - (i-1))$

$i = 3$ ($p_3 = 5$): $\pi(100/5) - 2 = \pi(20) - 2 = 8 - 2 = 6$. $i = 4$ ($p_4 = 7$): $\pi(100/7) - 3 = \pi(14) - 3 = 6 - 3 = 3$.

(Нет $i = 5$: $p_5 = 11 > \sqrt{100} = 10$

​=10.)

$P_2(100, 2) = 6 + 3 = 9$.

Итоговая формула: $\pi(100) = \phi(100, 2) + 2 - 1 - P_2(100, 2) = 33 + 1 - 9 = 25$

Верно: $\pi(100) = 25$.

Проверка $P_2$

$P_2(100, 2)$ считает числа $n = p \cdot q \le 100$ с $p, q > p_2 = 3$ ($p \le q$ простые $> 3$).

$p = 5$: $q$ простое, $5 \le q \le 20$: $\{5, 7, 11, 13, 17, 19\}$ — 6 значений. ✓ $p = 7$: $7 \le q \le 14$: $\{7, 11, 13\}$ — 3 значения. ✓ $p = 11$: $11 \le q \le 9.09$... нет (т.к. $11 > \sqrt{100} = 10$

​=10, значит $p \le \sqrt{100}$).

Итого $P_2 = 9$. Числа: $25, 35, 55, 65, 85, 95, 49, 77, 91$ — девять чисел. ✓

Реализация с мемоизацией $\phi$

Для ускорения кэшируем значения $\phi(v, a)$ при малых $a$ (типично $a \le 5$):

// Мемоизация через плоский хэш
// Ключ: (v * 8 + a) для a <= 7
// Для v <= sqrtx: используем 2D массив small_phi[v][a]
// Для v > sqrtx: используем unordered_map

struct PhiCache {
    int sqx;
    vector<array<ll, 8>> small_phi; // small_phi[v][a] для v <= sqx, a <= 7
    unordered_map<ll, ll> large_phi;

    PhiCache(int sq) : sqx(sq), small_phi(sq + 1) {
        for(auto& arr : small_phi) arr.fill(-1);
    }

    ll& get_small(int v, int a){ return small_phi[v][a]; }

    ll& get_large(ll v, int a){
        ll key = v * 8 + a;
        return large_phi[key];
    }
};

Для $a \le 7$ и $v \le \sqrt{x}$

​ это даёт $O(\sqrt{x})$ кэшированных значений. Для $a \le 7$ и $v > \sqrt{x}$ используем хэш: $v$ принимает $O(\sqrt{x}/p_7) = O(\sqrt{x}/17)$ значений, итого $O(\sqrt{x})$ записей.

Суммарная память кэша: $O(\sqrt{x} \cdot 7)$

​⋅7) — пренебрежимо мало.

---

Итог

Метод Мейсселя–Лемера позволяет вычислять $\pi(x)$ точно за $O(x^{2/3} / \log^2 x)$ времени с памятью $O(x^{1/3})$. Три ключевых компонента:

1. Функция $\phi(x, a)$: рекурсивно считает числа с наименьшим делителем $> p_a$; остановка при попадании в предподсчитанный диапазон.
2. Сумма $P_2(x, a)$: подсчёт чисел с ровно двумя простыми делителями $> p_a$; вычисляется через $\pi(x/p_i)$ в цикле.
3. Формула: $\pi(x) = \phi(x, a) + a - 1 - P_2(x, a)$ при $a = \pi(x^{1/4})$.

Сравнение с Min-25 Sieve показывает, что каждый метод имеет свою нишу: Лемер лучше для больших $N$ и ограниченной памяти, Min-25 — для задач с суммами мультипликативных функций.