BZOJ 4916 神犇和蒟蒻(杜教筛)

Description

很久很久以前，有一只神犇叫yzy;
很久很久之后，有一只蒟蒻叫lty;

---

Input

请你读入一个整数N;1<=N<=1E9,A、B模1E9+7;

---

Output

请你输出一个整数A=∑Ni=1μ(i2);
请你输出一个整数B=∑Ni=1φ(i2);

---

Sample Input

1

---

Sample Output

1
1

---

Solution

A=∑Ni=1μ(i2)是拿来搞笑的。
B=∑Ni=1φ(i2)化简一下，用杜教筛即可。

因为

φ(n)=n⋅Πpi−1pi，φ(n2)=n2⋅Πpi−1pi

其中

n=pk11pk22..pkmm,n2=p2k11p2k22..p2kmm

所以

φ(n2)=n⋅φ(n)

题目就是求

∑i=1nφ(i2)=∑i=1ni⋅φ(i)

这个利用狄利克雷卷积知识就可以知道，(id⋅φ)∗id=id2，
即

∑d|nd⋅φ(d)⋅nd=n∑d|nφ(d)=n2

我们设ϕ′(n)=∑ni=1i⋅φ(i)，利用平方前缀和公式，我们得到

n(n+1)(2n+1)6=∑i=1ni2=∑i=1n∑d|id⋅φ(d)⋅id=∑id=1id⋅∑d=1⎢⎣⎢⎢nid⎥⎦⎥⎥d⋅φ(d)=∑i=1ni⋅ϕ′(⌊ni⌋)

因此ϕ′(n)=n(n+1)(2n+1)6−∑ni=2i⋅ϕ′(⌊ni⌋)，预处理一下，用hash记忆化递归，时间就可以保证在O(n23)了。

---

Code

#include <cstdio>#include <cstdlib>#include <cstring>#include <iostream>#include <cmath>#include <algorithm>#define N 1000010#define M 2333333#define Mod 1000000007using namespace std;typedef long long LL;int n, cnt;bool vis[N];int prime[N];LL phi[N], hash_v[M], hash_n[M];void Da(){    vis[1] = true;    phi[1] = 1ll;    for(int i = 2; i < N; i++){      if(!vis[i]){        prime[++cnt] = i;        phi[i] = (LL)(i-1);      }      for(int j = 1; j <= cnt && i * prime[j] < N; j++){        vis[i * prime[j]] = true;        if(i % prime[j] == 0){          phi[i * prime[j]] = phi[i] * (LL)prime[j];          break;        }        else  phi[i * prime[j]] = phi[i] * (LL)(prime[j]-1);      }    }    for(int i = 1; i < N; i++){        phi[i] = phi[i] % Mod * i % Mod;      phi[i] = (phi[i] + phi[i-1]) % Mod;    }}LL Find(LL x){    int p = x % M;    while(hash_n[p] && hash_n[p] != x)  p = (p+1) % M;    if(!hash_n[p])  return -1;    return hash_v[p];}void Push(LL v, LL x){    int p = x % M;    while(hash_n[p])  p = (p+1) % M;    hash_v[p] = v;    hash_n[p] = x;}LL Sum(LL x){    if(x < N)  return phi[x];    LL res = Find(x);    if(~ res)  return res;    else{      LL temp1 = x * (x + 1), temp2 = 2 * x + 1;      if(!(temp1 % 6))  temp1 /= 6;      else if(!(temp2 % 6))  temp2 /= 6;      else if(!(temp1 % 2) && !(temp2 % 3))  temp1 /= 2, temp2 /= 3;      else if(!(temp2 % 2) && !(temp1 % 3))  temp1 /= 3, temp2 /= 2;      res = ((temp1 % Mod) * (temp2 % Mod)) % Mod;    }    LL last;    for(LL i = 2; i <= x; i = last+1){      last = x/(x/i);      res = (res - (((last+i) * (last-i+1) / 2) % Mod * Sum(x/i)) % Mod + Mod) % Mod;    }    Push(res, x);    return res;}int main(){    freopen("bzoj4916.in", "r", stdin);    freopen("bzoj4916.out", "w", stdout);    Da();    scanf("%lld", &n);    printf("1\n%lld\n", Sum(n));    return 0;}

---