[NOI2016] 循环之美

强烈建议先按照 $yxc$ 老师的进度学掉数学部分，不听课你会一脸懵的。

题意简述：给你 $n$, $m$ 其中 $n,m<=10^8$，要你求由 $x<=n$ 作为分子并且 $y<=m$ 作为分母的值不相等的的 $K$ 进制下的 能达成循环小数或者整数的总数。

如果你按照后面的提示做题可能会想偏这一点因人而异吧。

对于一个循环小数 $x$ ，我们把它进个 $l$ 位循环节的长度，它的值是不会改变的。

那么对于一个特定的循环节长度 $l$，不难发现：
$$\left[\dfrac{x\times K^l}{y}\right] =\ \left[\dfrac{x}{y}\right]$$

然后怎么办呢？？把小数部分表示成为原数字减掉下取整的形式看一看：

$$\frac{x\times K^l}{y} -\ \left\lfloor\dfrac{x\times K^l}{y}\right\rfloor=\frac{x}{y} -\ \left\lfloor\dfrac{x}{y}\right\rfloor$$

$${x\times K^l} -\ \left\lfloor\dfrac{x\times K^l}{y}\right\rfloor\times y={x} -\ \left\lfloor\dfrac{x}{y}\right\rfloor\times y$$

$$x(K^l-1)=y(\left\lfloor\dfrac{x\times K^l}{y}\right\rfloor-\left\lfloor\dfrac{x}{y}\right\rfloor)$$

然后是一个很惊悚的转化，让这道题突然变了质。

$$x(K^l-1)\equiv0\pmod{y}$$
$$x\times K^l\equiv x\pmod{y}$$
$$K^l\equiv 1\pmod{y}$$
不管怎么说，反正只要 $K$ 和 $y$ 互质就一定能找到一个 $l$ 使得 $K^l\equiv 1\pmod{y}$。

所以我们只需要求满足 $K$ 与 $y$ 互质的答案总数就可以了。
由于有些分数可能相等会导致算重复，因此我们钦定 $x$ 和 $y$ 互质。

因此 $Ans$ 不难用和式表示出来。

$$Ans=\sum_{i=1}^n\sum_{j=1}^m\ [i\perp j]\ [j\perp K]$$

$$Ans=\sum_{i=1}^n\sum_{j=1}^m\ [gcd(i,j)=1]\ [gcd(j,k)=1]$$

发现 $gcd$ 于是大力莫比乌斯反演。

把谁拆了呢？如果你把 $[gcd(i,j)=1]$ 你会像我一样推了一张纸的式子然后手足无措，
但如果你把 $[gcd(j,k)=1]$ 拆了的话，事情变得简单起来：
$$\sum_{i=1}^n\sum_{j=1}^m\ [gcd(i,j)=1]\ [gcd(j,k)=1]$$
$$\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1]\ \sum_{d \mid j}^{} \sum_{d \mid k}\mu(d)$$
把 $d$ 放在前面去。

$$\sum_{d \mid k}\ \mu(d) \sum_{i=1}^n \sum_{d \mid j}^{}\ [gcd(i,j)=1]$$
$$\sum_{d \mid k}\ \mu(d) \sum_{i=1}^n \sum_{d=1}^{\left\lfloor\dfrac{n}{d}\right\rfloor}\ [gcd(i,jd)=1]$$
$$\sum_{d \mid k}\ \mu(d) \sum_{i=1}^n \sum_{d=1}^{\left\lfloor\dfrac{n}{d}\right\rfloor}\ [gcd(i,j)=1]\ [gcd(i,d)=1]$$

于是你发现这里面右侧的和式出现了一个很像原答案的式子，
定义：
$$f(m,k)=\sum_{i=1}^n\mu(d) \sum_{j=1}^m\ [i\perp j]\ [j\perp K]$$

那么不难得到转移式子：
$$f(m,k)=\sum_{d \mid k}\mu(d)\ f(\left\lfloor\dfrac{n}{d}\right\rfloor,d)$$

关于 $k=0$ 的时候，$[gcd(j,k)=1]$ 的限制可以直接删掉，这变成又一个经典莫比乌斯反演问题，直接整除分块。

一定要开 $MAP$ 存已经算出来的值，否则慢的要命。。。

中间用到杜教筛筛 $\mu$

复杂度 $O(logn\times logm\times \sqrt{n}+n^{\frac{2}{3}})$

注意卡常，你看我从 $13.7s$ 卡到 $2s$ 以内

若有错误欢迎指出。

Code:

#include<map>
#include<cstdio>
#define min(a,b) ((a)<(b)?(a):(b))
typedef long long LL;

const int N=1e7+5,maxn=1e7;
int n,m,K;
int prime[700000],mu[N],sum[N],tot;
bool flag[N];
struct Fac{
    int p[2010];
}g[2010];
std::map<int,int> MP;

inline int gcd(int a,int b) { return b==0?a:gcd(b,a%b);}
inline void pre()
{
    mu[1]=1;
    for(register int i=2;i<=maxn;++i){
        if(!flag[i]){
            prime[++tot]=i;
            mu[i]=-1;
        }
        for(register int j=1;j<=tot;++j){
            if(i*prime[j]>maxn) break;
            flag[i*prime[j]]=1;
            if(i%prime[j]) mu[i*prime[j]]=-mu[i];
            else{
                mu[i*prime[j]]=0;
                break;
            }
        }
    }
    for(register int i=1;i<=maxn;++i) sum[i]=sum[i-1]+mu[i];
    for(register int k=2;k<=2000;++k){
        for(register int i=1;i<=k;++i)
            if(k%i==0)
                g[k].p[++g[k].p[0]]=i;
    }
}
inline LL Sum_mu(int n)
{
    if(MP[n]) return MP[n];
    //fprintf(stderr,"MU(%d)\n",n);
    LL ans=1;
    int n0;
    register int l=2,r;
    for(;l*l<=n;++l){
        n0=n/l;
        ans-=1llu*(n0<=maxn?sum[n0]:Sum_mu(n0));
    }
    for(;l<=n;l=r+1){
        n0=n/l;
        r=n/n0;
        ans-=1llu*(n0<=maxn?sum[n0]:Sum_mu(n0))*(r-l+1);
    }
    return MP[n]=ans;
}
inline LL solve(int n,int m)
{
    LL ans=0;
    int sx=min(n,m),n0,m0;
    register int l=1,r;
    for(;l<=sx;l=r+1){
        n0=n/l,m0=m/l;
        r=min(n/n0,m/m0);
        ans+=((r<=maxn?sum[r]:Sum_mu(r))-1llu*(l-1<=maxn?sum[l-1]:Sum_mu(l-1)))*n0*m0;
    }
    return ans;
}
inline LL S(int n,int m,int K)
{
    if(K==1) return solve(n,m);
    LL ans=0;
    for(register int i=1;i<=g[K].p[0];++i){
        int d=g[K].p[i];
        ans=ans+1llu*mu[d]*(mu[d]==0?0:(m<d?0:S(m/d,n,d)));
    }
    return ans;
}
signed main()
{
    pre();
    scanf("%d%d%d",&n,&m,&K);
    printf("%lld\n",S(n,m,K));
    return 0;
}