复算一条原初谱，并检验求解器

先预测，再运行

先把预期写下来。$R^2$ 高原取 $N_*=55$，最低阶公式给出 $n_s\approx0.9636$、$r\approx0.004$；$\tfrac12m^2\phi^2$ 给出 $n_s\approx0.9640$、$r\approx0.144$。

台阶的影响也先估一下：台阶放在结束前 $59$ 个 e 倍数，受影响的应该是比主尺度早约四个 e 倍数离开视界的模式，也就是 $k\approx10^{-3}\,\mathrm{Mpc}^{-1}$ 附近，主尺度本身不受影响。

对应阅读：Planck Collaboration, Planck 2018 results. X. Constraints on inflation, arXiv:1807.06211v2 · Original background and mode calculations with an independent Python reference

取得文件，固定运行环境

下载下方的 ZIP，解压到一个新目录。里面有参考求解器、基准脚本、比较脚本、五个输入模型和对应的预期输出，不含浏览器 Engine。

在虚拟环境里安装固定版本。预期结果由 Python 3.9.6、NumPy 2.0.2 和 SciPy 1.13.1 生成，每个输出文件都记下了产生它的软件版本。

python3 -m venv venv
. venv/bin/activate
pip install -r requirements.txt

对应阅读：Original background and mode calculations with an independent Python reference · CLASS explanatory.ini, lesgourg/class_public commit c5d93d01990116162b408b789b62ec5669fc4ea4

先用精确解检验求解器

幂律暴胀 $V=V_0e^{-\lambda\phi}$ 的 $\epsilon=\lambda^2/2$ 恒定，标量谱和张量谱都有闭式解，见 Stewart 与 Lyth 的式31、41、43。基准脚本拿数值谱和这个闭式解比较，而不是和慢滚估计比较。

我们的运行覆盖 $\lambda$ 从 $0.02$ 到 $0.3$。功率谱最大误差为 $2.0\times10^{-6}$，Wronskian 最大偏差为 $4.0\times10^{-7}$。最低阶慢滚公式在 $\epsilon=2\times10^{-4}$ 时偏差 $1.1\times10^{-4}$，在 $\epsilon=0.045$ 时偏差 $2.6\%$。检验因此要用闭式解，而不是慢滚近似。

只要有一项检查不通过，脚本就以错误退出。在你的机器上通过之前，不要往下做。

python verify_benchmarks.py my-benchmarks.json

对应阅读：Stewart & Lyth, A more accurate analytic calculation of the spectrum of cosmological perturbations produced during inflation (1993), gr-qc/9302019v1 · Original background and mode calculations with an independent Python reference

复算一个固定模型

每个输入文件写明势能、以 $M_{\mathrm{pl}}^4$ 为单位的振幅 $V_0$、初始场值、$N_*$、台阶参数，以及 $27$ 个以 $\mathrm{Mpc}^{-1}$ 为单位的波数。台阶采用 Adams、Cresswell 与 Easther 式1的 tanh 形式，位置和宽度以场值为单位。

输出先给背景：总 e 倍数、主尺度处的场值和 $\epsilon$；然后每个波数一行，含标量与张量功率、慢滚估计、离开视界的时刻和当时的场值。每个预设约需一秒。

在相同软件版本下，我们全新运行的结果与预期文件完全一致。换一个 SciPy 版本，最后几位数字可能变化；比较脚本的门槛是 $10^{-4}$，远大于求解器误差，也远小于本课程讨论的任何物理效应。

python recompute_inflation_spectrum.py inputs/step-ripples.json my-step.json
python compare_with_expected.py my-step.json expected/step-ripples.json

对应阅读：Adams, Cresswell & Easther, Inflationary perturbations from a potential with a step (2001), astro-ph/0102236v2 · Original background and mode calculations with an independent Python reference

亲手测出谱指数和张量比

在主尺度附近取三个点作输入，对功率的对数做差分。下面的代码以高原为例。

我们的结果：高原 $n_s=0.96490$、$r=0.00355$；加上默认台阶，小数点后四位不变；二次势 $n_s=0.96338$、$r=0.14403$；更短的再加热（$N_*=46$）$n_s=0.95826$、$r=0.00498$。

振幅要小心处理。高原的 $V_0$ 是按 $A_s=2.0989\times10^{-9}$ 选的，对应 Planck 的 $\ln(10^{10}A_s)=3.044$。同样的 $V_0$ 下，$N_*=46$ 的模型给出 $1.4880\times10^{-9}$。由于固定 $N_*$ 时 $\mathcal{P}_{\mathcal{R}}\propto V_0$，把 $V_0$ 乘以 $1.4105$ 就能恢复，谱指数和 $r$ 都不变。

import json, math, subprocess
model = json.load(open('inputs/smooth-plateau.json'))
model['kGridMpc'] = [0.049, 0.05, 0.051]
json.dump(model, open('tilt-input.json', 'w'))
subprocess.run(['python', 'recompute_inflation_spectrum.py', 'tilt-input.json', 'tilt-output.json'], check=True)
low, mid, high = json.load(open('tilt-output.json'))['spectrum']
print('n_s =', 1 + math.log(high['scalar'] / low['scalar']) / math.log(0.051 / 0.049))
print('r   =', mid['tensor'] / mid['scalar'])

对应阅读：Original background and mode calculations with an independent Python reference · Planck Collaboration, Planck 2018 results. X. Constraints on inflation, arXiv:1807.06211v2

与已发表的限制对照

Planck 2018 X 给出 $n_s=0.9649\pm0.0042$、$\ln(10^{10}A_s)=3.044\pm0.014$；PDG 2025 给出 $r<0.036$（$95\%$ 置信度）。高原正好落在谱指数色带中央。二次势的谱指数符合，$r$ 却是上限的四倍。$N_*=46$ 的高原谱指数偏低 $1.6\sigma$。

这些都是原初参数之间的比较，而这些参数本身是在幂律谱假设下拟合出来的。带台阶的模型必须经过 Boltzmann 程序和似然函数才能比较。CLASS 可以通过 external_Pk 选项接收你自己的谱表，这是自然的下一步；本课程停在原初谱。

对应阅读：Planck Collaboration, Planck 2018 results. X. Constraints on inflation, arXiv:1807.06211v2 · Ellis, Vennin & Wands, Review 23: Inflation, Review of Particle Physics (PDG 2025, revised August 2025) · CLASS explanatory.ini, lesgourg/class_public commit c5d93d01990116162b408b789b62ec5669fc4ea4

我们真正遇到过的失败

绝对容差。早期版本用了固定的 $\mathrm{atol}=10^{-30}$。$N_*=46$ 时，在这套单位里各波数的起始 $|\mathcal{R}|$ 约在 $2\times10^{-29}$ 到 $10^{-24}$ 之间，这个容差其实很松。在 $27$ 个波数上，功率最多错了 $0.55\%$，Wronskian 最多偏离 $0.88\%$，积分器却报告成功。$N_*=55$ 的高原误差仍低于 $10^{-4}$。修正办法是让容差随各分量的起始大小缩放。

台阶位置。台阶放在结束前 $58$ 个 e 倍数时，振荡一直延伸到主尺度，改变了报告的谱指数。默认台阶因此移到 $59$ 个 e 倍数，每次修改后都把主尺度谱指数和平滑模型比一比。

$r$ 的主尺度。默认台阶不改变 $0.05\,\mathrm{Mpc}^{-1}$ 处的 $r$，仍为 $0.00355$；但把 $0.002\,\mathrm{Mpc}^{-1}$ 处的值从 $0.00318$ 抬高到 $0.00441$。Lab 因此报告 $0.05\,\mathrm{Mpc}^{-1}$ 处的 $r$。

对应阅读：Original background and mode calculations with an independent Python reference · Planck Collaboration, Planck 2018 results. X. Constraints on inflation, arXiv:1807.06211v2

这次复算不能说明什么

模型只有一个标准标量场，引力是广义相对论，初态是 Bunch–Davies 真空，$\epsilon=1$ 之后没有动力学。它只给出原初谱，没有经过转移函数和似然函数，所以不声称拟合了数据。

Lab 积分器与这个参考程序在五个预设上相差不超过 $2.5\times10^{-4}$，说明两套代码解的是同一组方程。

对应阅读：Original background and mode calculations with an independent Python reference · Ellis, Vennin & Wands, Review 23: Inflation, Review of Particle Physics (PDG 2025, revised August 2025)