test/pcounter.test.cpp

• View this file on GitHub
• Last update: 2025-06-14 01:07:36+09:00
• Problem: https://judge.yosupo.jp/problem/counting_primes

Depends on

• k乗根(整数) (C++/math/kthrooti.hpp)
• 素数の個数 (C++/math/primecounter.hpp)
• Sieve of Eratosthenes (C++/math/primetable.hpp)

Code

#define PROBLEM "https://judge.yosupo.jp/problem/counting_primes"
#include <iostream>
#include "C++/math/primecounter.hpp"
int main() {
    int64_t n;
    std::cin >> n;
    std::cout << man::p_count(int64_t(1e11)).pi(n) << '\n';
}

#line 1 "test/pcounter.test.cpp"
#define PROBLEM "https://judge.yosupo.jp/problem/counting_primes"
#include <iostream>
#line 2 "C++/math/primecounter.hpp"

#include <vector>
#line 2 "C++/math/kthrooti.hpp"

#include <cstdint>
#include <limits>
#include <ranges>
#ifndef TEMPLATE
namespace man {
template <std::integral T, std::integral U> constexpr inline bool overflow_if_mul(const T a, const U b) noexcept { return (std::numeric_limits<T>::max()/a)<b; }
}
#endif
namespace man {
inline unsigned long long kthrooti(const unsigned long long n, const int k) {
    if(k == 1) {
		return n;
	}
	const auto chk = [&](const unsigned x) -> bool {
		unsigned long long mul = 1;
		for([[maybe_unused]] const auto _: std::views::iota(0, k)) {
            if(man::overflow_if_mul(mul, x)) {
                return false;
            }
            mul *= x;
        }
		return mul <= n;
	};
	unsigned long long ret = 0;
	for(const auto i: std::views::iota(0, 32) | std::views::reverse) {
		if(chk(ret | (1U << i))) {
			ret |= 1U << i;
		}
	}
	return ret;
}
}

/**
 * @brief k乗根(整数)
 */
#line 2 "C++/math/primetable.hpp"

#line 5 "C++/math/primetable.hpp"
namespace man {
struct p_table {
    std::vector<int> SoE;
    p_table(const int n): SoE(n + 1, 1) {
        SoE[0] = SoE[1] = 0;
        for(const long long i: std::views::iota(2, n + 1)) {
            if(!SoE[i]) {
                continue;
            }
            for(long long j = i * i; j <= n; j += i) {
                SoE[j] = 0;
            }
        }
    }
    std::vector<int> get() {
        std::vector<int> p;
        for(const auto i: std::views::iota(2, std::ssize(SoE))) {
            if(SoE[i]) {
                p.emplace_back(i);
            }
        }
        return p;
    }
};
}

/**
 * @brief Sieve of Eratosthenes
 */
#line 6 "C++/math/primecounter.hpp"
#ifndef TEMPLATE
namespace man {
template <class T> constexpr inline T sqr(const T x) noexcept { return x * x; }
template <class T> constexpr inline T cub(const T x) noexcept { return x * x * x; }
}
#endif
namespace man {
struct p_count {
private:
    long long sq;
    std::vector<int> prime;
    std::vector<long long> prime_sum, primes;
    constexpr inline long long p2(const long long x, const long long y) noexcept {
        if(x < 4) {
            return 0;
        }
        const long long a = pi(y), b = pi(kthrooti(x, 2));
        if(a >= b) {
            return 0;
        }
        long long sum = (a - 2) * (a + 1) / 2 - (b - 2) * (b + 1) / 2;
        for(const auto i: std::views::iota(a, b)) {
            sum += pi(x / primes[i]);
        }
        return sum;
    }
    constexpr inline long long phi(const long long m, const long long n) noexcept {
        if(m < 1) {
            return 0;
        }
        if(n > m) {
            return 1;
        }
        if(n < 1) {
            return m;
        }
        if(m <= sqr(primes[n - 1])) {
            return pi(m) - n + 1;
        }
        if(m <= cub(primes[n - 1]) && m <= sq) {
            const long long sx = pi(kthrooti(m, 2));
            long long ans = pi(m) - (sx + n - 2) * (sx - n + 1) / 2;
            for(const auto i: std::views::iota(n, sx)) {
                ans += pi(m / primes[i]);
            }
            return ans;
        }
        return phi(m, n - 1) - phi(m / primes[n - 1], n - 1);
    }
public:
    p_count(const long long lim): sq(kthrooti(lim, 2)), prime_sum(sq + 1) {
        prime = p_table(sq).SoE;
        for(const auto i: std::views::iota(1) | std::views::take(sq)) {
            prime_sum[i] = prime_sum[i - 1] + prime[i];
        }
        primes.reserve(prime_sum[sq]);
        for(const auto i: std::views::iota(1) | std::views::take(sq)) {
            if(prime[i]) {
                primes.emplace_back(i);
            }
        }
    }
    inline long long pi(const long long n) noexcept {
        if(n <= sq) {
            return prime_sum[n];
        }
        const long long m = kthrooti(n, 3);
        const long long a = pi(m);
        return phi(n, a) + a - 1 - p2(n, m);
    }
};
}

/**
 * @brief 素数の個数
 */
#line 4 "test/pcounter.test.cpp"
int main() {
    int64_t n;
    std::cin >> n;
    std::cout << man::p_count(int64_t(1e11)).pi(n) << '\n';
}

Back to top page