题解「Luogu1587 [NOI2016]循环之美」

题意

求满足 (1 leq x leq n,1 leq y leq m) 的在 (k) 进制下能写成纯循环小数的最简分数 (frac{x}{y}) 的个数。

题解

---

证明：

[ ext{Ans}=sum_{x=1}^{n}sum_{y=1}^{m}[x ot y][y ot k] ]

首先要知道， (k) 进制下 (frac{x}{y}) 小数点右移一位等于 (frac{x}{y} imes k) 。

假设 (frac{x}{y}) 小数部分的循环节长度为 (l) ，由于纯循环小数的特殊性质有：

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

这里的方括号代表取小数部分。

那么有：

[frac{xk^l}{y}-lfloorfrac{xk^l}{y} floor=frac{x}{y}-lfloorfrac{x}{y} floor ]

注意这里是下取整。

两边同时乘上 (y) ：

[xk^l-lfloorfrac{xk^l}{y} floor imes y=x-lfloorfrac{x}{y} floor imes y ]

变化成同余式：

[xk^l equiv x ({ m{mod}} y) ]

因为 (x ot y) ，则有：

[k^l equiv k equiv 1 ({ m{mod}} y) implies y ot k ]

那么分数 (frac{x}{y}) 当且仅当 (x ot y,y ot k) 才会对 ( ext{Ans}) 产生贡献。

证毕。

---

接下来就开始推式子。

将式子中 ([x ot y]) 用莫比乌斯函数替换：

[egin{align} sum_{x=1}^{n}sum_{y=1}^{m}[x ot y][y ot k] & = sum_{x=1}^{n}sum_{y=1}^{m}[y ot k]sum_{d|x,d|y} mu(d)\ &=sum_{d=1}^{{ m{min}}(n,m)}mu(d) lfloorfrac{n}{d} floorsum_{y=1}^{lfloorfrac{m}{d} floor} [yd ot k]=sum_{d=1}^{{ m{min}}(n,m)}mu(d) lfloorfrac{n}{d} floorsum_{y=1}^{lfloorfrac{m}{d} floor} [y ot k][d ot k]\ &=sum_{d=1}^{{ m{min}}(n,m)}mu(d)[d ot k] lfloorfrac{n}{d} floorsum_{y=1}^{lfloorfrac{m}{d} floor} [y ot k] end{align} ]

如果我们能快速求出 (sum_{y=1}^{lfloorfrac{m}{d} floor} [y ot k]) 和 (mu(d)[d ot k]) 的前缀和，那么用整除分块就容易得到 ( ext{Ans}) 。

首先不难看出 (k) 中多个相同的质因子对 (k) 与其他数是否互质没有影响，所以不妨假设 (k) 的所有质因子次数均为1。

令：

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

那么 (f(n,1)=sum_{i=1}^n1=n) 能直接得出， (g(n,1)=sum_{i=1}^nmu(i)) 能通过杜教筛得出。

考虑 (k ot = 1) 的情况：

对于 (k) 的一个质因子 (p) ，考虑 (f(n,k)) 如何由 (f(n,frac{k}{p})) 得来。

由于多了一个质因子，则与 (frac{k}{p}) 互质的数中有一些不与 (k) 互质，即能被 (p) 整除的数。所以(f(n,k)) 要从 (f(n,frac{k}{p})) 中去掉一部分：

[egin{align} f(n,k)&=sum_{i=1}^{n}[i ot k]=sum_{i=1}^n[i ot frac{k}{p}]-sum_{i=1}^n[iot frac{k}{p}][p|i]\ &=f(n,frac{k}{p})-sum_{i=1}^{lfloorfrac{n}{p} floor}[ipot frac{k}{p}] end{align} ]

因为 (p) 不存在于 (frac{k}{p}) 中，则有 ([ip ot frac{k}{p}]=[i ot frac{k}{p}]) ，则：

[f(n,k)=f(n,frac{k}{p})-sum_{i=1}^{lfloorfrac{n}{p} floor}[iot frac{k}{p}]=f(n,frac{k}{p})-f(lfloorfrac{n}{p} floor,frac{k}{p}) ]

类似的，有：

[egin{align} g(n,k)&=sum_{i=1}^nmu(i)[iot k]=sum_{i=1}^nmu(i)[iot frac{k}{p}]-sum_{i=1}^nmu(i)[iot frac{k}{p}][p|i]\ &=g(n,frac{k}{p})-sum_{i=1}^{lfloorfrac{n}{p} floor}mu(ip)[iot frac{k}{p}] end{align} ]

对于莫比乌斯函数，有：

[mu(ij)=mu(i)mu(j)[iot j] ]

代入上式得：

[egin{align} g(n,k)&=g(n,frac{k}{p})-sum_{i=1}^{lfloorfrac{n}{p} floor}mu(ip)[iot frac{k}{p}]\ &=g(n,frac{k}{p})-sum_{i=1}^{lfloorfrac{n}{p} floor}mu(i)mu(p)[iot p][iot frac{k}{p}]\ &=g(n,frac{k}{p})+sum_{i=1}^{lfloorfrac{n}{p} floor}mu(i)[iot k]\ &=g(n,frac{k}{p})+g(lfloorfrac{n}{p} floor,k) end{align} ]

递归求一下 (f,g) ，这道题就被完美解决了。

---

( ext{Code}:)

#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <unordered_map>
#define maxn 10000005
#define Rint register int
#define INF 0x3f3f3f3f
using namespace std;
typedef long long lxl;

template <typename T>
inline T read()
{
	T x=0,f=1;char ch=getchar();
	while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
	return x*f;
}

int prime[maxn],cnt;
bool flag[maxn];
lxl mu[maxn];
lxl N,M,K,pri[maxn],prc;

inline void sieve()
{
	mu[1]=1;
	for(int i=2;i<maxn;++i)
	{
		if(!flag[i]) prime[++cnt]=i,mu[i]=-1;
		for(int j=1;j<=cnt&&i*prime[j]<maxn;++j)
		{
			flag[i*prime[j]]=true;
			if(!(i%prime[j])) break;
			mu[i*prime[j]]=-mu[i];
		}
	}
	for(int i=1;i<maxn;++i)
		mu[i]+=mu[i-1];
	for(int i=1;i<=cnt;++i)
	{
		while(!(K%(prime[i]*prime[i]))) K/=prime[i];
		if(!(K%prime[i])) pri[++prc]=prime[i];
	}
}

unordered_map<int,lxl> mp;

inline lxl G1(int n)
{
	if(n<maxn) return mu[n];
	if(mp[n]) return mp[n];
	lxl res=1;
	for(lxl l=2,r=0;l<=n;l=r+1)
	{
		r=n/(n/l);
		res-=(r-l+1)*G1(n/l);
	}
	return mp[n]=res;
}

unordered_map<int,unordered_map<int,lxl> > F,G;

inline lxl f(int n,int k,int p=1)
{
	if(!n) return 0;
	if(k==1) return 1ll*n;
	if(F[n][k]) return F[n][k];
	return F[n][k]=f(n,k/pri[p],p+1)-f(n/pri[p],k/pri[p],p+1);
}

inline lxl g(int n,int k,int p=1)
{
	if(!n) return 0;
	if(k==1) return G1(n);
	if(G[n][k]) return G[n][k];
	return G[n][k]=g(n,k/pri[p],p+1)+g(n/pri[p],k,p);
}

int main()
{
#ifndef ONLINE_JUDGE
	freopen("P1587.in","r",stdin);
#endif
	N=read<int >(),M=read<int >(),K=read<int >();
	sieve();
	lxl res=0;
	int up=min(N,M);
	for(lxl l=1,r=0;l<=up;l=r+1)
	{
		r=min(N/(N/l),M/(M/l));
		res+=1ll*(N/l)*f(M/l,K)*(g(r,K)-g(l-1,K));
	}
	printf("%lld
",res);
	return 0;
}

---

参考文献：

神犇柳苏明的题解。

11Dimensions的课件。