题目链接
https://lydsy.com/JudgeOnline/problem.php?id=4652
题解
容易发现,若ij\frac{i}{j}ji在kkk进制下为纯循环小数,那么必定有
ikl=imod  jkl=1mod  j ik^l=i\mod{j}\\ k^l=1\mod{j} ikl=imodjkl=1modj
因此gcd(j,k)=1\gcd(j,k)=1gcd(j,k)=1。
那么题目要求的就是
F(n,m,k)=∑i=1n∑j=1m[gcd(i,j)=1][gcd(j,k)=1] F(n,m,k)=\sum_{i=1}^n \sum_{j=1}^m [\gcd(i,j)=1][\gcd(j,k)=1] F(n,m,k)=i=1∑nj=1∑m[gcd(i,j)=1][gcd(j,k)=1]
反演得到
F(n,m,k)=∑d∣kF(md,n,d) F(n,m,k)=\sum_{d|k}F(\frac{m}{d},n,d) F(n,m,k)=d∣k∑F(dm,n,d)
边界条件
F(0,m,k)=F(n,0,k)=0F(n,m,1)=∑i=1min(n,m)μ(i)⌊ni⌋⌊mi⌋ F(0,m,k)=F(n,0,k)=0\\ F(n,m,1)=\sum_{i=1}^{\min(n,m)}\mu(i)\lfloor\frac{n}{i}\rfloor\lfloor\frac{m}{i}\rfloor F(0,m,k)=F(n,0,k)=0F(n,m,1)=i=1∑min(n,m)μ(i)⌊in⌋⌊im⌋
注意需要记忆化。
代码
#include <map>
#include <cstdio>
#include <algorithm>
int read()
{
int x=0,f=1;
char ch=getchar();
while((ch<'0')||(ch>'9'))
{
if(ch=='-')
{
f=-f;
}
ch=getchar();
}
while((ch>='0')&&(ch<='9'))
{
x=x*10+ch-'0';
ch=getchar();
}
return x*f;
}
const int maxk=2000;
const int maxn=1000000;
int p[maxn+10],prime[maxn+10],cnt,mu[maxn+10],sum[maxn+10],d[maxk+10][maxk+10];
int getprime()
{
p[1]=mu[1]=1;
for(int i=2; i<=maxn; ++i)
{
if(!p[i])
{
prime[++cnt]=i;
mu[i]=-1;
}
for(int j=1; (j<=cnt)&&(i*prime[j]<=maxn); ++j)
{
int x=i*prime[j];
p[x]=1;
if(i%prime[j]==0)
{
mu[x]=0;
break;
}
mu[x]=-mu[i];
}
}
for(int i=1; i<=maxn; ++i)
{
sum[i]=sum[i-1]+mu[i];
}
for(int i=1; i<=maxk; ++i)
{
if(mu[i])
{
for(int j=i; j<=maxk; j+=i)
{
d[j][++d[j][0]]=i;
}
}
}
return 0;
}
std::map<int,int> mp;
int getsum(int n)
{
if(n<=maxn)
{
return sum[n];
}
if(mp.count(n))
{
return mp[n];
}
int ans=1;
for(int l=2,r; l<=n; l=r+1)
{
r=n/(n/l);
ans-=(r-l+1)*getsum(n/l);
}
return mp[n]=ans;
}
struct data
{
int n,m,k;
data(int _n=0,int _m=0,int _k=0):n(_n),m(_m),k(_k){}
bool operator <(const data &oth) const
{
return (n==oth.n)?((m==oth.m)?(k<oth.k):(m<oth.m)):(n<oth.n);
}
};
std::map<data,long long> mf;
long long F(int n,int m,int k)
{
if((n==0)||(m==0))
{
return 0;
}
data da(n,m,k);
if(mf.count(da))
{
return mf[da];
}
long long res=0;
if(k==1)
{
for(int l=1,r,s,t=0; l<=std::min(n,m); l=r+1,t=s)
{
r=std::min(n/(n/l),m/(m/l));
s=getsum(r);
res+=1ll*(s-t)*(n/l)*(m/l);
}
mf[data(m,n,k)]=res;
}
else
{
for(int i=1; (i<=d[k][0])&&(d[k][i]<=m); ++i)
{
res+=mu[d[k][i]]*F(m/d[k][i],n,d[k][i]);
}
}
return mf[da]=res;
}
int n,m,k;
int main()
{
getprime();
n=read();
m=read();
k=read();
printf("%lld\n",F(n,m,k));
return 0;
}