思路:
对于
xy
,我们肯定让
gcd(x,y)=1
,因为这样计算不会重复,而且比不互质时更优,且只要余数出现重复就是有循环节了,根据十进制小数转k进制的方法,我们可以得到
x⋅ks=x(mody)⇒ks=1(mody)
或
y|x
因为
gcd(x,y)=1
,所以
gcd(ks,y)=1
,也就是
gcd(k,y)=1
我们要求的式子就是
∑i=1n∑j=1m[gcd(i,j)=1][gcd(j,k)=1]
把
[gcd(i,j)=1]
移到后面,反演一下,再在前面更换指标
∑dμ(d)⌊nd⌋∑j=1⌊md⌋[gcd(dj,k)=1]
(然后我就反演了半天,最后被reflash和mrazer一眼秒了)
因为
gcd(dj,k)=1⇔gcd(d,k)=1,gcd(j,k)=1
,所以
∑gcd(d,k)=1μ(d)⌊nd⌋F(⌊md⌋)
,其中
F(m)=∑i=1m[gcd(i,k)=1]
关于F函数,因为
gcd(i,k)=gcd(i+k,k)
,所以很多都是重复的,可以预处理
k
以内的F,然后就能
(我原本没预处理,暴力算的,在uoj跑76分,问了reflash才知道可以
O(1)
算的)
前面部分根据d的取值分块算,现在问题就在于如何快速计算
∑gcd(d,k)=1μ(d)
%一发reflash
(这个我真没想到,卡在这里挺久的,因为没有把
μ(di)
拆开看)
设
g(n)=∑i=1n[gcd(i,k)=1]μ(i),M(n)=∑i=1nμ(i)
则
g(n)=M(n)−∑d=2n∑i=1n[gcd(i,k)=d]μ(i)=M(n)−∑d=2,d|kn∑i=1⌊nd⌋[gcd(i,kd)=1]μ(di)=M(n)−∑d=2,d|kn∑i=1⌊nd⌋[gcd(i,kd)=1]μ(d)μ(i)[gcd(i,d)=1]
两个互质条件合并一下,就是
[gcd(i,k)=1]
g(n)=M(n)−∑d=2,d|knμ(d)g(⌊nd⌋)
非常像杜教筛的处理方法,记忆化即可
小范围的时候可以暴力
代码:
#include<cstdio>
#include<iostream>
#define LL long long
#define ui unsigned int
using namespace std;
int gcd(int x,int y){return y?gcd(y,x%y):x;}
int n,m,k;
const int lim=2000000;
int prime[lim/10+5],f[2005],py[2005];
ui pr[2005];
short mu[lim+5],sum[lim+5];
char vis[lim+5];
void init()
{
mu[1]=1;
int t,i,j;
for (i=2;i<=lim;++i)
{
if (!vis[i])
prime[++prime[0]]=i,
mu[i]=-1;
for (j=1;j<=prime[0]&&(t=i*prime[j])<=lim;++j)
{
vis[t]=1;
if (i%prime[j]) mu[t]=-mu[i];
else break;
}
}
for (int i=1;i<=lim;++i) sum[i]=sum[i-1]+mu[i];
for (int i=2;i<=k;++i)
if (k%i==0&&mu[i]) pr[++pr[0]]=i;
for (int i=1;i<=k;++i)
{
f[i]=f[i-1];py[i]=py[i-1];
if (gcd(i,k)==1) ++f[i],py[i]+=mu[i];
}
}
LL S(ui m){return 1LL*(m/k)*f[k]+f[m%k];}
const int P=999983;
struct node{
int tot,first[P+5],next[200005],ha[200005];
ui data[200005];
void push(ui n,int t)
{
int p=n%P;
ha[++tot]=t;
data[tot]=n;
next[tot]=first[p];
first[p]=tot;
}
int get(ui n)
{
for (int i=first[n%P];i;i=next[i])
if (data[i]==n) return i;
return 0;
}
}mp_mu,mp;
int cal_mu(ui n)
{
if (n<=lim) return sum[n];
int id=mp_mu.get(n);
if (id) return mp_mu.ha[id];
int t=1;
for (ui last,p,i=2;i<=n;i=last+1)
last=n/(p=n/i),
t-=(last-i+1)*cal_mu(p);
mp_mu.push(n,t);
return t;
}
int cal(ui n)
{
if (n<=k) return py[n];
int id=mp.get(n);
if (id) return mp.ha[id];
int t=cal_mu(n);
for (int i=1;i<=pr[0]&&pr[i]<=n;++i)
t-=mu[pr[i]]*cal(n/pr[i]);
mp.push(n,t);
return t;
}
main()
{
scanf("%d%d%d",&n,&m,&k);
init();
LL ans=0;
int ls=0,rs,ll=min(n,m);
for (ui p,q,last,i=1;i<=ll;i=last+1)
last=min((ui)n/(p=n/i),(ui)m/(q=m/i)),
ans+=((rs=cal(last))-ls)*S(q)*p,
ls=rs;
printf("%lld\n",ans);
}