用Python和PyQUBO搞定整数分割问题：从Ising模型到QUBO矩阵的保姆级实战

用Python和PyQUBO搞定整数分割问题：从Ising模型到QUBO矩阵的保姆级实战

整数分割问题看似简单，却能在金融资产分配、物流装载优化等场景中发挥巨大价值。传统方法往往陷入组合爆炸的困境，而借助PyQUBO库，我们能够将这类离散优化问题转化为QUBO矩阵，利用模拟退火等启发式算法高效求解。本文将手把手带你完成从问题定义到结果分析的全流程实战。

1. 环境准备与工具链配置

在开始之前，确保你的Python环境已安装以下关键库：

pip install pyqubo neal numpy

• PyQUBO：核心建模工具，支持从高阶表达式自动生成QUBO矩阵
• neal：模拟退火求解器，适用于中小规模问题
• NumPy：基础数值计算支持

提示：建议使用Python 3.8+环境以避免潜在的库兼容性问题。若需要GPU加速，可额外安装dwave-neal替代标准neal库。

2. 问题建模：从整数分割到QUBO表达式

整数分割问题的典型形式是：给定一组正整数，如何将其划分为两个子集，使得两个子集的和尽可能接近。用数学语言描述：

设整数集合为S = {s₁, s₂, ..., sₙ}，寻找二进制变量xᵢ ∈ {0,1}，使得目标函数：

H = (∑ sᵢxᵢ - ∑ sᵢ(1-xᵢ))²

达到最小值。展开后可以发现，这等价于最小化：

H = 4(∑ sᵢxᵢ)² - 4(∑ sᵢ)(∑ sᵢxᵢ) + const.

在PyQUBO中，我们可以直接表达这个目标：

from pyqubo import Binary def build_model(numbers): x = [Binary(f'x_{i}') for i in range(len(numbers))] sum_A = sum(s * x_i for s, x_i in zip(numbers, x)) sum_B = sum(numbers) - sum_A H = (sum_A - sum_B)**2 model = H.compile() return model

注意：这里使用Binary而非Spin变量，因为金融场景中更习惯用0/1表示是否选择某资产。

3. 模型编译与参数调优

PyQUBO的compile()方法会将高阶表达式转换为QUBO矩阵。关键参数包括：

实际编译时可能需要调整strength：

model = H.compile(strength=5.0) # 增强约束条件的影响 qubo, offset = model.to_qubo() # 获取QUBO矩阵和常数偏移量

调试技巧：

• 打印qubo.items()检查非零元素分布
• 使用model.decode_sample()验证样本有效性

4. 模拟退火求解与结果分析

neal库提供了简单易用的模拟退火接口：

import neal def solve_qubo(qubo, num_reads=100): sampler = neal.SimulatedAnnealingSampler() response = sampler.sample_qubo(qubo, num_reads=num_reads) return response

关键结果分析方法：

1. 能量值对比：response.record['energy']展示各样本对应的目标函数值
2. 最优解提取：

   best_sample = response.first.sample best_energy = response.first.energy

3. 解的可视化：

   import matplotlib.pyplot as plt energies = [r.energy for r in response.record] plt.hist(energies, bins=20) plt.xlabel('Energy') plt.ylabel('Frequency')

典型问题排查：

• 如果所有解能量值相同，可能是strength设置不当
• 出现违反约束的解，需要检查目标函数构建逻辑

5. 实战案例：投资组合优化

假设需要将以下股票投资金额分配到两个策略中（单位：万元）：

stocks = [50, 75, 120, 60, 95, 35]

完整求解流程：

# 构建模型 model = build_model(stocks) qubo, offset = model.to_qubo() # 求解 response = solve_qubo(qubo, num_reads=500) # 解析结果 solution = model.decode_sample(response.first.sample, vartype='BINARY') selected = [stocks[i] for i, val in enumerate(solution['x']) if val == 1] print(f"策略A: {selected}, 总额: {sum(selected)}") print(f"策略B: {list(set(stocks)-set(selected))}, 总额: {sum(stocks)-sum(selected)}")

输出示例：

策略A: [75, 95, 35], 总额: 205 策略B: [50, 120, 60], 总额: 230

6. 性能优化技巧

当处理大规模问题时（n>50），需要考虑以下优化策略：

1. 变量缩减：

   • 提前排除明显应该/不应该被选中的元素
   • 合并相似大小的数字

2. 分层求解：

   def hierarchical_solve(numbers, threshold=10): if len(numbers) <= threshold: return direct_solve(numbers) else: groups = partition_numbers(numbers) representatives = [sum(group) for group in groups] return hierarchical_solve(representatives)

3. 混合求解策略：

   • 先用贪心算法获得初始解
   • 再以该解为起点进行退火优化

7. 进阶应用：多目标优化扩展

实际业务中往往需要平衡多个目标，例如：

• 分割金额接近
• 风险水平相当
• 行业分布均衡

可以通过加权组合构建复合目标函数：

H_total = α*H_amount + β*H_risk + γ*H_industry

权重系数建议通过敏感性分析确定：