LibreOJ 2085.「NOI2016」循环之美

首先,为了去重,我们可以增加一个限制:必须是最简分数。
然后我们看看什么情况会产生纯循环小数。

1010 进制下,当 gcd(x,y)=1,gcd(x,10)=1\gcd(x,y) = 1,\gcd(x,10) = 1 时,yx\dfrac y x 是纯循环小数。
于是我们猜测当 gcd(x,y)=1,gcd(x,k)=1\gcd(x,y) = 1,\gcd(x,k) = 1 时,yx\dfrac y xkk 进制下是纯循环小数。

证明:
yx\dfrac y xkk 进制下是纯循环小数,当且仅当 l0,yxyx=yklxyklx\exists l \neq 0,\dfrac y x - \left\lfloor\dfrac y x\right\rfloor = \dfrac {yk^l} x - \left\lfloor \dfrac {yk^l} x \right\rfloor
这玩意看着很像取模的形式,所以两边同乘 xx,变成

由于前面说了 gcd(x,y)=1\gcd(x,y) = 1,所以
于是由欧拉定理得 gcd(k,x)=1\gcd(k,x) = 1,一个特解为 l=φ(x)l = \varphi(x)

证毕。

问题变成了 i=1nj=1m[gcd(i,j)=1][gcd(j,k)=1]\sum\limits_{i = 1}^n\sum\limits_{j = 1}^m [\gcd(i,j) = 1][\gcd(j,k) = 1]

f(n)=i=1n[gcd(i,k)=1]f(n) = \sum\limits_{i = 1}^n [\gcd(i,k) = 1]
由于 gcd(i,k)=gcd(i+k,k)\gcd(i,k) = \gcd(i + k,k),有
发现 kk 非常小,暴力求出 x[0,k]x \in [0,k]f(x)f(x) 值即可。

于是有

如果数论分块的话,就剩下后面的两项。
考虑怎么算它们的前缀和。

S(n,k)=i=1n[gcd(i,k)=1]μ(i)S(n,k) = \sum\limits_{i = 1}^n [\gcd(i,k) = 1] \mu(i)

乍一看推不下去了。
但是我们反观 μ\mu 函数,根据定义容易发现当 gcd(x,y)1\gcd(x,y) \ne 1 时,μ(xy)=0\mu(xy) = 0

所以有

然后考虑递归求解 SS,容易发现 S(n,1)=i=1nμ(i)S(n,1) = \sum\limits_{i = 1}^n \mu(i),杜教筛即可。

于是最后对于 d=1min(n,m)ndf(md)[gcd(d,k)=1]μ(d)\sum\limits_{d = 1}^{\min(n,m)} \lfloor \dfrac n d \rfloor f(\lfloor \dfrac m d \rfloor) [\gcd(d,k) = 1] \mu(d) 数论分块即可。

略微卡常,少开 long long,算 SS 的时候如果系数是 00 就不用往下算。
unordered_map 居然不能用 pair,偷懒换成了 map(

代码:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
#include <cstdio>
#include <utility>
#include <algorithm>
#include <map>
#include <tr1/unordered_map>
using namespace std;
using namespace tr1;
const int N = 1e9;
const int K = 2e3;
const int CNT = 1e7;
int n,m,k;
int vis[CNT + 5],prime[CNT + 5],cnt,mu[CNT + 5];
int f[K + 5];
unordered_map<int,int> w;
map<pair<int,int>,long long> mem;
long long ans;
int calc(int n)
{
if(n <= CNT)
return mu[n];
if(w.count(n))
return w[n];
long long ret = 1;
for(register int l = 2,r;l <= n;l = r + 1)
{
r = n / (n / l);
ret -= (long long)(r - l + 1) * calc(n / l);
}
return w[n] = ret;
}
int F(int n)
{
return f[k] * (n / k) + f[n % k];
}
long long S(int n,int k)
{
if(n == 0 || n == 1)
return n;
if(k == 1)
return calc(n);
if(mem.count(make_pair(n,k)))
return mem[make_pair(n,k)];
long long ret = 0;
for(int i = 1;i * i <= k;++i)
if(!(k % i))
{
if(mu[i] - mu[i - 1])
ret += (mu[i] - mu[i - 1]) * (mu[i] - mu[i - 1]) * S(n / i,i);
if(i ^ (k / i) && mu[k / i] - mu[k / i - 1])
ret += (mu[k / i] - mu[k / i - 1]) * (mu[k / i] - mu[k / i - 1]) * S(n / (k / i),k / i);
}
return mem[make_pair(n,k)] = ret;
}
int main()
{
mu[1] = 1;
for(register int i = 2;i <= CNT;++i)
{
if(!vis[i])
mu[prime[++cnt] = i] = -1;
for(register int j = 1;j <= cnt && i * prime[j] <= CNT;++j)
{
vis[i * prime[j]] = 1;
if(!(i % prime[j]))
break;
else
mu[i * prime[j]] = -mu[i];
}
mu[i] += mu[i - 1];
}
scanf("%d%d%d",&n,&m,&k);
for(register int i = 1;i <= k;++i)
f[i] = f[i - 1] + (__gcd(i,k) == 1);
for(register int l = 1,r;l <= min(n,m);l = r + 1)
{
r = min(n / (n / l),m / (m / l));
ans += (S(r,k) - S(l - 1,k)) * (n / l) * F(m / l);
}
printf("%lld\n",ans);
}