经典柿子
[[gcd(i,j)=1]=sum_{wmid gcd(i,j)}mu(w)=sum_{wmid i,wmid j}mu(w)
]
POJ3904 Sky Code
description
给定(n) 个不超过(10000) 的正整数(a_1,a_2,cdots,a_n) ，求
[sum_{1le i<j<k<lle n}[gcd(a_i,a_j,a_k,a_l)=1]
]
(nle 10000,Tle10000)
solution
大力推柿子
[egin{aligned}
ans&=sum_{1le i<j<k<lle n}[gcd(a_i,a_j,a_k,a_l)=1]\
&=sum_{1le i<j<k<lle n}sum_{wmid a_i,wmid a_j,wmid a_k,wmid a_l}mu(w)\
&=sum_{w=1}^{10000}mu(w)sum_{1le i<j<k<lle n}[wmid a_i][wmid a_j][wmid a_k][wmid a_l]
end{aligned}
]
我们不妨假设有(c(w)) 个数是(w) 的倍数，那么后面那些东西就等于(dbinom{c(w)}{4}) 。
复杂度(mathcal O(nT)) 。
code
#include<cstdio>
#include<vector>
using namespace std;
const int N=1e4+5;
typedef long long ll;ll C[N];
bool flag[N];int mu[N];
vector<int>pr;
inline void pre()
{
int n=N-5;
flag[1]=1;mu[1]=1;
for(int i=2;i<=n;++i)
{
if(!flag[i]){pr.push_back(i);mu[i]=-1;}
for(int j=0;j<pr.size();++j)
{
int num=i*pr[j];if(num>n)break;
flag[num]=1;
if(i%pr[j])mu[num]=-mu[i];
else{mu[num]=0;break;}
}
}
for(int i=4;i<=n;++i)
C[i]=1ll*i*(i-1)*(i-2)*(i-3)/24;
}
int n,cnt[N];
int main()
{
pre();
while(scanf("%d",&n)==1)
{
int mx=0;
for(int i=1,d;i<=n;++i)
{
scanf("%d",&d);
for(int j=1;j*j<=d;++j)
if(d%j==0)
{
int w1=d/j,w2=j;
++cnt[w1];
if(w1^w2)++cnt[w2];
}
mx=max(mx,d);
}ll ans=0;
for(int i=1;i<=mx;++i)
{
if(!mu[i])continue;
ans+=mu[i]*C[cnt[i]];
}
printf("%lld
",ans);
for(int i=1;i<=mx;++i)cnt[i]=0;
}
return 0;
}
可见格点
description
三维直角坐标系中三个坐标均非负且不超过(n) 的整点中有多少点可以从原点看到。
(nle 10^6)
solution
原题相当于是求
[sum_{xge 0,yge0,zge0}[gcd(x,y,z)=1]
]
这是经典的莫反柿子。不过需要考虑(x,y,z) 中部分为零的情况，于是就有
[egin{aligned}
&=3+3sum_{x>0,y>0}[gcd(x,y)=1]+sum_{x>0,y>0,z>0}[gcd(x,y,z)=1]\
&=3+3sum_{x>0,y>0}sum_{wmid x,wmid y}mu(w)+sum_{x>0,y>0,z>0}sum_{wmid x,wmid y,wmid z}mu(w)\
&=3+3sum_{w=1}^nmu(w)lfloorfrac nw
floor^2+sum_{w=1}^nmu(w)lfloorfrac nw
floor^3\
&=3+sum_{w=1}^nmu(w)lfloorfrac nw
floor^2(lfloorfrac nw
floor+3)
end{aligned}
]
复杂度(mathcal O(nT)) 。
code
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1E+6+5,V=1e+6;
bool flag[N];vector<int>pr;
int mu[N];
inline void pre()
{
mu[1]=1,flag[1]=1;
for(int i=2;i<=V;++i)
{
if(!flag[i])pr.push_back(i),mu[i]=-1;
for(int c:pr)
{
int num=c*i;if(num>V)break;
flag[num]=1;
if(i%c)mu[num]=-mu[i];
else{mu[num]=0;break;}
}
}
}
int main()
{
pre();
int T;scanf("%d",&T);
while(T-->0)
{
int n;scanf("%d",&n);
ll ans=3;
for(int i=1;i<=n;++i)
{
if(!mu[i])continue;
int k=n/i;
ans+=1ll*mu[i]*k*k*(k+3);
}
printf("%lld
",ans);
}
return 0;
}
数一数a x b
description
令
[f(n)=sum_{i=1}^nsum_{j=1}^n[n
mid ij]
]
求
[g(n)=sum_{dmid n}f(d)
]
(nle 10^9,Tle 20000)
solution
不整除显然是难以计算的，我们考虑简单容斥，即
[f(n)=n^2-h(n)=n^2-sum_{i=1}^nsum_{j=1}^n[nmid ij]
]
因此就有
[g(n)=sum_{dmid n}d^2-sum_{dmid n}sum_{i=1}^dsum_{j=1}^d[dmid ij]
]
分开考虑。我们假设(n) 的质因数分解为(n=prod_limits{i=1}^kp_i^{a_i}) ，那么就有
[sum_{dmid n}d^2=prod_{i=1}^ksum_{j=0}^{a_i}(p_i^j)^2
]
因为次数之和大概是(mathcal O(log n)) 级别的，因此这部分复杂度为(mathcal O(log n)) 。
对于后面部分，我们枚举(w=gcd(d,i)) ，则有
[egin{aligned}
&=sum_{dmid n}sum_{i=1}^dsum_{j=1}^d[dmid ij]\
&=sum_{w=1}^nsum_{dmid n,wmid d}sum_{wmid i}^dsum_{j=1}^d[frac dwmid frac iwcdot j][gcd(frac dw,frac iw)=1]\
&=sum_{wmid n}sum_{dmid frac nw}sum_{j=1}^{wd}[dmid j]sum_{i=1}^{d}[gcd(d,i)=1]\
&=sum_{wmid n}sum_{dmid frac nw}sum_{j=1}^{wd}[dmid j]varphi(d)\
&=sum_{wmid n}sum_{dmid frac nw}sum_{j=1}^{w}varphi(d)\
&=sum_{wmid n}wsum_{dmid frac nw}varphi(d)\
&=sum_{wmid n}wfrac nw\
&=ncdot sigma_0(n)
end{aligned}
]
可以通过预处理质数来优化后续复杂度。至此，本题完结。
总复杂度为(mathcal O(sqrt n+Tlog n)) 。
code
#include<bits/stdc++.h>
using namespace std;
const int N=1e5+5;
typedef unsigned long long ll;
bool flag[N];vector<int>pr;
inline void pre(int n)
{
flag[1]=1;
for(int i=2;i<=n;++i)
{
if(!flag[i])pr.push_back(i);
for(int j=0;j<pr.size();++j)
{
int num=i*pr[j];if(num>n)break;
flag[num]=1;
if(i%pr[j]==0)break;
}
}
}
int main()
{
pre(4e4);int T;scanf("%d",&T);
while(T-->0)
{
int n;scanf("%d",&n);ll g=1,h=n;
for(int i=0;i<pr.size();++i)
{
int p=pr[i];ll ret=0;
if(n<p)break;if(n%p)continue;
int k=0;while(n%p==0)n/=p,++k;
for(int j=0,cur=1;j<=k;++j,cur*=p)ret+=1ull*cur*cur;
g*=ret;h*=k+1;
}
if(n>1)h<<=1,g*=(1+1ull*n*n);
printf("%llu
",g-h);
}
return 0;
}
墨菲斯
description
给定(n,m,p) ，求
[sum_{i=1}^nsum_{j=1}^m[h(gcd(i,j))le p]
]
其中(h(x)) 表示(x) 的唯一分解中质数个数。
solution
首先套路地推柿子，枚举(gcd)
[egin{aligned}
ans&=sum_{d=1}^n[h(d)le p]sum_{dmid i}^nsum_{dmid j}^m[gcd(frac id,frac jd)=1]\
&=sum_{d=1}^n[h(d)le p]sum_isum_jsum_{wmid i,wmid j}mu(w)\
&=sum_{d=1}^n[h(d)le p]sum_{w=1}^{lfloor frac nd
floor}mu(w)lfloor frac n{dw}
floorlfloor frac m{dw}
floor\
&=sum_{T=1}^nlfloor frac nT
floorlfloor frac mT
floorsum_{dmid T}[h(d)le p]mu(frac Td)
end{aligned}
]
似乎遇到了障碍。不过我们注意到(2^{19}>5 imes 10^5) ，因此当(pge 18) 时答案就是(nm) 。我们只用考虑(p<18) 的情况。
令
[f(p,d)=sum_{dmid T}[h(d)le p]mu(frac Td)
]
我们可以在(mathcal O(18nlog n)) 的时间内求出所有的(f(p,d)) ，然后对其求前缀和以便在整除分块时可以快速计算区间和。
复杂度(mathcal O(18nlog n+qsqrt n)) 。
code
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=5e5+5,MX=18;
bool flag[N];int mu[N],h[N];ll f[MX+5][N];
vector<int>pr;
inline void pre()
{
int n=N-5;
flag[1]=1;mu[1]=1;h[1]=0;
for(int i=2;i<=n;++i)
{
if(!flag[i])pr.push_back(i),mu[i]=-1,h[i]=1;
for(int j=0;j<pr.size();++j)
{
int num=i*pr[j];if(num>n)break;
flag[num]=1;h[num]=h[i]+1;
if(i%pr[j])mu[num]=-mu[i];
else{mu[num]=0;break;}
}
}
for(int p=0;p<=MX;++p)
{
for(int d=1;d<=n;++d)
{
if(h[d]>p)continue;
for(int t=d,c=1;t<=n;t+=d,++c)f[p][t]+=mu[c];
}
for(int i=1;i<=n;++i)f[p][i]+=f[p][i-1];
}
}
int main()
{
pre();int T;scanf("%d",&T);
while(T-->0)
{
int n,m,p;scanf("%d%d%d",&n,&m,&p);
if(p>MX){printf("%lld
",1ll*n*m);continue;}
if(n>m)swap(n,m);ll ans=0;
for(int l=1,r;l<=n;l=r+1)
{
r=min(n/(n/l),m/(m/l));
ans+=1ll*(n/l)*(m/l)*(f[p][r]-f[p][l-1]);
}
printf("%lld
",ans);
}
return 0;
}