【模板】min25 筛

使用条件:

  • 求积性函数单点前缀和。
  • \(p\) 为质数,则 \(f(p)\) 是低次多项式(或者可以拆成常数个完全积性函数,同时它们的前缀和都可以计算)
  • \(p\) 为质数,则要求 \(f(p^c)\) 能在比较低的时间计算。
  • \(n\leq 10^{10}\)(?)

算法流程建议直接看 OI-wiki 而不是我在这里复读一遍。你就找个时间认真地逐行读一下,注意一定是认真地读,搞清楚每个等号的由来。

代码中除了 \(idl, idg\) 记录了 \(x, n/x\) 对应的 id 之外还记了一个 \(idp_k=id(pri_k)\)。筛质数点值前缀和对我的代码负优化,没加,或者可能我写的有问题。\(id\) 的写法也很牛,你需要提前处理根号,不能做太多无用的除法。

点击查看代码

本题 \(f(p^k)=p^k(p^k-1)\)

#include <bits/stdc++.h>
using namespace std;
#ifdef LOCAL
#define debug(...) fprintf(stderr, ##__VA_ARGS__)
#else
#define endl "\n"
#define debug(...) void(0)
#endif
using LL = long long;
template <class T>
using must_int = enable_if_t<is_integral<T>::value, int>;
template <unsigned umod>
struct modint {/*{{{*/
  static constexpr int mod = umod;
  unsigned v;
  modint() : v(0) {}
  template <class T, must_int<T> = 0>
  modint(T x) {
    x %= mod;
    v = x < 0 ? x + mod : x;
  }
  modint operator+() const { return *this; }
  modint operator-() const { return modint() - *this; }
  friend int raw(const modint &self) { return self.v; }
  friend ostream& operator<<(ostream& os, const modint &self) { 
    return os << raw(self); 
  }
  modint &operator+=(const modint &rhs) {
    v += rhs.v;
    if (v >= umod) v -= umod;
    return *this;
  }
  modint &operator-=(const modint &rhs) {
    v -= rhs.v;
    if (v >= umod) v += umod;
    return *this;
  }
  modint &operator*=(const modint &rhs) {
    v = 1ull * v * rhs.v % umod;
    return *this;
  }
  modint &operator/=(const modint &rhs) {
    assert(rhs.v);
    return *this *= qpow(rhs, mod - 2);
  }
  template <class T, must_int<T> = 0>
  friend modint qpow(modint a, T b) {
    modint r = 1;
    for (; b; b >>= 1, a *= a)
      if (b & 1) r *= a;
    return r;
  }
  friend modint operator+(modint lhs, const modint &rhs) { return lhs += rhs; }
  friend modint operator-(modint lhs, const modint &rhs) { return lhs -= rhs; }
  friend modint operator*(modint lhs, const modint &rhs) { return lhs *= rhs; }
  friend modint operator/(modint lhs, const modint &rhs) { return lhs /= rhs; }
  bool operator==(const modint &rhs) const { return v == rhs.v; }
  bool operator!=(const modint &rhs) const { return v != rhs.v; }
};/*}}}*/
using mint = modint<int(1e9 + 7)>;
LL sqr(LL x) { return x * x; }
template <int N>
struct siever {
  int pri[N + 10], cnt;
  siever() : cnt(0) {
    fill(pri + 1, pri + N + 1, true);
    pri[0] = 1, pri[1] = false;
    for (int i = 1; i <= N; i++) {
      if (pri[i]) pri[++cnt] = i;
      for (int j = 1; j <= cnt && i * pri[j] <= N; j++) {
        pri[i * pri[j]] = false;
        if (i % pri[j] == 0) break;
      }
    }
  }
};
mint ID1(mint x) { return x * (x + 1) * 500000004; }
mint ID2(mint x) { return x * (x + 1) * (2 * x + 1) * 166666668; }
LL n, v[200010], lim;
mint f(mint pc) { return pc * (pc - 1); }
int m, idl[100010], idg[100010], idp[100010];
int id(LL x) { return x <= lim ? idl[x] : idg[n / x]; }
mint g1[200010], g2[200010], *g = g2;
siever<100010> sie;
mint dfs(LL n, int k) {
  if (n < sie.pri[k]) return 0;
  mint ans = g[id(n)] - g[idp[k - 1]];
  for (int i = k; sqr(sie.pri[i]) <= n; i++) {
    int pi = sie.pri[i];
    LL now = pi;
    for (int c = 1; now * pi <= n; c++, now *= pi) {
      ans += f(now) * dfs(n / now, i + 1) + f(now * pi);
    }
  }
  return ans;
}
int main() {
#ifndef LOCAL
  cin.tie(nullptr)->sync_with_stdio(false);  
#endif
  cin >> n;
  lim = sqrtl(n);
  for (LL l = 1, r; l <= n; l = r + 1) {
    v[++m] = n / l, r = n / v[m];
    (v[m] <= r ? idl[v[m]] : idg[r]) = m;
    g1[m] = ID1(v[m]) - 1, g2[m] = ID2(v[m]) - 1;
  }
  idp[0] = id(1);
  for (int k = 1; sqr(sie.pri[k]) <= n; k++) {
    idp[k] = id(sie.pri[k]);
    int p = sie.pri[k];
    for (int i = 1; i <= m && v[i] >= sqr(p); i++) {
      g1[i] -= p * (g1[id(v[i] / p)] - g1[idp[k - 1]]);
      g2[i] -= sqr(p) * (g2[id(v[i] / p)] - g2[idp[k - 1]]);
    }
  }
  for (int i = 1; i <= m; i++) g[i] -= g1[i];
  cout << dfs(n, 1) + 1 << endl;
  return 0;
}

LOJ6053【简单的函数】:这道题就是说 \(p\) 为质数时 \(f(p)=p\oplus 1=p+(-1)^{[p\neq 2]}\),就是说 \(f(2)=3\) 而其它 \(f(p)=p-1\)\(f(2)\) 要特判。这个特判要判的东西其实很少,首先 \(g\) 的部分根本不需要管,然后最后用 \(g\) 的时候 mint ans = g[id(n)] - g[idp[k - 1]]; 再把 \(2\) 的多出来一点点的答案加进去,实际上根本没啥要改的。

可以认为就是:\(g\) 部分就是算 lambda x: x ** k 这个函数的质数点处前缀和,然后你把它组合一下让它摇身一变变成原函数的质数点处前缀和。然后再暴力 dfs,去把每个数搜出来,按照质因子从小到大的顺序,对于质数用刚才的 \(g\) 直接前缀和求出,剩下的合数让它递归去。

posted @ 2024-09-10 21:09  caijianhong  阅读(24)  评论(0编辑  收藏  举报