BZOJ 4916: 神犇和蒟蒻 杜教筛 数学

我来写个题解造（骗）福（访）世（问）人（量）

第一问不会的出门左转百度μ是啥去

第二问的话显然答案是等于∑ni=1iφ(i)的，不知道的出门左转百度φ计算公式去……

然后考虑那个东西我们用杜教筛搞一下

考虑f(x)=x和g(x)=xφ(x)的狄利克雷卷积G(x)=∑d|xf(xd)g(d)

首先我们可以直接带入一下

G(t)=∑tx=1∑d|xxddφ(d)=∑d|xxφ(d)

G(t)=∑tx=1x∑d|xφ(d)

注意到后面的东西其实是x

所以其实G(t)=∑tx=1x2=t(t+1)(2t+1)6

然后考虑我们不带入然后算贡献

G(t)=∑tx=1∑d|xf(d)g(xd)

=∑tx=1f(x)∑nx|dg(dx)这步的话其实就是相当于我枚举了一个数和乘积然后算另外一个的贡献

=∑tx=1∑ndk=1g(k)

考虑当x=1的时候后面的东西就是g(n)即为所求，剩下的根号分块递归处理即可

最后令S(x)=∑ni=1iφ(i)

那么会有S(n)=G(i)−∑ni=2iS(ni)

手写了hash表并没有变快qwq

另附n=109时的答案为631887325

// 631887325
#include <map>
#include <stdio.h>
#include <string.h>
#include <iostream>
#include <algorithm>
using namespace std;
typedef long long LL;
const int N = 5e6+5;
const int Max = 5e6;
const int inv2 = 500000004;
const int inv6 = 166666668;
const int mod = 1e9+7;
const int mod1 = 5e6+5;
const int M = 5 * mod1;
int prime
, cnt, phi
, head[mod1+3];
bool F
;
map <int,int> mp;

struct graph {
int next, to, val;
graph () {}
graph (int _next,int _to,int _val)
:next(_next),to(_to), val(_val){}
}edge[M];

inline void add(int x,int y) {
int temp = x % mod1;
edge[++cnt] = graph(head[temp], x, y);
head[temp] = cnt;
}

inline int find(int x) {
int temp = x % mod1;
for(int i=head[temp];i;i=edge[i].next)
if(edge[i].to == x) return edge[i].val;
return -1;
}

inline void init() {
phi[1] = 1;
register int i, j;
for(i=2;i<=Max;++i) {
if(!F[i]) {
prime[++cnt] = i;
phi[i] = i-1;
}
for(j=1;prime[j]*i<=Max;++j) {
F[i*prime[j]] = 1;
if(i%prime[j]==0) {
phi[i*prime[j]] = phi[i] * prime[j];
break;
}
phi[i*prime[j]] = phi[i] * (prime[j]-1);
}
}
for(i=1;i<=Max;++i) phi[i] = (long long) i*phi[i]%mod;
for(i=1;i<=Max;++i) (phi[i] += phi[i-1]) %= mod;
}

inline int calc_3(int x) {
return (long long)x * (x+1) % mod * (2*x+1) % mod * inv6 % mod;
}

inline int calc_2(int x,int y) {
return (LL)(x+y)%mod*(y-x+1)%mod*inv2%mod;
}

inline int calc(int x) {
if(x <= Max) return phi[x];
int tt = find(x);
if(tt != -1) return tt;
int last = 0;
long long ans = 0;
(ans += calc_3(x)) %= mod;
for(int i=2;i<=x;i=last+1) {
last = x/(x/i);
(ans += mod-(LL)calc(x/i)*calc_2(i,last)%mod) %= mod;
}
add(x, ans);
return ans;
}

int main() {
init();
int n;
cin >> n;
puts("1");
cout << calc(n) << endl;
}