首先要知道什么樣的數才是"純循環數".打表可以發現,這樣的數當且僅當分母\(y\)和\(k\)互質,這是因為,首先考慮除法過程,每次先給當前余數\(*k\),然后對分母做帶余除法,如果把每次的余數寫成一個序列,那么"純循環數"就要使的第一個出現過兩次的余數正好為序列中第一個余數.設第一個余數為\(x\),循環節長度為\(a\),那么會有\(xk^a\equiv x \pmod y\),即\(k^a\equiv 1 \pmod y\),注意到\(i\ge 1\)時得到的\(k^i \bmod y\)都是\(gcd(k,y)\)的倍數,所以必須要有\(gcd(k,y)=1\)
然后要使得數值不同,所以其實要求的是這個東西
\(\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=1][\gcd(j,k)=1]\)
先把只有\(j\)和\(k\)的提前
\(\sum_{j=1}^{m}[\gcd(j,k)=1]\sum_{i=1}^{n}[\gcd(i,j)=1]\)
然后把\([\gcd(j,k)=1]\)轉化一下
\(\sum_{j=1}^{m}\sum_{d|j,d|k}\mu(d)\sum_{i=1}^{n}[\gcd(i,j)=1]\)
把\(d\)提前
\(\sum_{d=1}^{m}\mu(d)\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}\sum_{i=1}^{n}[\gcd(i,jd)=1]\)
\(\sum_{d=1}^{m}\mu(d)\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}\sum_{i=1}^{n}[\gcd(i,j)=1][\gcd(i,d)=1]\)
后面那個是不是有點眼熟?我們如果記\(s(i,j,k)=\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=1][\gcd(j,k)=1]\),那么就能得到
\(s(i,j,k)=\sum_{d|k}\mu(d)s(\lfloor\frac{m}{d}\rfloor,n,d)\)
這個遞歸處理就好了,邊界就是\(m=0\)或者\(n=0\)時為\(0\),\(k=1\)時為\(\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=1]\),數論分塊+杜教篩即可.復雜度大概是\(O(\log n\log m \sqrt{n}+n^{\frac{2}{3}})\),這個還是比較慢的,加了記憶化都要1800多ms 也可能是我寫丑了
// luogu-judger-enable-o2
#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdlib>
#include<cstdio>
#include<vector>
#include<cmath>
#include<ctime>
#include<queue>
#include<map>
#include<set>
#define LL long long
#define db double
using namespace std;
const int N=5e6+10,M=2000+10;
int rd()
{
int x=0,w=1;char ch=0;
while(ch<'0'||ch>'9'){if(ch=='-') w=-1;ch=getchar();}
while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+(ch^48);ch=getchar();}
return x*w;
}
int n,m,k,prm[N],tt;
bool pp[N];
LL mu[N];
map<int,LL> f;
int lm;
LL siv(int nn)
{
if(nn<=N-10) return mu[nn];
if(f.count(nn)) return f[nn];
LL &an=f[nn];
an=1;
for(int i=2,j;i<=nn;i=j+1)
{
j=nn/(nn/i);
an-=1ll*siv(nn/i)*(j-i+1);
}
return an;
}
LL ff(int nn,int mm)
{
lm=min(nn,mm);
LL an=0;
for(int i=1,j;i<=lm;i=j+1)
{
j=min(nn/(nn/i),mm/(mm/i));
an+=1ll*(siv(j)-siv(i-1))*(nn/i)*(mm/i);
}
return an;
}
vector<int> dd[M];
struct node
{
int n,m,k;
bool operator < (const node &bb) const {return n!=bb.n?n<bb.n:(m!=bb.m?m<bb.m:k<bb.k);}
};
map<node,LL> g;
LL sov(int n,int m,int k)
{
if(n<=0||m<=0) return 0;
node nw=(node){n,m,k};
if(g.count(nw)) return g[nw];
if(k==1) return g[nw]=ff(n,m);
LL an=0;
vector<int>::iterator it;
for(it=dd[k].begin();it!=dd[k].end();++it)
{
int i=*it;
an+=1ll*sov(m/i,n,i)*(mu[i]-mu[i-1]);
}
return g[nw]=an;
}
int main()
{
mu[1]=1;
for(int i=2;i<=N-10;++i)
{
if(!pp[i]) prm[++tt]=i,mu[i]=-1;
for(int j=1;j<=tt&&i*prm[j]<=N-10;++j)
{
pp[i*prm[j]]=1;
if(i%prm[j]==0) break;
mu[i*prm[j]]=-mu[i];
}
}
for(int i=2;i<=N-10;++i) mu[i]+=mu[i-1];
n=rd(),m=rd(),k=rd();
for(int i=1;i<=k;++i)
if(k%i==0) dd[k].push_back(i);
int nn=dd[k].size();
for(int i=0;i<nn-1;++i)
{
int x=dd[k][i];
vector<int>::iterator it;
for(it=dd[k].begin();it!=dd[k].end();++it)
{
int y=*it;
if(x<y) break;
if(x%y==0) dd[x].push_back(y);
}
}
printf("%lld\n",sov(n,m,k));
return 0;
}