BZOJ4652：循环之美（数学+杜教筛）

题面
这题天秀
题意：求有多少个互质的x，y，x≤n，y≤m
使得$\frac{x}{y}$的小数在k进制下从第一位开始无限循环

我看到这题，突然发现和我小学5年级奥数学的内容很像（我小学好厉害）

当时讲了怎么把一个循环小数化成分数
从第一位开始循环的小数
循环节长度就是分母中9的个数
分子就是循环节
(证明就是等比数列求和，然后求极限，好像当时也讲了）
比如

$$0.123456123456123456...=\frac{123456}{999999}$$

若不从第一位开始的就各种解方程的骚方法求出来
比如

$$0.123456666666...=\frac{6}{9}-0.54321$$
（小学奥数普及）

推导到k进制
分母应该只能是若干个(k-1)
约分后也是若干个(k-1)的约数
就必定与k互质了

正经的证明及以下的题解

设

$$f(n,m,k)=\sum_{i=1}^n\sum_{j=1}^m[(i,j)=1][(j,k)=1]$$
$$=\sum_{i=1}^n \sum_{jd=1}^m\sum_{d|jd,d|k}\mu(d)$$
一轮画柿子
$$=\sum_{d|k}\mu(d)f(\frac{m}{d},n,d)$$
k=1时为杜教筛
然后玄学45000ms卡过bzoj
洛谷T两个点

其实根据套路，玄学都要记忆化才能过
正解也是玄学
所以玄学的维数不要太多

有点T的玄学

#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
#include <ctime>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <map>

using namespace std;
#define mmst(a, b) memset(a, b, sizeof(a))
#define mmcp(a, b) memcpy(a, b, sizeof(b))

typedef long long LL;

const int N=6006000;
LL n,m,k;
int miu[N],smiu[N];
int prime[N],num;
bool b[N];

map<int,LL> mp;

int gm(int x)
{
    if(x<N)
    return smiu[x];
    if(mp.count(x))
    return mp[x];
    int res=1;
    for(int i=2,last;i<=x;i=last+1)
    {
        last=x/(x/i);
        res-=gm(x/i)*(last-i+1);
    }
    mp[x]=res;
    return res;
}

LL G(int x,int y)
{
    if(x>y)
    swap(x,y);
    LL res=0;
    for(int i=1,last;i<=x;i=last+1)
    {
        last=min(x/(x/i),y/(y/i));
        res+=(LL)(gm(last)-gm(i-1))*(LL)(x/i)*(y/i);
    }
    return res;
}

LL F(int x,int y,int d)
{
    if(!x||!y)
    return 0;
    if(d==1)
    return G(x,y);
    LL res=0;
    for(int i=1;i*i<=d;i++)
    if(d%i==0)
    {
        res+=(LL)miu[i]*F(y/i,x,i);
        if(i*i<d)
        res+=(LL)miu[d/i]*F(y/(d/i),x,d/i);
    }
    return res;
}

int main()
{
    miu[1]=1;
    for(int i=2;i<N;i++)
    {
        if(!b[i])
        prime[++num]=i,miu[i]=-1;
        for(int j=1;j<=num&&i*prime[j]<N;j++)
        {
            miu[i*prime[j]]=-miu[i];
            b[i*prime[j]]=1;
            if(i%prime[j]==0)
            {
                miu[i*prime[j]]=0;
                break;
            }
        }
    }
    for(int i=1;i<N;i++)
    smiu[i]=smiu[i-1]+miu[i];

    cin>>n>>m>>k;
    cout<<F(n,m,k)<<endl;

    return 0;
}

传说中的正解

#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
#include <ctime>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <map>

using namespace std;
#define mmst(a, b) memset(a, b, sizeof(a))
#define mmcp(a, b) memcpy(a, b, sizeof(b))

typedef long long LL;

const int N=6006000;
LL n,m,k,ans;
int u[N],miu[N],smiu[N],phi[N],f[N];
int prime[N],num;
bool b[N];

map<int,LL> mp;
map<LL,LL> mmp;

int gm(int x)
{
    if(x<N)
    return smiu[x];
    if(mp.count(x))
    return mp[x];
    int res=1;
    for(int i=2,last;i<=x;i=last+1)
    {
        last=x/(x/i);
        res-=gm(x/i)*(last-i+1);
    }
    mp[x]=res;
    return res;
}

int gcd(int x,int y)
{
    if(!y)
    return x;
    return gcd(y,x%y);
}

int lj(int x)
{
    return (x/k)*phi[k]+f[x%k];
}

LL F(int x,int y)
{
    if(mmp.count(2001ll*x+y))
    return mmp[2001ll*x+y];
    if(!x)
    return 0;
    if(y==1)
    return gm(x);
    LL res=F(x,y/u[y])+F(x/u[y],y);
    mmp[2001ll*x+y]=res;
    return res;
}

int main()
{
    phi[1]=miu[1]=1;
    for(int i=2;i<N;i++)
    {
        if(!b[i])
        prime[++num]=i,phi[i]=i-1,miu[i]=-1,u[i]=i;
        for(int j=1;j<=num&&i*prime[j]<N;j++)
        {
            miu[i*prime[j]]=-miu[i];
            b[i*prime[j]]=1;
            u[i*prime[j]]=prime[j];
            phi[i*prime[j]]=phi[i]*(prime[j]-1);
            if(i%prime[j]==0)
            {
                phi[i*prime[j]]=phi[i]*prime[j];
                miu[i*prime[j]]=0;
                break;
            }
        }
    }

    for(int i=1;i<N;i++)
    smiu[i]=smiu[i-1]+miu[i];

    cin>>n>>m>>k;
    int tu=1;
    while(k>1)
    {
        tu*=u[k];
        while(u[k/u[k]]==u[k])
        k/=u[k];
        k/=u[k];
    }
    k=tu;

    for(int i=1;i<=k;i++)
    if(gcd(i,k)==1)
    f[i]++;
    for(int i=1;i<=k;i++)
    f[i]+=f[i-1];

    for(int i=1,last;i<=min(n,m);i=last+1)
    {
        last=min(n/(n/i),m/(m/i));
        ans+=(F(last,k)-F(i-1,k))*(LL)(n/i)*lj(m/i);
    }

    cout<<ans<<endl;

    return 0;
}