简介一套专注于量子信息计算的MATLAB工具库适合科研人员、工程师及高年级学生用于量子算法仿真、量子态分析与量子测量实验。资源包内含279个文件以276个m格式脚本/函数为主体覆盖量子比特操作、量子门实现、量子线路模拟、密度矩阵求解、量子测量以及熵与纠缠计算等常用功能另有3个txt文档分别说明许可证、版本更新与快速上手指引压缩包整体仅163KB轻量且易集成到科研或教学流程。目前已有545人学习下载。库内目录模块划分明确如测量、纠缠与熵、工具函数等配合核心qlib库与startup.m启动配置可快速搭建仿真环境。借助量子熵、纠缠度等物理量计算用户能深入分析量子系统性质并基于MATLAB的数值计算与可视化能力直接用于量子信息教学演示、算法验证和论文复现。1. 量子信息计算库MATLAB不写 Python 也能正经做量子计算仿真先说一个反直觉的结论量子计算仿真这件事并不是 Python 生态的专利。很多搞通信、信号处理、电子信息方向的人平时主力工具就是 MATLAB遇到量子信息方向的需求——比如验证一个量子纠错码、仿真一个含噪信道下的量子态演化、跑一遍 Grover 搜索算法——第一反应是去装 QuTiP 或 Cirq结果被 Python 环境和依赖折腾一晚上。实际上MATLAB 上有一批开源的量子信息计算库能覆盖从量子态表示、矩阵运算、信道仿真到算法验证的完整链路而且和 MATLAB 原生的数值计算、绘图、并行工具箱无缝衔接。这类库解决的核心痛点有三个第一不用切语言在熟悉的 MATLAB 环境里就能复现量子算法论文里的实验第二量子态和算子的运算本质是线性代数MATLAB 的矩阵运算效率比 Python 原生循环高一个量级中小规模比特数15 比特以内的仿真完全跑得动第三绘图和数据分析是一条龙布洛赫球、保真度曲线、噪声扫描图直接在同一个脚本里出图。适合的人群很明确有 MATLAB 基础、需要把量子信息算法跑出数值结果来验证想法或完成课题的学生和工程师以及想快速评估某个量子方案可行性、不想被 Python 版本地狱拖累的研究人员。本篇就沿着「库怎么选 → 最小例子怎么跑 → 参数怎么设 → 噪声和常见坑在哪 → 怎么验证没算错」这条线展开。2. 都有哪些现成的库选型依据与边界2.1 开源社区的一块主力量子信息工具箱类库MATLAB 生态里做量子信息计算最常见的选择是某几个流传较广的开源工具箱。这类工具箱一般以 包Package的形式组织代码模仿 MATLAB 自带的工具箱风格安装后你在命令行敲quantum.state([1 0])或qinfo.state(...)就会返回一个量子态对象而不是一个冷冰冰的复矩阵。我一般会优先看这类库的三件事第一它是否实现了密度矩阵的数据结构而不是只有态矢量。因为真实噪声仿真必须用密度矩阵态矢量只能描述纯态遇到退极化信道或振幅阻尼信道就无能为力。第二它是否内置了信道仿真函数比如damping、depolarizing这些函数返回的是一个超级算子Superoperator或 Kraus 算子集可以直接作用在密度矩阵上。第三它是否提供测量和保真度相关的函数因为最终你要算的不是那个矩阵长什么样而是某个算法跑完之后的成功概率或保真度。选型时还要注意维护活跃度。有些库是好几年以前的毕业设计只支持到某个旧版本 MATLAB在新版本里会因为函数改名比如strsplit的行为变化而报错。我通常的验证方式是直接把库下载下来跑它自带的test或example脚本如果能一口气跑完再往自己的目标上靠如果连自带的示例都跑不通果断换下一个。2.2 官方与商业路线MATLAB 支持包的取舍除了开源库MATLAB 官方也有一个量子计算相关的支持包这可能是很多人不知道的。它在 Add-On Explorer 里可以直接搜到安装后提供量子比特qubit对象、量子电路circuit对象、测量操作和几个内置算法示例。这条路线的最大好处是和 MATLAB 版本绑定不会出现兼容性问题而且底层用了 MATLAB 的matlab.lang.makeValidName等一系列现代语言特性运行效率有保障。缺点是扩展性有限。官方支持包更像是「帮你快速上手量子电路可视化」的演示工具遇到论文里自定义的算符、非标准信道模型、或者需要对 Hamiltonian 做时间演化的场景它的抽象层次反而碍事。我的经验是如果你的目标是把一个已知算法比如 Grover、Shor 的入门版本跑通官方支持包够用如果你的目标是做研究、改信道模型、对比不同噪声参数下的性能那就得用开源工具箱这类更底层的东西。两者并不冲突甚至可以在同一个项目里共存——用开源库做核心矩阵运算用官方支持包的画图函数来出布洛赫球。2.3 怎么装让库在 MATLAB 里「永久可见」无论选哪条路线安装的关键动作只有一个把库所在的文件夹加入 MATLAB 路径。这一步看似简单实际是很多初学者的第一个坑——有人直接把整个库文件夹复制到toolbox目录下也不管结果每次启动 MATLAB 都提示找不到函数。% 将下载好的库文件夹路径加入 MATLAB 路径永久生效 libraryPath D:\Workspace\QuantumInfoLib; addpath(genpath(libraryPath)); % genpath 会递归包含所有子文件夹 savepath; % 保存路径下次启动 MATLAB 仍有效这段代码的逻辑是genpath(libraryPath)生成一个包含所有子文件夹的路径列表addpath把这个列表加入当前会话的搜索路径savepath把当前的路径设置写入 MATLAB 的路径配置文件pathdef.m。如果不加genpath只写addpath(D:\Workspace\QuantumInfoLib)那库目录下的二级子文件夹里的函数全部都无法访问这是最常见的「明明装了却找不到函数」的原因。参数说明libraryPath是库在你本机上的绝对路径注意不要带中文和空格某些版本的genpath对含空格路径的处理有玄学问题。如果你用的是官方支持包不需要手动addpath在 Add-On Explorer 里点安装后 MATLAB 自己会配好路径。2.4 装完先跑通一个内置示例再谈其他装完库的第一件事不是立刻写自己的代码而是跑一个库自带的示例脚本。这一步目的是确认三件事库函数能调、绘图没问题、矩阵运算核心没坏。以某开源库为例我一般先跑它的 GHZ 态制备示例% 运行库自带的 GHZ 态示例脚本 run example_ghz.m如果这个脚本能出图、能打印出态矢量的系数就说明库基本可用。如果报错先看错误类型Undefined function说明路径没加对Invalid syntax说明你的 MATLAB 版本太老库用了新语法特性。这里有个小习惯跑通后立刻执行which quantum.state或你用的库的核心构造函数看它返回的路径是不是在库目录下。如果which返回的是其他工具箱的同名函数说明函数名冲突了这是下一个大坑。3. 量子态初体验从 |0⟩ 到叠加态的完整计算3.1 构建量子态的最小命令态矢量与密度矩阵两种方式在 MATLAB 的量子信息库中创建量子态有两条路态矢量State Vector和密度矩阵Density Matrix。两者的区别很本质态矢量只能描述纯态密度矩阵既能描述纯态也能描述混合态。初学者刚上手时一定要同时学会这两种创建方式因为后续的噪声仿真全部要落在密度矩阵上。% 创建一个单量子比特的 |0 态态矢量方式 psi0 quantum.state([1; 0]); % 创建一个单量子比特的 |0 态密度矩阵方式 rho0 quantum.density([1 0; 0 0]);第一行代码中[1; 0]是列向量表示代表基态 |0⟩ 的概率幅为 1、|1⟩ 的概率幅为 0。第二行代码中[1 0; 0 0]是 2×2 矩阵表示密度矩阵只有一个非零对角元对应 |0⟩ 态的纯态密度矩阵。这里的核心逻辑是库内部会用quantum.state或quantum.density构造函数对输入数据做合法性检查——比如态矢量的范数是否为 1、密度矩阵是否满足迹为 1 且半正定。这也意味着你不小心传入一个未归一化的向量库会直接报错或警告而不是闷头算出一个错误结果。参数说明两种构造函数的入参都支持复数构造叠加态时可以直接传[1/sqrt(2); 1/sqrt(2)]。库的构造函数通常还接受第二参数字符串类型来标记这个态的名字方便在绘图时作为图例显示。3.2 单比特门操作用矩阵乘法实现量子逻辑门量子门在数学上就是酉矩阵作用在态上就是矩阵乘法。这是 MATLAB 做量子计算最舒服的一点——你不用写任何循环矩阵乘法的效率极高。以最常用的 Hadamard 门H 门为例它的作用是把 |0⟩ 变成叠加态 (|0⟩|1⟩)/√2。% 定义 Hadamard 门 H (1/sqrt(2)) * [1 1; 1 -1]; % 将 H 门作用在 |0 态上 psi_plus H * psi0;这里的H是 2×2 矩阵psi0是之前创建的态矢量列向量两者的矩阵维度自动匹配。计算结果是 psi_plus [1/√2; 1/√2]这正是等概率叠加态。如果你用的是密度矩阵方式量子门的操作不是直接左乘而是rho_new H * rho * H其中H是 H 的共轭转置。这个区别经常有人搞混——态矢量是 U|ψ⟩密度矩阵是 UρU†漏掉一个共轭转置结果就是完全错误的。参数说明实际项目中不要手写H (1/sqrt(2)) * [1 1; 1 -1]更稳妥的做法是调用库内置的quantum.gate.H或类似属性这样可以避免手算四舍五入误差。而且内置门定义更严格比如quantum.gate.X、quantum.gate.Z、quantum.gate.CNOT都是经过预定义验证的。3.3 布洛赫球可视化让量子态「看得到」做量子信息仿真和做通信仿真最大的不同是你看不到误码率曲线只能看到矩阵。布洛赫球是单量子比特态的可视化标准工具球面上的每一个点对应一个量子态。在 MATLAB 中绘图函数通常会读取态对象或密度矩阵计算其对应的布洛赫向量分量然后在三维球体上画出矢量箭头。% 将叠加态绘制到布洛赫球上 figure; quantum.plot.bloch(psi_plus); title(Bloch Sphere Representation of | State);代码执行后MATLAB 会弹出一个 3D 球体图球面上有一条从原点指向 (1,0,0) 方向的红色箭头这个坐标对应态 (|0⟩|1⟩)/√2。如果你绘制的是 |0⟩箭头指向 (0,0,1) 北极方向|1⟩ 指向 (0,0,-1) 南极方向。布洛赫球的一个重要用途是验证门操作是否正确比如对 |0⟩ 连续执行两次 H 门箭头应该回到北极如果没回来大概率是 H 门矩阵写错了或矩阵乘法方向错了。参数说明bloch类的绘图函数一般接受一个态对象或密度矩阵作为入参如果传入的是密度矩阵库会自动提取其 Bloch 向量。部分库还支持在同一张图上绘制多个态用元胞数组传入这在对比不同噪声强度下的态偏移时非常有用。3.4 测量与概率从量子态到经典比特的关键一步量子测量的结果是概率性的这是量子信息和经典信息最大的分水岭。在仿真中测量函数做的事很简单按概率分布随机采样一个结果然后返回测量后的态。% 对叠加态进行 Z 基测量返回测量结果和测量后状态 [measureOutcome, postMeasState] quantum.measure(psi_plus, Z); % 重复测量 1000 次统计 |0 的出现次数 count0 0; for idx 1:1000 [outcome, ~] quantum.measure(psi_plus, Z); if outcome 0 count0 count0 1; end end fprintf(测量到 |0 的次数%d频率%.3f\n, count0, count0/1000);这段代码的核心逻辑是quantum.measure接受态矢量和测量基字符串Z 表示计算基返回两个输出——第一个是测量结果0 或 1第二个是测量后坍缩的量子态。注意循环里每次测量都会重新随机采样所以 1000 次的频率会接近理论概率 50%。在实际算法验证中跑单次测量没有意义必须重复多次统计频率来验证概率分布。参数说明quantum.measure的第二个参数可选项包括 Z、X、Y分别对应泡利 Z、X、Y 基测量。如果你的仿真需要固定随机种子以便复现可以在循环前调用rng(42)这样每次跑出来的分布完全一致对调试非常关键。4. 多比特体系与纠缠态从单比特走向实用算法4.1 张量积的正确姿势把两个单比特合成一个双比特态真实量子算法不可能只跑单比特多比特体系的构建靠的是张量积运算。MATLAB 中kron函数可以直接做克罗内克积量子信息库通常会封装一个更语义化的函数比如quantum.tensor。% 构建双比特态第一个比特 |0第二个比特 |0 psiAB quantum.tensor(psi0, psi0); % 构建双比特 Bell 态 |Φ (|00 |11)/√2 bellState quantum.state((1/sqrt(2)) * [1; 0; 0; 1]);第一行代码中quantum.tensor(psi0, psi0)将两个 2×1 向量合成 4×1 向量对应四个计算基 |00⟩、|01⟩、|10⟩、|11⟩ 的概率幅分别是 1、0、0、0。第二行代码手动构造了一个 Bell 态4×1 列向量的第 1 个元素对应 |00⟩和第 4 个元素对应 |11⟩都是 1/√2这正是最大纠缠态。这里有个新手高频犯的错用kron时把两个向量的顺序搞反。kron(a, b)的结果是 a 的每个元素乘以整个 b对应量子信息中的习惯是第一个比特变化更慢。如果你的库里的约定是反过来的整个态就变了。所以用quantum.tensor而不是裸kron是更稳的做法。参数说明tensor函数接受任意数量的态矢量多个参数时从左到右依次执行张量积。4.2 双比特门的实现CNOT 是纠缠的发动机CNOT受控非门是双比特量子计算的核心它在经典对应的基础上增加了量子叠加的特性。它的矩阵是一个 4×4 酉矩阵当控制比特为 |1⟩ 时目标比特翻转控制比特为 |0⟩ 时目标比特不动。% 使用内置的 CNOT 门 CNOT quantum.gate.CNOT; % 将 CNOT 作用在 Bell 态的第一个比特控制和第二个比特目标上 psiCNOT CNOT * bellState;代码逻辑很简单CNOT是一个 4×4 矩阵bellState是 4×1 向量矩阵乘法直接得到新态。你可以验证如果把 CNOT 作用在 |11⟩ 态上结果变成 |10⟩作用在 |⟩|0⟩ 上则产生纠缠态。多比特门操作中最容易翻车的点不是门本身而是比特顺序。双比特态的向量中第一分量是 |00⟩第二分量是 |01⟩第三是 |10⟩第四是 |11⟩。CNOT 的控制比特可以是第一个比特也可以是第二个比特两者对应的 4×4 矩阵是不同的。库内置的CNOT通常默认控制比特为第一个、目标为第二个。如果你需要控制比特在第二个位置需要自己构造置换矩阵或直接用quantum.gate.CNOT的参数版本部分库支持传入控制比特和目标比特的索引。如果不支持就得手动重排矩阵行和列。4.3 用 CHSH 不等式验证纠缠仿真里的「贝尔测试」很多人做了 Bell 态仿真却不知道怎么验证这个态真的是纠缠的。一个可靠的量化指标是 CHSH 不等式对于局域隐变量模型CHSH 参数 S 的理论上限是 2而量子力学允许 S 最大达到 2√2 ≈ 2.828。在 MATLAB 里可以快速验证。% 定义两组测量角度单位弧度 thetaA [0, pi/2]; thetaB [pi/4, -pi/4]; % 计算关联函数 E(thetaA, thetaB) sigma_a ⊗ sigma_b E zeros(2, 2); for ii 1:2 for jj 1:2 % 构造测量算符单位向量方向上的泡利算符 sigmaA cos(thetaA(ii)) * quantum.gate.Z sin(thetaA(ii)) * quantum.gate.X; sigmaB cos(thetaB(jj)) * quantum.gate.Z sin(thetaB(jj)) * quantum.gate.X; % 计算期望值 measOp kron(sigmaA, sigmaB); E(ii, jj) real(bellState * measOp * bellState); end end % CHSH 参数 S E(0,π/4) - E(0,-π/4) E(π/2,π/4) E(π/2,-π/4) S E(1,1) - E(1,2) E(2,1) E(2,2); fprintf(CHSH 参数 S %.4f量子上限 2√2 ≈ 2.828\n, S);这段代码做的是对 Bell 态分别沿两组角度测量计算联合期望值最后组合成 CHSH 参数。核心逻辑是把测量角度映射为测量算符——角度 θ 对应的算符是 cosθ·Z sinθ·X——然后用kron构造联合测量算符计算其在 Bell 态下的期望值。如果你构造的不是最大纠缠态S 会明显低于 2.828如果是可分离态S 不会超过 2。参数说明这里的测量角度选取直接决定 S 值大小角度偏移一点S 就会明显下降这不是 bug是量子力学本身的性质。仿真中如果跑出 S 在 2.5~2.8 之间说明纠缠态制备基本正确。4.4 多比特体系的维度诅咒与内存边界MATLAB 仿真多比特体系最现实的约束是内存。一个 n 比特的态矢量有 2^n 个复数分量每个复数 double 占 16 字节。所以 20 比特需要 1,048,576 × 16 16 MB看起来不多但 28 比特需要 4 GB30 比特需要 16 GB这还不算门矩阵。在实际项目中我会在代码里加一个维度检查超过 25 比特直接报错。% 检查态矢量维度是否超出内存安全范围 nQubits 25; if nQubits 25 error(当前计算机内存无法安全仿真 %d 比特量子系统需要 %.1f GB, ... nQubits, 2^nQubits * 16 / 1024^3); end这个检查的核心逻辑是2^nQubits × 16 是态矢量所需的最小字节数在实际中还要考虑门矩阵、中间变量、绘图函数的额外开销所以我个人把安全线划在 25 比特。参数说明如果你用的是密度矩阵内存需求是态矢量的平方n 比特密度矩阵有 4^n 个元素所以密度矩阵仿真只能做到 14~15 比特。能用态矢量绝不用密度矩阵除非你在做噪声仿真。5. 含噪声信道仿真从理想量子态到贴近硬件的真实计算5.1 三种基础噪声模型振幅阻尼、退极化、相位阻尼真实的量子硬件和理想仿真最大的区别就是噪声。量子信息库中通常内置了几种标准噪声信道它们对密度矩阵的作用方式各不相同。最常用的三种是振幅阻尼对应能量耗散、退极化对应完全随机化、相位阻尼对应退相干。% 创建三种噪声信道的 Kraus 算符damping 参数是噪声强度 ampDamp quantum.channel.amplitudeDamping(0.1); depol quantum.channel.depolarizing(0.15); phaseDamp quantum.channel.phaseDamping(0.2); % 分别作用于同一个初始态上 rhoInitial quantum.density([1 0; 0 0]); rhoAfterAmp ampDamp(rhoInitial); rhoAfterDepol depol(rhoInitial); rhoAfterPhase phaseDamp(rhoInitial);这段代码的逻辑是信道对象本质上是一个「可以作用于密度矩阵的函数」quantum.channel.amplitudeDamping(0.1)返回一个 Kraus 算子集合用函数句柄封装起来。噪声强度参数的含义因信道而异振幅阻尼的参数 γ 表示能量耗散率0.1 意味着约 10% 的概率从 |1⟩ 掉到 |0⟩退极化的参数 p 表示完全随机化的概率0.15 时态有 15% 的概率变成完全混合态相位阻尼的参数 λ 表示退相干的强度。参数说明噪声强度的取值范围是 [0,1]0 表示无噪声1 表示完全噪声。这里有一个常见的理解误区把三类噪声参数当成同一类东西来调。实际仿真中振幅阻尼和相位阻尼是硬件里最主要的两种噪声退极化信道更多是一种理论简化模型。在做硬件相关仿真时优先调前两种。5.2 含噪态的可视化保真度曲线怎么画有了噪声信道下一步就是量化评估噪声对量子态的影响。最常用的指标是保真度 F ⟨ψ|ρ|ψ⟩它衡量的是实际态与理想态的重叠程度。在实际仿真中你会扫描噪声参数从 0 到 1画出保真度随噪声强度变化的曲线。% 扫描振幅阻尼噪声强度计算保真度变化 noiseStrength 0:0.02:1; fidelityValues zeros(size(noiseStrength)); rhoIdeal quantum.density([1; 0]); % 理想 |0 态 rhoInput quantum.density([0; 1]); % 输入 |1 态会受振幅阻尼影响 for idx 1:length(noiseStrength) channel quantum.channel.amplitudeDamping(noiseStrength(idx)); rhoNoisy channel(rhoInput); fidelityValues(idx) real(sqrt(rhoIdeal * rhoNoisy * rhoIdeal)); end % 绘制保真度曲线 figure; plot(noiseStrength, fidelityValues, LineWidth, 2); xlabel(振幅阻尼强度 γ); ylabel(保真度 F); grid on;这段代码的核心逻辑是准备一个理想目标态 |0⟩把 |1⟩ 态逐个通过不同强度的振幅阻尼信道然后计算输出态与 |0⟩ 的保真度。由于振幅阻尼的物理意义是 |1⟩ 向 |0⟩ 的衰减所以噪声越强输出态越靠近 |0⟩保真度反而越高——这个反直觉的结果经常让初学者困惑。参数说明这里的保真度用的是平方根形式取值范围 0~1。如果你的库返回负值或大于 1 的数值几乎可以确定是密度矩阵在数值上不再合法比如出现负特征值需要检查是否手动操作密度矩阵时漏掉了共轭转置。5.3 噪声仿真的避坑清单密度矩阵的物理合法性检查含噪声仿真中一个反复出现的坑是密度矩阵「物理合法性」检查。密度矩阵必须满足三个条件厄米性等于其共轭转置、半正定性所有特征值非负、迹为 1所有对角元之和等于 1。在含噪声环境下浮点运算误差和手动矩阵操作很容易破坏这三个条件。% 手动计算噪声后的密度矩阵这是错误的示范会破坏物理合法性 H quantum.gate.H; rhoWrong H * rhoInitial; % 漏掉了共轭转置结果不再是合法的密度矩阵 % 检查合法性 fprintf(迹 %.4f\n, trace(rhoWrong)); fprintf(厄米性误差 %.4e\n, norm(rhoWrong - rhoWrong)); % 正确方式酉变换是 rho - U * rho * U rhoCorrect H * rhoInitial * H;代码中第一个操作H * rhoInitial是典型的错误示范它只在态矢量表示时才正确对密度矩阵做变换必须左右两边各乘一次右边乘的是共轭转置。norm(rhoWrong - rhoWrong)的值如果远大于 1e-12说明矩阵不再厄米。使用库提供的信道函数时一般不会出这个问题但只要你手动用矩阵乘法去变换密度矩阵就必须执行合法性检查。注意很多库提供了quantum.validate或isValidState这类函数可以直接判断一个矩阵是否合法的密度矩阵。我的习惯是在任何手动构造或手动变换密度矩阵之后立刻调用一次这个检查比起后面拿着错误数据跑半天再回头找 bug多写这一行的成本几乎为零。5.4 一个完整的含噪 GHZ 态制备仿真流程把前面的概念串起来展示一个完整的小项目制备三比特 GHZ 态 (|000⟩|111⟩)/√2然后分别让它通过退极化信道最后测量保真度。% 创建三比特 GHZ 态H 门 两个 CNOT 门 psiInit quantum.tensor(quantum.state([1;0]), quantum.state([1;0]), quantum.state([1;0])); H quantum.gate.H; CNOT quantum.gate.CNOT; % 第一步对第一个比特施加 H 门作用于整个三比特空间H ⊗ I ⊗ I I eye(2); H3 kron(kron(H, I), I); psiAfterH H3 * psiInit; % 第二步CNOT(1,2)控制比特 1目标比特 2 CNOT12 kron(CNOT, I); psiAfterCNOT12 CNOT12 * psiAfterH; % 第三步CNOT(2,3)——这里需要自行构造控制第二比特翻转第三比特 % 先用置换矩阵将比特顺序交换再做 CNOT具体细节略这里展示了一个实际仿真会遇到的坎多比特门操作不是简单调用内置函数就行因为内置矩阵默认作用在「相邻的固定比特位」上。当你想让 CNOT 作用在第 2 和第 3 个比特上需要构造带单位阵的张量积或使用置换矩阵。参数说明kron(kron(H, I), I)的含义是 H 门作用在第一个比特上第二、第三比特不受影响。这个构造方式在多比特仿真中会反复使用。如果你的库提供了quantum.gate.applyTo或类似操作优先用库的方式来做——它内部会自动处理比特顺序和矩阵维度避免手写置换矩阵时出错。比如库的applyTo(psiAfterH, H, [1])明确表示把 H 门作用在第 1 个比特上。执行完三步操作后得到的态理论上就是 GHZ 态。要验证制备是否正确可以用 Fidelity 函数算一下它与理论 GHZ 态的保真度如果大于 0.999 可以认为制备成功。接着把密度矩阵版本通过退极化信道扫描不同噪声强度画出保真度衰减曲线——这就是你论文里的一张标准图。6. 避坑指南量子信息库换用了 5 个最值得记录的踩坑点6.1 函数名冲突which出来的路径不对现象调用quantum.state时 MATLAB 报错或返回的结果完全不符合预期比如创建出的态维数不对。用which quantum.state -all查看发现命中的是 MATLAB 自带工具箱或另一个无关库的同名函数。原因不同工具箱和库之间函数名撞车。量子信息库里的state、measure、gate这些名字太通用很容易与 MATLAB 其他工具箱的函数冲突。解决在addpath时把量子库的路径放在路径列表的最前面这样 MATLAB 搜索时会优先命中它。做法是addpath(genpath(libraryPath), -begin)其中-begin表示把新路径添加到搜索路径开头。如果还是不行就用which定位到具体冲突的函数手动改名后调用比如写成quantum.state(...)的包名调用形式。这是最稳妥的。6.2 矩阵维度对不上kron和tensor的区别没搞清楚现象用kron(psiA, psiB)构造双比特态后后续计算保真度时维度对不上或者算出来的结果与论文不一致。原因kron的结果是「数学上的克罗内克积」量子信息中的张量积在比特排序约定上可能与之相反。不同库对第一个参数对应「第一个比特」还是「最后一个比特」的定义不同。解决不要直接用kron去构造量子态改用库封装的quantum.tensor。如果必须用kron先在文档中查清楚库的比特排序约定然后用一个已知结果做验证——比如计算 |01⟩ 的 tensor 表示看第二分量是否为 1。这一步验证只需要几秒钟但能省几个小时的排查时间。6.3 密度矩阵变换漏了共轭转置结果看起来合理其实全错现象对密度矩阵做 H 门变换后保真度和布洛赫球箭头方向都不对但程序不报错、矩阵维度也对看起来像「能跑出结果」的错误。原因态矢量的酉变换是 U|ψ⟩密度矩阵的酉变换是 UρU†U† 是 U 的共轭转置。很多人从态矢量的写法直接迁移到密度矩阵只做了左乘漏了右乘。解决实在怕漏就用信道函数来做变换channel quantum.channel.unitary(H)然后rhoNew channel(rho)。这样本质上把 UρU† 的变换封装在库函数内部你想漏都漏不掉。如果一定要手写写完后立刻检查trace(rhoNew)是否仍等于 1近似以及rhoNew是否厄米。6.4 内存爆炸试图用密度矩阵仿真 20 比特以上系统现象跑一个 18 或 20 比特的密度矩阵仿真MATLAB 直接内存溢出或电脑风扇狂转、响应极慢。原因密度矩阵的元素数量是 4^nn20 时是 1.1 万亿个复数需要 PB 级内存这显然不可能。这不是库的问题是物理规律限制了经典计算机仿真量子系统的能力。解决能改用态矢量就绝不用密度矩阵。先确认你的仿真目标是否真的需要密度矩阵——只有涉及噪声、混合态的时候才需要纯态算法仿真Grover、Shor 等理想情况全程用态矢量就够了。如果确需密度矩阵且比特数超过 14考虑改为采样法Monte Carlo 波函数法多次运行态矢量仿真并对结果取平均这种方式能推高到 25 比特左右。6.5 随机测量结果不可复现没有固定随机种子现象同一个脚本跑两次测量统计结果不同用if分支判断测量结果时两次运行走进了不同的代码路径导致后续调试困难。原因量子测量函数内部调用了随机数生成器默认每次启动 MATLAB 随机种子都不同。这不是 bug是概率仿真的正常行为。解决在脚本开头或测量循环前调用rng(42)42 是自己选的种子任意整数都行这样整个会话的随机数序列固定测量结果序列也就固定了。注意rng是全局的如果在同一个脚本里同时做噪声采样和测量采样两者共享随机数序列要复现就必须保证整个脚本的调用顺序不变。7. 验证方法怎么确定你算出来的量子态没算错7.1 第一性验证与理论值直接对比最直接的验证方式是把仿真结果与已知的解析结果对比。比如 H 门作用在 |0⟩ 上理论结果是 (|0⟩|1⟩)/√2对应向量 [1/√2; 1/√2]你可以直接用norm(psiResult - expectedPsi)来检查。这个差值的范数如果小于 1e-10说明门操作完全正确。% 验证 H 门结果的正确性 H quantum.gate.H; psi0 quantum.state([1;0]); psiComputed H * psi0; psiExpected [1/sqrt(2); 1/sqrt(2)]; errorNorm norm(psiComputed - psiExpected); if errorNorm 1e-10 disp(H 门结果验证通过); else fprintf(验证失败误差 %.4e\n, errorNorm); end这段代码的逻辑很简单算完一个操作后立即与理论值做差看误差范数。关键在于选择合适的误差阈值——浮点运算下误差不可能为 01e-10 是一个合理的安全线。如果你的结果与理论值差 0.1 这种量级不要犹豫一定是某个门矩阵写错或比特顺序搞反了。参数说明norm默认计算 2-范数对量子态误差来说这是最常用的度量。对于密度矩阵可以换成trace(rhoComputed - rhoExpected)的绝对值或者用保真度函数来度量两个态的距离。7.2 第二验证利用酉矩阵的数学性质检查门操作量子门必须满足酉性U†U I。这个性质可以用来检查你手动构造的门矩阵是否正确。比如你从论文里抄了一个看似正确的 4×4 矩阵第一件事就是检查酉性和行列式。% 检查手动构造的门矩阵是否是合法的酉矩阵 U [1 0 0 0; 0 1 0 0; 0 0 0 1; 0 0 1 0]; % 这是一个 CNOT 门 unitarityError norm(U * U - eye(4)); fprintf(酉性误差 %.4e\n, unitarityError); if unitarityError 1e-12 error(矩阵不是酉矩阵检查构造是否正确); end这段代码的作用是U是 U 的共轭转置U * U如果等于单位矩阵则 U 是酉矩阵。norm计算差矩阵的 2-范数如果大于 1e-12说明矩阵构造有问题。这个检查之所以重要是因为很多量子算法论文里的门矩阵可能印刷有误或抄录时手误在跑完整算法之前先排除这种低级错误能省大量时间。7.3 第三验证用已知量子算法做端到端测试验单门操作正确不算数真正稳妥的验证是跑一个已知输出的完整算法。最经典、也最简单的是 Deutsch-Jozsa 算法——它只需要几个门操作输出是确定性的适合做端到端验证。另一个推荐的是 Grover 搜索的两比特版本因为它的成功概率理论值精确可算。% 用两比特 Grover 算法验证库的正确性搜索 |11 态 % 初始化两个比特都为 |0 psi quantum.tensor(quantum.state([1;0]), quantum.state([1;0])); % 第一步对两个比特分别施加 H 门 psi kron(H, H) * psi; % 第二步应用 Oracle将 |11 的相位翻转其他不变 Oracle diag([1 1 1 -1]); psi Oracle * psi; % 第三步应用扩散算子Grover 扩散算符 D 2 * (psi * psi) / (psi * psi) - eye(4); % 注意这里应用的是投影算子 % 但 Grover 扩散算子的正确构造是 D 2|ss| - I|s 为均匀叠加态 s 0.5 * ones(4,1); GroverDiffuser 2 * (s * s) - eye(4); psi GroverDiffuser * psi; % 测量 prob11 abs(psi(4))^2; fprintf(搜索到 |11 态的概率 %.4f理论值应为 1.0\n, prob11);这段代码有一个重点扩散算子的构造用的是均匀叠加态 |s⟩ (|00⟩|01⟩|10⟩|11⟩)/2 对应的投影算子 2|s⟩⟨s| - I。有些初学者会混淆地使用当前态的投影算子——那是错的。正确逻辑是Grover 扩散算子是固定算符与当前态无关。跑完这个算法后测量概率如果接近 1.0说明库的矩阵运算、张量积操作、门定义全部正确。如果概率偏离很大就能定位到库的某个函数行为与理论预期不符合排查时先回到单门验证那一步。参数说明这里只用了一次 Grover 迭代。两比特搜索需要刚好一次迭代就能以概率 1 找到目标态这是教科书标准结论。如果你的仿真结果明显小于 1优先检查 Oracle 的矩阵构造——diag([1 1 1 -1])矩阵对角线上第 4 个元素是 -1对应 |11⟩ 相位翻转如果写成diag([-1 1 1 1])那就变成翻转 |00⟩ 了结果完全不同。我个人的经验教训是每写一段量子仿真代码至少配一种验证手段。验证的次序也很讲究——先验证门矩阵酉性再验证单比特态和理论值做差最后才验证多比特算法跑端到端。如果跳过了前两步直接跑 Grover一旦结果不对排查范围是整个流程非常痛苦。养成每步都检查的习惯之后量子仿真其实比经典仿真更好调试因为每一步的数学性质都非常明确不存在「跑通了但不知道为什么对」的模糊地带。希望这些方法和踩坑记录能帮你在 MATLAB 上顺利跑起自己的量子信息计算实验。本文还有配套的精品资源点击获取
简介一套专注于量子信息计算的MATLAB工具库适合科研人员、工程师及高年级学生用于量子算法仿真、量子态分析与量子测量实验。资源包内含279个文件以276个m格式脚本/函数为主体覆盖量子比特操作、量子门实现、量子线路模拟、密度矩阵求解、量子测量以及熵与纠缠计算等常用功能另有3个txt文档分别说明许可证、版本更新与快速上手指引压缩包整体仅163KB轻量且易集成到科研或教学流程。目前已有545人学习下载。库内目录模块划分明确如测量、纠缠与熵、工具函数等配合核心qlib库与startup.m启动配置可快速搭建仿真环境。借助量子熵、纠缠度等物理量计算用户能深入分析量子系统性质并基于MATLAB的数值计算与可视化能力直接用于量子信息教学演示、算法验证和论文复现。1. 量子信息计算库MATLAB不写 Python 也能正经做量子计算仿真先说一个反直觉的结论量子计算仿真这件事并不是 Python 生态的专利。很多搞通信、信号处理、电子信息方向的人平时主力工具就是 MATLAB遇到量子信息方向的需求——比如验证一个量子纠错码、仿真一个含噪信道下的量子态演化、跑一遍 Grover 搜索算法——第一反应是去装 QuTiP 或 Cirq结果被 Python 环境和依赖折腾一晚上。实际上MATLAB 上有一批开源的量子信息计算库能覆盖从量子态表示、矩阵运算、信道仿真到算法验证的完整链路而且和 MATLAB 原生的数值计算、绘图、并行工具箱无缝衔接。这类库解决的核心痛点有三个第一不用切语言在熟悉的 MATLAB 环境里就能复现量子算法论文里的实验第二量子态和算子的运算本质是线性代数MATLAB 的矩阵运算效率比 Python 原生循环高一个量级中小规模比特数15 比特以内的仿真完全跑得动第三绘图和数据分析是一条龙布洛赫球、保真度曲线、噪声扫描图直接在同一个脚本里出图。适合的人群很明确有 MATLAB 基础、需要把量子信息算法跑出数值结果来验证想法或完成课题的学生和工程师以及想快速评估某个量子方案可行性、不想被 Python 版本地狱拖累的研究人员。本篇就沿着「库怎么选 → 最小例子怎么跑 → 参数怎么设 → 噪声和常见坑在哪 → 怎么验证没算错」这条线展开。2. 都有哪些现成的库选型依据与边界2.1 开源社区的一块主力量子信息工具箱类库MATLAB 生态里做量子信息计算最常见的选择是某几个流传较广的开源工具箱。这类工具箱一般以 包Package的形式组织代码模仿 MATLAB 自带的工具箱风格安装后你在命令行敲quantum.state([1 0])或qinfo.state(...)就会返回一个量子态对象而不是一个冷冰冰的复矩阵。我一般会优先看这类库的三件事第一它是否实现了密度矩阵的数据结构而不是只有态矢量。因为真实噪声仿真必须用密度矩阵态矢量只能描述纯态遇到退极化信道或振幅阻尼信道就无能为力。第二它是否内置了信道仿真函数比如damping、depolarizing这些函数返回的是一个超级算子Superoperator或 Kraus 算子集可以直接作用在密度矩阵上。第三它是否提供测量和保真度相关的函数因为最终你要算的不是那个矩阵长什么样而是某个算法跑完之后的成功概率或保真度。选型时还要注意维护活跃度。有些库是好几年以前的毕业设计只支持到某个旧版本 MATLAB在新版本里会因为函数改名比如strsplit的行为变化而报错。我通常的验证方式是直接把库下载下来跑它自带的test或example脚本如果能一口气跑完再往自己的目标上靠如果连自带的示例都跑不通果断换下一个。2.2 官方与商业路线MATLAB 支持包的取舍除了开源库MATLAB 官方也有一个量子计算相关的支持包这可能是很多人不知道的。它在 Add-On Explorer 里可以直接搜到安装后提供量子比特qubit对象、量子电路circuit对象、测量操作和几个内置算法示例。这条路线的最大好处是和 MATLAB 版本绑定不会出现兼容性问题而且底层用了 MATLAB 的matlab.lang.makeValidName等一系列现代语言特性运行效率有保障。缺点是扩展性有限。官方支持包更像是「帮你快速上手量子电路可视化」的演示工具遇到论文里自定义的算符、非标准信道模型、或者需要对 Hamiltonian 做时间演化的场景它的抽象层次反而碍事。我的经验是如果你的目标是把一个已知算法比如 Grover、Shor 的入门版本跑通官方支持包够用如果你的目标是做研究、改信道模型、对比不同噪声参数下的性能那就得用开源工具箱这类更底层的东西。两者并不冲突甚至可以在同一个项目里共存——用开源库做核心矩阵运算用官方支持包的画图函数来出布洛赫球。2.3 怎么装让库在 MATLAB 里「永久可见」无论选哪条路线安装的关键动作只有一个把库所在的文件夹加入 MATLAB 路径。这一步看似简单实际是很多初学者的第一个坑——有人直接把整个库文件夹复制到toolbox目录下也不管结果每次启动 MATLAB 都提示找不到函数。% 将下载好的库文件夹路径加入 MATLAB 路径永久生效 libraryPath D:\Workspace\QuantumInfoLib; addpath(genpath(libraryPath)); % genpath 会递归包含所有子文件夹 savepath; % 保存路径下次启动 MATLAB 仍有效这段代码的逻辑是genpath(libraryPath)生成一个包含所有子文件夹的路径列表addpath把这个列表加入当前会话的搜索路径savepath把当前的路径设置写入 MATLAB 的路径配置文件pathdef.m。如果不加genpath只写addpath(D:\Workspace\QuantumInfoLib)那库目录下的二级子文件夹里的函数全部都无法访问这是最常见的「明明装了却找不到函数」的原因。参数说明libraryPath是库在你本机上的绝对路径注意不要带中文和空格某些版本的genpath对含空格路径的处理有玄学问题。如果你用的是官方支持包不需要手动addpath在 Add-On Explorer 里点安装后 MATLAB 自己会配好路径。2.4 装完先跑通一个内置示例再谈其他装完库的第一件事不是立刻写自己的代码而是跑一个库自带的示例脚本。这一步目的是确认三件事库函数能调、绘图没问题、矩阵运算核心没坏。以某开源库为例我一般先跑它的 GHZ 态制备示例% 运行库自带的 GHZ 态示例脚本 run example_ghz.m如果这个脚本能出图、能打印出态矢量的系数就说明库基本可用。如果报错先看错误类型Undefined function说明路径没加对Invalid syntax说明你的 MATLAB 版本太老库用了新语法特性。这里有个小习惯跑通后立刻执行which quantum.state或你用的库的核心构造函数看它返回的路径是不是在库目录下。如果which返回的是其他工具箱的同名函数说明函数名冲突了这是下一个大坑。3. 量子态初体验从 |0⟩ 到叠加态的完整计算3.1 构建量子态的最小命令态矢量与密度矩阵两种方式在 MATLAB 的量子信息库中创建量子态有两条路态矢量State Vector和密度矩阵Density Matrix。两者的区别很本质态矢量只能描述纯态密度矩阵既能描述纯态也能描述混合态。初学者刚上手时一定要同时学会这两种创建方式因为后续的噪声仿真全部要落在密度矩阵上。% 创建一个单量子比特的 |0 态态矢量方式 psi0 quantum.state([1; 0]); % 创建一个单量子比特的 |0 态密度矩阵方式 rho0 quantum.density([1 0; 0 0]);第一行代码中[1; 0]是列向量表示代表基态 |0⟩ 的概率幅为 1、|1⟩ 的概率幅为 0。第二行代码中[1 0; 0 0]是 2×2 矩阵表示密度矩阵只有一个非零对角元对应 |0⟩ 态的纯态密度矩阵。这里的核心逻辑是库内部会用quantum.state或quantum.density构造函数对输入数据做合法性检查——比如态矢量的范数是否为 1、密度矩阵是否满足迹为 1 且半正定。这也意味着你不小心传入一个未归一化的向量库会直接报错或警告而不是闷头算出一个错误结果。参数说明两种构造函数的入参都支持复数构造叠加态时可以直接传[1/sqrt(2); 1/sqrt(2)]。库的构造函数通常还接受第二参数字符串类型来标记这个态的名字方便在绘图时作为图例显示。3.2 单比特门操作用矩阵乘法实现量子逻辑门量子门在数学上就是酉矩阵作用在态上就是矩阵乘法。这是 MATLAB 做量子计算最舒服的一点——你不用写任何循环矩阵乘法的效率极高。以最常用的 Hadamard 门H 门为例它的作用是把 |0⟩ 变成叠加态 (|0⟩|1⟩)/√2。% 定义 Hadamard 门 H (1/sqrt(2)) * [1 1; 1 -1]; % 将 H 门作用在 |0 态上 psi_plus H * psi0;这里的H是 2×2 矩阵psi0是之前创建的态矢量列向量两者的矩阵维度自动匹配。计算结果是 psi_plus [1/√2; 1/√2]这正是等概率叠加态。如果你用的是密度矩阵方式量子门的操作不是直接左乘而是rho_new H * rho * H其中H是 H 的共轭转置。这个区别经常有人搞混——态矢量是 U|ψ⟩密度矩阵是 UρU†漏掉一个共轭转置结果就是完全错误的。参数说明实际项目中不要手写H (1/sqrt(2)) * [1 1; 1 -1]更稳妥的做法是调用库内置的quantum.gate.H或类似属性这样可以避免手算四舍五入误差。而且内置门定义更严格比如quantum.gate.X、quantum.gate.Z、quantum.gate.CNOT都是经过预定义验证的。3.3 布洛赫球可视化让量子态「看得到」做量子信息仿真和做通信仿真最大的不同是你看不到误码率曲线只能看到矩阵。布洛赫球是单量子比特态的可视化标准工具球面上的每一个点对应一个量子态。在 MATLAB 中绘图函数通常会读取态对象或密度矩阵计算其对应的布洛赫向量分量然后在三维球体上画出矢量箭头。% 将叠加态绘制到布洛赫球上 figure; quantum.plot.bloch(psi_plus); title(Bloch Sphere Representation of | State);代码执行后MATLAB 会弹出一个 3D 球体图球面上有一条从原点指向 (1,0,0) 方向的红色箭头这个坐标对应态 (|0⟩|1⟩)/√2。如果你绘制的是 |0⟩箭头指向 (0,0,1) 北极方向|1⟩ 指向 (0,0,-1) 南极方向。布洛赫球的一个重要用途是验证门操作是否正确比如对 |0⟩ 连续执行两次 H 门箭头应该回到北极如果没回来大概率是 H 门矩阵写错了或矩阵乘法方向错了。参数说明bloch类的绘图函数一般接受一个态对象或密度矩阵作为入参如果传入的是密度矩阵库会自动提取其 Bloch 向量。部分库还支持在同一张图上绘制多个态用元胞数组传入这在对比不同噪声强度下的态偏移时非常有用。3.4 测量与概率从量子态到经典比特的关键一步量子测量的结果是概率性的这是量子信息和经典信息最大的分水岭。在仿真中测量函数做的事很简单按概率分布随机采样一个结果然后返回测量后的态。% 对叠加态进行 Z 基测量返回测量结果和测量后状态 [measureOutcome, postMeasState] quantum.measure(psi_plus, Z); % 重复测量 1000 次统计 |0 的出现次数 count0 0; for idx 1:1000 [outcome, ~] quantum.measure(psi_plus, Z); if outcome 0 count0 count0 1; end end fprintf(测量到 |0 的次数%d频率%.3f\n, count0, count0/1000);这段代码的核心逻辑是quantum.measure接受态矢量和测量基字符串Z 表示计算基返回两个输出——第一个是测量结果0 或 1第二个是测量后坍缩的量子态。注意循环里每次测量都会重新随机采样所以 1000 次的频率会接近理论概率 50%。在实际算法验证中跑单次测量没有意义必须重复多次统计频率来验证概率分布。参数说明quantum.measure的第二个参数可选项包括 Z、X、Y分别对应泡利 Z、X、Y 基测量。如果你的仿真需要固定随机种子以便复现可以在循环前调用rng(42)这样每次跑出来的分布完全一致对调试非常关键。4. 多比特体系与纠缠态从单比特走向实用算法4.1 张量积的正确姿势把两个单比特合成一个双比特态真实量子算法不可能只跑单比特多比特体系的构建靠的是张量积运算。MATLAB 中kron函数可以直接做克罗内克积量子信息库通常会封装一个更语义化的函数比如quantum.tensor。% 构建双比特态第一个比特 |0第二个比特 |0 psiAB quantum.tensor(psi0, psi0); % 构建双比特 Bell 态 |Φ (|00 |11)/√2 bellState quantum.state((1/sqrt(2)) * [1; 0; 0; 1]);第一行代码中quantum.tensor(psi0, psi0)将两个 2×1 向量合成 4×1 向量对应四个计算基 |00⟩、|01⟩、|10⟩、|11⟩ 的概率幅分别是 1、0、0、0。第二行代码手动构造了一个 Bell 态4×1 列向量的第 1 个元素对应 |00⟩和第 4 个元素对应 |11⟩都是 1/√2这正是最大纠缠态。这里有个新手高频犯的错用kron时把两个向量的顺序搞反。kron(a, b)的结果是 a 的每个元素乘以整个 b对应量子信息中的习惯是第一个比特变化更慢。如果你的库里的约定是反过来的整个态就变了。所以用quantum.tensor而不是裸kron是更稳的做法。参数说明tensor函数接受任意数量的态矢量多个参数时从左到右依次执行张量积。4.2 双比特门的实现CNOT 是纠缠的发动机CNOT受控非门是双比特量子计算的核心它在经典对应的基础上增加了量子叠加的特性。它的矩阵是一个 4×4 酉矩阵当控制比特为 |1⟩ 时目标比特翻转控制比特为 |0⟩ 时目标比特不动。% 使用内置的 CNOT 门 CNOT quantum.gate.CNOT; % 将 CNOT 作用在 Bell 态的第一个比特控制和第二个比特目标上 psiCNOT CNOT * bellState;代码逻辑很简单CNOT是一个 4×4 矩阵bellState是 4×1 向量矩阵乘法直接得到新态。你可以验证如果把 CNOT 作用在 |11⟩ 态上结果变成 |10⟩作用在 |⟩|0⟩ 上则产生纠缠态。多比特门操作中最容易翻车的点不是门本身而是比特顺序。双比特态的向量中第一分量是 |00⟩第二分量是 |01⟩第三是 |10⟩第四是 |11⟩。CNOT 的控制比特可以是第一个比特也可以是第二个比特两者对应的 4×4 矩阵是不同的。库内置的CNOT通常默认控制比特为第一个、目标为第二个。如果你需要控制比特在第二个位置需要自己构造置换矩阵或直接用quantum.gate.CNOT的参数版本部分库支持传入控制比特和目标比特的索引。如果不支持就得手动重排矩阵行和列。4.3 用 CHSH 不等式验证纠缠仿真里的「贝尔测试」很多人做了 Bell 态仿真却不知道怎么验证这个态真的是纠缠的。一个可靠的量化指标是 CHSH 不等式对于局域隐变量模型CHSH 参数 S 的理论上限是 2而量子力学允许 S 最大达到 2√2 ≈ 2.828。在 MATLAB 里可以快速验证。% 定义两组测量角度单位弧度 thetaA [0, pi/2]; thetaB [pi/4, -pi/4]; % 计算关联函数 E(thetaA, thetaB) sigma_a ⊗ sigma_b E zeros(2, 2); for ii 1:2 for jj 1:2 % 构造测量算符单位向量方向上的泡利算符 sigmaA cos(thetaA(ii)) * quantum.gate.Z sin(thetaA(ii)) * quantum.gate.X; sigmaB cos(thetaB(jj)) * quantum.gate.Z sin(thetaB(jj)) * quantum.gate.X; % 计算期望值 measOp kron(sigmaA, sigmaB); E(ii, jj) real(bellState * measOp * bellState); end end % CHSH 参数 S E(0,π/4) - E(0,-π/4) E(π/2,π/4) E(π/2,-π/4) S E(1,1) - E(1,2) E(2,1) E(2,2); fprintf(CHSH 参数 S %.4f量子上限 2√2 ≈ 2.828\n, S);这段代码做的是对 Bell 态分别沿两组角度测量计算联合期望值最后组合成 CHSH 参数。核心逻辑是把测量角度映射为测量算符——角度 θ 对应的算符是 cosθ·Z sinθ·X——然后用kron构造联合测量算符计算其在 Bell 态下的期望值。如果你构造的不是最大纠缠态S 会明显低于 2.828如果是可分离态S 不会超过 2。参数说明这里的测量角度选取直接决定 S 值大小角度偏移一点S 就会明显下降这不是 bug是量子力学本身的性质。仿真中如果跑出 S 在 2.5~2.8 之间说明纠缠态制备基本正确。4.4 多比特体系的维度诅咒与内存边界MATLAB 仿真多比特体系最现实的约束是内存。一个 n 比特的态矢量有 2^n 个复数分量每个复数 double 占 16 字节。所以 20 比特需要 1,048,576 × 16 16 MB看起来不多但 28 比特需要 4 GB30 比特需要 16 GB这还不算门矩阵。在实际项目中我会在代码里加一个维度检查超过 25 比特直接报错。% 检查态矢量维度是否超出内存安全范围 nQubits 25; if nQubits 25 error(当前计算机内存无法安全仿真 %d 比特量子系统需要 %.1f GB, ... nQubits, 2^nQubits * 16 / 1024^3); end这个检查的核心逻辑是2^nQubits × 16 是态矢量所需的最小字节数在实际中还要考虑门矩阵、中间变量、绘图函数的额外开销所以我个人把安全线划在 25 比特。参数说明如果你用的是密度矩阵内存需求是态矢量的平方n 比特密度矩阵有 4^n 个元素所以密度矩阵仿真只能做到 14~15 比特。能用态矢量绝不用密度矩阵除非你在做噪声仿真。5. 含噪声信道仿真从理想量子态到贴近硬件的真实计算5.1 三种基础噪声模型振幅阻尼、退极化、相位阻尼真实的量子硬件和理想仿真最大的区别就是噪声。量子信息库中通常内置了几种标准噪声信道它们对密度矩阵的作用方式各不相同。最常用的三种是振幅阻尼对应能量耗散、退极化对应完全随机化、相位阻尼对应退相干。% 创建三种噪声信道的 Kraus 算符damping 参数是噪声强度 ampDamp quantum.channel.amplitudeDamping(0.1); depol quantum.channel.depolarizing(0.15); phaseDamp quantum.channel.phaseDamping(0.2); % 分别作用于同一个初始态上 rhoInitial quantum.density([1 0; 0 0]); rhoAfterAmp ampDamp(rhoInitial); rhoAfterDepol depol(rhoInitial); rhoAfterPhase phaseDamp(rhoInitial);这段代码的逻辑是信道对象本质上是一个「可以作用于密度矩阵的函数」quantum.channel.amplitudeDamping(0.1)返回一个 Kraus 算子集合用函数句柄封装起来。噪声强度参数的含义因信道而异振幅阻尼的参数 γ 表示能量耗散率0.1 意味着约 10% 的概率从 |1⟩ 掉到 |0⟩退极化的参数 p 表示完全随机化的概率0.15 时态有 15% 的概率变成完全混合态相位阻尼的参数 λ 表示退相干的强度。参数说明噪声强度的取值范围是 [0,1]0 表示无噪声1 表示完全噪声。这里有一个常见的理解误区把三类噪声参数当成同一类东西来调。实际仿真中振幅阻尼和相位阻尼是硬件里最主要的两种噪声退极化信道更多是一种理论简化模型。在做硬件相关仿真时优先调前两种。5.2 含噪态的可视化保真度曲线怎么画有了噪声信道下一步就是量化评估噪声对量子态的影响。最常用的指标是保真度 F ⟨ψ|ρ|ψ⟩它衡量的是实际态与理想态的重叠程度。在实际仿真中你会扫描噪声参数从 0 到 1画出保真度随噪声强度变化的曲线。% 扫描振幅阻尼噪声强度计算保真度变化 noiseStrength 0:0.02:1; fidelityValues zeros(size(noiseStrength)); rhoIdeal quantum.density([1; 0]); % 理想 |0 态 rhoInput quantum.density([0; 1]); % 输入 |1 态会受振幅阻尼影响 for idx 1:length(noiseStrength) channel quantum.channel.amplitudeDamping(noiseStrength(idx)); rhoNoisy channel(rhoInput); fidelityValues(idx) real(sqrt(rhoIdeal * rhoNoisy * rhoIdeal)); end % 绘制保真度曲线 figure; plot(noiseStrength, fidelityValues, LineWidth, 2); xlabel(振幅阻尼强度 γ); ylabel(保真度 F); grid on;这段代码的核心逻辑是准备一个理想目标态 |0⟩把 |1⟩ 态逐个通过不同强度的振幅阻尼信道然后计算输出态与 |0⟩ 的保真度。由于振幅阻尼的物理意义是 |1⟩ 向 |0⟩ 的衰减所以噪声越强输出态越靠近 |0⟩保真度反而越高——这个反直觉的结果经常让初学者困惑。参数说明这里的保真度用的是平方根形式取值范围 0~1。如果你的库返回负值或大于 1 的数值几乎可以确定是密度矩阵在数值上不再合法比如出现负特征值需要检查是否手动操作密度矩阵时漏掉了共轭转置。5.3 噪声仿真的避坑清单密度矩阵的物理合法性检查含噪声仿真中一个反复出现的坑是密度矩阵「物理合法性」检查。密度矩阵必须满足三个条件厄米性等于其共轭转置、半正定性所有特征值非负、迹为 1所有对角元之和等于 1。在含噪声环境下浮点运算误差和手动矩阵操作很容易破坏这三个条件。% 手动计算噪声后的密度矩阵这是错误的示范会破坏物理合法性 H quantum.gate.H; rhoWrong H * rhoInitial; % 漏掉了共轭转置结果不再是合法的密度矩阵 % 检查合法性 fprintf(迹 %.4f\n, trace(rhoWrong)); fprintf(厄米性误差 %.4e\n, norm(rhoWrong - rhoWrong)); % 正确方式酉变换是 rho - U * rho * U rhoCorrect H * rhoInitial * H;代码中第一个操作H * rhoInitial是典型的错误示范它只在态矢量表示时才正确对密度矩阵做变换必须左右两边各乘一次右边乘的是共轭转置。norm(rhoWrong - rhoWrong)的值如果远大于 1e-12说明矩阵不再厄米。使用库提供的信道函数时一般不会出这个问题但只要你手动用矩阵乘法去变换密度矩阵就必须执行合法性检查。注意很多库提供了quantum.validate或isValidState这类函数可以直接判断一个矩阵是否合法的密度矩阵。我的习惯是在任何手动构造或手动变换密度矩阵之后立刻调用一次这个检查比起后面拿着错误数据跑半天再回头找 bug多写这一行的成本几乎为零。5.4 一个完整的含噪 GHZ 态制备仿真流程把前面的概念串起来展示一个完整的小项目制备三比特 GHZ 态 (|000⟩|111⟩)/√2然后分别让它通过退极化信道最后测量保真度。% 创建三比特 GHZ 态H 门 两个 CNOT 门 psiInit quantum.tensor(quantum.state([1;0]), quantum.state([1;0]), quantum.state([1;0])); H quantum.gate.H; CNOT quantum.gate.CNOT; % 第一步对第一个比特施加 H 门作用于整个三比特空间H ⊗ I ⊗ I I eye(2); H3 kron(kron(H, I), I); psiAfterH H3 * psiInit; % 第二步CNOT(1,2)控制比特 1目标比特 2 CNOT12 kron(CNOT, I); psiAfterCNOT12 CNOT12 * psiAfterH; % 第三步CNOT(2,3)——这里需要自行构造控制第二比特翻转第三比特 % 先用置换矩阵将比特顺序交换再做 CNOT具体细节略这里展示了一个实际仿真会遇到的坎多比特门操作不是简单调用内置函数就行因为内置矩阵默认作用在「相邻的固定比特位」上。当你想让 CNOT 作用在第 2 和第 3 个比特上需要构造带单位阵的张量积或使用置换矩阵。参数说明kron(kron(H, I), I)的含义是 H 门作用在第一个比特上第二、第三比特不受影响。这个构造方式在多比特仿真中会反复使用。如果你的库提供了quantum.gate.applyTo或类似操作优先用库的方式来做——它内部会自动处理比特顺序和矩阵维度避免手写置换矩阵时出错。比如库的applyTo(psiAfterH, H, [1])明确表示把 H 门作用在第 1 个比特上。执行完三步操作后得到的态理论上就是 GHZ 态。要验证制备是否正确可以用 Fidelity 函数算一下它与理论 GHZ 态的保真度如果大于 0.999 可以认为制备成功。接着把密度矩阵版本通过退极化信道扫描不同噪声强度画出保真度衰减曲线——这就是你论文里的一张标准图。6. 避坑指南量子信息库换用了 5 个最值得记录的踩坑点6.1 函数名冲突which出来的路径不对现象调用quantum.state时 MATLAB 报错或返回的结果完全不符合预期比如创建出的态维数不对。用which quantum.state -all查看发现命中的是 MATLAB 自带工具箱或另一个无关库的同名函数。原因不同工具箱和库之间函数名撞车。量子信息库里的state、measure、gate这些名字太通用很容易与 MATLAB 其他工具箱的函数冲突。解决在addpath时把量子库的路径放在路径列表的最前面这样 MATLAB 搜索时会优先命中它。做法是addpath(genpath(libraryPath), -begin)其中-begin表示把新路径添加到搜索路径开头。如果还是不行就用which定位到具体冲突的函数手动改名后调用比如写成quantum.state(...)的包名调用形式。这是最稳妥的。6.2 矩阵维度对不上kron和tensor的区别没搞清楚现象用kron(psiA, psiB)构造双比特态后后续计算保真度时维度对不上或者算出来的结果与论文不一致。原因kron的结果是「数学上的克罗内克积」量子信息中的张量积在比特排序约定上可能与之相反。不同库对第一个参数对应「第一个比特」还是「最后一个比特」的定义不同。解决不要直接用kron去构造量子态改用库封装的quantum.tensor。如果必须用kron先在文档中查清楚库的比特排序约定然后用一个已知结果做验证——比如计算 |01⟩ 的 tensor 表示看第二分量是否为 1。这一步验证只需要几秒钟但能省几个小时的排查时间。6.3 密度矩阵变换漏了共轭转置结果看起来合理其实全错现象对密度矩阵做 H 门变换后保真度和布洛赫球箭头方向都不对但程序不报错、矩阵维度也对看起来像「能跑出结果」的错误。原因态矢量的酉变换是 U|ψ⟩密度矩阵的酉变换是 UρU†U† 是 U 的共轭转置。很多人从态矢量的写法直接迁移到密度矩阵只做了左乘漏了右乘。解决实在怕漏就用信道函数来做变换channel quantum.channel.unitary(H)然后rhoNew channel(rho)。这样本质上把 UρU† 的变换封装在库函数内部你想漏都漏不掉。如果一定要手写写完后立刻检查trace(rhoNew)是否仍等于 1近似以及rhoNew是否厄米。6.4 内存爆炸试图用密度矩阵仿真 20 比特以上系统现象跑一个 18 或 20 比特的密度矩阵仿真MATLAB 直接内存溢出或电脑风扇狂转、响应极慢。原因密度矩阵的元素数量是 4^nn20 时是 1.1 万亿个复数需要 PB 级内存这显然不可能。这不是库的问题是物理规律限制了经典计算机仿真量子系统的能力。解决能改用态矢量就绝不用密度矩阵。先确认你的仿真目标是否真的需要密度矩阵——只有涉及噪声、混合态的时候才需要纯态算法仿真Grover、Shor 等理想情况全程用态矢量就够了。如果确需密度矩阵且比特数超过 14考虑改为采样法Monte Carlo 波函数法多次运行态矢量仿真并对结果取平均这种方式能推高到 25 比特左右。6.5 随机测量结果不可复现没有固定随机种子现象同一个脚本跑两次测量统计结果不同用if分支判断测量结果时两次运行走进了不同的代码路径导致后续调试困难。原因量子测量函数内部调用了随机数生成器默认每次启动 MATLAB 随机种子都不同。这不是 bug是概率仿真的正常行为。解决在脚本开头或测量循环前调用rng(42)42 是自己选的种子任意整数都行这样整个会话的随机数序列固定测量结果序列也就固定了。注意rng是全局的如果在同一个脚本里同时做噪声采样和测量采样两者共享随机数序列要复现就必须保证整个脚本的调用顺序不变。7. 验证方法怎么确定你算出来的量子态没算错7.1 第一性验证与理论值直接对比最直接的验证方式是把仿真结果与已知的解析结果对比。比如 H 门作用在 |0⟩ 上理论结果是 (|0⟩|1⟩)/√2对应向量 [1/√2; 1/√2]你可以直接用norm(psiResult - expectedPsi)来检查。这个差值的范数如果小于 1e-10说明门操作完全正确。% 验证 H 门结果的正确性 H quantum.gate.H; psi0 quantum.state([1;0]); psiComputed H * psi0; psiExpected [1/sqrt(2); 1/sqrt(2)]; errorNorm norm(psiComputed - psiExpected); if errorNorm 1e-10 disp(H 门结果验证通过); else fprintf(验证失败误差 %.4e\n, errorNorm); end这段代码的逻辑很简单算完一个操作后立即与理论值做差看误差范数。关键在于选择合适的误差阈值——浮点运算下误差不可能为 01e-10 是一个合理的安全线。如果你的结果与理论值差 0.1 这种量级不要犹豫一定是某个门矩阵写错或比特顺序搞反了。参数说明norm默认计算 2-范数对量子态误差来说这是最常用的度量。对于密度矩阵可以换成trace(rhoComputed - rhoExpected)的绝对值或者用保真度函数来度量两个态的距离。7.2 第二验证利用酉矩阵的数学性质检查门操作量子门必须满足酉性U†U I。这个性质可以用来检查你手动构造的门矩阵是否正确。比如你从论文里抄了一个看似正确的 4×4 矩阵第一件事就是检查酉性和行列式。% 检查手动构造的门矩阵是否是合法的酉矩阵 U [1 0 0 0; 0 1 0 0; 0 0 0 1; 0 0 1 0]; % 这是一个 CNOT 门 unitarityError norm(U * U - eye(4)); fprintf(酉性误差 %.4e\n, unitarityError); if unitarityError 1e-12 error(矩阵不是酉矩阵检查构造是否正确); end这段代码的作用是U是 U 的共轭转置U * U如果等于单位矩阵则 U 是酉矩阵。norm计算差矩阵的 2-范数如果大于 1e-12说明矩阵构造有问题。这个检查之所以重要是因为很多量子算法论文里的门矩阵可能印刷有误或抄录时手误在跑完整算法之前先排除这种低级错误能省大量时间。7.3 第三验证用已知量子算法做端到端测试验单门操作正确不算数真正稳妥的验证是跑一个已知输出的完整算法。最经典、也最简单的是 Deutsch-Jozsa 算法——它只需要几个门操作输出是确定性的适合做端到端验证。另一个推荐的是 Grover 搜索的两比特版本因为它的成功概率理论值精确可算。% 用两比特 Grover 算法验证库的正确性搜索 |11 态 % 初始化两个比特都为 |0 psi quantum.tensor(quantum.state([1;0]), quantum.state([1;0])); % 第一步对两个比特分别施加 H 门 psi kron(H, H) * psi; % 第二步应用 Oracle将 |11 的相位翻转其他不变 Oracle diag([1 1 1 -1]); psi Oracle * psi; % 第三步应用扩散算子Grover 扩散算符 D 2 * (psi * psi) / (psi * psi) - eye(4); % 注意这里应用的是投影算子 % 但 Grover 扩散算子的正确构造是 D 2|ss| - I|s 为均匀叠加态 s 0.5 * ones(4,1); GroverDiffuser 2 * (s * s) - eye(4); psi GroverDiffuser * psi; % 测量 prob11 abs(psi(4))^2; fprintf(搜索到 |11 态的概率 %.4f理论值应为 1.0\n, prob11);这段代码有一个重点扩散算子的构造用的是均匀叠加态 |s⟩ (|00⟩|01⟩|10⟩|11⟩)/2 对应的投影算子 2|s⟩⟨s| - I。有些初学者会混淆地使用当前态的投影算子——那是错的。正确逻辑是Grover 扩散算子是固定算符与当前态无关。跑完这个算法后测量概率如果接近 1.0说明库的矩阵运算、张量积操作、门定义全部正确。如果概率偏离很大就能定位到库的某个函数行为与理论预期不符合排查时先回到单门验证那一步。参数说明这里只用了一次 Grover 迭代。两比特搜索需要刚好一次迭代就能以概率 1 找到目标态这是教科书标准结论。如果你的仿真结果明显小于 1优先检查 Oracle 的矩阵构造——diag([1 1 1 -1])矩阵对角线上第 4 个元素是 -1对应 |11⟩ 相位翻转如果写成diag([-1 1 1 1])那就变成翻转 |00⟩ 了结果完全不同。我个人的经验教训是每写一段量子仿真代码至少配一种验证手段。验证的次序也很讲究——先验证门矩阵酉性再验证单比特态和理论值做差最后才验证多比特算法跑端到端。如果跳过了前两步直接跑 Grover一旦结果不对排查范围是整个流程非常痛苦。养成每步都检查的习惯之后量子仿真其实比经典仿真更好调试因为每一步的数学性质都非常明确不存在「跑通了但不知道为什么对」的模糊地带。希望这些方法和踩坑记录能帮你在 MATLAB 上顺利跑起自己的量子信息计算实验。本文还有配套的精品资源点击获取
