求解 H₂ 基态能量

适用版本 · ArcQML 0.1.0

本篇展示关键代码片段,完整程序和数据随 ArcQML 发行包提供。配对代码框内的 Rust / Python 按钮会同步切换全篇示例,并在 VQE、QNN 之间记住选择;单段代码和命令保持显示。

本页目录

ArcQML 是一个用 Rust 语言原生实现的量子机器学习框架,提供了量子电路构建、状态向量模拟、可观测量期望值计算、自动微分和参数优化等核心能力。框架既支持单个量子态的计算,也支持批量状态模拟,并可使用 Pauli 算符及其线性组合表示可观测量。量子电路计算得到的期望值以可微分张量形式输出,可以继续参与损失函数及其他经典可微运算;执行反向传播时,梯度能够沿完整计算图传递至量子电路参数。对于含参量子电路,ArcQML 可以利用伴随法计算电路参数梯度,再与 Adam、SGD 等经典优化器配合完成混合量子—经典算法的训练。

这些能力可以用于变分量子本征求解、量子神经网络、量子态分析和小规模酉矩阵综合等任务。由于 Rust 语言具备原生编译、内存安全和高效并行等特性,ArcQML 能在保持可靠性的同时提供良好的计算性能。

本教程使用 ArcQML 实现一个完整的变分量子本征求解器 (Variational Quantum Eigensolver,VQE),估计氢分子 H₂ 在给定几何结构、基组和活性空间下的基态能量。教程将依次说明电子结构问题如何转换为量子比特哈密顿量 (Hamiltonian),如何构造量子电路与制备量子态,以及如何计算能量、获得梯度并更新参数。完整程序分别位于 examples/rust/h2_vqe.rsexamples/python/h2_vqe.py

1. 准备运行环境

请先按安装与环境配置准备发行包和运行环境。Rust 用户完成 Rust 与 Runtime 配置;Python 用户完成 wheel 安装即可,无需安装 Rust 或设置 Runtime 环境变量。

2. VQE 拟解决的问题

2.1 完整计算流程

一个分子可以处于许多能量不同的状态,其中能量最低的状态称为基态。分子的基态能量是预测分子稳定性、化学反应和材料性质的重要基础量。随着体系规模增大,直接在经典计算机上求大分子的最低能量所需的存储空间和计算时间呈指数级上涨。变分量子本征求解器 (Variational Quantum Eigensolver,VQE) 将量子设备用于制备参数化量子态并测量能量期望值,将经典计算机用于更新变分电路参数。通过这种量子和经典混合的方式,VQE 有效避免了在经典计算机中存储和对角化完整高维矩阵所带来的巨大开销。

VQE 是一种用于估计量子系统基态能量的量子—经典混合算法。设系统哈密顿量 HH 的最小本征值为 E0E_0,则对于任意归一化量子态 ψ\lvert\psi\rangle,变分原理给出:

ψHψE0.\langle\psi\rvert H\lvert\psi\rangle\ge E_0.

因此,可以构造包含参数 θ\boldsymbol{\theta} 的量子电路 U(θ)U(\boldsymbol{\theta}) 。选取易于制备的量子参考态 ϕ0\lvert\phi_{\mathrm{0}}\rangle,利用参数化量子电路将其演化为参数化量子态 ψ(θ)\lvert\psi(\boldsymbol{\theta})\rangle

ψ(θ)=U(θ)ϕ0.\lvert\psi(\boldsymbol{\theta})\rangle =U(\boldsymbol{\theta})\lvert\phi_{\mathrm{0}}\rangle.

量子化学中的电子哈密顿量通常先写成由费米子产生算符和湮灭算符组成的二次量子化形式。选定活性电子和活性轨道后,再通过 Jordan–Wigner 等费米子到量子比特的映射,将其转换为量子计算可以处理的量子比特哈密顿量:

H=jcjPj,H=\sum_j c_jP_j,

其中 cjc_j 是实数系数,PjP_j 是由 Pauli X、Y、Z 和单位算符组成的 Pauli 字符串。最终对应的目标函数是该参数化量子态对哈密顿量的能量期望值:

E(θ)=ψ(θ)Hψ(θ)=jcjψ(θ)Pjψ(θ).E(\boldsymbol{\theta})= \langle\psi(\boldsymbol{\theta})\rvert H \lvert\psi(\boldsymbol{\theta})\rangle =\sum_j c_j \langle\psi(\boldsymbol{\theta})\rvert P_j \lvert\psi(\boldsymbol{\theta})\rangle.

经典优化器依据测量结果迭代更新参数,使 E(θ)E(\boldsymbol{\theta}) 逐步降低。最终得到的最小能量是基态能量 E0E_0 的变分上界;当参数化电路能够表示基态且优化充分收敛时,该结果可逼近真实基态能量。

2.2 化学设置

本例采用以下量子化学设置:

配置 数值
分子几何 H (0, 0, -0.35 Å),H (0, 0, 0.35 Å)
原子间距 0.70 Å
基组 STO-3G
活性电子数 2
活性轨道数 2
费米子映射 Jordan–Wigner
量子比特数 4
Hamiltonian 项数 15
固定粒子数子空间维数 6
Hartree–Fock 能量 -1.1173490350562805 Ha
精确基态能量 -1.1361894542078266 Ha

3. 初始化设置

首先导入 ArcQML 和标准库中的必要接口:

Rust / Python
use arcqml::prelude::*;

// `AppResult<T>` 表示函数成功时返回 T,失败时返回可打印的错误
// 后续代码中的 `?` 会把错误自动交给调用者处理
type AppResult<T> = std::result::Result<T, Box<dyn std::error::Error>>;
from __future__ import annotations

import json
import math
from pathlib import Path

import arcqml

示例中将主要数值集中声明为常量:

Rust / Python
const NUM_QUBITS: usize = 4;       // H₂ 映射后需要 4 个量子比特
const ACTIVE_ELECTRONS: usize = 2; // Hartree–Fock 初态占据前 2 个位置
const LAYERS: usize = 6;           // 可训练电路重复 6 层
const STEPS: usize = 100;          // 共更新参数 100 次

const EXACT_GROUND_ENERGY: f64 = -1.136_189_454_207_826_6;    // 精确基态能量
const HARTREE_FOCK_ENERGY: f64 = -1.117_349_035_056_280_5;    // Hartree–Fock 能量(初态能量)

const CHEMICAL_ACCURACY: f64 = 1.6e-3;    // 化学精度误差阈值
NUM_QUBITS = 4       # H₂ 映射后需要 4 个量子比特
ACTIVE_ELECTRONS = 2 # Hartree–Fock 初态占据前 2 个位置
LAYERS = 6           # 可训练电路重复 6 层
STEPS = 100          # 共更新参数 100 次

EXACT_GROUND_ENERGY = -1.136_189_454_207_826_6    # 精确基态能量
HARTREE_FOCK_ENERGY = -1.117_349_035_056_280_5    # Hartree–Fock 能量(初态能量)

CHEMICAL_ACCURACY = 1.6e-3    # 化学精度误差阈值

其中层数和训练步数是本示例选择的超参数,可以在后续实验中调整。

4. 哈密顿量构造

本示例的哈密顿量由第三方量子计算工具 PennyLane 生成。Jordan–Wigner 变换后,H₂ Hamiltonian 是 15 个 Pauli 字符串的实系数线性组合。ArcQML 的 Rust API 使用 SparsePauliOp,Python API 使用 PauliSum 保存若干带实数系数的 Pauli 项。Pauli X、Y、Z 是描述单个量子比特变化或测量的基本算符;把它们及其系数相加,就能表达本例的分子能量。Rust 通过 include_str! 在编译时读取 JSON,Python 则在运行时读取并逐项构造对象:

Rust / Python
let hamiltonian = SparsePauliOp::from_json(include_str!(
    "../data/h2_hamiltonian.json"
))?;

// 确保导入的哈密顿量是 4 量子比特以及 15 项
assert_eq!(hamiltonian.num_qubits(), 4);
assert_eq!(hamiltonian.len(), 15);
DATA_PATH = Path(__file__).resolve().parents[1] / "data" / "h2_hamiltonian.json"


def load_hamiltonian(path: Path) -> arcqml.PauliSum:
    # Python 公开 API 目前暂未实现直接接收 JSON 文本,因此逐项构造 PauliSum
    payload = json.loads(path.read_text(encoding="utf-8"))
    hamiltonian = arcqml.PauliSum(num_qubits=payload["num_qubits"])

    for term in payload["terms"]:
        coefficient = term["coefficient"]
        operations = term["paulis"]
        if not operations:
            # 单位算符
            hamiltonian.add_identity(coefficient=coefficient)
            continue
        hamiltonian.add_term(
            paulis="".join(operation["pauli"] for operation in operations),
            qubits=[operation["qubit"] for operation in operations],
            coefficient=coefficient,
        )
    return hamiltonian


hamiltonian = load_hamiltonian(path=DATA_PATH)
# 确保导入的哈密顿量是 4 量子比特以及 15 项
assert hamiltonian.num_qubits == 4
assert hamiltonian.num_terms == 15

ArcQML 的 JSON 格式由顶层 num_qubitsterms 组成;每个 Pauli 操作显式记录 qubitpauli。例如常数项和 X(q0)X(q1)Y(q2)Y(q3) 项写作:

json
{
  "num_qubits": 4,
  "terms": [
    {
      "coefficient": -0.042078985845795724,
      "paulis": []
    },
    {
      "coefficient": -0.04475014386992153,
      "paulis": [
        { "qubit": 0, "pauli": "X" },
        { "qubit": 1, "pauli": "X" },
        { "qubit": 2, "pauli": "Y" },
        { "qubit": 3, "pauli": "Y" }
      ]
    }
  ]
}

paulis 数组表示单位算符。

5. 初始量子态制备

Circuit::new(4) 创建一条 4 量子比特电路。模拟器的默认初始状态是 0000\lvert0000\rangle,即四个位置都未被占据。Jordan–Wigner 映射下,本例的两个活性电子占据编号最小的两个位置,因此对量子比特 0 和 1 施加 X 门,把对应位从 0 翻转为 1:

Rust / Python
let mut circuit = Circuit::new(NUM_QUBITS)?;

// 依次对 qubit 0 和 qubit 1 施加 X 门
for qubit in 0..ACTIVE_ELECTRONS {
    circuit.x(qubit)?;
}
circuit = arcqml.Circuit(num_qubits=NUM_QUBITS)

# 依次对 qubit 0 和 qubit 1 施加 X 门
for qubit in range(ACTIVE_ELECTRONS):
    circuit.x(qubit=qubit)

ArcQML 显示状态字符串时,量子比特编号从大到小排列,即 q3q2q1q0q_3q_2q_1q_0。所以量子比特 0 和 1 被置为 1 后,状态显示为 0011\lvert0011\rangle

6. 变分量子电路构造

6.1 角度初始化

代码中的初始化公式是:

θk=0.04sin(0.37k),k=1,2,,48.\theta_k=0.04\sin(0.37k),\qquad k=1,2,\ldots,48.

其中 0.04 是角度幅度,单位为弧度。它将所有旋转门的初始角限制在 [0.04,0.04][-0.04,0.04],避免一开始就施加很大的单量子比特旋转。而 0.37 是相邻参数的固定相位步长,用于让相邻参数得到不同的非对称初值,避免 48 个角度全部相同。正弦函数让角度随参数编号平滑变化,并保证数值有正有负。

6.2 电路构造

本例使用 6 层硬件有效参数化量子线路。每层首先在每个量子比特上施加 RY、RZ,随后使用环形 CNOT 纠缠。对应代码如下:

Rust / Python
let mut parameter_index = 0usize;

for _layer in 0..LAYERS {
    for qubit in 0..NUM_QUBITS {
        // 初始化角度
        let ry = 0.04 * (0.37 * (parameter_index + 1) as f64).sin();
        parameter_index += 1;

        let rz = 0.04 * (0.37 * (parameter_index + 1) as f64).sin();
        parameter_index += 1;

        // ry 和 rz 函数会把传入角度登记为可训练参数
        circuit.ry(ry, qubit)?;
        circuit.rz(rz, qubit)?;
    }

    for control in 0..NUM_QUBITS {
        // 环形 CNOT 纠缠
        let target = (control + 1) % NUM_QUBITS;
        circuit.cnot(control, target)?;
    }
}
parameter_index = 0

for _layer in range(LAYERS):
    for qubit in range(NUM_QUBITS):
        # 初始化角度
        ry = 0.04 * math.sin(0.37 * (parameter_index + 1))
        parameter_index += 1

        rz = 0.04 * math.sin(0.37 * (parameter_index + 1))
        parameter_index += 1

        # ry 和 rz 函数会把传入角度登记为可训练参数
        circuit.ry(angle=ry, qubit=qubit)
        circuit.rz(angle=rz, qubit=qubit)

    for control in range(NUM_QUBITS):
        # 环形 CNOT 纠缠
        circuit.cnot(
            control=control,
            target=(control + 1) % NUM_QUBITS,
        )

电路结构如下:

电路结构

该硬件有效参数化量子线路不严格保持粒子数。它适合演示通用 VQE 和 ArcQML 自动微分。在化学计算中也可以改用保持粒子数的激发电路。

7. 模拟器创建和优化器选择

7.1 状态向量模拟器

Rust / Python
let simulator = StateVectorSimulator::new(NUM_QUBITS)?;
simulator = arcqml.StateVectorSimulator(num_qubits=NUM_QUBITS)

模拟器从 0000\lvert0000\rangle 开始执行整条电路。Hartree–Fock 制备门已经放在线路最前面,所以每次计算都会先得到 0011\lvert0011\rangle,再执行后续可训练量子电路。

Rust 的 simulator.run(&circuit, &hamiltonian) 与 Python 的 simulator.run(circuit=circuit, observable=hamiltonian) 都不会永久改变模拟器保存的初始状态,因此同一个模拟器可以在全部训练步骤中重复使用。

7.2 Adam 优化器

普通梯度下降对所有参数使用同一个固定更新尺度。Adam 会为每个参数分别记录近期梯度的平均趋势和梯度平方的平均趋势,据此自动调节每个参数的实际更新幅度。

Rust / Python
let mut optimizer = Adam::new(
    0.05,  // 学习率
    0.9,   // beta1
    0.999, // beta2
    1e-8,  // epsilon
    0.0,   // weight_decay
)?;
# 其余可选参数采用与 Rust 代码相同的默认值
optimizer = arcqml.Adam(learning_rate=0.05)

该超参数组合不适用于所有分子和电路。如果训练效果不好,该部分超参数也可能需要重新调整。

8. 执行 VQE 优化

一次训练迭代包含四个步骤:

  1. 前向计算当前能量。
  2. 从能量反向计算 48 个参数的梯度。
  3. Adam 根据梯度更新参数。
  4. 清空本轮梯度,避免它们与下一轮结果累加。
Rust / Python
for step in 1..=STEPS {
    // 前向计算当前能量
    let energy = simulator.run(&circuit, &hamiltonian)?;

    // 反向计算所有可训练角度的梯度
    energy.backward()?;

    // 优化器读取梯度并更新参数
    optimizer.step(circuit.parameters())?;

    // 清空本轮梯度
    optimizer.zero_grad(circuit.parameters());

    // 第 1 步和之后每 10 步打印一次更新后的能量
    if step == 1 || step % 10 == 0 {
        let current = simulator.run(&circuit, &hamiltonian)?.value()?;
        println!("step {step:>3}: energy = {current:.12} Ha");
    }
}
for step in range(1, STEPS + 1):
    # 前向计算当前能量
    energy_tensor = simulator.run(circuit=circuit, observable=hamiltonian)

    # 反向计算所有可训练角度的梯度
    energy_tensor.backward()

    # 优化器读取梯度并更新参数
    optimizer.step(circuit=circuit)

    # 清空本轮梯度
    optimizer.zero_grad(circuit=circuit)

    # 第 1 步和之后每 10 步打印一次更新后的能量
    if step == 1 or step % 10 == 0:
        with arcqml.no_grad():
            current = simulator.run(
                circuit=circuit,
                observable=hamiltonian,
            ).item()
        print(f"step {step:>3}: energy = {current:.12f} Ha")

run 返回的是 Tensor,调用 backward() 时,ArcQML 根据这份记录自动计算梯度。第一次迭代前参数还没有梯度,所以不需要预先调用 zero_grad;第一次反向之后,每轮都必须清零。

9. 运行完整示例

9.1 运行 Rust 版本

安装页的 Rust 配置说明已经设置 ARCQML_RUNTIME_LIB_DIR 的同一个 Windows PowerShell 或 Linux 终端中执行:

text
cargo run --release -p arcqml --example h2_vqe

9.2 运行 Python 版本

先按照安装页的 Python 配置说明安装 wheel 并激活 Conda 环境,再从 ArcQML 根目录执行:

text
conda activate arcqml-example
python examples/python/h2_vqe.py

9.3 输出示例

程序每 10 步输出一次能量以及相对精确基态能量的绝对误差,最后判断是否达到化学精度:

text
H2/STO-3G, Jordan-Wigner, 4 qubits
Hamiltonian terms : 15
Trainable params  : 48
Hartree-Fock      : -1.117349035056 Ha
Exact ground      : -1.136189454208 Ha
step   1: energy = -0.401257301811 Ha, |error| = 7.349e-1 Ha
...
step 100: energy = -1.136145412038 Ha, |error| = 4.404e-5 Ha

VQE energy        : -1.136145412038 Ha
Absolute error    : 4.404217e-5 Ha
Chemical accuracy: reached (threshold 1.6e-3 Ha)

10. 如何理解结果

11. 下一步

可以在此示例上继续尝试: