「NOI 2016」循环之美

Posted on | Edited on | In , ,

Problem

Description

给定三个十进制数 n, m, k ，求在 k 进制下，有多少个数值互不相等的纯循环小数可以用分数 \displaystyle\frac{x}{y} \left(x, y \in \mathbb{N^{+}}, x \le n, y \le m\right) 表示。
一个数是纯循环的，当且仅当其可以写成以下形式：

a.\dot{c_1} c_2 c_3 \ldots c_{p - 1} \dot{c_p}

其中， a 是一个整数， p \ge 1 ；对于 1 \le i \le p ， c_i 是 k 进制下的一位数字。

需要特别注意的是，我们认为一个整数是纯循环的，因为它的小数部分可以表示成 0 的循环或是 k-1 的循环。

Constraints

1 \le n, m \le 10^{9}, 2 \le k \le 2000

Solution

Analysis

对于 x \perp y ，观察 \displaystyle\frac{x}{y} 进制转换的过程，若其为纯循环小数，设其循环节长度为 l ，则有：

\begin{aligned} z &= \frac{x}{y} k^{l} - \frac{x}{y} \in \mathbb{Z} \\ xk^{l} - x &= yz \\ xk^{l} &\equiv x \pmod{y} \\ k^{l} &\equiv 1 \pmod{y} \end{aligned}

那么 k^{l - 1} 即为 k 在模 y 意义下的逆元，而这是 k \perp y 的充要条件。

那么答案即为：

\begin{aligned} \sum_{i = 1}^{n} \sum_{j = 1}^{m} [i \perp j] [j \perp k] &= \sum_{i = 1}^{m} [i \perp k] \sum_{j = 1}^{n} [i \perp j] \\ &= \sum_{d = 1}^{\min(n, m)} [d \perp k] \mu(d) (n / d) \sum_{id \le m} [i \perp k] \end{aligned}

现在问题转化为求 [n \perp k] \mu(n) 和 [n \perp k] 的前缀和。设：

\begin{aligned} f(n) &= [n \perp k] \\ F(n) &= \sum_{i = 1}^{n} f(i) \\ g(n) &= [n \perp k] \mu(n) \\ G(n) &= \sum_{i = 1}^{n} g(i) \end{aligned}

对于 F(n) ，由于 \gcd(x, y) = \gcd(y, x \bmod{y}) ，有：

\begin{aligned} F(n) &= \sum_{i = 1}^{n} [i \perp k] \\ &= \sum_{i = 1}^{n} [(i \bmod{k}) \perp k] \\ &= (n / k) F(k) + F(n \bmod{k}) \end{aligned}

预处理 F(1), F(2), \dots, F(k) 即可 O(1) 回答，时间复杂度 O(k) - O(1) 。

对于 G(n) ，发现 [n \perp k] 一项难以处理，考虑将其消掉。
观察到 [x \perp z] [y \perp z] = [xy \perp z] ，故可以将 g 卷上 f 来杜教筛：

\begin{aligned} (g \ast f)(n) &= \sum_{d | n} [d \perp k] \mu(d) \left[\frac{n}{d} \perp k\right] \\ &= \sum_{d | n} [n \perp k] \mu(d) \\ &= [n \perp k] \epsilon(n) \\ &= \epsilon(n) \\ f(1)G(n) &= \sum_{i = 1}^{n} (g \ast f)(i) - \sum_{i = 2}^{n} f(i) G(n / i) \\ G(n) &= 1 - \sum_{i = 2}^{n} f(i) G(n / i) \end{aligned}

时间复杂度 O\left(n^{\frac{2}{3}} + k\right) 。

Code

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
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
#include <cstdio>

using i64 = long long;

constexpr int maxk = 2000;
constexpr int maxs = 1000000; // n^{2 / 3}

template <typename _Tp>
inline _Tp Min(const _Tp &x, const _Tp &y) {
	return x < y ? x : y;
}

inline int io() {
	static int v;
	return scanf("%d", &v), v;
}

int pd[6], dcnt;

inline void divide(int k) {
	for (int i = 2; i <= k; ++i)
		if (k % i == 0) {
			pd[++dcnt] = i;
			do k /= i; while (k % i == 0);
		}
}

inline bool perp(const int &n) {
	for (int i = 1; i <= dcnt; ++i)
		if (n % pd[i] == 0)
			return false;
	return true;
}

int pri[maxs / 7], pcnt,
	lpf[maxs + 1],
	f[maxk + 1],
	g[maxs + 1];

inline void precalcF(const int &k) {
	f[1] = 1;
	for (int i = 2; i <= k; ++i)
		f[i] = f[i - 1] + perp(i);
}

inline void sieve(const int &n) {
	g[1] = 1; i64 v;
	for (int i = 2; i <= n; ++i) {
		if (lpf[i] == 0)
			pri[lpf[i] = ++pcnt] = i, g[i] = -1;
		for (int j = 1; j < lpf[i] && (v = (i64)i * pri[j]) <= n; ++j)
			lpf[v] = j, g[v] = -g[i];
		if ((v = (i64)i * pri[lpf[i]]) <= n)
			lpf[v] = lpf[i];
		perp(i) ? g[i] += g[i - 1] : (g[i] = g[i - 1]);
	}
}

const int global_n = io(), m = io(), k = io();

inline int F(const int &n) {
	return (n / k) * f[k] + f[n % k];
}

namespace hash {

constexpr int hashsz = 1048576;

struct hashlnk {
	int key, val; hashlnk *las;
	inline hashlnk* init(const int &ik, const int &iv, hashlnk *const &ls) {
		return key = ik, val = iv, las = ls, this;
	}
} *las[hashsz];

inline int& find(const int &key) {
	static hashlnk pool[maxs], *alloc_p = pool - 1;
	const int pos = (key & 1048575);
	for (hashlnk *o = las[pos]; o; o = o->las)
		if (o->key == key)
			return o->val;
	las[pos] = (++alloc_p)->init(key, 0, las[pos]);
	return las[pos]->val;
}

}

inline int G(const int &n) {
	if (n <= maxs) return g[n];
	int &mem = hash::find(n);
	if (mem) return mem;
	int ans = 1;
	for (int i = 2, j; i <= n; i = j + 1) {
		j = n / (n / i);
		ans -= (F(j) - F(i - 1)) * G(n / i);
	} return mem = ans;
}

int main() {
	const int lim = Min(global_n, m);

	divide(k);
	precalcF(k);
	sieve(Min(lim, maxs));

	i64 ans = 0;
	for (int i = 1, j; i <= lim; i = j + 1) {
		j = Min(global_n / (global_n / i), m / (m / i));
		ans += (i64)(G(j) - G(i - 1)) * (global_n / i) * F(m / i);
	} printf("%lld", ans);

	return 0;
}