密度矩阵在薛定谔绘景下的演化

上述结果是在一个特殊的旋转框架（相互作用绘景）下得到的。为了理解系统密度矩阵在薛定谔绘景下的演化，我们先来回顾薛定谔绘景和相互作用绘景之间的关系。

薛定谔绘景和相互作用绘景之间的关系：相会作用绘景由一个旋转算符定义

$R=e^{i(H_s+H_b)t}\quad\qquad\quad R^\dagger=e^{-i(H_s+H_b)t}$

根据 $|\phi(t)\rangle=R(t)|\psi(t)\rangle$ 和 $\rho^{\mathrm{Sch}}=|\psi(t)\rangle\langle\psi(t)|$ ，得到

$\rho^{\mathrm{Int}}=|\phi(t)\rangle\langle\phi(t)|=R\rho^{\mathrm{Sch}}R^{\dagger}$

系统的约化密度矩阵 $\rho_s^{\text{Sch}}(t)$ ，是通过对总密度矩阵 $\rho^{\text{Sch}}(t)$ 中热库自由度进行偏迹运算而得到的

$\begin{aligned} \rho_{s}^{\mathrm{Sch}}(t) & =\mathrm{Tr}_{b}\left(\rho^{\mathrm{Sch}}(t)\right) =\mathrm{Tr}_{b}\left[R^{\dagger}(t)\rho^{\mathrm{Int}}(t)R(t)\right] \\[1em] & =\mathrm{Tr}_{b}\left[e^{-i(H_{s}+H_{b})t}\rho^{\mathrm{Int}}(t)e^{i(H_{s}+H_{b})t}\right] \\[1em] & =\mathrm{Tr}_{b}\left[e^{-iH_{s}t}e^{-iH_{b}t}\rho^{\mathrm{Int}}(t)e^{iH_{b}t}e^{iH_{s}t}\right] \\[1em] & =e^{-iH_{s}t}\mathrm{Tr}_{b}\left[e^{-iH_{b}t}\rho^{\mathrm{Int}}(t)e^{iH_{b}t}\right]e^{iH_{s}t} \\[1em] & =e^{-iH_{s}t}\rho_{s}^{\mathrm{Int}}(t)e^{iH_{s}t} \end{aligned}$

这表明，系统密度矩阵在薛定谔绘景与相互作用绘景之间，是通过一个仅涉及系统自由哈密顿量 $H_s$ 的幺正变换联系起来的。

对上式取时间导数

$\begin{aligned} \frac{d}{dt}\rho_{s}^{\mathrm{Sch}} & =\frac{d}{dt}e^{-iH_{s}t}\rho_{s}^{\mathrm{Int}}e^{iH_{s}t} \\[1em] & =-iH_{s}e^{-iH_{s}t}\rho_{s}^{\mathrm{Int}}e^{iH_{s}t}+e^{-iH_{s}t}\frac{d}{dt}(\rho_{s}^{\mathrm{Int}})e^{iH_{s}t}+e^{-iH_{s}t}\rho_{s}^{\mathrm{Int}}e^{iH_{s}t}(iH_{s}) \\[1em] & =-iH_{s}\rho_{s}^{\mathrm{Sch}}+e^{-iH_{s}t}\frac{d}{dt}(\rho_{s}^{\mathrm{Int}})e^{iH_{s}t}+\rho_{s}^{\mathrm{Sch}}iH_{s} \\[1em] & =-i[H_{s},\rho_{s}^{\mathrm{Sch}}]+e^{-iH_{s}t}\frac{d}{dt}(\rho_{s}^{\mathrm{Int}})e^{iH_{s}t} \end{aligned}$

我们已经在相互作用绘景中建立了 $\rho^{\text{Int}}_s (t)$ 的运动方程，该方程包含了源于热库相互作用的耗散项和噪声项

$\frac{d}{dt}\rho_s^{\text{Int}}(t) = \frac{\gamma}{2}N \cdot \left( 2a^{\dagger\text{Int}}\rho_s^{\text{Int}}(t)a^{\text{Int}} - a^{\dagger\text{Int}}a^{\text{Int}}\rho_s^{\text{Int}}(t) - \rho_s^{\text{Int}}(t)a^{\dagger\text{Int}}a^{\text{Int}} \right) + \dots$

将上式代入，我们可以得到薛定谔绘景下的主方程

$\begin{align*} e^{-iH_s t} \frac{d}{dt}\rho_s^{\text{Int}} e^{iH_s t} &= \cdots e^{-iH_s t} \left( 2a^{\dagger\text{Int}}\rho_s^{\text{Int}}a^{\text{Int}} + \cdots \right) e^{iH_s t} \\[1em] &= \cdots \left( 2e^{-iH_s t}a^{\dagger\text{Int}}\textcolor{blue}{e^{iH_s t} e^{-iH_s t}}\rho_s^{\text{Int}}\textcolor{blue}{e^{iH_s t} e^{-iH_s t}}a^{\text{Int}}e^{iH_s t} + \cdots \right. \\[1em] &= \cdots \left( 2\left(e^{-iH_s t}a^{\dagger\text{Int}}e^{iH_s t}\right) \left(e^{-iH_s t}\rho_s^{\text{Int}}e^{iH_s t}\right) \left(e^{-iH_s t}a^{\text{Int}}e^{iH_s t}\right) + \cdots \right. \\[1em] &= \frac{\gamma}{2}N \cdot \underbrace{\left( 2a^{\dagger}\rho_s a - a^{\dagger}a\rho_s - \rho_s a^{\dagger}a + \cdots \right.}_{\text{Sch picture}} \end{align*}$

在对热库自由度求迹并应用马尔可夫近似之后，我们得到了薛定谔绘景下系统密度矩阵的最终的、时间局域的运动方程，这就是著名阻尼谐振子的林德布拉德 $\text{Lindblad}$ 主方程：

$\boxed{ \begin{aligned} \frac{d}{dt}\rho_s &= -i[H_s, \rho_s] \\[1em] &\quad + \frac{\gamma}{2}N(\omega_s) \left( 2a^{\dagger}\rho_s a - aa^{\dagger}\rho_s - \rho_s aa^{\dagger} \right) \\[1em] &\quad + \frac{\gamma}{2}[N(\omega_s) + 1] \left( 2a\rho_s a^{\dagger} - a^{\dagger}a\rho_s - \rho_s a^{\dagger}a \right) \end{aligned} }$

林德布拉德主方程可以写成适用于多个耗散通道的通用形式：

$\frac{d\rho}{dt}=-i[H,\rho]+\sum_m\frac{\gamma_m}{2}\left(2L_m\rho L_m^\dagger-L_m^\dagger L_m\rho-\rho L_m^\dagger L_m\right)$

其中 $L_m$ 被称为林德布拉德算符（或“跳跃算符”），用于描述不同的衰减过程，而 $\gamma_m$ 则是它们对应的速率。

林德布拉德方程的广义形式与有效性：

该推导假设了 $H_s$ 包含一个自由部分 $H_{s0} = \omega_s a^\dagger a$ 。但该结果对于更一般的系统哈密顿量 $H_s = H_{s0} + H_{s1}$ （例如：一个由激光驱动的原子）同样有效。为了验证这一点，可以使用旋转参考系变换 $R = e^{i(H_{s0} + H_b)t}$ 来重复该推导过程。

该主方程的推导依赖于几个关键假设，这些假设界定了它的适用范围。

马尔可夫近似

相比于系统的演化时间尺度 ( $\gamma^{-1}$ )，热库的关联时间 $\tau_b$ 非常短。这意味着热库几乎瞬间“遗忘”了系统传递给它的任何信息，使得系统的未来状态仅取决于其当前状态，而与其过去的历史无关。这一点隐含在取 $\kappa(\omega) = \kappa$ 的近似中。

何时失效： 如果系统与热库的耦合非常强，热库的关联时间可能就变得不可忽略。在这种情况下，必须采用非马尔可夫方法，例如布朗运动主方程。

热库关联性质

我们假设了 $\langle b^\dagger b^\dagger \rangle = 0$ 且 $\langle bb \rangle = 0$ 。对于处于热平衡或真空平衡态的热库，这通常是成立的。

何时失效： 对于压缩热库，这些反常关联是非零的 ( $\langle bb \rangle \neq 0$ )。这会导致主方程中出现额外的、不同的项，并且这是量子光学中用于产生特殊量子态的一个活跃研究领域。

至此，我们完成了林德布拉德主方程的推导，这是模拟耗散量子系统的核心工具。接下来，我们将利用它来求解特定量子光学系统的动力学。

腔光子模式的衰减

我们已经反复断言，将腔场耦合到热库会导致其衰减。让我们针对低温和高温热库两种情形，明确地验证这一说法。

低温下的腔场损失

假设腔在初始时刻 $t=0$ 包含 $n$ 个光子，由密度矩阵描述

$\rho(0)=|n\rangle\langle n|$

如果热库是冻结的（实际上处于零温， $N=0$ ），林德布拉德主方程简化为：

$\begin{aligned} \frac{d}{dt}\rho_s &= \frac{\gamma}{2}N \cdot (2a^\dagger \rho_s a - a a^\dagger \rho_s - \rho_s a a^\dagger) + \frac{\gamma}{2}(N+1) \cdot (2a\rho_s a^\dagger - a^\dagger a \rho_s - \rho_s a^\dagger a) \\[1em] &= \frac{\gamma}{2} \cdot (2a\rho_s a^\dagger - a^\dagger a \rho_s - \rho_s a^\dagger a) \end{aligned}$

其中 $\gamma$ 是腔场的衰减速率。我们知道

$\begin{gather*} a|n\rangle\langle n|a^\dagger = \sqrt{n}|n-1\rangle\langle n-1|\sqrt{n} = n|n-1\rangle\langle n-1| \\[1em] a^\dagger a|n\rangle\langle n| = n|n\rangle\langle n| \\[1em] |n\rangle\langle n|a^\dagger a = |n\rangle\langle n|n \end{gather*}$

将这些结果代入主方程，得到微分增量

$\frac{\rho_s(\Delta t) - \rho_s(0)}{\Delta t} = \frac{\gamma}{2} \cdot 2n|n-1\rangle\langle n-1| - \frac{\gamma}{2} \cdot 2n|n\rangle\langle n|$

在微小时间间隔 $\Delta t$ 后，腔场密度矩阵为

$\rho_s(\Delta t) = n\gamma\Delta t|n-1\rangle\langle n-1| + (1 - n\gamma\Delta t)|n\rangle\langle n|$

第一项代表光子以正比于 $\gamma n$ 的速率流失到热库，而第二项描述了初始 $|n\rangle$ 态的存活概率。

初始时刻 $t=0$ 时，腔场处于纯态 $|n\rangle$ 。经过时间 $\Delta t$ 后，腔场变成了混合态，包含 $|n\rangle$ 态和 $|n-1\rangle$ 态的概率分布。这反映了光子从腔场流失到热库的过程。

高温下的腔场损失

现在，考虑一个处于有限温度且具有非零热光子数 $N > 0$ 的热库，这时必须使用完整的主方程

$\frac{d}{dt}\rho_s=\frac{\gamma}{2}N\cdot(2a^\dagger\rho_sa-aa^\dagger\rho_s-\rho_saa^\dagger)+\frac{\gamma}{2}(N+1)\cdot(2a\rho_sa^\dagger-a^\dagger a\rho_s-\rho_sa^\dagger a)$

我们再次计算初始态 $\rho_s(0) = |n\rangle\langle n|$ 的变化。我们需要这些新项的作用

$\begin{aligned} a^\dagger \rho_s a &= a^\dagger |n\rangle\langle n|a = (n+1)|n+1\rangle\langle n+1| \\[1em] a a^\dagger \rho_s &= (n+1)|n\rangle\langle n| \\[1em] \rho_s a a^\dagger &= (n+1)|n\rangle\langle n| \end{aligned}$

因此，经过短时间 $\Delta t$ 后的密度矩阵为

$\begin{aligned} \rho_s(\Delta t) &= N(n+1)\gamma\Delta t|n+1\rangle\langle n+1| + (N+1)n\gamma\Delta t|n-1\rangle\langle n-1| \\[1em] &\quad+ (1 - N(n+1)\gamma\Delta t - (N+1)n\gamma\Delta t)|n\rangle\langle n| \end{aligned}$

这一结果揭示了丰富且相互竞争的动力学

光子吸收：与 $N(n + 1)$ 成正比的项代表从热库吸收一个光子，从而将腔光子数增加到 $n + 1$ 。该速率正比于热光子数 $N$ 以及末态的“可用性”因子 $(n + 1)$ 。
光子发射：与 $(N + 1)n$ 成正比的项代表发射，它将光子数减少到 $n-1$ 。该速率包含一个自发部分 ( $1\cdot n$ ) 和一个受激部分 ( $N\cdot n$ )。
平衡： 系统将演化趋向热平衡，在此状态下吸收速率与发射速率相等。稳态将是一个热态，其平均光子数由热库温度决定。

这展示了主方程的强大之处：不仅能描述衰减，还能描述热化以及自发与受激过程之间精妙的平衡。

这还导致了依赖于温度的趋向平衡过程：

对于 $N > n$ ，即热库比腔更热，跃迁至 $|n+1\rangle \langle n+1|$ 的概率更高
对于 $N < n$ ，即热库比腔更冷，跃迁至 $|n-1\rangle \langle n-1|$ 的概率更高

原子的自发辐射

我们已经成功模拟了腔衰减，现在将林德布拉德主方程应用于量子光学中最基础的过程之一：二能级原子的自发辐射。这将会体现其普适性，并提供对光与物质相互作用的深入洞察。

$\text{Lindblad}$ 主方程是一个固定的结构

$\frac{d\rho}{dt}=-i[H,\rho]+\sum_m\frac{\gamma_m}{2}\left(2L_m\rho L_m^\dagger-L_m^\dagger L_m\rho-\rho L_m^\dagger L_m\right)$

这里的 $L_m$ 是跳跃算符

如果系统是光，衰减就是光子消失，所以 $L=a$
如果系统是原子，衰减就是原子从激发态掉下来，所以 $L=\sigma_{12}$ (即 $\sigma_-$ )。

考虑一个二能级原子模型，其基态为 $|1\rangle$ ，激发态为 $|2\rangle$ ，并将其耦合到电磁场模式的连续谱（热库），原子的林德布拉德主方程为

$\frac{d}{dt}\rho_{s}=\frac{\gamma}{2}N\cdot(2\sigma_{21}\rho_{s}\sigma_{12}-\sigma_{12}\sigma_{21}\rho_{s}-\rho_{s}\sigma_{12}\sigma_{21})+\frac{\gamma}{2}(N+1)\cdot(2\sigma_{12}\rho_{s}\sigma_{21}-\sigma_{21}\sigma_{12}\rho_{s}-\rho_{s}\sigma_{21}\sigma_{12})$

为了求解该方程，我们在原子基底 $|1\rangle, |2\rangle$ 下将其表示出来，密度矩阵和相关算符为

$\rho_s=\begin{pmatrix}\rho_{11}&\rho_{12}\\[1em]\rho_{21}&\rho_{22}\end{pmatrix},\quad\sigma_{12}=\begin{pmatrix}0&1\\[1em]0&0\end{pmatrix}=|1\rangle\langle2|,\quad\sigma_{21}=\begin{pmatrix}0&0\\[1em]1&0\end{pmatrix}=|2\rangle\langle1|$

接下来以矩阵方式计算每一项林德布拉德项对 $\rho_s$ 的作用。经过简单的计算，我们得到

$\begin{aligned} \frac{d}{dt}\begin{bmatrix} \rho_{11} & \rho_{12} \\ \rho_{21} & \rho_{22} \end{bmatrix} &= \begin{bmatrix} -\gamma N \cdot \rho_{11} & -\frac{\gamma}{2}N \cdot \rho_{12} \\[1em] -\frac{\gamma}{2}N \cdot \rho_{21} & \gamma N \cdot \rho_{11} \end{bmatrix} + \begin{bmatrix} \frac{\gamma}{2}(N+1) \cdot \rho_{22} & -\frac{\gamma}{2}(N+1) \cdot \rho_{12} \\[1em] -\frac{\gamma}{2}(N+1) \cdot \rho_{21} & -\gamma(N+1) \cdot \rho_{22} \end{bmatrix} \\[2em] &= \begin{bmatrix} -\gamma N \cdot \rho_{11} + \gamma(N+1) \cdot \rho_{22} & -\frac{\gamma}{2}(2N+1) \cdot \rho_{12} \\[1em] -\frac{\gamma}{2}(2N+1) \cdot \rho_{21} & \gamma N \cdot \rho_{11} - \gamma(N+1) \cdot \rho_{22} \end{bmatrix} \end{aligned}$

令 $d\rho_s/dt = 0$ 求解稳态，给出了两个独立的方程

$\begin{aligned} (\dot{\rho}_{11} = 0): &\quad -\gamma N\rho_{11} + \gamma(N+1)\rho_{22} = 0 \\[1em] (\dot{\rho}_{22} = 0): &\quad \gamma N\rho_{11} - \gamma(N+1)\rho_{22} = 0 \end{aligned}$

这两式都得到了相同的关系

$\frac{\rho_{22}}{\rho_{11}}=\frac{N}{N+1}$

该结果表明，原子激发态的相对布居数完全由热库在跃迁频率处的热光子数 $N$ 决定。在高温下（ $N \to \infty$ ）， $\rho_{22}/\rho_{11} \to 1$ ，这与经典均分定理相符。在 $T \to 0$ 时， $N \to 0$ ， $\rho_{22} \to 0$ ，激发态为空。

对于处于热平衡且温度为 $T$ 的热库，光子数服从玻色-爱因斯坦分布

$N(\omega_2)=\frac{1}{e^{\hbar\omega_2/k_BT}-1}$

将其代入我们的稳态比率中

$\frac{\rho_{22}}{\rho_{11}}=\frac{\dfrac{1}{e^{\hbar\omega_{2}/k_{B}T}-1}}{\dfrac{1}{e^{\hbar\omega_{2}/k_{B}T}-1}+1}=\dfrac{1}{e^{\hbar\omega_{2}/k_{B}T}}=e^{-\hbar\omega_{2}/k_{B}T}$