BZOJ 4314 倍数？倍数！

发表于 | 分类于

Description

要求你从\(Z_n\)（模\(n\)的剩余系）里选出\(k\)个不重复元素，使得他们的和模\(n\)为0。求方案数。

\(1\leq n\leq 10^9, 1\leq k\leq 10^3\)。

Solution

看上去就很像单位根反演的题？

先放一下众所周知的单位根反演的式子： \[ \frac{1}{n}\sum_{i = 0}^{n - 1}(\xi_{n}^i)^k=[n | k] \] 证明很容易：假如有\(n | k\)，那么和式中每一项都是1，加起来除\(n\)后就事1了（这种情况下等比数列公比为1，所以不可以采用等比数列求和公式）；反之，我们采用等比数列求和公式可知原式中的和式为\(\frac{1 - \xi_n^{kn}}{1 - \xi_n^k}\)，这个东西的分子为0，那么自然和式左边为0。

那么考虑构造二元答案多项式\(f(x, y)\)： \[ f(x, y) = \prod_{i = 0}^{n - 1}(1 + x^iy) \] 这个多项式的每一项中\(x\)的次数表示和，\(y\)表示选数的数量。我们显然要求所有含\(y^k\)且\(x\)的次数为\(n\)的倍数的项的系数之和，那么考虑将\(y\)视为常数，对\(x\)进行单位根反演： \[ f(\xi_n^t, y)=\prod_{i = 0}^{n - 1}(1 + (\xi_n^{t})^iy) \] 观察到如果\(t\)和\(n\)不互质，那么两者可以同时约去一约数（反之，若有\(t\perp n\)，那么称\(\xi_n^t\)为一个\(n\)次本原单位根）。假设约完之后的单位根为\(\xi_a^b\)，那么我们发现当积式中枚举的\(i\)大于\(a\)时，会出现循环。所以我们只需要取循环节（显然长度为\(a\)，因为我们知道所有\(a\)阶单位根构成一个循环群，而\(a\)阶本原单位根一定是该群的生成元。原因事首先显然用\(\xi_a^1\)可以生成所有\(a\)阶单位根，所以所有\(a\)阶单位根构成一个循环群；然后你根据数论里的欧拉定理可以知道\(x^{\phi(p)}\equiv 1\pmod{p}(x\perp p)\)，那么对于任意本原单位根\(\xi_a^b\)有\(\xi_a^{b^{\phi(a)}} = \xi_a^1\)，式子左边的东西是\(\xi_a^b\)的若干次方，所以用\(\xi_a^b\)可以生成\(\xi_a^1\)，自然就可以生成整个循环群）的若干次方就行了，那么有： \[ f(\xi_a^b, y) = (\prod_{i = 0}^{a - 1}(1 + (\xi_a^b)^iy))^\frac{n}{a} \]

根据上面的说法，考虑里面的积式一定可以枚举到所有\(a\)阶单位根，所以这个积式也可以写成： \[ f(\xi_a^b, y) = (\prod_{i = 0}^{a - 1} (1 + \xi_a^iy))^\frac{n}{a} \] 这么一来，我们发现对于任意\(a\)阶本原单位根\(\xi_a^b\)，对答案的贡献都是一样的，而我们知道\(a\)阶本原单位根有\(\phi(a)\)个。只要我们对所有\(n\)的约数求出其所有本原单位根的贡献的和，那么这题就做完了。

那么我们考虑，如果我们知道了\(a\)，那么怎么快速求一个本原单位根的贡献呢？

考虑那个多项式（里面的积式），首先很容易发现其0次项为1。然后我们从零点入手……

然后我们发现这个多项式几乎没法求零点……那么我们做一步小转换吧： \[ \prod_{i = 0}^{a - 1}\xi_a^i(\xi_a^{a - i} + y) \]

然后我们发现\(\xi_a^{a - i}\)就枚举了所有单位根，所以零点集合就是所有单位根相反数的集合。

如果\(a\)为偶数，那么把所有单位根取相反数之后得到的集合其实和原集合事一样的……这个大致可以理解为所有单位根构成的图形关于原点中心对称。而我们知道\(x^a - 1 = 0\)的解就是全体\(a\)阶单位根，因此答案多项式和\(x^a - 1\)只有常数项的不同，而取反之后得到了\(1 - x^a\)，其零点、常数项都符合我们的要求。

如果\(a\)为奇数，那么考虑\(x^a\)，往里面带一个单位根的相反数的话，因为\(a\)为奇数所以会得到\(-1\)，所以该多项式的零点就是\(x^a + 1 = 0\)的全体解。而左边多项式的常数项也符合我们的要求，因此此时多项式为\(1 + x^a\)。

综上可得： \[ \prod_{i = 0}^{a - 1}(1 + \xi_a^iy) = 1 - (-y)^a \] 而最后的答案多项式为： \[ \frac{1}{n}\sum_{d | n}\phi(d)(1 - (-y)^d)^\frac{n}{d} \] 这个东西的\(y^k\)项的系数就是答案了……

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
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <functional>
#include <utility>
typedef long long ll;
const ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
  ll ans = 1LL, res = a;
  while(b) {
    if(1LL & b) ans = (ans * res) % ha;
    res = (res * res) % ha; b >>= 1;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}
 
ll sqrt(ll n) {
  ll ans = 1;
  while((ans + 1LL) * (ans + 1LL) <= n) ans ++;
  return ans;
}
int phi(int n) {
  int ans = n, m = sqrt(n);
  for(int i = 2; i <= m; i ++) {
    if(n % i == 0) {
      ans = (ans / i) * (i - 1);
      while(n % i == 0) n /= i;
    }
  }
  if(n > 1) ans = (ans / n) * (n - 1);
  return ans;
}
ll C(int n, int m) {
  if(n < m) return 0;
  ll ans = 1;
  for(int i = 1; i <= m; i ++) {
    ans = (ans * (ll(n - i + 1))) % ha;
    ans = (ans * inv(i)) % ha;
  }
  return ans;
}
ll query(int d, int nd, int k) {
  if(k % d != 0) return 0;
  ll ret = phi(d);
  if(d % 2 == 0 && (k / d) % 2 != 0) {
    ret = (ha - ret) % ha;
  }
  ret = (ret * C(nd, k / d)) % ha;
  return ret;
}
ll calc(int n, int k) {
  int m = sqrt(n); ll ans = 0;
  for(int i = 1; i <= m; i ++) {
    if(n % i != 0) continue;
    ans = (ans + query(i, n / i, k)) % ha;
    if(i * i != n) ans = (ans + query(n / i, i, k)) % ha;
  }
  ans = (ans * inv(n)) % ha;
  return ans;
}
 
int main() {
  int n, k; scanf("%d%d", &n, &k);
  printf("%lld\n", calc(n, k));
  return 0;
}