问题重述

我们有一个未知的正整数集合 S ⊆ {1,2,...,n}，对于任意 x ≤ n，f(x) 表示用 S 中的元素（可重复使用）组成 x 的方案数。现在已知所有 f(i) mod p 的值，要恢复出 S。

数学基础

1. 生成函数表示

设集合 S 的生成函数为：

$$F(x) = \sum_{i=0}^\infty f(i)x^i$$

其中 f(0) = 1（空集表示0，有1种方式）。

关键定理：如果 S = {s₁, s₂, ..., sₘ}，那么：

$$F(x) = \prod_{s \in S} \frac{1}{1 - x^s}$$

证明：

• 每个元素 s ∈ S 可以出现 0,1,2,... 次
• 对应的生成函数是 $\frac{1}{1-x^s} = 1 + x^s + x^{2s} + \cdots$
• 所有元素的生成函数相乘就得到总的生成函数

2. 取对数变换

对两边取对数：

$$\ln F(x) = \ln\left(\prod_{s \in S} \frac{1}{1 - x^s}\right) = -\sum_{s \in S} \ln(1 - x^s) $$

利用展开式 $\ln(1 - y) = -\sum_{k=1}^\infty \frac{y^k}{k}$：

$$\ln F(x) = \sum_{s \in S} \sum_{k=1}^\infty \frac{x^{ks}}{k} $$

交换求和顺序：

$$\ln F(x) = \sum_{m=1}^\infty \left(\sum_{s \in S, s|m} 1\right) \frac{x^m}{m} $$

令 $c_m = |\{s \in S : s \mid m\}|$，即能被 S 中某个元素整除的 m 的计数，那么：

$$\ln F(x) = \sum_{m=1}^\infty c_m \frac{x^m}{m}$$

3. 莫比乌斯反演

现在我们已知 F(x)（即 f(i)），通过多项式对数可以求出 $\ln F(x)$ 的系数。设：

$$g_m = [x^m] \ln F(x)$$

那么有：

$$g_m = \frac{c_m}{m} \quad \Rightarrow \quad c_m = m \cdot g_m $$

现在我们知道每个 m 的 c_m 值，需要恢复出 S。

关键观察：c_m 是 S 中能整除 m 的元素个数。我们想要知道的是：哪些数在 S 中？

使用莫比乌斯反演：

$$a_m = \sum_{d|m} \mu(d) c_{m/d}$$

其中 μ 是莫比乌斯函数。可以证明：a_m > 0 当且仅当 m ∈ S。

算法实现细节

1. 多项式对数计算

给定 f[0..n]，要计算 g = ln(f)。

公式：如果 G(x) = ln F(x)，那么：

$$F'(x) = G'(x) F(x) \quad \Rightarrow \quad G'(x) = \frac{F'(x)}{F(x)} $$

算法步骤：

1. 计算 f 的导数：d[i] = (i+1) * f[i+1]
2. 计算 f 的逆：inv_f = 1/f
3. 计算 g' = d * inv_f
4. 积分得到 g：g[i] = g'[i-1] / i

il void ln(int *x, int n, int *res) {
    int inv[maxn], d[maxn];
    // 计算导数
    For(i, 1, n) d[i - 1] = 1LL * i * x[i] % p;
    memset(inv, 0, sizeof inv);
    // 计算逆
    inverse(x, n + 1, inv);
    // 相乘得到g'
    mul(d, n, inv, n, res, n);
    // 积分得到g
    Rep(i, n, 1) res[i] = 1LL * res[i - 1] * power(i, p - 2) % p;
}

2. 多项式求逆

使用牛顿迭代法：

• 如果已知 $F(x)G(x) = 1 \pmod{x^m}$
• 那么 $G_{2m}(x) = 2G_m(x) - F(x)G_m^2(x) \pmod{x^{2m}}$

il void inverse(int *x, int y, int *res) {
    if (y == 1)
        return res[0] = power(x[0], p - 2), void();
    
    inverse(x, (y + 1) >> 1, res);
    mul(x, y - 1, res, y - 1, tmp, y - 1);
    tmp[0] = 2 - tmp[0] + p;
    For(i, 1, y - 1) tmp[i] = p - tmp[i];
    mul(tmp, y - 1, res, y - 1, res, y - 1);
}

3. 快速傅里叶变换（FFT）

用于加速多项式乘法。由于在模 p 下运算，使用三次FFT技巧：

将系数 a 拆分为 $a = a_1 \cdot M + a_2$，其中 $M = \lfloor \sqrt{p} \rfloor$

这样可以将大数乘法转化为多个小数乘法，避免精度问题。

4. 莫比乌斯反演

For(i, 1, n) for (int j = i; j <= n; j += i)
    a[j] += mu[i] * g[j / i];

这里 g 已经乘以了对应的系数，所以直接使用。

复杂度分析

• 多项式求逆：O(n log n)
• 多项式对数：O(n log n)
• 多项式乘法：O(n log n)
• 莫比乌斯反演：O(n log n)

总复杂度：O(n log n)

正确性证明

为什么这样能恢复 S？

设 S 是我们要找的集合。定义：

• $c_m = |\{s \in S : s \mid m\}|$
• 那么 $\ln F(x) = \sum_{m=1}^\infty c_m \frac{x^m}{m}$

我们计算出的 g_m 应该满足：$m \cdot g_m = c_m$

现在考虑莫比乌斯反演：

$$a_m = \sum_{d|m} \mu(d) c_{m/d}$$

代入 $c_{m/d} = \sum_{s \in S, s|(m/d)} 1$：

$$a_m = \sum_{d|m} \mu(d) \sum_{s \in S, s|(m/d)} 1 = \sum_{s \in S} \sum_{d|(m/s)} \mu(d) $$

根据莫比乌斯函数性质：$\sum_{d|k} \mu(d) = 1$ 当 k=1，否则为0。

所以：

$$a_m = \sum_{s \in S} [m/s = 1] = \sum_{s \in S} [s = m] $$

即：$a_m = 1$ 当 m ∈ S，否则为 0。

字典序最小的保证

由于我们按从小到大的顺序遍历 1 到 n，当发现 a[i] ≠ 0 时就输出 i，这自然保证了字典序最小。

样例详细验证

样例1：S = {1,2,3,4,5}

生成函数：

$$F(x) = \frac{1}{(1-x)(1-x^2)(1-x^3)(1-x^4)(1-x^5)} $$

计算 f(1)=1, f(2)=2, f(3)=3, f(4)=5, f(5)=7 与题目一致。

样例2：S = {1,2,5}

生成函数：

$$F(x) = \frac{1}{(1-x)(1-x^2)(1-x^5)}$$

计算 f(4)=3, f(5)=4, f(6)=5, f(7)=6 与题目一致。

这个算法通过严谨的数学推导，高效地解决了集合恢复问题，保证了正确性和最优性