自定义势能的数值求解器

上一节我们对着谐振子把有限差分法走通了一遍。但每次换一个势都重新摆弄网格、单位和检查项，太不成体统。这一节做两件事：把方法打包成一个随取随用的求解器；再用它啃一块解析方法啃不动的硬骨头——双势阱——并从数值结果里读出真实的物理：氨分子的微波激射器，以及化学键的雏形。

求解器本体

代码几乎就是上一节的翻版，唯一的新意是把「势」作为参数传进去——任何能对数组求值的函数都行：

import numpy as np
from scipy.linalg import eigh_tridiagonal

def solve(V, L, N, k=6):
    """解 -psi''/2 + V(x) psi = E psi（无量纲单位 hbar = m = 1）。

    V : 势能函数，接受 numpy 数组
    L : 求解区间 [-L/2, L/2] 的宽度（两端视为无限高墙）
    N : 内部格点数
    k : 要求的最低能级个数
    返回：格点 x、能级 E[0..k-1]、波函数 psi[:, i]（已归一化）
    """
    x = np.linspace(-L/2, L/2, N + 2)[1:-1]
    h = x[1] - x[0]
    diag = 1.0 / h**2 + V(x)
    off  = -0.5 / h**2 * np.ones(N - 1)
    E, psi = eigh_tridiagonal(diag, off, select='i', select_range=(0, k - 1))
    return x, E, psi / np.sqrt(h)

先跑一遍已知答案，确认没抄错：

x, E, psi = solve(lambda x: 0.5 * x**2, L=16, N=1000)
print(np.round(E, 4))    # [0.5  1.5  2.4999  3.4998  4.4997  5.4995]

谐振子能级 $n+\frac12$ 如期而至。

**量纲怎么办？**代码用的是 $\hbar=m=1$ 的无量纲单位——这不是偷懒，是数值计算的标准操作：避免让 $10^{-34}$ 量级的数在浮点运算里乱窜。换算规则一句话：选定长度单位 $L_0$，能量单位就自动是 $\varepsilon_0=\hbar^2/(mL_0^2)$。对电子取 $L_0=1$ nm：

$$\varepsilon_0=\frac{(\hbar c)^2}{mc^2\,L_0^2}=\frac{197.3^2}{511000\times1^2}\approx0.0762\ \text{eV}\tag{4.5.1}$$

把真实势按 $\varepsilon_0$ 缩成无量纲数喂进去，算出的 $E$ 再乘回 $\varepsilon_0$，就是 eV。

验收实战：双势阱

现在上主菜。取

$$V(x)=V_0\,(x^2-1)^2\tag{4.5.2}$$

两个极小值在 $x=\pm1$，中间隔着高 $V_0$ 的驼峰。这是解析方法的坟场——没有封闭解——却是物理的富矿：氨分子里氮原子在氢平面两侧的翻转、化学键中共享的电子、超导量子比特的双阱回路，全是它的变奏。

对求解器来说它毫无特殊之处，一行换势：

x, E, psi = solve(lambda x: 20.0 * (x**2 - 1)**2, L=8, N=2000)
print(np.round(E, 3))    # [6.036  6.056  16.23  17.111  23.379  27.933]

盯住这串数字，它在讲一个故事：

• $E_0=6.036$ 与 $E_1=6.056$ 挤成一对，只差 $\Delta E=0.019$；而 $E_1$ 到 $E_2$ 隔着足足 $10.2$。能级不再等距铺开，而是成对出现、对内几乎简并。
• 检查波函数：$\psi_0$ 是偶函数——左右两阱上各一个鼓包，同号相连；$\psi_1$ 是奇函数——同样的两个鼓包，右边反号，中间穿零。清单第 3 条通过：势对称，偶奇交替。
• 谐振子极限也对得上（清单第 5 条）：单个阱底附近 $V''(\pm1)=8V_0=160$，等效频率 $\omega=\sqrt{160}\approx12.6$，零点能约 $\omega/2\approx6.3$——与 $E_0\approx6.0$ 相符，差的那点正是非谐修正。

为什么挤成一对？想象驼峰无限高：左阱和右阱是两个隔绝的世界，各有一个基态，能量完全相同——严格二重简并。驼峰有限时，两阱的波函数尾巴在峰下隧穿相遇，简并被打破，一个能级劈成两个。$\Delta E$ 就叫隧穿劈裂，它直接量度两阱间「串门」的难易。

数值可以当场验证这层解释：把驼峰从 $V_0=20$ 降到 $10$，劈裂从 $0.019$ 跳到 $0.125$——驼峰高度减半，劈裂涨了六倍多。这种不成比例的剧烈响应正是 2.11 节隧穿指数 $\ee^{-2\kappa a}$ 的指纹：劈裂对势垒参数指数敏感。

从一对能级到来回振荡

近简并的一对能级不只是能谱上的装饰——它藏着一台钟。

推导：粒子在两阱间的来回振荡基础~8 min展开

第一步：造出「在左阱」的态。$\psi_0$（偶）和 $\psi_1$（奇）在左阱同号、在右阱反号，于是叠加

$$\psi_L=\frac{\psi_0+\psi_1}{\sqrt2},\qquad \psi_R=\frac{\psi_0-\psi_1}{\sqrt2}\tag{4.5.3}$$

$\psi_L$ 在右阱处两项相消、左阱处相长——它局域在左阱；$\psi_R$ 反之。（系数 $1/\sqrt2$ 保持归一化；两式正交。）注意：$\psi_L$ 是两个不同能量本征态的叠加，不是定态——好戏由此开场。

**第二步：让它演化。**初态取 $\Psi(0)=\psi_L$。按第 2 章的老规矩给每个本征成分挂上相位钟：

$$\Psi(t)=\frac{1}{\sqrt2}\left[\psi_0\,\ee^{-\ii E_0t/\hbar}+\psi_1\,\ee^{-\ii E_1t/\hbar}\right]\tag{4.5.4}$$

**第三步：问它还在不在左边。**投影回 $\psi_L$（利用 $\psi_0,\psi_1$ 正交归一）：

$$\braket{\psi_L}{\Psi(t)}=\frac12\left[\ee^{-\ii E_0t/\hbar}+\ee^{-\ii E_1t/\hbar}\right] =\ee^{-\ii(E_0+E_1)t/2\hbar}\cos\frac{\Delta E\,t}{2\hbar}\tag{4.5.5}$$

第二步提出了公共相位，剩下实实在在的余弦。概率为

$$\boxed{\ P_L(t)=\cos^2\frac{\Delta E\,t}{2\hbar}\ }\qquad P_R(t)=\sin^2\frac{\Delta E\,t}{2\hbar}\tag{4.5.6}$$

读结果。粒子以周期 $T_{\text{osc}}=2\pi\hbar/\Delta E$ 在两阱间整体来回搬家：$t=\pi\hbar/\Delta E$ 时它百分之百在右阱。驱动这台钟的唯一参数是隧穿劈裂 $\Delta E$——劈裂越小（垒越高），搬家越慢。驼峰无限高的极限下 $\Delta E\to0$，周期无穷长：粒子被永远困在一边，回到「两个隔绝世界」的图像，自洽。

把它变成你的工具

这个求解器从此是你的实验台。几个值得亲手跑的势（每个都记得过验收清单）：

• 不对称双阱 $V_0(x^2-1)^2+\lambda x$：加一点点倾斜 $\lambda$，看近简并的一对如何被拆散、波函数如何各自缩回一个阱——4.1 节「势不对称则宇称失效」的数值版。
• 斜底阱 $V=\alpha x$（$x>0$，左端为墙）：三角形阱，半导体异质结界面处电子的真实处境。
• 截断谐振子：把 $x^2/2$ 在某高度削平，数一数还剩几个束缚态，对照 4.1 节的 $z_0$ 计数经验。

写论文级的代码还会再加两样：Numerov 方法把精度从 $O(h^2)$ 提到 $O(h^4)$，虚时间演化处理二维三维——但骨架不变：离散化、对角化、验收。

∎本节关键公式

• 无量纲化

  $$\varepsilon_0=\frac{\hbar^2}{mL_0^2}=\frac{(\hbar c)^2}{mc^2L_0^2}\ \xrightarrow{\text{电子},\,L_0=1\,\text{nm}}\ 0.0762\ \text{eV}$$

  长度单位选定后能量单位自动确定；算完乘回去

• 隧穿劈裂

  $$E_1-E_0=\Delta E\ \propto\ \ee^{-2\int\kappa\,\dd x}$$

  对势垒高度与宽度指数敏感；垒无限高时回到二重简并

• 左右局域态

  $$\psi_{L,R}=\frac{\psi_0\pm\psi_1}{\sqrt2}$$

  定态的叠加、本身非定态；偶奇宇称在单侧相长/相消

• 两阱振荡

  $$P_L(t)=\cos^2\frac{\Delta E\,t}{2\hbar},\qquad T_{\text{osc}}=\frac{2\pi\hbar}{\Delta E}$$

  氨分子 ΔE ≈ 9.8×10⁻⁵ eV → 24 GHz，微波激射器的工作频率

?自测共 4 题

1. 1.

   双势阱最低两个能级挤成近简并的一对（ΔE 远小于到下一能级的间距）。这对「挤在一起」的物理根源是：

   • 势垒无限高时两阱严格简并，有限的隧穿把简并微微劈开
   • 有限差分法的数值误差把一个能级算成了两个
   • 两个能级分属两个不同的阱
   • 宇称守恒要求能级成对

   提交
2. 2.

   关于局域态 ψ_L = (ψ₀ + ψ₁)/√2，正确的说法有：（多选）

   多选题

   • 它不是能量本征态，会随时间演化
   • 以它为初态，概率按 cos²(ΔE·t/2ħ) 留在左阱，周期性地整体搬到右阱
   • 它没有确定宇称
   • 它的平均能量高于 E₁

   提交
3. 3.

   数值解一个新势后，下列哪项「验收」最能揭穿「盒子 L 取小了」的错误？

   • 加大 L 重算，看低能级是否移动
   • 检查本征矢量相互正交
   • 数波函数的节点数
   • 检查能量是实数

   提交
4. 4.

   把本节求解器用于谐振子 V(x) = x²/2（无量纲单位），网格取得足够密。输出的第 5 个能级 E₄ 应该非常接近多少？

   允许 5% 相对误差

   提交

本章结束

盘点一下这一章添置的家当：宇称用一个对易子把束缚态整理成偶奇交替，还免费划掉一批积分；相移把散射浓缩成一个角度，共振与束缚态在 Levinson 定理里握手；传递矩阵把多层结构变成 $2\times2$ 乘法，交出了共振隧穿和能带；有限差分把「算符即矩阵」落成代码，任意一维势从此都是十行程序的事。一维世界，至此可以说尽在掌握。

但真实的原子不是一维的。氢原子里的电子活在三维库仑势中，而三维带来一个一维根本没有的角色：转动。粒子可以绕着核转，转动有快慢、有方向——角动量登场。好消息是你的工具全都用得上：中心势下的三维方程可以分离出一个径向方程，它恰恰活在半轴上，正是 4.2 节那个「墙加势」的舞台。坏消息（其实是最精彩的消息）是角向部分自成天地：三个角动量分量两两不对易，$[\hat L_x,\hat L_y]=\ii\hbar\hat L_z$——第 3 章的对易子代数将不再只是判断「能否同时测准」的工具，而是亲手把角动量的取值一格一格造出来。

第 5 章，从转动开始。