bzoj 4176 Lucas的数论 莫比乌斯反演（杜教筛）

Posted 2020-10-07 Troy Ricardo

Description

去年的Lucas非常喜欢数论题，但是一年以后的Lucas却不那么喜欢了。

在整理以前的试题时，发现了这样一道题目“求Sigma(f(i)),其中1<=i<=N”，其中 表示i的约数个数。他现在长大了，题目也变难了。

求如下表达式的值：

 

 

一行一个整数ans，表示答案模1000000007的值。

Sample Input

2

Sample Output

8

HINT

 

 对于100%的数据n <= 10^9。

题解：

　　解锁新技能：杜教筛。

　　再复习一下：

　　若$F(n)=\sum_{i=1}^{n}f(i)，g(i)=\sum_{j|i}f(i),G(n)=\sum_{i=1}^{n}g(i)$，

　　则有：$G(n)=\sum_{i=1}^{n}F(\lfloor\frac{n}{i}\rfloor)$。

　　即：$F(n)=G(n)-\sum_{i=2}^{n}F(\lfloor\frac{n}{i}\rfloor)$

　　然后就可以上杜教筛了。

　　话归正题：
　　对于本题，其实和bzoj3994式子一样……

$\sum_{i=1}^{n}\sum_{j=1}^{n}f(ij)$

$=\sum_{i=1}^{n}\sum_{j=1}^{n}\sum_{p=1}^{n^{2}}[p|ij]$

$=\sum_{i=1}^{n}\sum_{j=1}^{n}\sum_{p=1}^{n^{2}}[\frac{p}{(p,i)}|j]$

$=\sum_{t=1}^{n}\sum_{i=1}^{\lfloor\frac{n}{t}\rfloor}\sum_{j=1}^{n}\sum_{p=1}^{\lfloor\frac{n^{2}}{t}\rfloor}[p|j]*[(i,p)==0]$

$=\sum_{t=1}^{n}\sum_{i=1}^{\lfloor\frac{n}{t}\rfloor}\sum_{p=1}^{\lfloor\frac{n^{2}}{t}\rfloor}\sum_{d|(i,j)}\mu(d)\lfloor \frac{n}{p} \rfloor$

$=\sum_{d=1}^{n}\mu(d)\sum_{t=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{i=1}^{\lfloor\frac{n}{td}\rfloor}\sum_{p=1}^{\lfloor\frac{n^{2}}{td}\rfloor}\lfloor\frac{n}{pd}\rfloor$

$n>=pd,n^{2}>=tdn$

$\therefore Ans=\sum_{d=1}^{n}\mu(d)\sum_{t=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{i=1}^{\lfloor\frac{n}{td}\rfloor}\sum_{p=1}^{\lfloor\frac{n}{td}\rfloor}\lfloor\frac{n}{pd}\rfloor$

$Ans=\sum_{d=1}^{n}\mu(d)\sum_{t=1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{n}{td}\rfloor\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{n}{pd}\rfloor$

$设f(i)=\sum_{x=1}^{i}\lfloor\frac{i}{x} \rfloor$

$\therefore Ans=\sum_{d=1}^{n}\mu(d)f(\lfloor\frac{n}{d}\rfloor)$

　　然后分块+杜教筛即可。

代码：

　　

 1 #include<cstdio>
 2 using namespace std;
 3 typedef long long ll;
 4 const ll N=(ll)3e6+1;
 5 const ll mod=1000000007ll;
 6 short miu[N];
 7 int prim[N/10],num;
 8 bool vis[N];
 9 int sum[N];
10 inline void init(){
11     miu[1]=sum[1]=1;
12     for(int i=2;i<N;i++){
13         if(!vis[i]){
14             prim[++num]=i;
15             miu[i]=-1;
16         }for(int j=1;j<=num&&prim[j]*i<N;j++){
17             vis[i*prim[j]]=1;
18             if(i%prim[j]==0){
19                 miu[i*prim[j]]=0;
20                 break;
21             }
22             else
23                 miu[i*prim[j]]=-miu[i];
24         }
25       //  printf("miu[%d]=%d\n",i,miu[i]);
26     }
27     for(int i=2;i<N;i++)
28         sum[i]=sum[i-1]+miu[i];
29 }
30 struct edges{
31     ll v;int w;edges *last;
32 }edge[N/10],*head[76545];int cnt;
33 const int limit=76543;
34 inline void push(int u,ll v,int w){
35     edge[++cnt]=(edges){v,w,head[u]};head[u]=edge+cnt;
36 }
37 inline ll Get_sum(ll x){
38     if(x<N) return sum[x];
39     int t=x%limit;
40     for(edges *i=head[t];i;i=i->last)
41         if(i->v==x) return i->w;
42     ll ans=1;
43     for(ll i=2,pos;i<=x;i=pos+1){
44         pos=x/i;
45         pos=x/pos;
46         ans-=(pos-i+1)*Get_sum(x/pos);
47     }
48     push(t,x,ans);
49     return ans;
50 }
51 inline ll Get_F(ll x){
52     ll ans=0;
53     for(ll i=1,pos;i<=x;i=pos+1){
54         pos=x/i;
55         pos=x/pos;
56         ll t=x/i;
57         ans+=t*(pos-i+1);
58         ans%=mod;
59     }return ans;
60 }
61 inline ll solve(ll x){
62     ll ans=0;
63     for(ll i=1,pos;i<=x;i=pos+1){
64         pos=x/i;
65         pos=x/pos;
66         ll t=Get_F(x/i);
67         ans+=(Get_sum(pos)-Get_sum(i-1))%mod*t%mod*t%mod;
68         ans%=mod;
69     }
70     return (ans%mod+mod)%mod;
71 }
72 int main(){
73     init();
74     ll a;
75     scanf("%lld",&a);
76     printf("%lld\n",solve(a));
77 }