BZOJ4916: 神犇和蒟蒻【杜教篩】

Description

好久好久之前，有一隻神犇叫yzy;

好久好久以後，有一隻蒟蒻叫lty;

Input

請你讀入一個整數N;1<=N<=1E9,A、B模1E9+7;

Output

請你輸出一個整數A=\sum_{i=1}^N{\mu (i^2)};

請你輸出一個整數B=\sum_{i=1}^N{\varphi (i^2)};

Sample Input

1

Sample Output

1
1

---

思路

首先發現第一個必定是1.。。。

而後發現第二個其實能夠表示成
\[ \sum_{i = 1}^n\phi(i)*i \]
而後咱們令
\[ f(i)=\phi(i)*i \\ g(i)=i \]
那麼能夠獲得
\[ ans=Sum(n)=\sum_{i=1}^nf(i) \]
又由於
\[ \sum_{i = 1}^n\sum_{d|i}f(d)g(\frac{i}{d})=\sum_{i=1}^n i^2=\frac{n *(n + 1)*(2n+1)}{6} \]
且
\[ \sum_{i = 1}^n\sum_{d|i}f(d)g(\frac{i}{d})=\sum_{k = 1}^ng(k)\sum_{d = 1}^{\lfloor\frac{n}{k}\rfloor}f(d)=\sum_{k = 1}^ng(k)Sum(\lfloor\frac{n}{k}\rfloor) \]
因此有
\[ Sum(n)=\frac{n *(n + 1)*(2n+1)}{6}-\sum_{k = 2}^ng(k)Sum(\lfloor\frac{n}{k}\rfloor) \]
而後上杜教篩板子。。

---

#include<bits/stdc++.h>

using namespace std;

typedef long long ll;

const ll Mod = 1e9 + 7;
const ll N = 1e7 + 10;
const ll inv6 = 166666668;
const ll inv2 = 500000004;

ll prime[N], cnt = 0;
ll phi[N], sum[N], vis[N];

map<ll, ll> mp;

ll add(ll a, ll b) {
  return (a += b) >= Mod ? a - Mod : a;
}

ll sub(ll a, ll b) {
  return (a -= b) < 0 ? a + Mod : a;
}

ll mul(ll a, ll b) {
  return a * b % Mod;
}

void get_prime() {
  phi[1] = 1;
  for (ll i = 2; i < N; i++) {
    if (!vis[i]) {
      phi[i] = i - 1;
      prime[++cnt] = i;
    }
    for (ll j = 1; j <= cnt && i * prime[j] < N; j++) {
      vis[i * prime[j]] = 1;
      if (i % prime[j] == 0) {
        phi[i * prime[j]] = phi[i] * prime[j];
        break;
      } else {
        phi[i * prime[j]] = phi[i] * (prime[j] - 1);
      }
    }
  }
  for (ll i = 1; i < N; i++)
    sum[i] = add(sum[i - 1], mul(i, phi[i]));
}

ll solve(ll n) {
  if (n < N) return sum[n];
  if (mp.count(n)) return mp[n];
  ll res = mul(mul(n, n + 1), mul(2 * n + 1, inv6));
  for (ll i = 2; i <= n; i++) {
    ll j = n / (n / i);
    res = sub(res, mul(solve(n / i), mul(inv2, mul(i + j, j - i + 1))));
    i = j;
  }
  return mp[n] = res;
}

int main() {
  get_prime();
  ll n; cin >> n;
  cout << 1 << "\n" << solve(n);
  return 0;
}