几个有用的结论：
记狄利克雷卷积： $(f*g)(n)=sumlimits_{d|n} f(d) * g(frac{n}{d})$
则有$mu * 1 = [n=1]$与$phi *1 = id$
若$F(n)=sumlimits_{d|n} f(d)$
则$f(n)=sumlimits_{d|n} mu(d) * F(frac{n}{d})$
bzoj 2190 仪仗队
题目大意：
求$n*n$的矩形中能从左下角被直接看到的点的个数
思路：
设左下角$(0,0)$ 则相当于求$sumlimits_{i=1}^{n-1} sumlimits_{j=1}^{n-1} [gcd(i,j)==1]+2$（$(1,0),(0,1)$
将式子转化为$1+ 2* sumlimits_{i=1}^{n-1} sumlimits_{j=1}^{i} [gcd(i,j)==1] $（由于有$(1,0),(0,1),(1,1)$三个特殊点存在
然后发现满足欧拉函数的定义，则所求为$1+2* sumlimits_{i=1}^{n-1} phi(i) $ 筛出欧拉函数即可
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 50010 15 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 16 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 17 #define ren for(register int i=fst[x];i;i=nxt[i]) 18 #define pb(i,x) vec[i].push_back(x) 19 #define pls(a,b) (a+b)%MOD 20 #define mns(a,b) (a-b+MOD)%MOD 21 #define mul(a,b) (1LL*(a)*(b))%MOD 22 using namespace std; 23 inline int read() 24 { 25 int x=0,f=1;char ch=getchar(); 26 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 27 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 28 return x*f; 29 } 30 int n,tot,p[MAXN+100],ntp[MAXN+100],phi[MAXN+100]; 31 void mem() 32 { 33 phi[1]=1; 34 rep(i,2,MAXN) 35 { 36 if(!ntp[i]) p[++tot]=i,phi[i]=i-1; 37 rep(j,1,tot) if(i*p[j]>MAXN) break; 38 else {ntp[i*p[j]]=1;if(i%p[j]) phi[i*p[j]]=phi[i]*phi[p[j]];else {phi[i*p[j]]=phi[i]*p[j];break;}} 39 } 40 rep(i,2,n-1) phi[i]+=phi[i-1]; 41 } 42 int main() 43 { 44 n=read();mem();printf("%d ",n>1?phi[n-1]<<1|1:0); 45 }
bzoj 2818 Gcd
题目大意：
求有序数对$(i,j),i,j leq n$满足$gcd(i,j)$为素数的数对个数
思路：
法1：和上一题很相似，设$p$为素数则$sumlimits_{i=1}^n sumlimits_{j=1}^n [gcd(i,j)==p] =sumlimits_{i=1}^{n/p} sumlimits_{j=1}^{n/p} [gcd(i,j)==1]$
这样我们枚举每一个素数就可以解决了
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 10001000 15 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 16 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 17 #define ren for(register int i=fst[x];i;i=nxt[i]) 18 #define pb(i,x) vec[i].push_back(x) 19 #define pls(a,b) (a+b)%MOD 20 #define mns(a,b) (a-b+MOD)%MOD 21 #define mul(a,b) (1LL*(a)*(b))%MOD 22 using namespace std; 23 inline int read() 24 { 25 int x=0,f=1;char ch=getchar(); 26 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 27 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 28 return x*f; 29 } 30 int n,tot,p[MAXN],ntp[MAXN+5],phi[MAXN]; 31 ll ans,sum[MAXN]; 32 void mem() 33 { 34 phi[1]=1; 35 rep(i,2,n) 36 { 37 if(!ntp[i]) p[++tot]=i,phi[i]=i-1; 38 rep(j,1,tot) if(i*p[j]>n) break; 39 else {ntp[i*p[j]]=1;if(i%p[j]) phi[i*p[j]]=phi[i]*phi[p[j]];else {phi[i*p[j]]=phi[i]*p[j];break;}} 40 } 41 rep(i,1,n) sum[i]=sum[i-1]+phi[i]; 42 } 43 int main() 44 { 45 n=read();mem(); 46 rep(i,1,tot) ans+=sum[n/p[i]]*2-1;printf("%lld ",ans); 47 }
法2：因为$gcd(i,j)==1$的形式可以转化为$sumlimits_{d|gcd(i,j)} mu(d)$
由于对于每一个$d$，相对应的$i,j$有$n/i,n/j$个。所以前两个枚举是没有必要的，设$m$为$n/p$，式子为$sumlimits_{d=1}^m mu(d)*frac{m}{d}*frac{m}{d}$
枚举素数后数论分块
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 10001000 15 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 16 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 17 #define ren for(register int i=fst[x];i;i=nxt[i]) 18 #define pb(i,x) vec[i].push_back(x) 19 #define pls(a,b) (a+b)%MOD 20 #define mns(a,b) (a-b+MOD)%MOD 21 #define mul(a,b) (1LL*(a)*(b))%MOD 22 using namespace std; 23 inline int read() 24 { 25 int x=0,f=1;char ch=getchar(); 26 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 27 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 28 return x*f; 29 } 30 int n,tot,p[MAXN],ntp[MAXN+5],mu[MAXN]; 31 ll ans,sum[MAXN]; 32 void mem() 33 { 34 mu[1]=1; 35 rep(i,2,n) 36 { 37 if(!ntp[i]) p[++tot]=i,mu[i]=-1; 38 rep(j,1,tot) if(i*p[j]>n) break; 39 else {ntp[i*p[j]]=1;if(i%p[j]) mu[i*p[j]]=-mu[i];else {mu[i*p[j]]=0;break;}} 40 } 41 rep(i,1,n) sum[i]=sum[i-1]+mu[i]; 42 } 43 ll solve(int x,ll res=0) 44 { 45 int pos,lmt=n/x;rep(i,1,lmt) {pos=lmt/(lmt/i);res+=1LL*(sum[pos]-sum[i-1])*(lmt/i)*(lmt/i);i=pos;} 46 return res; 47 } 48 int main() 49 { 50 n=read();mem(); 51 rep(i,1,tot) ans+=solve(p[i]);printf("%lld ",ans); 52 }
bzoj 1101 Zap
题目大意：
求有序数对$(i,j),ileq a,j leq b$满足$gcd(i,j)==d$的数对个数
思路：
和上一题基本一样
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 60100 15 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 16 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 17 #define ren for(register int i=fst[x];i;i=nxt[i]) 18 #define pb(i,x) vec[i].push_back(x) 19 #define pls(a,b) (a+b)%MOD 20 #define mns(a,b) (a-b+MOD)%MOD 21 #define mul(a,b) (1LL*(a)*(b))%MOD 22 using namespace std; 23 inline int read() 24 { 25 int x=0,f=1;char ch=getchar(); 26 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 27 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 28 return x*f; 29 } 30 int n,tot,p[MAXN],ntp[MAXN+5],mu[MAXN]; 31 ll ans,sum[MAXN]; 32 void mem(int n) 33 { 34 mu[1]=1; 35 rep(i,2,n) 36 { 37 if(!ntp[i]) p[++tot]=i,mu[i]=-1; 38 rep(j,1,tot) if(i*p[j]>n) break; 39 else {ntp[i*p[j]]=1;if(i%p[j]) mu[i*p[j]]=-mu[i];else {mu[i*p[j]]=0;break;}} 40 } 41 rep(i,1,n) sum[i]=sum[i-1]+mu[i]; 42 } 43 ll solve(int x,int y,ll res=0) 44 { 45 int pos,lmt=min(x,y);rep(i,1,lmt) 46 {pos=min(x/(x/i),(y/(y/i)));res+=1LL*(sum[pos]-sum[i-1])*(x/i)*(y/i);i=pos;} 47 return res; 48 } 49 int main() 50 { 51 int T=read(),a,b,d;mem(50010);while(T--) 52 {a=read(),b=read(),d=read();printf("%lld ",solve(a/d,b/d));} 53 }
bzoj 4804 欧拉心算
题目大意：
求$sumlimits_i^n sumlimits_j^n phi(gcd(i,j))$
思路：
枚举$gcd$ 得到$sumlimits_{k=1}^n phi(k) sumlimits_{i=1}^n sumlimits_{j=1}^n [gcd(i,j)==k]$
后面的式子转化为第一题的样子即原式为$sumlimits_{k=1}^n phi(k) (-1+2*sumlimits_{i=1}^{n/k}phi(i))$
设$f(k)=-1+2*sumlimits_{i=1}^k phi(i)$ 则原式为$sumlimits_{k=1}^n phi(k)*f(frac{n}{k})$
使用数论分块即可
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 10000001 15 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 16 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 17 #define ren for(register int i=fst[x];i;i=nxt[i]) 18 #define pb(i,x) vec[i].push_back(x) 19 #define pls(a,b) (a+b)%MOD 20 #define mns(a,b) (a-b+MOD)%MOD 21 #define mul(a,b) (1LL*(a)*(b))%MOD 22 using namespace std; 23 inline int read() 24 { 25 int x=0,f=1;char ch=getchar(); 26 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 27 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 28 return x*f; 29 } 30 int n,tot,p[MAXN],ntp[MAXN]; 31 ll ans,phi[MAXN],f[MAXN]; 32 void mem(int n) 33 { 34 phi[1]=1; 35 rep(i,2,n) 36 { 37 if(!ntp[i]) p[++tot]=i,phi[i]=i-1; 38 rep(j,1,tot) if(i*p[j]>n) break; 39 else {ntp[i*p[j]]=1;if(i%p[j]) phi[i*p[j]]=phi[i]*phi[p[j]];else {phi[i*p[j]]=phi[i]*p[j];break;}} 40 } 41 rep(i,1,n) phi[i]+=phi[i-1],f[i]=-1+(phi[i]<<1); 42 } 43 ll solve(int n,ll res=0) 44 { 45 int pos;rep(i,1,n) {pos=n/(n/i);res+=1LL*(phi[pos]-phi[i-1])*f[n/i];i=pos;}return res; 46 } 47 int main() 48 { 49 int T=read(),a,b,d;mem(10000000);while(T--) 50 {a=read();printf("%lld ",solve(a));} 51 }
（卡我一个数组的空间常数怕不是有毒
bzoj 3529 数表
题目大意：
有一张 n×m 的数表，其第 i 行第 j 列（1 <= i <= n, 1 <= j <= m）的数值为能同时整除 i 和 j 的所有自然数之和
q次询问 每次给出n m a 求数表中不大于 a 的数之和
思路：
首先设$f(i)=sumlimits_{d|i} d$
则所求为$sumlimits_{i=1}^n sumlimits_{j=1}^m f(gcd(i,j)) *[f(gcd(i,j))<=a]$
先不管最后那一个限制，继续推式子，我们枚举$gcd$：$sumlimits_{k=1}^n f(k) sumlimits_{i=1}^{n/k} sumlimits_{j=1}^{m/k} [gcd(i,j)==1]$
后面的式子好熟啊，变成莫比乌斯函数形式：$sumlimits_{k=1}^n f(k) sumlimits_d^{min(n/k,m/k)} mu(d) frac{n/k}{d} frac{m/k}{d}$
如果后面把$k$写在分数线下面，发现我们可以枚举$k*d$设为$t$ 得到：$sumlimits_{t=1}^{min(n,m)} frac{n}{t} frac{m}{t} sumlimits_{i|t} f(i)*mu(frac{t}{i})$
把后面的拿出来设$h(t)=sumlimits_{i|t} f(i)*mu(frac{t}{i})$ 原式为$sumlimits_{t=1}^{min(n,m)} h(t) frac{n}{t} frac{m}{t}$我们可以数论分块
对于$h(i)$，我们可以把询问的$a$限制从小到大排序，将所有$i$按照$f(i)$排序后一次加入
每一个$i$，枚举其倍数$t=k *i $即对$h(t)$造成的影响，其复杂度为调和级数，暴力使用树状数组维护对每个倍数的位置$pos(i*j)+=f(i)*mu(j)$
数论分块的时候前缀和用树状数组查询即可
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 100100 15 #define MOD 2147483647 16 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 17 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 18 #define ren for(register int i=fst[x];i;i=nxt[i]) 19 #define pb(i,x) vec[i].push_back(x) 20 #define pls(a,b) (a+b)%MOD 21 #define mns(a,b) (a-b+MOD)%MOD 22 #define mul(a,b) (1LL*(a)*(b))%MOD 23 using namespace std; 24 inline int read() 25 { 26 int x=0,f=1;char ch=getchar(); 27 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 28 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 29 return x*f; 30 } 31 int n,tot,p[MAXN],ntp[MAXN],mu[MAXN],lmt; 32 ll phi[MAXN],f[MAXN],g[MAXN],ans[MAXN],af[MAXN]; 33 struct Ask{int id,n,m,a;}q[MAXN]; 34 bool cmpa(Ask x,Ask y) {return x.a<y.a;} 35 bool cmp(int a,int b) {return f[a]<f[b];} 36 ll c[MAXN]; 37 int lowbit(int x) {return x&(-x);} 38 void mdf(int x,int val) {for(;x<=lmt;x+=lowbit(x)) (c[x]+=val)&=MOD;} 39 ll query(int x,ll res=0) {for(;x;x-=lowbit(x)) (res+=c[x])&=MOD;return res;} 40 void mem(int n) 41 { 42 mu[1]=1;rep(i,2,n) 43 { 44 if(!ntp[i]) p[++tot]=i,mu[i]=-1; 45 rep(j,1,tot) if(i*p[j]>n) break; 46 else {ntp[i*p[j]]=1;if(i%p[j]) mu[i*p[j]]=-mu[i];else {mu[i*p[j]]=0;break;}} 47 } 48 rep(i,1,n) for(int j=af[i]=i;j<=n;j+=i) f[j]=(f[j]+i)&MOD; 49 } 50 ll solve(int n,int m,ll res=0) 51 { 52 int pos;rep(i,1,min(n,m)) 53 {pos=min(n/(n/i),m/(m/i));(res+=1LL*(query(pos)-query(i-1))*(((n/i)*(m/i))&MOD))&=MOD;i=pos;} 54 return res; 55 } 56 int main() 57 { 58 int num=read(),a,b,d;mem(lmt=100000); 59 rep(i,1,num) q[i].n=read(),q[i].m=read(),q[i].a=read(),q[i].id=i; 60 sort(q+1,q+num+1,cmpa);sort(af+1,af+lmt+1,cmp);int pos=1; 61 rep(i,1,num) 62 { 63 while(pos<lmt&&f[af[pos]]<=q[i].a) 64 {rep(j,1,lmt/af[pos]) mdf(af[pos]*j,f[af[pos]]*mu[j]);pos++;} 65 ans[q[i].id]=solve(q[i].n,q[i].m); 66 } 67 rep(i,1,num) printf("%d ",ans[i]); 68 }
bzoj 2154 Crash的数字表格
题目大意：
求$sumlimits_{i=1}^n sumlimits_{j=1}^m lcm(i,j)$
思路：
先将$lcm$转化为$gcd$：$sumlimits_{i=1}^n sumlimits_{j=1}^m frac{i*j}{gcd(i,j)}$
令$n$为较小的那个数，然后枚举一波gcd：$sumlimits_{d=1}^n sumlimits_{i=1}^{n/d} sumlimits_{j=1}^{m/d} frac{id*jd}{d} *[gcd(i,j)==1]$
把$d$提到前面去接化简后面那个很熟的式子：$sumlimits_{d=1}^n d* sumlimits_{i=1}^{n/d} sumlimits_{j=1}^{m/d} sumlimits_{t|gcd(i,j)} i*j*mu(t)$
对于第二个$gcd$继续枚举：$sumlimits_{d=1}^n d * sumlimits_{t=1}^{n/d} sumlimits_{i=1}^{n/(dt)} sumlimits_{j=1}^{m/(dt)} it *jt*mu(t)$
提出$t$： $sumlimits_{d=1}^n d *( sumlimits_{t=1}^{n/d} t^2 * mu(t)) sumlimits_{i=1}^{n/(dt)} i *sumlimits_{j=1}^{m/(dt)} j$
设$s(t)=sumlimits_{i=1}^t i$，则：$sumlimits_{d=1}^n d * sumlimits_{t=1}^{n/d} t^2 * mu(t) * s(frac{n}{dt}) *s(frac{m}{dt})$
与前一题类似的，我们设$g=d*t$，得到：$sumlimits_{g=1}^n sumlimits_{d|g} d * mu(frac{g}{d}) * (frac{g}{d})^2 * s(frac{n}{g}) * s(frac{m}{g})$
把$g$提出来：$sumlimits_{g=1}^n g * s(frac{n}{g}) * s(frac{m}{g})sumlimits_{d|g} frac{g}{d}* mu(frac{g}{d}) $
把后面的部分单独拿出来，设$f(g)=sumlimits_{d|g} d * mu(d)$，发现$f$是$d*mu(d)$与$1$的狄利克雷卷积结果
由于那两个都是积性函数，因此$f$也是积性函数，可以通过筛法来求
对于素数$p$，两个因数分别$1,p$，因此$f(p)=-p+1$
对于线筛中$i%p[j]==0$的情况，考虑新增添的因数一定含有平方因子，因此其莫比乌斯函数一定为0，不需要考虑
最终式子为$sumlimits_{g=1}^n g * s(frac{n}{g}) * s(frac{m}{g}) *f(g) $然后愉快的数论分块即可
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 10001000 15 #define MOD 20101009 16 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 17 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 18 #define ren for(register int i=fst[x];i;i=nxt[i]) 19 #define pb(i,x) vec[i].push_back(x) 20 #define pls(a,b) (a+b+MOD)%MOD 21 #define mns(a,b) (a-b+MOD)%MOD 22 #define mul(a,b) (1LL*(a)*(b))%MOD 23 using namespace std; 24 inline int read() 25 { 26 int x=0,f=1;char ch=getchar(); 27 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 28 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 29 return x*f; 30 } 31 int n,m,tot,p[MAXN],ntp[MAXN],lmt; 32 ll f[MAXN],s[MAXN]; 33 void mem(int n) 34 { 35 f[1]=1;rep(i,2,n) 36 { 37 if(!ntp[i]) p[++tot]=i,f[i]=-i+1; 38 rep(j,1,tot) if(i*p[j]>n) break; 39 else {ntp[i*p[j]]=1;if(i%p[j]) f[i*p[j]]=mul(f[i],f[p[j]]);else {f[i*p[j]]=f[i];break;}} 40 } 41 rep(i,1,n) s[i]=pls(s[i-1],i),f[i]=pls(f[i-1],mul(f[i],i)); 42 } 43 ll solve(int n,int m,ll res=0) 44 { 45 int pos;rep(i,1,min(n,m)) 46 {pos=min(n/(n/i),m/(m/i));res=pls(res,mul(f[pos]-f[i-1],mul(s[n/i],s[m/i])));i=pos;} 47 return res; 48 } 49 int main() 50 { 51 n=read(),m=read();mem(max(n,m));printf("%lld ",solve(n,m)); 52 }