[ASTRA] 合数模下的区间线性基：带时间标记的 p-进阶梯表

· · 算法·理论

本文介绍 GPT-6-ASTRA 完成的一项算法研究成果。

完整代码、证明、测试和实验结果已经整理在仓库：ASTRA4OI。

结论

可以做到接近普通区间线性基的复杂度。设已知

m=\prod_{s=1}^{w}p_s^{k_s},\qquad K=\sum_s k_s.
模式 预处理的主要模运算量 单次查询的主要模运算量 基结构空间
单个素数幂，按右端点离线 O(nd^2k) O(d^2) O(d^2k)
一般整数模数，按右端点离线 O(nd^2K) O(wd^2) O(d^2K)
保存所有前缀，在线回答任意历史区间 同上 同上 O(nd^2K)

特别地，m=p² 只需每列两个槽位，与普通域上线性基有相同的主要渐近复杂度。k=1 时算法退化成普通的带最新位置的线性基。

本文介绍一个合数模下的区间线性表示算法。

给定向量序列

a_i\in(\mathbb Z/m\mathbb Z)^d,

需要回答：

x\in\left\langle a_l,a_{l+1},\ldots,a_r\right\rangle_{\mathbb Z/m\mathbb Z}

是否成立，也就是是否存在模 m 的系数 \lambda_i，使得

x=\sum_{i=l}^{r}\lambda_i a_i.

本文对应的成果目录为：composite-modulus-basis。

其中包括：

一、为什么普通线性基不能直接推广

在域上，线性基的核心操作是用主元消去某一列。

但是 \mathbb Z/m\mathbb Z 在 m 为合数时通常不是域，存在零因子，非零元素也可能不可逆。

例如模 4 下，只有一个向量

a=(2,1).

普通的消元会把它作为第一列主元，但是

2a=(0,2)

又在第二列产生了新的信息。

如果只保留一个主元，就会漏掉向量 (0,2)。

另一个问题是，不能只保存每一列最新的原始向量。

考虑模 4 的一维序列：

a_1=1,\qquad a_2=2.

整个区间 [1,2] 可以表示 1，但区间 [2,2] 只能表示 0 和 2。

因此：

对于模 p^k，一维序列

1,p,p^2,\ldots,p^{k-1}

会产生 k 个不同的后缀层次。

所以，在合数模下，每一列不能只保存一个主元，而需要保存多个 p-进层次。

二、先考虑模 p^k

先令

R=\mathbb Z/p^k\mathbb Z,

其中 p 是素数。

所有坐标都规范化到

0,1,\ldots,p^k-1.

对于非零元素 x，定义

v_p(x)=\max\{v:p^v\mid x\}.

例如模 8 下：

v_2(1)=0,\quad v_2(2)=1,\quad v_2(4)=2.

对一个非零向量 b，设它第一个非零坐标的位置为 j，并令

v=v_p(b_j).

定义它的首项位置为

\operatorname{lead}(b)=(j,v).

首项位置按照先比较坐标 j，再比较 p-进赋值 v 的字典序排列。

如果一个向量的第一个非零坐标是 p^v u，其中 u 与 p 互质，那么可以乘以单位 u^{-1}，把主元归一化成 p^v。

注意，这里只需要对单位求逆，绝不会对非单位 p^v 求逆。

三、p-进阶梯表

维护一个二维表：

B[j][v],

其中：

因此总共有 dk 个槽位。

一个槽位中保存一行：

(b,t,e),

其中：

表中行需要满足：

  1. 前 j 个坐标为零；
  2. 第 j 个坐标为 p^v；
  3. 该行的首项位置为 (j,v)。

时间戳的含义是：这条行只使用了位置不小于 t 的原始向量。

因此，处理完右端点 r 后，区间 [l,r] 只需要使用时间戳满足

t\ge l

的槽位。

四、p-数字基

这里的“基”不是通常意义上的自由模基。

我们称一组行是一个 p-数字基，如果其中每个元素都可以唯一写成

\sum_b c_b b, \qquad c_b\in\{0,1,\ldots,p-1\}.

每个系数只取一个 p 进制数字。

引理一：数字组合的唯一性

设不同的行具有不同的首项位置。

考虑一个非零组合：

\sum_b c_b b=0, \qquad c_b\in\{0,1,\ldots,p-1\}.

从所有 c_b\ne0 的行中，选取首项位置最小的那一行，记它的首项为 (j,v)。

在第 j 列：

由于

1\le c_b\le p-1,

所以 c_b p^v 不被 p^{v+1} 整除，无法被其他项抵消。

矛盾。

因此，不同的 p-数字组合一定得到不同的向量。

引理二：使用阶梯表进行成员查询

假设当前表是某个子模 H 的 p-数字基。

给定目标向量 x，从左到右处理每个非零坐标。

若当前第 j 列的 p-进赋值为 v，就寻找槽位

B[j][v].

如果槽位不存在，则 x 不属于 H。

如果槽位存在，设主元行是 b，则令

c=\frac{x_j}{p^v}.

由于 x_j 和 p^v 都是普通整数，且 p^v\mid x_j，这里不需要求逆。

执行：

x\leftarrow x-cb\pmod {p^k}.

这样第 j 列会被完全清零。

每次操作都会清空当前最左侧非零列，因此最多进行 d 次消元。

如果最终全部清零，就得到了一份表示 x 的线性组合。

反之，如果某一步找不到对应槽位，则根据 p-数字基的唯一性，x 不可能属于当前子模。

所以成员查询是正确的。

五、加入一个新向量会产生什么

设当前已经有一个较新的子模 H，现在加入向量 a。

考虑商群：

(H+Ra)/H.

它由 a+H 生成，因此是一个循环 p-群。

设它的阶为

p^h.

那么：

a+H,\ pa+H,\ p^2a+H,\ldots,p^{h-1}a+H

全部非零，而

p^h a\in H.

对这些向量分别相对于 H 做消元。

得到的非零余行会具有严格递增的首项位置，并且可以写成：

b_e\equiv u_e p^e a\pmod H, \qquad 0\le e<h,

其中 u_e 是单位。

因此，加入一个向量时，恰好增加 h 个新的 p-数字层。

这 h 行的数字组合一共有

p^h

种，正好覆盖循环商群

(H+Ra)/H

的全部元素。

于是：

如果 H 已经有 p-数字基，那么加入一个向量 a 后，保留 a,pa,\ldots,p^{k-1}a 中相对于 H 非零的部分，就可以得到 H+Ra 的 p-数字基。

这就是整个算法的核心。

六、时间戳算法

设当前要加入的是第 r 个向量 a_r。

首先将以下至多 k 条任务加入优先队列：

(a_r,r,0), (pa_r,r,1), (pa_r,r,2), \ldots, (p^{k-1}a_r,r,k-1).

队列按照以下顺序处理：

  1. 时间戳较大的优先；
  2. 时间戳相同时，幂层较小的优先。

每条任务还带有一个起始列 start，表示从这一列开始继续寻找主元。

伪代码如下：

while 队列非空：
    取出时间戳最大的任务 (row, t, e, start)

    从 start 开始寻找 row 的第一个非零坐标 j
    如果不存在这样的 j：
        丢弃该任务
        continue

    v = valuation_p(row[j])
    old = B[j][v]

    如果 old 不存在：
        将 row 乘以单位，归一化为主元 p^v
        写入 B[j][v]
        结束当前任务

    如果 old 的时间戳大于 t：
        用 old 消去 row 的第 j 列
        从 j+1 继续处理 row

    如果 t 大于 old 的时间戳：
        将 row 归一化并替换 old
        用新 row 消去 old
        将 old 的余行以旧时间戳重新加入队列
        结束当前任务

同一时间戳的不同幂层不会占用同一个槽位，这是上一节循环商群结论保证的。

被替换掉的旧行必须重新进入队列，不能直接沿着当前流程继续处理。因为所有较新的时间戳必须先稳定下来，旧行才能相对于新的较新子模完成消元。

七、正确性证明

设

H_{t+1}^{(r)} = \left\langle a_{t+1},a_{t+2},\ldots,a_r \right\rangle_R.

处理完前缀 a_1,\ldots,a_r 后，维护以下不变量。

不变量一

对于任意左端点 l，所有时间戳满足 t\ge l 的表中行，构成区间

[a_l,a_{l+1},\ldots,a_r]

生成子模的 p-数字基。

不变量二

设 a_t 在商群

R^d/H_{t+1}^{(r)}

中的阶为

p^{h_t}.

那么时间戳恰好为 t 的行有且仅有 h_t 条，对应幂层

e=0,1,\ldots,h_t-1.

并且每一条都满足：

b_{t,e} \equiv u_{t,e}p^e a_t \pmod {H_{t+1}^{(r)}},

其中 u_{t,e} 是单位。

归纳证明

当表为空时，不变量显然成立。

现在假设插入新向量后，所有时间戳大于 t 的行已经稳定。

这些行生成了新的较新子模：

H=H_{t+1}^{(r)}.

考虑时间戳为 t 的行。

在插入新向量之前，时间 t 的行对应于旧的商群层次。

加入新的较新向量后，子模 H 变大，因此商群中的阶不会变大。

如果新的阶指数为 h_{\text{new}}，旧的阶指数为 h_{\text{old}}，则：

h_{\text{new}}\le h_{\text{old}}.

对于幂层 e\ge h_{\text{new}}，有：

p^e a_t\in H,

所以这些层最终会被较新的行完全消去。

对于 e<h_{\text{new}}，对应的陪集仍然非零，并且会在某个 H 缺少的首项位置停止。

如果该位置原本有更旧时间戳的行，那么新行替换旧行，旧行的余量重新进入队列。

由于队列始终优先处理较新的时间戳，时间 t 的所有任务都会在所有更旧任务之前稳定。

根据前面的循环商群引理，时间 t 最终恰好保留：

e=0,1,\ldots,h_{\text{new}}-1

这些层，并且它们的首项位置互不相同。

因此，时间 t 的不变量成立。

从较新的时间戳向较旧的时间戳逐层归纳，得到全部不变量。

于是，处理完右端点 r 后，筛选时间戳至少为 l 的行，正好得到区间 [l,r] 的 p-数字基。

根据引理二，查询算法返回 YES 当且仅当：

x\in \left\langle a_l,\ldots,a_r \right\rangle_R.

所以算法正确。

八、复杂度证明

一次插入最开始只产生至多 k 条幂链。

一条任务链每次消元都会清空一个坐标，并从下一列继续，因此一条链最多访问 d 列。

每次访问一列需要进行一次长度为 d 的向量操作，复杂度为 O(d)。

所以单次插入的主要复杂度为：

O(kd^2).

处理全部 n 个向量：

O(nkd^2).

一次查询最多清空 d 个坐标，每次消元需要 O(d) 操作，因此：

O(d^2).

所以对于模 p^k：

操作 复杂度
全部预处理 O(nkd^2)
单次查询 O(d^2)
当前表空间 O(kd^2)

此外还包括：

这些成本与 p^k 的二进制长度有关。

九、推广到任意合数模

设

m=\prod_{s=1}^{w}p_s^{k_s},

其中各个 p_s 互不相同。

根据中国剩余定理：

\mathbb Z/m\mathbb Z \cong \prod_{s=1}^{w} \mathbb Z/p_s^{k_s}\mathbb Z.

因此，一个向量在模 m 下可表示，当且仅当它在每个素数幂分量下都可表示。

具体做法是：

  1. 对每个 p_s^{k_s} 建立一套 p-进阶梯表；
  2. 分别进行区间查询；
  3. 只有所有分量都返回 YES 时，整体才返回 YES。

充分性需要注意一个细节。

不同素数幂分量中得到的系数可能不同：

\lambda_i^{(s)}.

但是对于每个固定下标 i，可以通过中国剩余定理，把这些系数组合成一个模 m 的系数 \lambda_i。

因此，不要求不同分量先得到完全相同的系数。

令

K=\sum_{s=1}^{w}k_s.

总复杂度为：

O(nd^2K)

的预处理，以及：

O(wd^2)

的单次查询。

在线保存所有前缀快照时，空间复杂度为：

O(nd^2K).

离线按右端点处理时，基结构空间为：

O(d^2K).

因数分解本身不包含在上述复杂度中。

代码默认使用试除法，因此大模数应通过 --factors 显式提供素因子分解。

十、在线版本和离线版本

离线版本

将所有查询按照右端点 r 分桶。

从左到右扫描向量：

for r = 1 .. n:
    插入 a_r
    回答所有右端点为 r 的查询

只需要保存当前的阶梯表。

在线版本

如果需要在不断追加向量的过程中查询任意历史区间，可以在每次追加后保存一个表快照。

快照只复制槽位指针，行本身是不可变对象，因此：

当前实现支持：

不支持修改中间向量或删除向量。

十一、实现与验证

仓库中的实现只使用 Python 标准库：

solver.py

测试采用了独立的穷举生成子模作为参考答案，没有使用待测算法生成预期结果。

当前测试包括：

实际测试结果包括：

运行方式：

cd results/composite-modulus-basis

python3 test_solver.py
python3 solver.py < example.in
python3 solver.py --online < example.in
python3 solver.py --factors 2:2 < example.in

十二、适用范围和限制

本文解决的主要问题是：

(\mathbb Z/m\mathbb Z)^d

上的区间生成子模成员查询。

需要注意：

  1. 合数 m 的因数分解成本需要单独计算；
  2. --factors 提供的底数必须确实是素数；
  3. 一般有限交换环部分给出的是理论归约，不是一个通用环输入接口；
  4. 在线版本支持追加和历史查询，不支持修改和删除；
  5. 复杂度中的模运算还需要乘上大整数位复杂度。

总结

合数模线性代数和域上的线性代数有本质区别：

解决方法是：

  1. 对模 p^k 建立每列多个 p-进槽位；
  2. 使用 p-数字基描述有限生成子模；
  3. 用时间戳同时维护所有后缀；
  4. 用优先队列保证新旧行按照正确顺序消元；
  5. 最后通过中国剩余定理处理任意整数模数。

完整成果见：

https://github.com/Zi-Gao/ASTRA4OI

本文及对应算法成果由 GPT-6-ASTRA 完成。