最近我一直在琢磨一件事量子计算模拟器到底难不难写如果抛开 Qiskit、Cirq 这类现成框架完全从零开始用纯 Python 搭建一个能跑量子线路、能测量、能复现经典实验的模拟器代码量会不会大得吓人带着这个疑问我把 DREAMVFIA 这个开源项目整理成了完整的实现笔记。这篇是上篇先把最核心的框架、数据结构和基本算子讲透——量子态怎么表示、量子门怎么作用、测量怎么坍缩、贝尔态怎么跑通。整个实现只用 NumPy 一个第三方库核心代码不足三百行非常适合两类读者一类是想搞懂量子计算底层数学逻辑的 Python 开发者另一类是正在学量子信息、想找一个能动手改的模拟器教学项目的学生。看完这篇你会对“量子比特就是一个复向量”、“量子门就是矩阵乘法”、“测量就是概率采样”这三句话有切身的体感。1. 项目定位与整体设计思路1.1 为什么非要自己写一个模拟器市面上的量子模拟器其实不少。Qiskit Aer 能模拟上百个量子比特Cirq 背靠 Google 的量子生态ProjectQ 性能也很强。那为什么还要自己做我的原因很简单也很“笨”不自己实现一遍量子计算的很多概念永远是黑箱。用 Qiskit 跑一个贝尔态三行代码出结果但你不会知道 1000 次测量里为什么恰好一半是 00、一半是 11。这个“为什么”不是靠文档能消化的必须亲手把量子态向量写出来、把 H 门乘进去、把 CNOT 门作用上去亲眼看到振幅是怎么从[1, 0, 0, 0]变成[0.707, 0, 0, 0.707]的才算真正理解。DREAMVFIA 的定位就是教学型模拟器它刻意不追求性能和功能堆砌而是追求可读性和可扩展性。每一个矩阵、每一步位运算都没有隐藏魔法读代码的过程就是复习量子计算基础的过程。对于给学生上课的场景这种“能拆开看内部结构”的代码比任何商业框架都有价值。1.2 分层架构与模块规划写模拟器之前我先把物理模型映射成软件模块。量子计算的运行过程可以分成四件事准备量子态、施加量子门、测量、统计结果。于是 DREAMVFIA 的核心也分成四层模块职责对应物理概念state.py量子态的存储与初始化态向量、振幅gates.py常见量子门的矩阵定义量子门、酉变换circuit.py线路的构建与量子门施加量子线路measure.py测量与坍缩逻辑测量、Born 规则这样的分层有一个明显的好处每一层都可以独立测试。我写代码的时候是严格按层推进的先完成量子态模块并跑通单元测试再写量子门再写线路调度最后做测量。任何一个环节出错都能快速定位到具体的层而不是在几百行代码里大海捞针。1.3 技术选型Python 与 NumPy 的组合逻辑选 Python 不是因为性能而是因为它能把线性代数写得像数学公式。量子计算本质上是复向量空间上的线性变换NumPy 提供的complex128复数类型、向量化运算和广播机制几乎是为此量身定制的。选 NumPy 而非纯 Python 列表也是同理。模拟器的核心操作是矩阵乘法和振幅遍历如果用嵌套列表手动实现不仅代码冗长还容易在“列表拷贝”这种细节上踩坑。NumPy 的ndarray切片和拷贝语义清晰配合dtypenp.complex128一个振幅向量的内存占用和读写效率都远优于原生列表。当然代价是有的。Python 的循环逐索引处理 2 的 n 次方个振幅时性能天花板很低。这个问题的应对策略我在第 5 节专门讲这里先卖个关子。总而言之教学项目的第一目标是“让人看懂”在这个前提下Python NumPy 就是最优解。2. 量子计算模拟的核心原理2.1 量子比特从比特到状态向量经典比特只有 0 和 1 两个状态量子比特则多了一个“叠加”的可能性。数学上一个量子比特的状态用一个二维复向量表示|ψ⟩ α|0⟩ β|1⟩其中 α 和 β 都是复数分别叫“|0⟩ 态振幅”和“|1⟩ 态振幅”。它们满足归一化条件|α|² |β|² 1这个条件的物理含义是对量子比特进行测量得到 0 的概率是|α|²得到 1 的概率是|β|²总概率必须等于 1。你可以把量子比特想象成一根箭头在二维复平面里的方向向量只不过这根箭头的“长度”固定为 1但“方向”可以连续变化。经典比特只能指向两个固定方向0 或 1而量子比特可以指向圆周上的任意一个位置。正是这种“连续自由度”让 n 个量子比特能编码 2 的 n 次方个复数振幅。n 个量子比特的联合状态是一个 2ⁿ 维的复向量|ψ⟩ Σ c_i |i⟩其中 i 从 0 到 2ⁿ-1每个下标 i 对应一个计算基态比如三个量子比特时下标 5二进制 101对应的基态是|101⟩。整个模拟器最核心的数据结构就是存放这 2ⁿ 个复数振幅的一维数组。2.2 量子门酉矩阵与矩阵乘法量子门是对量子态施加的变换数学上就是一个酉矩阵酉矩阵满足 U†U I即共轭转置乘自身等于单位阵。酉矩阵保证变换前后向量的长度不变也就是归一化条件不会因为施加量子门而被破坏。单比特量子门是一个 2×2 酉矩阵作用在量子比特上就是一次矩阵乘法。以最经典的 Hadamard 门简称 H 门为例H 1/√2 [[1, 1], [1, -1]]H 门的作用是把|0⟩变成(|0⟩ |1⟩)/√2也就是将确定态变成等概率叠加态。这个门的矩阵写法是H|0⟩ 1/√2 [[1, 1], [1, -1]] [1, 0]ᵀ [1/√2, 1/√2]ᵀ看到没有本质上就是高中数学里的矩阵乘向量。Pauli-X 门就是量子版 NOT 门的矩阵是[[0, 1], [1, 0]]作用在|0⟩上得到|1⟩作用在|1⟩上得到|0⟩和经典取反一模一样。两比特门中最重要的是 CNOT 门受控非门。它有两个输入控制比特和目标比特。如果控制比特是 1就对目标比特取反如果控制比特是 0目标比特保持不变。CNOT 是一个 4×4 矩阵具体形式在实现一节里给出。这里有一个初学者最容易忽略的点量子门必须保证可逆。经典与非门会丢弃信息但量子门不会因为酉矩阵一定可逆。这也是为什么模拟器的所有门操作都只是数值上的线性变换不存在“if 分支丢弃数据”这种写法测量除外。2.3 测量Born 规则与波函数坍缩测量是量子计算里最特殊的一环也是模拟器里唯一有“随机性”和“破坏性”的操作。对一个量子比特测量得到结果 k 的概率由振幅的模平方决定P(k) |c_k|²这个规则叫 Born 规则它把抽象的复数振幅和现实世界里的观测概率联系了起来。模拟器的测量模块要做两件事第一按概率随机采样。用np.random.random()生成一个 [0, 1) 的随机数再计算累积概率落在哪个区间就返回哪个测量结果。第二坍缩量子态。测量一旦发生量子态不再是叠加态而是“坍缩”到与测量结果对应的基态上。在数学上这相当于把所有不满足测量结果的振幅置为零然后重新归一化。测量是不可逆的破坏性操作这就是为什么模拟器里测量之后量子态会“变掉”。如果多次测量同一个量子比特第一次测量会决定后续所有结果——这一点在模拟器里天然成立因为坍缩后其他振幅已经归零了。3. DREAMVFIA 核心代码实现3.1 量子态的数据结构量子态是整个模拟器的心脏。我选择把 n 个量子比特的状态存成一个长度为 2ⁿ 的一维复数数组下标 i 的二进制表示对应计算基态|i⟩的振幅。初始化时全部振幅置 0仅下标 0 处置 1表示所有量子比特都处在|0...0⟩态。import numpy as np class QuantumState: def __init__(self, n_qubits): self.n_qubits n_qubits self.dim 1 n_qubits # 2^n self.amp np.zeros(self.dim, dtypenp.complex128) self.amp[0] 1.0 0.0j # |00...0 def probability(self, index): return abs(self.amp[index]) ** 2dtypenp.complex128是必须写清楚的。如果不指定NumPy 默认用float64复数振幅会被截断成实数整个模拟器直接废掉。这里还有一个隐含约定下标 0 对应所有比特为 0 的状态这也是绝大多数量子模拟器的默认约定。3.2 量子门定义与线性代数基础我在gates.py里把常用量子门定义成函数返回 NumPy 矩阵。这样做的可读性比直接写一坨带数字的数组好得多后续如果要加参数化门比如相位门只需要在这个函数里加一个角度参数。def hadamard(): return np.array([[1, 1], [1, -1]]) / np.sqrt(2) def pauli_x(): return np.array([[0, 1], [1, 0]]) def pauli_z(): return np.array([[1, 0], [0, -1]]) def cnot(): return np.array([[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 0, 1], [0, 0, 1, 0]])注意 CNOT 矩阵的行列顺序输入态按|00⟩, |01⟩, |10⟩, |11⟩排列矩阵把|10⟩映射到|11⟩把|11⟩映射到|10⟩这就是“控制比特为 1 时目标比特取反”的矩阵表达。这些矩阵都是酉矩阵。写完之后建议用np.allclose(U U.conj().T, np.eye(len(U)))自检一下这是我在开发中坚持做的第一层验证。3.3 单比特门的施加与比特序有了量子态和量子门最关键的问题来了怎么把一个 2×2 的矩阵“作用到”一个 2ⁿ 维向量中的某一个量子比特上我的做法是逐对索引处理。假设目标比特是 t那么任意下标 i 都有一个“孪生兄弟”j i ^ (1 t)它们的二进制位只在第 t 位不同。施加单比特门 U 时对每一对 (i, j) 做局部矩阵乘法def apply_single_qubit_gate(state, gate, target): new_amp state.amp.copy() u00, u01, u10, u11 gate[0, 0], gate[0, 1], gate[1, 0], gate[1, 1] for i in range(state.dim): j i ^ (1 target) if i j: # 每对只处理一次 a_i, a_j state.amp[i], state.amp[j] new_amp[i] u00 * a_i u01 * a_j new_amp[j] u10 * a_i u11 * a_j state.amp new_amp这里有一个非常重要的约定比特序是“小端”的即 qubit 0 对应二进制的最低位。比如 2 个量子比特时下标 1二进制 01表示 qubit 0 为 1、qubit 1 为 0。为什么i j就能保证每对只处理一次因为当第 t 位是 0 时j i 2^t i当第 t 位是 1 时j i - 2^t i。所以循环只在“第 t 位为 0 的那个下标”处完成这对振幅的更新不会重复计算也不会漏掉任何一对。3.4 CNOT 门的施加逻辑CNOT 是两比特门作用于控制比特 c 和目标比特 t。当控制位为 1 时交换目标位为 0 和 1 的振幅控制位为 0 时不做任何事。def apply_cnot(state, control, target): new_amp state.amp.copy() for i in range(state.dim): if ((i control) 1) 1: # 控制位为 1 j i ^ (1 target) if i j: new_amp[i] state.amp[j] new_amp[j] state.amp[i] state.amp new_amp同样用i j避免把同一对振幅交换两次。这里再强调一次交换操作必须基于原始state.amp读取写入到new_amp否则会在同一轮循环里把已经交换过的值又读出来造成“二次交换”结果完全错误。3.5 测量与坍缩的实现测量的实现分三步计算目标比特为 1 的总概率、按概率随机采样、坍缩并归一化。def measure(state, qubit_index): prob_one 0.0 for i in range(state.dim): if (i qubit_index) 1: prob_one abs(state.amp[i]) ** 2 outcome 1 if np.random.random() prob_one else 0 mask 1 qubit_index if outcome 1: for i in range(state.dim): if (i mask) 0: state.amp[i] 0.0 else: for i in range(state.dim): if (i mask) ! 0: state.amp[i] 0.0 norm np.sqrt(np.sum(np.abs(state.amp) ** 2)) state.amp / norm return outcome坍缩后必须要做归一化。因为置零之后向量的模长不再等于 1如果不除以新的模长后续所有的概率计算都会失真。state.amp / norm这一步是在模拟“测量后量子态重新成为合法量子态”的物理过程。4. 实操验证跑通你的第一个量子线路4.1 环境准备与项目结构代码只依赖 Python 3.9 和 NumPy。建议用虚拟环境隔离依赖python -m venv venv source venv/bin/activate pip install numpy pytest项目结构我按模块职责拆成下面这样每个文件只干一件事测试文件也一一对应dreamvfia/ ├── dreamvfia/ │ ├── __init__.py │ ├── state.py │ ├── gates.py │ ├── circuit.py │ └── measure.py ├── examples/ │ └── bell_state.py └── tests/ ├── test_state.py ├── test_gates.py └── test_bell_state.pycircuit.py我封装了一个简单的线路类把量子态、量子门和测量串起来from .state import QuantumState from . import gates from .gates import hadamard, pauli_x from .measure import measure as measure_qubit class QuantumCircuit: def __init__(self, n_qubits): self.n_qubits n_qubits self.state QuantumState(n_qubits) self.instructions [] def h(self, target): self.instructions.append((h, target)) apply_single_qubit_gate(self.state, hadamard(), target) def x(self, target): self.instructions.append((x, target)) apply_single_qubit_gate(self.state, pauli_x(), target) def cnot(self, control, target): self.instructions.append((cnot, control, target)) apply_cnot(self.state, control, target) def measure(self, target): return measure_qubit(self.state, target)线路类把每一次操作记录到instructions列表里这一步对调试验证很有用。跑完一个实验可以打印这个列表确认门是按预期顺序被施加的。4.2 构建贝尔态并测量贝尔态Bell State是量子纠缠里最经典的例子。线路只有两步先把 qubit 0 用 H 门变成叠加态再用 CNOT 门让 qubit 1 与 qubit 0 纠缠from dreamvfia import QuantumCircuit qc QuantumCircuit(2) qc.h(0) qc.cnot(control0, target1)运行之后直接看一眼振幅向量print(qc.state.amp.round(4)) # [0.70710.j 0.0.j 0.0.j 0.70710.j]这个输出信息量很大。下标 0 和下标 3 的振幅非零且都等于1/√2说明量子态是(|00⟩ |11⟩)/√2——两个量子比特处于纠缠态。下标 1 和下标 2 的振幅为 0意味着无论怎么测量都不可能得到 01 或 10 的结果。再对两个量子比特做 1000 次测量统计结果分布from collections import Counter results Counter() shots 1000 for _ in range(shots): qc QuantumCircuit(2) qc.h(0) qc.cnot(0, 1) b0 qc.measure(0) b1 qc.measure(1) results[(b0, b1)] 1 print(results) # 典型输出: Counter({(0, 0): 498, (1, 1): 502})完美复现了贝尔态的理论预言只可能得到 00 或 11各约 50% 概率永远不可能出现 01 或 10。这背后就是量子纠缠的数学本质——测量 qubit 0 会瞬间决定 qubit 1 的状态因为它们的振幅被 CNOT 门耦合在了一起。4.3 单元测试与自检方法我一个人写模拟器时最大的感受是量子计算太容易“看起来对”了。很多错误不会让程序崩溃只会让概率统计悄悄偏离一点点。所以自检必须做成自动化测试。下面是我在项目里保留的几组核心测试def test_hadamard_on_zero(): qc QuantumCircuit(1) qc.h(0) p0 abs(qc.state.amp[0]) ** 2 p1 abs(qc.state.amp[1]) ** 2 assert abs(p0 - 0.5) 1e-12 assert abs(p1 - 0.5) 1e-12 def test_x_gate_flips_state(): qc QuantumCircuit(1) qc.x(0) assert np.allclose(qc.state.amp, [0, 1]) def test_cnot_on_control_one(): qc QuantumCircuit(2) qc.x(0) # qubit0 置为 |1⟩ qc.cnot(0, 1) assert np.allclose(qc.state.amp, [0, 0, 0, 1]) # |11⟩ def test_h_twice_is_identity(): qc QuantumCircuit(1) qc.h(0) qc.h(0) assert np.allclose(qc.state.amp, [1, 0])这几个测试覆盖了绝大部分底层逻辑。尤其是H门连续作用两次等于单位变换这一点如果单比特门的比特序搞错了这个测试会立刻报错。我在开发过程中就是靠这组测试把 index 操作和矩阵方向上的各种小错误一网打尽的。5. 开发中踩过的坑与排查经验5.1 复数精度与归一化第一个坑是我自己踩得最深的state.amp.copy()拷贝出来的数组类型是complex128但如果在初始化时写了np.zeros(self.dim, dtypefloat)后面所有复数振幅都会被截断成实部H 门算出来的1/√2倒是没问题但相位门这类依赖虚部的门会全军覆没。还有一个隐藏坑概率计算必须用abs(x) ** 2而不是x ** 2。复数(abj)的平方是(abj)² a² - b² 2abj得到的结果还是一个复数取出来当概率用简直是在随机数生成器上跳舞。正确的算法是(abj) * (a-bj)也就是模平方abs()函数做的就是这件事。5.2 比特序混乱比特序是模拟器开发里最容易出问题、也最难排查的问题。同一个振幅数组如果 qubit 0 对应最低位下标 1 表示|01⟩如果对应最高位下标 1 表示|10⟩。两种约定本身没有对错但代码内部必须全程一致。我建议在项目 README 的第一行就写明“qubit i 对应二进制第 i 位最低位为 qubit 0”并且在所有变更振幅的代码里都用位运算(i target) 1提取比特状态不要靠手算数值。检查比特序问题的一个实用技巧构造一个|1000⟩态下标 8分别对 qubit 0 和 qubit 3 施加 X 门看结果下标是 9 还是 1——这一步就能快速确认你的比特序约定。5.3 动态范围与数值稳定性量子线路深了之后振幅会变成非常小的浮点数比如 30 层线路后振幅可能只有 1e-15 量级。这时候如果某个振幅应该是 0 但实际是 1e-16测量概率计算中abs()之后几乎无影响但归一化时可能会引入微小误差。我的处理方法是在关键测试里统一用np.allclose(actual, expected, atol1e-12)而不是精确相等给数值误差留出余地。同时在测量坍缩之后立刻重新归一化避免误差在多次测量中累积。5.4 性能瓶颈在哪里最后聊一下性能。这个实现每个门都是遍历 2ⁿ 个振幅时间复杂度 O(2ⁿ)12 量子比特以下体验流畅20 比特以上开始卡顿。内存方面complex128每个振幅占 16 字节可以估算一下量子比特数振幅数内存占用约101,02416 KB1665,5361 MB201,048,57616 MB2533,554,432512 MB301,073,741,82416 GB到 30 量子比特个人电脑基本就到极限了。教学场景里 10~15 比特完全够用这也是 DREAMVFIA 现阶段的目标。如果真需要模拟更大的规模优化方向是把循环改成 NumPy 向量化操作或者对稀疏状态用字典存储非零振幅再进一步可以上 GPU 或分布式计算——这些都是下篇会展开的内容。5.5 问题速查表症状可能原因排查/解决概率总和不等于 1忘记归一化振幅用了float类型检查 dtype测量后除以模长X 门作用后结果不对比特序约定混乱用 CNOT 后纠缠态消失交换时读到了已更新的值交换必须基于旧数组副本概率出现复数用x**2而不是abs(x)**2统一用abs()取模平方16 比特以上卡顿循环逐索引处理 2ⁿ 个振幅向量化降低规模换下篇的实现6. 扩展方向与下篇预告6.1 从“能跑”到“好用”的三步当前实现能跑通贝尔态但距离“好用”还有三步要走。第一步是补齐量子门库把 Toffoli 门、SWAP 门、相位门 S/T、参数化旋转门 Rx/Ry/Rz 都加进去这会覆盖大部分量子算法的教学需求。第二步是支持部分测量和复位操作很多算法比如量子隐形传态中间需要测量某个比特并据此决定后续操作这要求线路模型支持经典比特和控制流。第三步是测量统计的可视化用 Matplotlib 画概率直方图教学演示时直观很多。6.2 性能优化路线如果要往更大规模走我有几条明确的优化路线。最直接的是把单比特门施加从 for 循环改成向量化操作用数组切片和矩阵广播一次性更新所有振幅对实测能提速几十倍。其次是引入稀疏振幅表示很多算法比如 GHZ 态制备的中间态只有少部分非零振幅用字典或稀疏矩阵存储能大幅扩展规模。最后是并行化不同振幅对之间的更新完全独立天然适合多线程或 GPU 加速这也是 DREAMVFIA 下篇要展示的重头戏。6.3 下篇能学到的内容说了这么多这篇文章其实只覆盖了“静态量子线路”——量子态按固定的门序列演化最后一次性测量。下篇我们会做一个质的飞跃实现带经典控制流的动态线路跑通量子隐形传态和 Deutsch-Jozsa 算法再把性能优化做一个完整的基准测试对比。到那时候你会看到同一个模拟器从 12 比特模拟到 20 比特以上靠的完全是向量化技巧和数据结构设计而不是换语言。我从这个项目里最大的体会是模拟量子计算机这件事门槛远没有想象中高但细节的坑远比想象中多。一个看似简单的小问题——比如比特序约定——就能让你调试一晚上。好在只要把每个环节拆开、弄懂、测试一切都变得清晰可控。希望这篇上篇笔记能帮你迈过“量子计算很神秘”这道心理门槛如果你自己也跟着写了一遍可以在评论区聊聊尤其是你踩过的那些坑大概率我当年也踩过。