📚 课程目录第4章数值分析工具链搭建

4. 数值分析工具链搭建

超导量子比特的能级数值分析，本质上是在求解一个有限维哈密顿量的本征问题。无论是透射子比特的电荷色散关系、磁通子比特的隧穿劈裂，还是多能级耦合下的能级排斥，最终都归结为“构造矩阵—求本征值/本征矢—提取物理量”这条链路。本节的目标不是讲物理，而是把这条链路所需的计算环境搭起来，并用可复现的判据确认它确实可用。

一、工具链的分层与选型依据

整个环境可以分成三层。底层是 Python 科学计算栈，负责数组运算、稀疏矩阵存储与线性代数；中层是 QuTiP，它提供量子对象（Qobj）、张量积、常用算符与求解器封装，适合快速搭建少自由度模型；上层是自定义稀疏矩阵求解环境，用于处理 QuTiP 在高维截断、特殊边界条件或大规模参数扫描下不够灵活的场景。三层之间通过 NumPy 的 ndarray 与 scipy.sparse 的稀疏矩阵格式互通，这是整个工具链能够拼接的关键接口。

选型上有一个常见误区：认为 QuTiP 能覆盖全部需求。实际上 QuTiP 的 Qobj 在维度达到数千以上时，张量积构造和默认稠密化操作会显著拖慢速度；而超导比特的能级分析经常需要扫描几十个参数点、每个点截断到 20～40 个能级，此时直接使用 scipy.sparse 构造哈密顿量、调用 eigsh 求部分本征值，效率会高出一个量级。因此合理的做法是：小模型用 QuTiP 验证物理正确性，大模型用自定义稀疏环境做批量计算，两者结果在重叠区间互相校验。

二、环境配置的关键步骤

建议使用独立的虚拟环境，避免与系统包冲突。核心依赖包括 numpy、scipy、qutip、matplotlib，其中 QuTiP 5.x 起对 NumPy 2.x 有版本要求，安装时应让 pip 自行解析依赖，不要手工锁定旧版 NumPy。安装完成后，用一段最小脚本确认版本与后端：

import numpy as np, scipy, qutip
print(np.__version__, scipy.__version__, qutip.__version__)
print(qutip.settings.core)   # 查看默认数据层与稀疏后端

QuTiP 5 默认使用 qutip.core.data 的稀疏数据层，若输出显示 Dense 或出现警告，说明某些操作触发了稠密化，需要在后续建模中留意。另一个容易被忽略的配置是 BLAS 线程数：在参数扫描时，多线程 BLAS 与多进程并行会互相争抢核心，通常建议在脚本开头设置 OMP_NUM_THREADS=1，把并行度交给外层的进程池控制。

三、用透射子比特做首次验证

验证环境是否真正可用，不能只看“能 import”，而要跑一个物理结果已知的最小模型。以电荷基底下的透射子比特为例，哈密顿量写作

H = 4E_C(n − n_g)² − E_Jcos φ

在电荷基底中，n 为对角矩阵，cos φ 只连接相邻电荷态。取 E_C = 0.3 GHz、E_J = 20 GHz、n_g = 0.5，截断到 −20～+20 共 41 个电荷态，用 QuTiP 构造并求前几个本征能量，应得到近似等间距的能级，且基态与第一激发态之差约在 4～5 GHz 量级（具体数值随截断略有变化，属于示例性合理预期，不是实测值）。若结果出现能级顺序错乱或出现明显负值，通常说明 cos φ 的构造用了错误的基底索引。

随后把同一哈密顿量用 scipy.sparse 重新构造：对角项用 sparse.diags，耦合项用 sparse.diags 的偏移参数生成次对角与超对角，再调用 eigsh(H, k=6, which='SA')。两条路径的前六个本征值之差应小于 1e-8（相对误差），这就是工具链互通的判据。若差异偏大，优先检查 QuTiP 的 Qobj 转稀疏矩阵时是否发生了维度重排。

四、验证判据与常见故障定位

把两条路径的结果对齐后，需要看三类判据。第一是数值一致性：前六个本征值的相对误差应低于 1e-8，若在 1e-4 量级徘徊，多半是截断维度不同或 cos φ 的系数差了 1/2。第二是物理合理性：基态与第一激发态之差应随 E_J/E_C 增大而增大，且在 n_g=0.5 附近对 n_g 不敏感（电荷噪声保护点），这是透射子比特的标志性特征。第三是性能基线：41 维模型单次求解应在毫秒量级，若耗时达到秒级，说明某处发生了稠密化，可用 H.nnz 检查非零元数量是否与预期相符（对角 41 个、次对角各 40 个，共约 121 个）。

建议把上述验证写成一个 test_environment.py，用 pytest 组织成三个断言：版本可导入、两路径本征值一致、稀疏矩阵非零元数量正确。每次更换环境或升级 QuTiP 后先跑这个测试，比事后排查物理结果异常要省时得多。

五、面向后续章节的目录约定

为了让后面各节的脚本能直接复用，建议固定一套目录与命名约定：models/ 存放哈密顿量构造函数，统一签名 build_H(params, dim) -> sparse matrix；solvers/ 存放求解封装，返回本征值与按能级排序的本征矢；validation/ 存放交叉校验脚本；data/ 存放参数扫描结果，用 .npz 格式保存以便复现。所有构造函数只依赖 NumPy 与 SciPy，QuTiP 仅出现在验证脚本中，这样上层求解环境不会被 QuTiP 的版本变动牵连。

完成本节后，应当能够独立回答三个问题：当前环境里 QuTiP 与 SciPy 的版本是否匹配、透射子比特的前六个能级在两条路径下是否一致、稀疏矩阵的非零元数量是否符合预期。这三个问题都通过，工具链才算真正搭好，后续的能级色散分析、耦合比特能级排斥与参数扫描才有可靠的计算基础。