NicheLibrary

This documentation is automatically generated by NotLeonian/competitive-verifier (forked from competitive-verifier/competitive-verifier)

View the Project on GitHub NotLeonian/NicheLibrary

:heavy_check_mark: $m$ で割った余りごとの $N$ 以下の素数の個数の列挙 (math/multiplicative-function/prime-counting-modulo.hpp)

概要

  • $N$ 以下の素数を、$m$ で割った余りごとに数える。
  • Lucy DP のテーブルを余りごとに持つ。
  • prime_counting_modulo_mf_prefix_sum_table<T> は、各余りについて Black Algorithm 用の Fprime を返す。
  • $m$ が合成数でも特別扱いしない。

使い方

  • prime_counting_modulo_table(N, m)
    • pair<vector<long long>, vector<vector<long long>>> を返す。
    • first[i] はテーブルの値 $x$ である。
    • second[k][i] は、$x$ 以下の素数で $m$ で割った余りが $k$ であるものの個数を返す。
    • 前提: $N\ge 0,\;m>0$。
    • 備考: second の 1 つ目の添字は余りである。
  • prime_counting_modulo(N, m)
    • 長さ $m$ の vector<long long> を返す。
    • 返り値の $k$ 番目は、$N$ 以下の素数で $m$ で割った余りが $k$ であるものの個数である。
    • 前提: $N\ge 0,\;m>0$。
  • prime_counting_modulo_mf_prefix_sum_table<T>(N, m)
    • 長さ $m$ の vector<vector<T>> を返す。
    • 返り値の $k$ 番目は、$m$ で割った余りが $k$ である素数だけを対象にする Fprime である。
    • 各要素は long long から T に変換される。
    • 前提: $N\ge 0,\;m>0$。
    • 備考: 複数の余りをまとめる場合は、返り値をユーザー側で足す。
    • 備考: $N=0$ では各行は空である。

計算量

  • 時間: $O(m N^{3/4}/\log N)$
  • 空間: $O(m\sqrt{N})$

Verified with

Code

#ifndef MATH_MULTIPLICATIVE_FUNCTION_PRIME_COUNTING_MODULO_HPP
#define MATH_MULTIPLICATIVE_FUNCTION_PRIME_COUNTING_MODULO_HPP

// N 以下の素数を、m で割った余りごとに数える。
// Lucy DP のテーブルを余りごとに持ち、素数 x によるふるいを同時に行う。
// N は非負、m は正を仮定する。
// 戻り値のテーブルの 1 つ目の添字は m で割った余りである。
// 計算量 O(m N^{3/4} / log N)、空間 O(m sqrt(N))。

#include <cassert>
#include <utility>
#include <vector>

namespace prime_counting_modulo_internal {
inline long long integer_sqrt(long long n) {
    assert(n >= 0);
    long long ok = 0, ng = 1;
    while (ng <= n / ng) {
        ng <<= 1;
    }
    while (ng - ok > 1) {
        const long long mid = ok + (ng - ok) / 2;
        if (mid <= n / mid) {
            ok = mid;
        } else {
            ng = mid;
        }
    }
    return ok;
}

inline long long count_residue_2_to_n(long long n, long long m, long long r) {
    assert(n >= 0);
    assert(m > 0);
    assert(0 <= r && r < m);
    long long res = 0;
    if (r <= n) {
        res = (n - r) / m + 1;
    }
    if (r == 0) {
        --res;
    }
    if (n >= 1 && r == 1 % m) {
        --res;
    }
    return res;
}

inline void add_mod(long long &x, long long a, long long m) {
    assert(0 <= x && x < m);
    assert(0 <= a && a < m);
    if (a == 0) {
        return;
    }
    if (x >= m - a) {
        x -= m - a;
    } else {
        x += a;
    }
}
} // namespace prime_counting_modulo_internal

inline std::pair<std::vector<long long>, std::vector<std::vector<long long>>>
prime_counting_modulo_table(long long N, long long m) {
    assert(N >= 0);
    assert(m > 0);
    using i64 = long long;
    std::vector<i64> ns{0};
    for (i64 i = N; i > 0;) {
        ns.push_back(i);
        const i64 q = N / i;
        if (q == N) {
            break;
        }
        i = N / (q + 1);
    }
    const i64 sq = prime_counting_modulo_internal::integer_sqrt(N);
    const i64 nsz = static_cast<i64>(ns.size());
    std::vector<std::vector<i64>> h(m, std::vector<i64>(nsz));
    for (i64 r = 0; r < m; ++r) {
        for (i64 i = 0; i < nsz; ++i) {
            h[r][i] = prime_counting_modulo_internal::count_residue_2_to_n(
                ns[i], m, r);
        }
    }
    for (i64 x = 2; x <= sq; ++x) {
        const i64 x_mod = x % m;
        const i64 x_idx = nsz - x;
        const i64 prev_idx = nsz - x + 1;
        if (h[x_mod][x_idx] == h[x_mod][prev_idx]) {
            continue;
        }
        const i64 x2 = x * x;
        for (i64 i = 1; i < nsz && ns[i] >= x2; ++i) {
            const i64 n = ns[i];
            const i64 q = n / x;
            const i64 q_idx = i <= sq / x ? i * x : nsz - q;
            i64 to = 0;
            for (i64 r = 0; r < m; ++r) {
                h[to][i] -= h[r][q_idx] - h[r][prev_idx];
                prime_counting_modulo_internal::add_mod(to, x_mod, m);
            }
        }
    }
    return {ns, h};
}

inline std::vector<long long> prime_counting_modulo(long long N, long long m) {
    assert(N >= 0);
    assert(m > 0);
    std::vector<long long> res(m);
    if (N == 0) {
        return res;
    }
    const auto table = prime_counting_modulo_table(N, m).second;
    for (long long r = 0; r < m; ++r) {
        res[r] = table[r][1];
    }
    return res;
}

template <class T>
std::vector<std::vector<T>>
prime_counting_modulo_mf_prefix_sum_table(long long N, long long m) {
    assert(N >= 0);
    assert(m > 0);
    std::vector<std::vector<T>> res(m);
    if (N == 0) {
        return res;
    }
    const auto table = prime_counting_modulo_table(N, m).second;
    for (long long r = 0; r < m; ++r) {
        res[r].resize(table[r].size());
        for (long long i = 0; i < static_cast<long long>(table[r].size());
             ++i) {
            res[r][i] = static_cast<T>(table[r][i]);
        }
    }
    return res;
}

#endif
#line 1 "math/multiplicative-function/prime-counting-modulo.hpp"



// N 以下の素数を、m で割った余りごとに数える。
// Lucy DP のテーブルを余りごとに持ち、素数 x によるふるいを同時に行う。
// N は非負、m は正を仮定する。
// 戻り値のテーブルの 1 つ目の添字は m で割った余りである。
// 計算量 O(m N^{3/4} / log N)、空間 O(m sqrt(N))。

#include <cassert>
#include <utility>
#include <vector>

namespace prime_counting_modulo_internal {
inline long long integer_sqrt(long long n) {
    assert(n >= 0);
    long long ok = 0, ng = 1;
    while (ng <= n / ng) {
        ng <<= 1;
    }
    while (ng - ok > 1) {
        const long long mid = ok + (ng - ok) / 2;
        if (mid <= n / mid) {
            ok = mid;
        } else {
            ng = mid;
        }
    }
    return ok;
}

inline long long count_residue_2_to_n(long long n, long long m, long long r) {
    assert(n >= 0);
    assert(m > 0);
    assert(0 <= r && r < m);
    long long res = 0;
    if (r <= n) {
        res = (n - r) / m + 1;
    }
    if (r == 0) {
        --res;
    }
    if (n >= 1 && r == 1 % m) {
        --res;
    }
    return res;
}

inline void add_mod(long long &x, long long a, long long m) {
    assert(0 <= x && x < m);
    assert(0 <= a && a < m);
    if (a == 0) {
        return;
    }
    if (x >= m - a) {
        x -= m - a;
    } else {
        x += a;
    }
}
} // namespace prime_counting_modulo_internal

inline std::pair<std::vector<long long>, std::vector<std::vector<long long>>>
prime_counting_modulo_table(long long N, long long m) {
    assert(N >= 0);
    assert(m > 0);
    using i64 = long long;
    std::vector<i64> ns{0};
    for (i64 i = N; i > 0;) {
        ns.push_back(i);
        const i64 q = N / i;
        if (q == N) {
            break;
        }
        i = N / (q + 1);
    }
    const i64 sq = prime_counting_modulo_internal::integer_sqrt(N);
    const i64 nsz = static_cast<i64>(ns.size());
    std::vector<std::vector<i64>> h(m, std::vector<i64>(nsz));
    for (i64 r = 0; r < m; ++r) {
        for (i64 i = 0; i < nsz; ++i) {
            h[r][i] = prime_counting_modulo_internal::count_residue_2_to_n(
                ns[i], m, r);
        }
    }
    for (i64 x = 2; x <= sq; ++x) {
        const i64 x_mod = x % m;
        const i64 x_idx = nsz - x;
        const i64 prev_idx = nsz - x + 1;
        if (h[x_mod][x_idx] == h[x_mod][prev_idx]) {
            continue;
        }
        const i64 x2 = x * x;
        for (i64 i = 1; i < nsz && ns[i] >= x2; ++i) {
            const i64 n = ns[i];
            const i64 q = n / x;
            const i64 q_idx = i <= sq / x ? i * x : nsz - q;
            i64 to = 0;
            for (i64 r = 0; r < m; ++r) {
                h[to][i] -= h[r][q_idx] - h[r][prev_idx];
                prime_counting_modulo_internal::add_mod(to, x_mod, m);
            }
        }
    }
    return {ns, h};
}

inline std::vector<long long> prime_counting_modulo(long long N, long long m) {
    assert(N >= 0);
    assert(m > 0);
    std::vector<long long> res(m);
    if (N == 0) {
        return res;
    }
    const auto table = prime_counting_modulo_table(N, m).second;
    for (long long r = 0; r < m; ++r) {
        res[r] = table[r][1];
    }
    return res;
}

template <class T>
std::vector<std::vector<T>>
prime_counting_modulo_mf_prefix_sum_table(long long N, long long m) {
    assert(N >= 0);
    assert(m > 0);
    std::vector<std::vector<T>> res(m);
    if (N == 0) {
        return res;
    }
    const auto table = prime_counting_modulo_table(N, m).second;
    for (long long r = 0; r < m; ++r) {
        res[r].resize(table[r].size());
        for (long long i = 0; i < static_cast<long long>(table[r].size());
             ++i) {
            res[r][i] = static_cast<T>(table[r][i]);
        }
    }
    return res;
}


Back to top page