企业官网建设流程全解析

钣金下料这事儿，几乎是把集合划分问题原封不动搬到台面上的典型场景。一张标准板材上要切出几十个零件，每种零件有好几种摆法，每种摆法对应一块子区域和一份废料率，你要从这堆摆法里挑一组出来，让每个零件不多不少正好被切一次，同时废料最少。这个“不多不少正好一次”的硬性要求，就是集合划分问题的核心骨架，也是它比一般的背包、覆盖问题更难缠的地方。我这两年在做组合优化建模时反复回到这个结构上，原因很朴素：它的约束形式特别干净，干净到可以被 QUBO（二次无约束二值优化）完整吃下去，不需要引入任何辅助变量。这篇文章就把这条链路完整走一遍，从“问题长什么样”到“矩阵怎么填”，再到“惩罚系数怎么定”和“代码怎么落地”。适合已经了解 0/1 规划基础、想把手里的组合优化问题搬进量子退火或 QAOA 这类求解框架的读者，也适合只想搞清楚 QUBO 到底怎么建、坑在哪里的工程同学。

1. 从一张钢板的切割清单说起：集合划分到底在选什么

先把业务语言翻译成标准形式，这一步不翻译清楚，后面所有矩阵都会填错。给定一个元素集合 E，比如所有待切零件的编号；再给定一个候选子集族 S，每个子集 S_j 是 E 的一个子集，代表一种可行的摆放方案；每个子集带一个代价 c_j，可以是废料面积、可以是耗时、可以是成本。你要做的是从 S 里挑出若干个子集，使得 E 中每个元素恰好出现在一个被选中的子集里，并且总代价最小。

注意这里的关键词是“恰好”，不是“至少”。至少一次是集合覆盖问题，恰好一次是集合划分问题。这两个问题在建模上只差一个字，在解空间结构上却差得很远：覆盖问题允许冗余，解空间更宽松，往往更容易找到可行解；划分问题的可行域是一个个孤立的组合，稍不留神就会滑到不可行区域里，而且可行解之间在二值超立方体上的距离往往很远。这一点直接决定了后面惩罚系数不能给得太抠。

1.1 用“恰好覆盖一次”把业务约束翻译成数学

把上面的描述写成标准整数规划，形式是这样的：

变量 x_j ∈ {0, 1}，表示第 j 个候选子集是否被选中。

目标函数是 min Σ_j c_j x_j。

约束是对每一个元素 e ∈ E，都有 Σ_{j: e ∈ S_j} x_j = 1。

看起来平平无奇，但这个约束矩阵有一个非常有用的性质：它的每一列（对应一个子集）里的元素全是 0 和 1，没有其他系数。这个性质在后面的 QUBO 展开里会带来极大便利，因为我们不需要处理任意实数系数，只需要处理“属于”和“不属于”两种状态。很多其他组合优化问题一旦写进 QUBO，系数会变成五花八门的浮点数，展开时容易出错，而集合划分的展开结果可以用一句非常直观的话概括：对角项扣掉“覆盖了几个元素”，交叉项加上“两个子集重叠了几个元素”。

我给这个性质起了个外号叫“交集计价”。你选了方案 A 和方案 B，如果它们抢了同一个零件，就要为这个零件付一次罚款；抢得越多，罚得越狠。这个直觉比记公式可靠得多，后面写代码时我基本是拿这句话做自检的。

1.2 三个变体：等和划分、装箱、加权覆盖

集合划分不是一个孤立的问题，它有三个常见的变体，理解它们之间的差异能帮你在建模时少走弯路。

第一个是等和划分，也叫数划分。把一组数字分成两堆，让两堆的和尽量接近。这个问题的特殊之处在于它没有显式的元素约束，约束是隐含在“所有数都必须被分配”里的。它的 QUBO 形式非常经典，也是我建议所有人入门时先手推一遍的案例，因为它的推导过程能让你彻底理解“平方展开 + 二值变量平方归零”这套动作。

第二个是装箱的变体，每个箱子有容量上限，每个物品必须且只能放进一个箱子。这里的约束从“恰好覆盖一次”变成了两部分：每个物品恰好被分配一次，以及每个箱子的总重量不超过容量。后一条是不等式，进 QUBO 需要引入松弛变量，这是后面第五章要专门讲的内容。

第三个是加权覆盖，某些元素允许被覆盖多次，但要付出额外代价。这种问题的约束是软的，处理起来反而简单，因为不需要严格的惩罚，直接折算进目标函数即可。

三个变体的共同点是：它们的决策变量都是二值的，约束都是线性的，代价都是可加的。这三点凑齐，就说明它可以被 QUBO 表达。

1.3 为什么这类问题天生适合二次型表达

QUBO 的能量函数形式是 x^T Q x，展开后每一项要么是 Q_ii x_i（一次项，因为 x_i² = x_i），要么是 Q_ij x_i x_j（二次项，i < j）。它最多只能表达二阶交互，这是硬性限制。那集合划分为什么刚好落在限制之内？

因为它的约束是线性的，而惩罚一个线性约束的标准手法是把它平方。线性表达式平方之后，最高次项恰好是二次。如果约束写成 Σ x_j - 1 = 0，平方后是 (Σ x_j)² - 2Σ x_j + 1，其中 (Σ x_j)² 展开只有 x_j² 和 x_j x_k 两类项，全部落在二次以内。这就是它“天生合适”的根本原因，不是什么巧合，是数学结构决定的。

反过来看，如果约束写成“最多选一个”这类涉及到“三个里面最多选一个”的三元关系，或者约束本身带三次项，那就超出了二次型的表达能力，必须引入辅助变量做降次。降次会膨胀变量数，也会引入新的对称性，是 QUBO 建模里最不愿意做的操作之一。所以在建模初期判断一个约束能不能保持二阶，比后面调参重要得多。

还有一个不那么显眼但很实用的点：蕴含约束，比如“选了 A 就必须选 B”，可以写成 x_A ≤ x_B。这个不等式直接转成 x_A(1 - x_B) = 0，惩罚项就是 P · x_A - P · x_A x_B，一次项和二次项各一个，完全不需要任何变换。这类小技巧在集合划分的边角约束里经常用到，比如“某个方案和另一个方案互斥”，直接写 P · x_A x_B 就行，一个交叉项解决。

2. QUBO 的骨架：只有对角线和交叉项的能量函数

很多人第一次接触 QUBO 会被符号吓到，其实拆开看就那么回事。整个模型就是一张表，对角线上的数管“单独选这个变量要付多少”，非对角线上的数管“同时选这两个变量要付多少”。求解器要做的就是在这张表上找一个 0/1 组合，让总和最小。没有任何其他东西。

真正需要动脑的地方在于：你的业务问题里那些“必须”“只能”“不超过”，怎么变成这张表上的数字。这是个翻译过程，翻译得好，求解器能自己找到合法解；翻译得不好，求解器会兴高采烈地给你一个目标值很低但完全非法的答案。

2.1 二值变量的“平方归零”性质

x_i ∈ {0,1}，所以 x_i² = x_i。这一条看似废话，却是整个 QUBO 推导的基石，几乎所有展开化简都靠它。它的实际含义是：二次项 x_i² 在二值世界里退化成了线性项，所以真正携带“交互信息”的只有 i ≠ j 的交叉项。

这个性质第一次让我吃亏是在手写展开的时候。我当时推导一个覆盖约束，展开了 (Σx)² 之后得到一堆 x_i²，下意识地当成二次项写进了矩阵的非对角位置，结果跑出来的解全乱。后来才反应过来，x_i² 必须并到对角线上去。这个错误在新手里极常见，而且因为程序不会报错，只会给出一个安静的错误答案，所以特别隐蔽。

正确的处理流程是：先把所有项展开成最原始的形式，然后逐项扫描，凡是 x_i² 一律改成 x_i 并累加到对角项 Q_ii；剩下形如 x_i x_j（i < j）的，累加到 Q_ij。这个流程走一遍，就不会出错。

2.2 约束靠罚不靠加

QUBO 的全称里有个词叫“无约束”，这是它的定义决定的：求解器只认能量函数，不认约束。所以所有约束都必须以惩罚项的形式揉进目标函数，变成“违反约束就付出代价”。

这个转换带来一个非常重要的后果：可行解不再是被强制保证的，而是被“经济激励”诱导出来的。如果你的惩罚系数不够大，违反约束省下来的代价超过了罚款，求解器就会理直气壮地违反约束。这不是求解器的 bug，是你建模的问题。

我见过太多人第一次跑 QUBO 时拿到一个明显违反约束的解，然后怀疑是求解器实现有问题。实际上九成以上的情况是惩罚系数给得太小。这一点在第 3 章会展开讲，这里先建立意识：QUBO 里的约束是“软”的，软到什么程度由你给的系数决定。

2.3 覆盖约束展开后的两个系数结论

现在把集合划分的约束完整展开一遍，这是全文最核心的一段推导，值得逐行看清楚。

对每个元素 e，定义 A_e = {j : e ∈ S_j}，即所有包含 e 的候选子集的编号集合。约束是 Σ_{j∈A_e} x_j = 1。惩罚项写成 P 乘以偏差的平方：

P · (Σ_{j∈A_e} x_j - 1)²

展开平方：

(Σ x_j)² - 2Σ x_j + 1

其中 (Σ x_j)² = Σ x_j² + 2Σ_{j<k} x_j x_k。用 x_j² = x_j 替换：

= Σ x_j + 2Σ_{j<k} x_j x_k - 2Σ x_j + 1 = -Σ_{j∈A_e} x_j + 2Σ_{j<k, j,k∈A_e} x_j x_k + 1

于是单个元素的约束贡献了三部分：每个 j ∈ A_e 的对角项减去 P；每对 j < k 都在 A_e 里的交叉项加上 2P；外加一个常数 P。常数对所有解都一样，直接丢掉不用管。

把所有元素累加起来，就得到两个可以直接写进代码的结论：

其中 deg(j) 就是子集 S_j 的大小，也就是它覆盖了多少个元素。

这个结论漂亮得让人想拍桌子：你完全不需要去遍历元素、展开平方、逐项化简，只要数两个数——每个子集有几个元素，每对子集有几个共同元素。这两个数都是集合运算，一行代码的事。我在实际项目里就是拿这个结论直接建矩阵的，不再走符号推导的流程，速度快而且不容易错。

不过这里有个细节要提醒：交叉项的系数是 2P 乘以交集大小，这个 2 是从平方展开的交叉项系数 2 来的。如果你用的是对称矩阵表示（即 Q_ij = Q_ji，能量里算 2Q_ij x_i x_j 的约定），那填进矩阵的值应该是 P 乘以交集大小，不是 2P。这个约定差异是第 7 章第一个坑的主题，后面细说。

3. 惩罚系数 P：建模成败的分水岭

如果只能给一条 QUBO 建模的建议，我会说：把 80% 的调试时间花在惩罚系数上。这个数字不是夸张，我在几个不同规模的问题上试过，P 取值区间往往只有一到两个数量级是可用的，调对了求解器像开了窍，调错了怎么改算法都白搭。

P 的作用是给“违反一次约束”标价。价格定得太低，违反约束变得划算；定得太高，整个能量地形被压成一片平原，求解器找不到方向。这两头的表现完全不同，需要分别识别。

3.1 P 偏小时解合法却次优的典型症状

P 偏小的表现非常具有欺骗性：求解器给出的解，目标函数值可能很不错，但一校验就发现某个元素被覆盖了两次，或者某个元素压根没被覆盖。更气人的是，这种情况往往不是偶发，而是稳定复现——因为从能量角度看，那个非法解确实比所有合法解都便宜。

举个具体数字感受一下。假设有 20 个候选方案，单个方案代价在 1 到 10 之间，总代价最多 200。如果 P 取 5，那么违反一次覆盖约束的罚款是 5，但省下来的代价可能是一个 10 块钱的方案不被选。求解器当然选择违规。这种情况下正确的 P 应该至少大于“所有方案代价之和”，也就是 200 左右，保险起见取 400 到 600。

这里有个容易忽略的点：判断 P 够不够，不能看单个元素的罚金，要看整体。因为一个非法解可能同时违反好几个元素的约束，罚金是累加的，但与此同时它也省下了好几个方案的代价。所以安全的做法是拿“目标函数可能的最大波动范围”去比较，而不是拿单次罚款去比较。

3.2 P 偏大时能隙塌缩与数值退化

另一个方向的问题没那么直观，但同样致命。P 取得过大，比如比目标项大出三四个数量级，会带来两个后果。

第一个后果是能隙塌缩。所有合法解的能量都变成了一个巨大的负数（因为对角项里减去了大量的 P），合法解之间的能量差异，相对于这个巨大的基数来说变得微乎其微。求解器在比较两个合法解时，需要分辨的差值可能只占总能量的百万分之一，这在浮点数精度和硬件噪声面前基本上等于无法分辨。表现为：解是合法的，但质量很差，而且多次运行结果高度一致地差，说明搜索根本没有有效区分优劣。

第二个后果是数值退化。很多求解器内部用定点数或者有限精度浮点表示系数，P 太大时，目标项的小系数会被直接截断成零或者被舍入掉。我遇到过一次，目标代价全在 0.1 到 1 之间，P 取了 10000，结果跑出来的解完全无视目标函数，只要合法就随便挑一个。查了半天才发现是系数精度被吃掉了。

判断方法很简单：把 Q 矩阵里所有非零元素的绝对值排序，看看最大值和最小值的比值。如果超过 10^6，就要警惕了。理想情况下这个比值控制在一万以内，最好在一千以内。

3.3 一个可以直接套用的 P 下界与归一化流程

经过多次试错，我固化下来一套流程，基本能保证第一次就把 P 定在合理区间。

第一步，先把所有代价归一化。把 c_j 除以 max|c_j|，让它们全部落在 [-1, 1] 区间内。这一步的意义是把问题尺度标准化，后面所有系数都基于这个尺度来定，换数据集时不用重新调参。代价如果是非负的，归一化后落在 [0, 1] 更好。

第二步，算一个下界。让 P 至少大于“任何可行解的代价上界”。可行解最多选 n 个方案，每个代价不超过 1，所以上界是 n。因此 P > n 是一个必要下界。但实际用的时候我一般取 P = 2n 到 4n，因为要留出与非法解竞争的余量。

第三步，做一次验证。建好矩阵后，构造一个故意违反约束的解（比如全选 1），算它的能量；再构造一个明显合法的解，算能量。前者必须显著大于后者。如果前者反而更小，说明 P 不够，加倍再来。

这套流程我在十几个不同规模的问题上用过，基本一次就能定下来，剩下的只是微调。唯一要注意的是，如果代价允许为负（比如某些方案带来收益），那 n 这个上界不成立，需要用“所有正代价之和”作为上界，或者把负代价单独处理。

还有一个提升稳定性的小技巧：在归一化之后，把所有系数统一放大到一个合适的整数区间，比如乘以 1000 后取整。这能让定点数求解器工作得更精确，也方便你把矩阵存成整数格式做快速验证。代价是矩阵里会出现一些大整数，但相对于浮点数的不确定性，这点代价是值得的。

4. 抛硬币式的翻译：QUBO 到 Ising 的系数对应

QUBO 是给优化问题用的语言，Ising 是给物理硬件用的语言，两者之间需要一次翻译。翻译的桥梁是一个极简单的代换：把 0/1 变量换成 ±1 变量。这个代换本身没什么难度，难的是翻译过程中常数项、一次项、二次项的归位，一旦搞错，硬件上跑出来的东西和你以为的完全不是一回事。

我之所以把这一章单独拿出来讲，是因为我见过太多人卡在这里。他们建好了 QUBO，矩阵看着也没问题，一送进硬件结果就不对了，查来查去发现是 Ising 转换时把常数项丢了，或者把交叉项系数算成了四分之一的两倍。

4.1 x=(1+s)/2 代换的三步推导

代换关系是 x_i = (1 + s_i) / 2，其中 s_i ∈ {-1, +1}。当 x_i = 1 时 s_i = 1，当 x_i = 0 时 s_i = -1。这个关系是一一对应的，没有信息损失。

把 QUBO 的两类项分别代换进去。

一次项：Q_ii x_i = Q_ii (1 + s_i)/2 = Q_ii/2 + Q_ii s_i/2。

二次项：Q_ij x_i x_j = Q_ij (1+s_i)(1+s_j)/4 = Q_ij/4 · (1 + s_i + s_j + s_i s_j)。

把二次项展开成四项：常数项 Q_ij/4，s_i 项 Q_ij s_i/4，s_j 项 Q_ij s_j/4，以及交互项 Q_ij s_i s_j / 4。

现在做归位。所有常数项收集起来，记作 C；所有 s_i 的一次项系数收集起来，记作 h_i，这就是 Ising 模型里的局部场；所有 s_i s_j 的系数收集起来，记作 J_ij，这就是耦合强度。

4.2 常数项、局部场、耦合矩阵的落位

按照上面的展开，对每个变量 i：

h_i = Q_ii / 2 + Σ_{j ≠ i} Q_ij / 4

这里要注意求和符号的范围。Q_ij 中 i 和 j 是对称出现的，对于固定的 i，它在所有包含 i 的交叉项里都会贡献一个 Q_ij/4。如果 Q 是上三角存储（只存 i < j 的部分），那对 i 来说，它既作为较小下标出现在 (i, j) 里，也作为较大下标出现在 (k, i) 里，两种情况都要算进去，所以是求和所有 j ≠ i。

J_ij = Q_ij / 4，只对 i < j 定义，且 J 矩阵对称。

C = Σ_i Q_ii/2 + Σ_{i<j} Q_ij/4。

最终的能量表达式是 E = Σ_i h_i s_i + Σ_{i<j} J_ij s_i s_j + C。

这里有个非常重要的实践点：常数项虽然不影响最优解的位置，但它会影响你读能量的绝对值。如果你拿硬件返回的能量和本地计算的 QUBO 值做对比，忘了常数项就会对不上，然后开始怀疑硬件有问题。我建议在代码里把常数项显式算出来、存下来，每次对比都带上它。

4.3 用两变量例子验证翻译是否掉项

光看公式容易糊，拿一个两变量的例子跑一遍就知道对不对了。

设 Q_11 = -3，Q_22 = -5，Q_12 = 4，没有其他项。QUBO 能量是 f(x) = -3x_1 - 5x_2 + 4x_1 x_2。

四个组合的能量分别是：f(0,0)=0，f(1,0)=-3，f(0,1)=-5，f(1,1)=-3-5+4=-4。最小值在 (0,1)，能量 -5。

现在转 Ising。h_1 = -3/2 + 4/4 = -1.5 + 1 = -0.5。h_2 = -5/2 + 4/4 = -2.5 + 1 = -1.5。J_12 = 4/4 = 1。C = -3/2 - 5/2 + 4/4 = -1.5 - 2.5 + 1 = -3。

验证 s = (-1, 1) 对应 x = (0, 1)：E = h_1(-1) + h_2(1) + J_12(-1)(1) + C = 0.5 - 1.5 - 1 - 3 = -5。和 QUBO 的 -5 完全一致。

再验证 s = (1, 1) 对应 x = (1, 1)：E = -0.5 - 1.5 + 1 - 3 = -4，也对上了。

这个验证步骤花不了两分钟，但能挡掉绝大部分低级错误。我的习惯是每次写完转换函数，先拿两三个变量的随机小例子跑一遍全枚举比对，确认无误再上规模。

5. 不等式约束怎么塞进二次型：松弛变量的二进制展开

集合划分的纯形式只有等式约束，但现实问题很少这么干净。容量限制、时间窗、最多选几个，这些都是不等式。不等式进 QUBO 比等式麻烦一个量级，因为它需要一个连续的松弛量，而 QUBO 里只允许 0/1 变量存在。

处理办法是把这个连续松弛量用二进制展开表示出来。听起来复杂，逻辑其实很直白：你要表示一个 0 到 C 之间的整数，就用若干个二进制位拼出来，每一位是一个 0/1 变量。

5.1 为什么必须引入 slack

拿装箱举例。约束是“第 b 个箱子里所有物品的总重量不超过容量 W_b”。写成式子是对每个箱子 b，Σ_{j ∈ 物品} w_j · y_{jb} ≤ W_b，其中 y_{jb} 表示物品 j 是否放进箱子 b。

这是个不等式，不能直接平方。要把它变成等式，就得补一个非负的松弛量 s_b，使得：

Σ w_j y_{jb} + s_b = W_b

然后对这个等式做和前面一样的平方惩罚。注意 s_b 的符号：如果约束是“不超过”，松弛量加在左边；如果是“不少于”，松弛量要减在左边。这个符号搞反了，整个惩罚项就变成了鼓励违规，问题会变得极其诡异。

我吃过一次这个亏。当时是“至少选 k 个”的约束，我照着“不超过”的模板写，松弛量的符号没改，结果求解器疯狂地选满所有变量。查了两个小时才定位到是一个正负号的问题。从那以后我养成了一个习惯：每次写完不等式约束，先手算一个明显违规的解和一个明显合法的解，看惩罚项是不是前者大，不验证不往下走。

5.2 位宽算错的两种表现

松弛量的二进制展开是这样的：s = Σ_{k=0}^{K-1} 2^k · z_k，其中 z_k ∈ {0, 1}。K 是位数，必须满足 2^K - 1 ≥ C，也就是 K ≥ ⌈log2(C + 1)⌉。

这里有两个容易出错的地方。

第一个是位宽算小。如果 C = 10，你需要能表示 0 到 10 共 11 个值，2^3 - 1 = 7 不够，2^4 - 1 = 15 才够，所以 K = 4。很多人图省事写成 K = ⌈log2(C)⌉ = 4，刚好蒙对；但如果 C = 8，⌈log2(8)⌉ = 3，2^3 - 1 = 7 < 8，就表示不了 8 这个值，导致某些可行解在模型里根本不存在。表现是求解器给出的解看起来总是差一点，怎么调都到不了理论最优。

第二个是位宽算大得离谱。C 只有 20，你给了 20 位，凭空多出 20 个变量。变量数每增加一个，解空间翻倍，搜索难度指数上升。而且多余的位会制造大量等价解（高位全 0 和某些低位组合表示同一个松弛值），让搜索在原地打转。

我的做法是写一个小函数专门算位宽，单元测试里固定几个边界值（C = 1, 2, 3, 7, 8, 15, 16），跑一遍确认不越界。这种小工具函数看着不起眼，但它是整个建模流程里最容易埋雷的地方。

5.3 变量膨胀后的取舍

引入松弛变量会让变量数从 n 涨到 n + Σ K_b。如果问题里有几十个箱子，每个箱子要 5 位，那就是一两百个额外变量，规模直接翻倍。这时候需要权衡。

一种取舍是放松约束的精度。比如容量是 1000，你不需要精确到个位，只需要精确到 10，那就把原问题除以 10，位宽从 10 位降到 7 位。代价是可能丢掉一些精细的最优解，但对于大规模问题，这点损失通常可以接受。

另一种取舍是把不等式约束改成软约束，用不对称的惩罚替代硬约束。具体做法是只在“超过容量”时施加惩罚，而不是强制等式成立。这在 QUBO 里可以用 ReLU 式的写法近似，但会引入非线性，需要额外的辅助变量，未必划算。

还有一种思路是把问题分解。如果箱子之间没有耦合（比如每个物品只能进一个箱子，但箱子之间完全独立），那其实可以拆成多个小问题分别求解再合并。这个思路在第七章会详细讲，它往往比死磕一个大 QUBO 有效得多。

6. 一套能跑通的实现：从建矩阵到校验闭环

前面都是纸面推导，这一章把它落地。我用 Python 写一套最小可用的实现，包含三个部分：建矩阵、求解、校验。代码不长，但每个部分都有它存在的理由，缺一个就会在某个环节翻车。

这套代码我刻意没有依赖任何特定的量子计算框架，因为建矩阵和校验这两步跟用什么求解器无关。你可以把建好的 Q 矩阵丢给模拟退火、丢给禁忌搜索、丢给量子退火硬件，或者丢给门模型上的 QAOA 编译流程，输入格式都是同一份。

6.1 用字典存稀疏 Q，别急着上大矩阵

小规模问题用二维数组存 Q 没问题，但一旦变量数上百，稠密矩阵就是灾难。集合划分的 Q 矩阵通常很稀疏，因为只有共享元素的子集之间才有交叉项，绝大多数子集对是不相交的，交叉项为零。

我用字典来存，键是 (i, j) 且保证 i ≤ j，值是对应的系数。这样内存只跟非零项数成正比，遍历也快。

from collections import defaultdict from itertools import combinations def build_set_partition_qubo(subsets, costs, penalty): """subsets: list[set]，每个元素是一个候选方案覆盖的元素集合 costs: list[float]，对应方案的代价 penalty: 惩罚系数 P 返回 (Q, n)，Q 是 dict[(i,j)] -> float，只存 i <= j 的项 """ n = len(subsets) Q = defaultdict(float) # 目标项先落到对角线上 for j in range(n): Q[(j, j)] += costs[j] # 建立元素 -> 覆盖它的方案列表 cover = defaultdict(list) for j, s in enumerate(subsets): for e in s: cover[e].append(j) # 逐元素施加覆盖约束的惩罚 for e, idxs in cover.items(): # 每个被覆盖的方案，对角项减 P for j in idxs: Q[(j, j)] -= penalty # 每一对共同覆盖该元素的方案，交叉项加 2P for j, k in combinations(sorted(idxs), 2): Q[(j, k)] += 2.0 * penalty return dict(Q), n

这段代码就是第 2.3 节那个结论的直译。Q[(j, j)] -= penalty对应“每个元素贡献 -P 到对角”，Q[(j, k)] += 2P对应“共同覆盖一个元素加 2P”。注意这里用的是上三角约定，即 Q 矩阵只存 i ≤ j 的部分，能量计算时交叉项只用乘一次系数，不是两次。如果你的求解器要求对称矩阵（能量里交叉项算两次），那这里要改成加 P 而不是 2P，具体看下一节的能量函数怎么写。

我建议在函数里加一行断言，检查所有被覆盖的元素确实出现在至少一个方案里，否则那个约束永远无法满足，问题无解。这种前置校验能省掉后面大量的困惑。

6.2 翻转增量的手写模拟退火

有了 Q，下一步是求解。我写一个最朴素的单比特翻转模拟退火，不追求性能，追求的是可读和可复现。它的作用是给你一个基线：如果你的量子方案跑不过这个几十行的东西，那说明还有很长的路要走。

import random, math def qubo_energy(Q, x, n): e = 0.0 for (i, j), v in Q.items(): if i == j: if x[i]: e += v else: if x[i] and x[j]: e += v return e def flip_delta(Q, x, i, n): """翻转第 i 位后能量的变化量（上三角约定）""" d = 0.0 xi = x[i] if (i, i) in Q: d += -Q[(i, i)] if xi else Q[(i, i)] for j in range(n): if j == i: continue key = (i, j) if i < j else (j, i) v = Q.get(key, 0.0) if v == 0.0: continue if x[j]: d += -v if xi else v return d def anneal(Q, n, sweeps=50000, t_hi=5.0, t_lo=0.01, seed=0): rng = random.Random(seed) x = [rng.randint(0, 1) for _ in range(n)] cur = qubo_energy(Q, x, n) best_x, best_e = x[:], cur for k in range(sweeps): t = t_hi * (t_lo / t_hi) ** (k / sweeps) i = rng.randrange(n) d = flip_delta(Q, x, i, n) if d <= 0 or rng.random() < math.exp(-d / max(t, 1e-12)): x[i] ^= 1 cur += d if cur < best_e - 1e-12: best_x, best_e = x[:], cur return best_x, best_e

flip_delta这段是性能关键。它只遍历跟 i 有关的项，而不是重算整个能量，复杂度从 O(nnz) 降到 O(deg(i))。对于稀疏的 Q 矩阵，这个优化非常明显。

温度调度用的是几何降温，从 t_hi 按指数衰减到 t_lo。t_hi 的选取有讲究：它应该和“典型翻转的能量差”同量级。我一般先用一次随机游走采样几百个翻转差值，取绝对值的 90 分位数作为 t_hi，这样比拍脑袋定一个数靠谱。

还有一点：best_e的更新用了容差1e-12，因为浮点累加会有微小误差，严格小于号可能导致错过等价的更好解。这类细节在长时间运行里会显现出来。

6.3 校验器与非法解拦截

求解器返回的解必须校验，这一步不能省。QUBO 的约束是软的，求解器没有任何义务给你一个合法解。

def check_partition(x, subsets, elements): cnt = {e: 0 for e in elements} for j, s in enumerate(subsets): if x[j]: for e in s: cnt[e] += 1 bad = {e: c for e, c in cnt.items() if c != 1} return len(bad) == 0, bad

返回一个布尔值和非法元素的分布。如果非法，先看非法元素的个数和具体分布：如果是零星一两个，说明惩罚系数差一点点，往上调 20% 再试；如果是大面积非法，那就不是系数问题，是矩阵建错了，回去检查。

我习惯把校验器做成流水线的一部分，每次求解都自动跑一遍，把非法率记下来。批量实验时看非法率随 P 变化的曲线，能很快定位到合适的区间。

6.4 一组实测数据与参数敏感性

拿一个中等规模的例子跑一下。20 个元素、60 个候选子集，每个子集随机包含 2 到 4 个元素，代价随机在 0 到 10 之间。这个规模的全枚举不可行，用我的模拟退火跑 50000 步，重复 20 次。

不同惩罚系数下的表现大致是这样的：

这个表很说明问题：可用区间是 200 到 2000 之间，两头都会退化。500 附近是甜点区，大概是总代价上界（60 × 10 = 600）的 0.8 倍。这个比例和我在其他问题上观察到的经验值基本吻合，一般取总代价上界的 0.5 到 1.5 倍之间。

还有个小细节：P 取 200 的时候非法率还有 5%，但那 5% 的解的代价明显低于平均，说明它们是通过违反约束换来的收益。这些都是很好的诊断信号，能帮你判断 P 还差多少。

7. 我在这个模型上踩过的四个坑

前面几章讲的都是“应该怎么做”，这一章讲讲“实际上会怎么错”。这些坑有的来自数学约定，有的来自问题结构，有的来自对工具的误解，共同点是：踩的时候都很自信，事后看都很低级。

7.1 上三角约定导致的系数差 2 倍

QUBO 的交叉项在能量函数里可能出现一次也可能出现两次，取决于你怎么定义。有两种常见写法：

写法 A：f(x) = Σ_i Q_ii x_i + Σ_{i<j} Q_ij x_i x_j，Q 只存上三角。 写法 B：f(x) = Σ_{i,j} Q_ij x_i x_j，Q 是完整对称矩阵，交叉项在 (i,j) 和 (j,i) 各出现一次。

两者的关系是：写法 B 里的 Q_ij 等于写法 A 里的 Q_ij 的一半。也就是说，如果你按写法 A 建好了矩阵，直接喂给按写法 B 实现的求解器，所有交叉项都会变成两倍。

这个错误的后果是所有重叠惩罚被放大一倍，模型倾向于选出交集更少的方案组合。解还是合法的（如果惩罚够大），但目标值明显偏离最优。更难查的是，当 P 本身就在临界值附近时，这个 2 倍误差可能刚好让模型从“勉强可行”变成“完全可行但次优”，看起来一切正常。

我的防御办法是在代码里显式声明约定，并且写一个双变量的手工测试用例做断言。比如固定 Q_12 = 4，构造 x = (1,1)，检查能量函数返回的是 4 还是 8。这个断言写一次，后面所有改动都能被挡住。

7.2 重复子集制造的等价解沼泽

真实数据里经常出现重复的候选方案：两个方案覆盖的元素完全一样，代价也完全一样，或者代价略有差异。如果不做预处理直接建矩阵，这两列会形成一个对称结构：选 A 不选 B 和选 B 不选 A 的能量完全相同。如果有 k 个这样的重复，就有 2^k 个等价解。

这会让搜索效率大打折扣。退火过程在等价解之间来回跳，浪费时间；更糟的是，如果求解器有某种“偏好”（比如总是从编号小的变量开始翻转），它可能长期卡在某个特定的等价类里出不来，误以为已经收敛。

处理办法是建矩阵之前先去重：把覆盖集合相同、代价差异在容差内的候选方案合并成一个，只保留代价最小的那个。如果覆盖集合相同但代价不同，保留代价小的，把代价大的直接删掉——因为它在任何情况下都不会是最优选择。这一步能把变量数削减 10% 到 30%，具体取决于数据质量。

还有一种不那么明显的等价：两个方案覆盖集合不同，但在当前的约束结构下效果一样。比如某个元素只被这两个方案覆盖，而这两个方案其他部分都不被任何约束涉及，那它们就是可以互换的。这种等价需要做连通性分析才能发现，一般是把问题拆成连通分量之后顺带解决的。

7.3 连通分量分解：被忽略的降维手段

这是我觉得最被低估的一个技巧。把候选子集建成一张图：两个子集如果共享至少一个元素，就连一条边。这张图的连通分量之间是完全独立的——因为没有任何约束同时涉及两个分量的变量，它们的能量贡献可以分开计算。

也就是说，如果整张图有 5 个连通分量，你就可以把一个大问题拆成 5 个小问题，分别求解，然后把结果拼起来。拼接的结果一定是全局最优，因为分量之间没有任何耦合。

我在一个 200 变量的实际数据上试过，图分解之后最大的分量只有 30 来个变量，剩下的都是零星几个的小分量。分解后求解时间从分钟级降到秒级，解的质量还更好，因为每个子问题都能充分搜索。

代码实现上，就是标准的并查集：遍历每个元素，把它覆盖的所有方案合并到一个集合里。最后按根节点分组，每组单独建 QUBO、单独求解。

要注意的是，这个方法只在约束是“恰好覆盖一次”时严格成立。如果目标函数里存在跨分量的项（比如“选了这个方案就不能选那个方案”，而这两个方案属于不同分量），那分量之间就有了耦合，不能独立求解。所以做分解之前要先检查目标函数的结构，确认它是可分离的。

7.4 别把 QUBO 当成能自动解决 NP-hard 的魔法

这是态度层面的一条，但我觉得比任何技术细节都重要。QUBO 是一套问题表示格式，不是算法，更不是加速器。它做的事情是把问题统一成一张矩阵，让不同的求解器可以在同一个输入上竞争。它不改变问题的复杂度，也不保证任何求解速度上的优势。

集合划分是 NP-hard 问题，这一点不会因为用了 QUBO 就改变。当前量子退火硬件能处理的变量数有限，连通性也有限制，一个中等规模的工业问题往往需要先做图嵌入，把逻辑变量映射到物理比特上，这个映射过程本身就可能消耗大量额外比特。门模型上的 QAOA 目前更多是原理验证阶段，在组合优化上的实际表现还没有稳定超过经典启发式。

所以务实的做法是：把 QUBO 当成一个统一的表达层，先用它把你的问题理清楚，然后用成熟的经典启发式（模拟退火、禁忌搜索、局部搜索加扰动）跑一个基线。只有当问题规模、结构或者求解时间确实成为瓶颈，并且你有明确的理由相信量子方法可能有优势时，再去考虑上硬件。营销材料里那些“量子加速”的说法，在真实的组合优化问题上，绝大多数情况下还只是潜力，不是现实。

我个人在实际操作中的体会是，把问题建成 QUBO 这件事本身的收益，往往比用它去跑量子硬件更大。因为建模过程逼着你把约束理清楚、把代价量化、把结构看清楚，这些工作做完之后，你就算用最朴素的模拟退火，也能拿到相当不错的结果。

最后再分享一个小技巧：如果你的集合划分问题规模大到 QUBO 都建不下，试试先做一次贪心构造，得到一个可行解作为初始点，然后把 QUBO 的搜索限制在这个解附近的邻域里。退火从全空间随机初始化，往往不如从一个不错的可行解出发，在小邻域里做精细调整。这个技巧在几百个变量的规模上效果特别明显，而且实现成本几乎为零，只需要把初始解改成贪心结果，把温度上限调低就行。