tests/Counting_Primes.test.cpp

View this file on GitHub
Last update: 2025-11-28 10:18:48+07:00
Problem: https://judge.yosupo.jp/problem/counting_primes

Depends on

Code

#define PROBLEM "https://judge.yosupo.jp/problem/counting_primes"

#include "../misc/macros.h"
#include "../math/Min25.h"

void solve() {
  Min25<i64> solver;
  i64 n;
  cin >> n;
  solver.init(n);
  cout << solver.g0[solver.id(n)];
}

signed main() {
  ios::sync_with_stdio(false);
  cin.tie(0);
  int tc = 1;
  // cin >> tc;
  while (tc--) solve();
}

#line 1 "tests/Counting_Primes.test.cpp"
#define PROBLEM "https://judge.yosupo.jp/problem/counting_primes"

#line 1 "misc/macros.h"
// #pragma GCC optimize("Ofast,unroll-loops")       // unroll long, simple loops
// #pragma GCC target("avx2,fma")                   // vectorizing code
// #pragma GCC target("lzcnt,popcnt,abm,bmi,bmi2")  // for fast bitset operation

#include <bits/extc++.h>
#include <tr2/dynamic_bitset>

using namespace std;
using namespace __gnu_pbds;  // ordered_set, gp_hash_table
// using namespace __gnu_cxx; // rope

// for templates to work
#define all(x) (x).begin(), (x).end()
#define sz(x) (int) (x).size()
#define pb push_back
#define eb emplace_back
using i32 = int32_t;
using u32 = uint32_t;
using i64 = int64_t;
using u64 = uint64_t;
using i128 = __int128_t;
using u128 = __uint128_t;
using ld = long double;
using pii = pair<i32, i32>;
using vi = vector<i32>;

// fast map
const int RANDOM = chrono::high_resolution_clock::now().time_since_epoch().count();
struct chash {  // customize hash function for gp_hash_table
  int operator()(int x) const { return x ^ RANDOM; }
};
gp_hash_table<int, int, chash> table;

/* ordered set
    find_by_order(k): returns an iterator to the k-th element (0-based)
    order_of_key(k): returns the number of elements in the set that are strictly less than k
*/
template <class T>
using ordered_set = tree<T, null_type, less<T>, rb_tree_tag, tree_order_statistics_node_update>;

/*  rope
    rope <int> cur = v.substr(l, r - l + 1);
    v.erase(l, r - l + 1);
    v.insert(v.mutable_begin(), cur);
*/
#line 1 "math/Min25.h"
template <class T>
struct Min25 {
  i64 n;
  int sq;
  vector<int> primes, id1, id2;
  vector<i64> vals;
  vector<T> g0, g1;  // g0: sum p^0, g1: sum p^1
  int id(i64 x) { return x <= sq ? id1[x] : id2[n / x]; }
  void init(i64 N) {
    n = N, sq = sqrt(n);
    primes.clear();
    vector<bool> is_p(sq + 1, true);
    for (int i = 2; i <= sq; ++i) {
      if (is_p[i]) {
        primes.pb(i);
        for (int j = i * 2; j <= sq; j += i) is_p[j] = false;
      }
    }
    vals.clear(), id1.assign(sq + 1, 0), id2.assign(sq + 1, 0);
    for (i64 l = 1, r; l <= n; l = r + 1) {
      i64 v = n / l;
      r = n / v;
      vals.pb(v);
      if (v <= sq) id1[v] = sz(vals) - 1;
      else id2[n / v] = sz(vals) - 1;
    }
    g0.resize(sz(vals)), g1.resize(sz(vals));
    T inv2 = T(1) / T(2);
    for (int i = 0; i < sz(vals); ++i) {
      T v = T(vals[i]);
      g0[i] = v - 1;
      g1[i] = v * (v + 1) * inv2 - 1;
    }
    for (int p : primes) {
      T sp0 = g0[id(p - 1)], sp1 = g1[id(p - 1)];
      i64 p2 = (i64) p * p;
      T tp = T(p);
      for (int i = 0; i < sz(vals); ++i) {
        if (vals[i] < p2) break;
        int k = id(vals[i] / p);
        g0[i] -= g0[k] - sp0;
        g1[i] -= tp * (g1[k] - sp1);
      }
    }
  }
  // A, B: f(p) = A*1 + B*p
  // func: (p, e) -> f(p^e) trả về T
  template <class Func>
  T solve(T A, T B, Func f_pe) {
    vector<T> s_fp(sz(primes) + 1);
    for (int i = 0; i < sz(primes); ++i)
      s_fp[i + 1] = s_fp[i] + A + B * T(primes[i]);

    auto S = [&](auto&& self, i64 x, int j) -> T {
      if (x <= 1 || (j < sz(primes) && primes[j] > x)) return 0;
      int k = id(x);
      T ans = A * g0[k] + B * g1[k];
      ans -= s_fp[j];
      for (int i = j; i < sz(primes); ++i) {
        i64 p = primes[i];
        if (p * p > x) break;
        i64 pe = p;
        for (int e = 1; pe * p <= x; ++e) {
          ans += f_pe(p, e) * self(self, x / pe, i + 1);
          ans += f_pe(p, e + 1);
          pe *= p;
        }
      }
      return ans;
    };
    return S(S, n, 0) + 1;
  }
};
#line 5 "tests/Counting_Primes.test.cpp"

void solve() {
  Min25<i64> solver;
  i64 n;
  cin >> n;
  solver.init(n);
  cout << solver.g0[solver.id(n)];
}

signed main() {
  ios::sync_with_stdio(false);
  cin.tie(0);
  int tc = 1;
  // cin >> tc;
  while (tc--) solve();
}