题目：QOJ #8690 Archaeological Recovery，The 2022 ICPC World Finals Luxor。

代码已通过 QOJ #3010917：1 ms，4068 KB，64 个正式测试和 2 个额外测试全部 Accepted。算法与实现由 ChatGPT 辅助推导及编写。

记 $N=3^k$。下面给出确定性的 $O((n+k)N)$ 时间、$O(N+n+k)$ 空间算法。所有计算均为整数运算，无浮点数或哈希碰撞。先恢复每个无向操作 $\{v,-v\}$ 的数量，再用逐因子反卷积贪心选择方向，无需方向子集和 DP。

与已有方法的关系

这不是首次对本题使用 Fourier 变换。需要区分以下公开基线：

因此，相对赛事解题分析，本文消除了 $N^2$ 项；相对已检查的快速变换实现，方向恢复消除了时间中的 $k$ 因子和空间中的 $n$ 因子，并给出不依赖模数碰撞假设的构造证明。Fourier、无向类可识别性及翻转与平移的关系已有理论基础，见 Glaudo–Kravitz 和 Ciprietti–Glaudo。本文不主张这些基础事实或该复杂度在所有未公开解法中首次出现。

1. 用整数表示三元 Fourier 变换

将 A,E,I 编码为 $0,1,2$，将操作 0,+,- 编码为 $0,1,2$。所有状态构成加法群 $G=\mathbb F_3^k$。输入 $C[x]$ 是 $n$ 个未知操作向量 $u_1,\ldots,u_n$ 的子集和重数。

取 $\omega^2+\omega+1=0$。定义

$$F[t]=\sum_{x\in G}C[x]\omega^{t\cdot x}=\prod_{j=1}^n(1+\omega^{t\cdot u_j}).$$

每个数用整数对 $(a,b)$ 表示 $a+b\omega$。乘以 $\omega$ 和 $\omega^2$ 分别为

$$ (a,b)\mapsto(-b,a-b),\qquad(a,b)\mapsto(b-a,-a). $$

因此每一维的蝶形

$$ (X,Y,Z)\mapsto(X+Y+Z,\ X+\omega Y+\omega^2Z,\ X+\omega^2Y+\omega Z) $$

只需常数次整数加减。对 $k$ 个三元位各做一层，耗时 $O(kN)$。

设 $z[t]$ 是满足 $t\cdot u_j=0$ 的操作数，$S=\sum_j u_j$。由于 $1+\omega^d=-\omega^{-d}$ 对 $d=1,2$ 成立，

$$F[t]=2^{z[t]}(-1)^{n-z[t]}\omega^{-t\cdot S}.$$

故每个 $F[t]$ 都是 $2^{z[t]}$ 乘以六次单位根。它的整数对一定是 $2^{z[t]}$ 乘以 $(\pm1,0),(0,\pm1),(1,1),(-1,-1)$ 之一。

所以 $z[t]$ 恰好等于 abs(a) | abs(b) 的二进制尾零数。abs(a) | abs(b) 非零，因此尾零计数有定义；不需要计算复数模长或浮点对数。

2. 第二次变换恢复无向操作

再对 $z$ 做同一个变换：

$$H[v]=\sum_{t\in G}z[t]\omega^{t\cdot v}.$$

交换求和次序，考察一个操作 $u$ 的贡献：

$$\sum_{t:t\cdot u=0}\omega^{t\cdot v}.$$

若 $u=0$，贡献仅在 $v=0$ 时为 $N$。若 $u\ne0$，由子空间上的特征标正交性，贡献在 $v\in\{0,u,-u\}$ 时为 $N/3$，其余为零。

设零操作数为 $m_0$，于是

$$m_0=\frac{H[0]/(N/3)-n}{2},$$

$$m(v)+m(-v)=\frac{H[v]}{N/3}\quad(v\ne0).$$

这些 $H[v]$ 都是普通整数，其 $\omega$ 系数为零。每对 $\{v,-v\}$ 只枚举一次，先任选一个方向，按恢复出的数量列出 $n$ 个操作。这一步只确定无向多重集。

3. 逐因子反卷积选择方向

定义平移 $(T_vD)[x]=D[x-v]$。增加操作 $v$ 等价于将剩余分布 $D$ 变成 $(I+T_v)D$。

因为 $3v=0$，有 $T_v^3=I$，从而

$$ (I+T_v)^{-1}=\frac{I-T_v+T_v^2}{2}. $$

对下一个无向操作 $\{v,-v\}$，直接计算

$$ D[x]=\frac{C[x]-C[x-v]+C[x+v]}2. $$

这是精确除法。假设当前某个合法表示包含 $v$，则 $D$ 是删除该操作后的分布；若它包含 $-v$，则 $D$ 是剩余分布的平移，因为

$$I+T_{-v}=T_{-v}(I+T_v).$$

因此 $D$ 始终是非负整数分布，但零点可能不在其支撑中。下面的引理解决这个问题。

引理： 一个合法子集和分布的平移，只要零点重数为正，就仍是同一无向操作多重集的某种合法定向产生的分布。

证明： 若 $(T_aF)[0]=F[-a]>0$，则存在操作子集 $B$ 的和为 $-a$。翻转 $B$ 中全部操作的方向，会使整个子集和多重集平移 $-\sum B=a$。这正好产生 $T_aF$。证毕。

每个合法分布都包含空集，所以当前 $C[0]>0$。又因

$$C[0]=D[0]+D[-v],$$

可作如下贪心选择：

• 若 $D[0]>0$，输出 $v$，令剩余分布为 $D$。
• 否则 $D[-v]>0$，输出 $-v$，令剩余分布为 $T_vD$。

两种选择都使新分布的零点重数为正。引理保证剩余无向操作仍存在合法定向，而且第二种选择满足

$$ (I+T_{-v})(T_vD)=(I+T_v)D=C. $$

由归纳法，每一步都精确去掉一个最终输出的因子。去掉全部 $n$ 个操作后，重数总和从 $2^n$ 变为 $1$，且零点重数为正，故剩下的正好是空集分布 $\delta_0$。把所有输出因子乘回去，就得到输入的完整分布，证明输出正确。零操作也适用：此时公式直接变成 $D=C/2$。

4. 实现与复杂度

不能在每个状态转移时重新拆 $k$ 个三元位，否则会多出一个 $k$ 因子。对固定 $v$，预先在 $O(N)$ 时间构造排列 $p[x]=x+v$：先令 $p[0]=v$，按位扩展已构造的前缀；若本位权重为 $h$、$v$ 的本位为 $d$，则

$$p[x+h]=p[x]+((d+1)\bmod3-d)h,$$ $$p[x+2h]=p[x]+((d+2)\bmod3-d)h.$$

扩展的总规模为 $1+3+\cdots+3^{k-1}=O(N)$。于是 $p[x]=x+v$，$p[p[x]]=x-v$；一轮反卷积、检查零点及必要平移共 $O(N)$。

两次变换为 $O(kN)$，全部 $n$ 轮反卷积为 $O(nN)$。输入、负向量计算与输出的 $O(kN+nk)$ 已被包含，故总时间为

$$\boxed{O((n+k)3^k)},\qquad\boxed{\text{空间 }O(3^k+n+k)}.$$

复杂度按能存储输入频数的机器字上的整数运算计。原题 $n\le40,k\le5$：第一次变换的每个部分和是互斥原始项的加权和，两个整数系数绝对值均不超过 $2^n$；第二次不超过 $nN$；反卷积临时分子绝对值不超过 $3\cdot2^n$。所有值都远小于有符号 64 位范围。

C++17 实现

以下代码与对应的 QOJ 提交一致。

#include <cassert>
#include <cstdlib>
#include <iostream>
#include <string>
#include <vector>
using namespace std;
using i64 = long long;

// Exact arithmetic in Z[w], w^2+w+1=0.
struct E {
    i64 a = 0, b = 0;
    E operator+(const E &o) const { return {a + o.a, b + o.b}; }
    E w() const { return {-b, a - b}; }
    E w2() const { return {b - a, -a}; }
};

void transform(vector<E> &a) {
    int n = (int)a.size();
    for (int h = 1; h < n; h *= 3)
        for (int b = 0; b < n; b += 3 * h)
            for (int i = b; i < b + h; ++i) {
                E x = a[i], y = a[i + h], z = a[i + 2 * h];
                a[i] = x + y + z;
                a[i + h] = x + y.w() + z.w2();
                a[i + 2 * h] = x + y.w2() + z.w();
            }
}

// All digitwise additions x+v in F_3^k, in O(3^k) time.
vector<int> translations(int v, int N) {
    vector<int> p(N);
    p[0] = v;
    for (int h = 1; h < N; h *= 3) {
        int d = v / h % 3;
        for (int x = 0; x < h; ++x) {
            p[x + h] = p[x] + ((d + 1) % 3 - d) * h;
            p[x + 2 * h] = p[x] + ((d + 2) % 3 - d) * h;
        }
    }
    return p;
}

vector<int> reconstruct(int n, int k, const vector<i64>& counts) {
    int N = 1;
    for (int i = 0; i < k; ++i) N *= 3;
    vector<E> f(N);
    for (int x = 0; x < N; ++x)
        f[x].a = counts[x];
    transform(f);
    for (E &x : f) {
        unsigned long long bits = (unsigned long long)llabs(x.a)
                                | (unsigned long long)llabs(x.b);
        assert(bits);
        x = {(i64)__builtin_ctzll(bits), 0};
    }
    transform(f);

    int zero = (f[0].a / (N / 3) - n) / 2;
    vector<int> levers(zero, 0), neg(N);
    for (int x = 1; x < N; ++x) {
        int y = x, p = 1;
        for (int j = 0; j < k; ++j, y /= 3, p *= 3)
            neg[x] += (3 - y % 3) % 3 * p;
        if (x < neg[x]) {
            assert(f[x].b == 0 && f[x].a % (N / 3) == 0);
            int count = f[x].a / (N / 3);
            for (int j = 0; j < count; ++j) levers.push_back(x);
        }
    }
    assert((int)levers.size() == n);

    // Divide out each unsigned factor and choose the sign that keeps zero
    // in the residual support. A supported translation of subset sums can
    // always be realized by flipping signs of a subset of the remaining levers.
    vector<i64> current = counts, quotient(N);
    for (int &v : levers) {
        vector<int> plus = translations(v, N);
        for (int x = 0; x < N; ++x) {
            int minus = plus[plus[x]];
            i64 twice = current[x] - current[minus] + current[plus[x]];
            assert(twice >= 0 && twice % 2 == 0);
            quotient[x] = twice / 2;
        }
        if (quotient[0] > 0) {
            current.swap(quotient);
        } else {
            // C=(1+T_v)D=(1+T_{-v})(T_v D), T_v D[x]=D[x-v].
            v = neg[v];
            for (int x = 0; x < N; ++x)
                current[x] = quotient[plus[plus[x]]];
        }
        assert(current[0] > 0);
    }
    assert(current[0] == 1);
    for (int x = 1; x < N; ++x) assert(current[x] == 0);
    return levers;
}

#ifndef ARCHAEOLOGY_LIBRARY
int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    int n, k, t;
    if (!(cin >> n >> k >> t)) return 0;
    int N = 1;
    for (int i = 0; i < k; ++i) N *= 3;
    vector<i64> counts(N);
    for (int i = 0; i < t; ++i) {
        string s;
        i64 count;
        cin >> s >> count;
        int x = 0, p = 1;
        for (int j = 0; j < k; ++j, p *= 3)
            x += (s[j] == 'A' ? 0 : s[j] == 'E' ? 1 : 2) * p;
        counts[x] = count;
    }
    for (int v : reconstruct(n, k, counts)) {
        for (int j = 0; j < k; ++j, v /= 3)
            cout << "0+-"[v % 3];
        cout << '\n';
    }
}
#endif