[bzoj4652] [NOI2016]循环之美

Description

牛牛是一个热爱算法设计的高中生。在他设计的算法中，常常会使用带小数的数进行计算。牛牛认为，如果在 k

进制下，一个数的小数部分是纯循环的，那么它就是美的。现在，牛牛想知道：对于已知的十进制数 n 和 m，在

kk 进制下，有多少个数值上互不相等的纯循环小数，可以用分数 xy 表示,其中 1≤x≤n,1≤y≤m，且 x,y是整数

。一个数是纯循环的，当且仅当其可以写成以下形式：a.c1˙c2c3…cp-1cp˙其中，a 是一个整数，p≥1；对于 1

≤i≤p，ci是 kk 进制下的一位数字。例如，在十进制下，0.45454545……=0.4˙5˙是纯循环的，它可以用 5/11

、10/22 等分数表示；在十进制下，0.1666666……=0.16˙则不是纯循环的，它可以用 1/6 等分数表示。需要特

别注意的是，我们认为一个整数是纯循环的，因为它的小数部分可以表示成 0 的循环或是 k?1 的循环；而一个小

数部分非 0 的有限小数不是纯循环的。

Input

只有一行，包含三个十进制数N,M,K意义如题所述,保证 1≤n≤10^9,1≤m≤109,2≤k≤2000

Output

一行一个整数，表示满足条件的美的数的个数。

Sample Input

2 6 10

Sample Output

4

Solution

\(k\)进制下\(x/y\)为循环小数，即：

\[\big[\frac{xk^l}{y}\big]=\big[\frac{x}{y}\big] \]

其中中括号表示小数部分，\(l\)是\(k\)进制下循环节长度。

中括号展开可得：

\[\frac{xk^l}{y}-\lfloor\frac{xk^l}{y}\rfloor=\frac{x}{y}-\lfloor\frac{x}{y}\rfloor\\ xk^l-y\lfloor\frac{xk^l}{y}\rfloor=x-y\lfloor\frac{x}{y}\rfloor \]

考虑取整那一项不好处理，可以玩一波骚操作，两边对\(y\)取模：

\[xk^l\equiv x\pmod{y} \]

因为\((x,y)=1\)，所以两边除掉：

\[k^l\equiv1\pmod{y} \]

可得：

\[(k,y)=1 \]

所以答案可以形式化的写成：

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

然后尝试着对\([(i,j)=1]\)进行莫比乌斯反演：

\[\begin{align} ans=&\sum_{i=1}^n\sum_{j=1}^m[(j,k)=1]\sum_{d|i\&d|j}\mu(d)\\ =&\sum_{d=1}^{\min(n,m)}\mu(d)[(d,k)=1]\lfloor\frac{n}{d}\rfloor\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[(j,k)=1] \end{align} \]

此处省略了中间交换求和符号之类的过程。

然后考虑下这个式子后面一块，设：

\[f(n)=\sum_{i=1}^n[(i,k)=1] \]

考虑分成几块，每块长度都是\(k\)，显然：

\[f(n)=\lfloor\frac{n}{k}\rfloor f(k)+f(n\bmod k) \]

再考虑下答案的前半段，由于后两项要数论分块，所以要求前两项的前缀和，设：

\[g(n,k)=\sum_{d=1}^n\mu(d)[(d,k)=1] \]

对后面莫比乌斯反演下：

\[\begin{align} g(n,k)=&\sum_{d=1}^n\mu(d)\sum_{t|d\&t|k}\mu(t)\\ =&\sum_{t=1}^n\mu(t)[t|k]\sum_{d=1}^{\lfloor\frac{n}{t}\rfloor}\mu(dt) \end{align} \]

然后最玄学的一步来了，，若想要使\(\mu(dt)\ne 0\)，显然\((d,t)=1\)，此时\(\mu(dt)=\mu(d)\mu(t)\)，所以后面一项可以写成：

\[\begin{align} g(n,k)=&\sum_{t=1}^n\mu(t)[t|k]\sum_{d=1}^{\lfloor\frac{n}{t}\rfloor}\mu(d)\mu(t)[(d,t)=1]\\ =&\sum_{t=1}^n\mu^2(t)[t|k]\sum_{d=1}^{\lfloor\frac{n}{t}\rfloor}\mu(d)[(d,t)=1]\\ =&\sum_{t=1}^n\mu^2(t)[t|k]g(\lfloor\frac{n}{d}\rfloor,t)\\ =&\sum_{t|k}\mu^2(t)g(\lfloor\frac{n}{d}\rfloor,t)\\ \end{align} \]

这是个递归的形式，预处理下\(k\)以内所有数的约数，直接记忆化暴力算就好了。

对于\(k=1\)的情况，\(g\)就变成了\(\sum_{i=1}^n \mu(i)\)，这个可以杜教筛出来。

答案可以写成：

\[ans=\sum_{d=1}^{\min(n,m)}\mu(d)[(d,k)=1]\lfloor\frac{n}{d}\rfloor f(\lfloor\frac{m}{d}\rfloor) \]

后面数论分块就行了。

#pragma GCC optimize(3)
#include<bits/stdc++.h>
using namespace std;

#ifdef ONLINE_JUDGE
#define getchar() ((p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin)),p1==p2)?EOF:*p1++)
#endif

namespace fast_IO {
	char buf[1<<21],*p1=buf,*p2=buf;

	template <typename T> inline void read(T &x) {
		x=0;T f=1;char ch=getchar();
		for(;!isdigit(ch);ch=getchar()) if(ch=='-') f=-f;
		for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';x*=f;
	}
	template <typename T,typename... Args> inline void read(T& x,Args& ...args) {
		read(x),read(args...);
	}

	char buf2[1<<21],a[80];int p,p3=-1;

	inline void flush() {fwrite(buf2,1,p3+1,stdout),p3=-1;}
	template <typename T> inline void write(T x) {
		if(p3>(1<<20)) flush();
		if(x<0) buf2[++p3]='-',x=-x;
		do {a[++p]=x%10+48;} while(x/=10);
		do {buf2[++p3]=a[p];} while(--p);
		buf2[++p3]='\n';
	}
	template <typename T,typename... Args> inline void write(T x,Args ...args) {
		write(x),write(args...);
	}
}

using fast_IO :: read;
using fast_IO :: write;
using fast_IO :: flush;

#define ll long long 

const int maxn = 5e3+10;
const int N = 5e6+10;

int pri[N],tot,mu[N],vis[N],F[maxn],Mu[N];

void sieve() {
	mu[1]=1;
	for(int i=2;i<N;i++) {
		if(!vis[i]) pri[++tot]=i,mu[i]=-1;
		for(int j=1;j<=tot&&i*pri[j]<N;j++) {
			vis[i*pri[j]]=1;
			if(i%pri[j]==0) break;
			mu[i*pri[j]]=-mu[i];
		}
	}
	for(int i=1;i<N;i++) Mu[i]=Mu[i-1]+mu[i];
}

map<int ,int > mp;

int sum_mu(int n) {
	if(n<N) return Mu[n];
	if(mp[n]) return mp[n];
	int res=1,T=2;
	while(T<=n) {
		int pre=T;T=n/(n/T);
		res=res-(T-pre+1)*sum_mu(n/T);T++;
	}
	return mp[n]=res;
}

int dv[2001][201];
map<int ,int > G[2001];

int g(int n,int k) {
	if(k==1) return sum_mu(n);int ans=0;
	if(G[k][n]) return G[k][n];
	if(!n) return 0;
	for(int i=1;i<=dv[k][0];i++)
		ans+=abs(mu[dv[k][i]])*g(n/dv[k][i],dv[k][i]);
	return G[k][n]=ans;
}

int n,m,k;

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

int main() {
	read(n,m,k);
	for(int i=1;i<=k;i++) F[i]=F[i-1]+(__gcd(i,k)==1);
	for(int i=1;i<=k;i++)
		for(int j=1;j<=i;j++)
			if(i%j==0) dv[i][++dv[i][0]]=j;
	sieve();int T=1;ll ans=0;
	while(T<=n&&T<=m) {
		int pre=T;T=min(n/(n/T),m/(m/T));
		ans=ans+1ll*(g(T,k)-g(pre-1,k))*(n/T)*f(m/T);T++;
	}write(ans);
	flush();
	return 0;
}