Introduction
An atom coupled to the electromagnetic field is an open quantum system: its state becomes mixed as it interacts with the environment. The density operator provides the appropriate formalism for describing such open systems, and the master equation governs their dynamics. These tools are essential for understanding spontaneous emission, decoherence, and the collective behavior of driven, dissipative atomic systems.
1. The Density Operator and Open Quantum Systems
1.1 Pure States, Mixed States, and the Density Matrix
A quantum system is said to be in a pure state if its state can be described by a single ket $| \psi \rangle$ in the Hilbert space $\mathcal{H}$. The expectation value of any observable $\hat{O}$ is $\langle \hat{O} \rangle = \langle \psi | \hat{O} | \psi \rangle$. Equivalently, we can define the density operator \begin{equation} \hat{\rho} = | \psi \rangle \langle \psi |, \end{equation} which satisfies $\hat{\rho}^2 = \hat{\rho}$ (idempotence) and $\operatorname{Tr}(\hat{\rho}) = 1$. Expectation values are computed as $\langle \hat{O} \rangle = \operatorname{Tr}(\hat{O}\hat{\rho})$.
A mixed state arises when the system is in a statistical ensemble of pure states $\{ | \psi_i \rangle \}$ with probabilities $p_i$ ($\sum_i p_i = 1$, $p_i \ge 0$). The density operator is then \begin{equation} \hat{\rho} = \sum_i p_i \, | \psi_i \rangle \langle \psi_i |. \label{eq:mixed_state} \end{equation} For a mixed state, $\hat{\rho}^2 \neq \hat{\rho}$ in general; the "purity" is quantified by $\operatorname{Tr}(\hat{\rho}^2) \le 1$, with equality holding only for pure states. The density operator is Hermitian ($\hat{\rho}^\dagger = \hat{\rho}$), positive semi-definite ($\langle \phi | \hat{\rho} | \phi \rangle \ge 0$ for all $| \phi \rangle$), and has unit trace.
1.2 Liouville–von Neumann Equation
For an isolated system governed by a Hamiltonian $\hat{H}$, each pure state in the ensemble evolves according to the Schrödinger equation $i\hbar\,\partial_t | \psi_i \rangle = \hat{H} | \psi_i \rangle$. The density operator then evolves as \begin{equation} \frac{\mathrm{d}\hat{\rho}}{\mathrm{d} t} = \sum_i p_i \left( \frac{\mathrm{d} | \psi_i \rangle}{\mathrm{d} t} \langle \psi_i | + | \psi_i \rangle \frac{\mathrm{d} \langle \psi_i |}{\mathrm{d} t} \right) = -\frac{i}{\hbar} [\hat{H}, \hat{\rho}]. \label{eq:Liouville} \end{equation} This is the Liouville–von Neumann equation. It is the quantum analog of the classical Liouville equation for the phase-space probability density. The commutator structure ensures that the eigenvalues of $\hat{\rho}$ (the probabilities $p_i$) are constant in time—unitary evolution preserves purity.
1.3 Open Quantum Systems and the Concept of a Reservoir
An open quantum system is a system $S$ of interest that interacts with an external environment (or reservoir, or bath) $R$. The total Hamiltonian is \begin{equation} \hat{H} = \hat{H}_S + \hat{H}_R + \hat{H}_{SR}, \end{equation} where $\hat{H}_S$ acts only on $S$, $\hat{H}_R$ acts only on $R$, and $\hat{H}_{SR}$ couples them. The total system $S+R$ is assumed to be closed and evolves unitarily according to (\ref{eq:Liouville}) with the full Hamiltonian. However, we are typically interested only in observables of the subsystem $S$. The reduced state of $S$ is obtained by tracing over the reservoir degrees of freedom: \begin{equation} \hat{\rho}_S(t) = \operatorname{Tr}_R \big[ \hat{\rho}_{SR}(t) \big]. \end{equation} The central goal of the theory of open quantum systems is to derive a closed, autonomous equation of motion for $\hat{\rho}_S(t)$—a master equation—by eliminating the reservoir degrees of freedom under suitable approximations. The resulting equation is, in general, non-unitary: it describes dissipation, decoherence, and the approach to thermal equilibrium.
2. The Born–Markov Master Equation
2.1 Interaction Picture and Formal Perturbation Theory
We work in the interaction picture with respect to $\hat{H}_0 = \hat{H}_S + \hat{H}_R$. The total density operator in the interaction picture is \begin{equation} \tilde{\rho}_{SR}(t) = e^{i\hat{H}_0 t/\hbar} \, \hat{\rho}_{SR}(t) \, e^{-i\hat{H}_0 t/\hbar}, \end{equation} and it satisfies \begin{equation} \frac{\mathrm{d}}{\mathrm{d} t} \tilde{\rho}_{SR}(t) = -\frac{i}{\hbar} [\tilde{H}_{SR}(t), \tilde{\rho}_{SR}(t)], \label{eq:IP_Liouville} \end{equation} where $\tilde{H}_{SR}(t) = e^{i\hat{H}_0 t/\hbar} \hat{H}_{SR} e^{-i\hat{H}_0 t/\hbar}$. Integrating (\ref{eq:IP_Liouville}) formally and substituting back into itself yields \begin{equation} \frac{\mathrm{d}}{\mathrm{d} t} \tilde{\rho}_{SR}(t) = -\frac{i}{\hbar} [\tilde{H}_{SR}(t), \tilde{\rho}_{SR}(0)] - \frac{1}{\hbar^2} \int_0^t \mathrm{d}\tau \, [\tilde{H}_{SR}(t), [\tilde{H}_{SR}(\tau), \tilde{\rho}_{SR}(\tau)]]. \label{eq:formal_expansion} \end{equation} Tracing over the reservoir, \begin{equation} \frac{\mathrm{d}}{\mathrm{d} t} \tilde{\rho}_S(t) = -\frac{i}{\hbar} \operatorname{Tr}_R \big\{ [\tilde{H}_{SR}(t), \tilde{\rho}_{SR}(0)] \big\} - \frac{1}{\hbar^2} \int_0^t \mathrm{d}\tau \, \operatorname{Tr}_R \big\{ [\tilde{H}_{SR}(t), [\tilde{H}_{SR}(\tau), \tilde{\rho}_{SR}(\tau)]] \big\}. \label{eq:exact_master} \end{equation} This equation is exact but not closed: the right-hand side depends on the full $\tilde{\rho}_{SR}(\tau)$, not just $\tilde{\rho}_S(\tau)$.
2.2 The Born Approximation
The Born approximation assumes that the coupling between system and reservoir is weak, and that the reservoir is sufficiently large that it is essentially unaffected by the system. We therefore write the total state in the factorized form \begin{equation} \tilde{\rho}_{SR}(t) \approx \tilde{\rho}_S(t) \otimes \hat{\rho}_R, \label{eq:Born} \end{equation} where $\hat{\rho}_R$ is the stationary state of the reservoir (e.g., the thermal equilibrium state or the vacuum state), which we take to be time-independent in the interaction picture: $[\hat{H}_R, \hat{\rho}_R] = 0$. We also assume that the first term in (\ref{eq:exact_master}) vanishes, which is guaranteed if $\operatorname{Tr}_R\{\tilde{H}_{SR}(t)\hat{\rho}_R\} = 0$ (this can always be arranged by absorbing any mean-field shift into $\hat{H}_S$). Under the Born approximation, the master equation becomes \begin{equation} \frac{\mathrm{d}}{\mathrm{d} t} \tilde{\rho}_S(t) = -\frac{1}{\hbar^2} \int_0^t \mathrm{d}\tau \, \operatorname{Tr}_R \big\{ [\tilde{H}_{SR}(t), [\tilde{H}_{SR}(\tau), \tilde{\rho}_S(\tau) \otimes \hat{\rho}_R]] \big\}. \label{eq:Born_master} \end{equation} This is a closed integro-differential equation for $\tilde{\rho}_S(t)$, but it is non-Markovian: the derivative at time $t$ depends on the history $\tilde{\rho}_S(\tau)$ for $0 \le \tau \le t$.
2.3 The Markov Approximation
The Markov approximation exploits the separation of timescales between the system and the reservoir. The reservoir correlation time $\tau_R$—the timescale over which the reservoir two-time correlation functions decay—is assumed to be much shorter than the characteristic timescale $\tau_S$ over which $\tilde{\rho}_S$ evolves appreciably ($\tau_S \sim 1/\Gamma$, where $\Gamma$ is the dissipation rate). In this regime, we can:
- Replace $\tilde{\rho}_S(\tau)$ by $\tilde{\rho}_S(t)$ inside the integral, since the memory kernel is sharply peaked around $\tau = t$.
- Extend the upper limit of the integration to $\infty$ (since the kernel vanishes for $\tau \gg \tau_R$).
- Change variable to $s = t - \tau$: \begin{equation} \int_0^t \mathrm{d}\tau \, \tilde{\rho}_S(\tau) \, K(t,\tau) \approx \int_0^\infty \mathrm{d}s \, \tilde{\rho}_S(t) \, K(t, t-s). \end{equation}
These steps constitute the Born–Markov approximation. The resulting master equation is local in time (Markovian) and is first-order in the system–reservoir coupling. The general form of the Born–Markov master equation in the Schrödinger picture is \begin{equation} \frac{\mathrm{d}}{\mathrm{d} t} \hat{\rho}_S(t) = -\frac{i}{\hbar} [\hat{H}_S, \hat{\rho}_S(t)] + \mathcal{L}[\hat{\rho}_S(t)], \label{eq:Born_Markov_general} \end{equation} where $\mathcal{L}$ is a superoperator (the Lindbladian) that encodes the dissipative dynamics.
2.4 The Lindblad Form
The most general Markovian master equation that preserves the trace, Hermiticity, and complete positivity of the density operator is the Gorini–Kossakowski–Sudarshan–Lindblad (GKSL) equation: \begin{equation} \boxed{\frac{\mathrm{d}}{\mathrm{d} t} \hat{\rho}_S = -\frac{i}{\hbar} [\hat{H}_{\text{eff}}, \hat{\rho}_S] + \sum_k \gamma_k \left( \hat{L}_k \hat{\rho}_S \hat{L}_k^\dagger - \frac{1}{2} \{ \hat{L}_k^\dagger \hat{L}_k, \hat{\rho}_S \} \right)}. \label{eq:Lindblad} \end{equation} Here, $\hat{H}_{\text{eff}}$ is the effective Hamiltonian of the system, which may include Lamb-shift terms induced by the reservoir. The $\hat{L}_k$ are Lindblad operators (or jump operators) that describe the various dissipative channels, and $\gamma_k \ge 0$ are the corresponding decay rates. The curly brackets denote the anti-commutator: $\{\hat{A},\hat{B}\} = \hat{A}\hat{B} + \hat{B}\hat{A}$.
The Lindblad form is the cornerstone of quantum optics and open quantum systems theory. It guarantees that $\hat{\rho}_S(t)$ remains a valid density operator for all times. The terms $\hat{L}_k \hat{\rho}_S \hat{L}_k^\dagger$ represent "quantum jumps"—the probabilistic, discontinuous change of the system state upon emission or absorption of a quantum of energy—while the anti-commutator term $-\frac{1}{2}\{\hat{L}_k^\dagger\hat{L}_k,\hat{\rho}_S\}$ ensures probability conservation and gives rise to the damping of coherences.
3. Master Equation for a Two-Level Atom
3.1 Derivation from the Wigner–Weisskopf Model
We now specialize the Born–Markov formalism to the case of a two-level atom coupled to the electromagnetic vacuum. This was the system studied in the Wigner–Weisskopf theory. The system Hamiltonian is \begin{equation} \hat{H}_S = \frac{\hbar\omega_0}{2} \hat{\sigma}_z \quad \text{(with } E_g = -\hbar\omega_0/2, \; E_e = +\hbar\omega_0/2\text{)}. \end{equation} The interaction Hamiltonian in the rotating-wave approximation is \begin{equation} \hat{H}_{SR} = \hbar \sum_{\mathbf{k},\lambda} \left( g_{\mathbf{k}\lambda}^* \hat{\sigma}_+ \hat{a}_{\mathbf{k}\lambda} + g_{\mathbf{k}\lambda} \hat{\sigma}_- \hat{a}_{\mathbf{k}\lambda}^\dagger \right), \end{equation} where we have redefined the coupling to absorb factors of $i$ for convenience.
The reservoir is the multimode electromagnetic field in the vacuum state $\hat{\rho}_R = | \{0\} \rangle \langle \{0\} |$. The reservoir correlation functions are \begin{align} \langle \hat{a}_{\mathbf{k}\lambda}(t) \hat{a}_{\mathbf{k}'\lambda'}^\dagger(t') \rangle_R &= \delta_{\mathbf{k}\mathbf{k}'} \delta_{\lambda\lambda'} e^{-i\omega_k (t-t')}, \\ \langle \hat{a}_{\mathbf{k}\lambda}^\dagger(t) \hat{a}_{\mathbf{k}'\lambda'}(t') \rangle_R &= 0 \quad (\text{vacuum}), \end{align} where $\langle \cdot \rangle_R = \operatorname{Tr}_R(\cdot \hat{\rho}_R)$. Substituting $\hat{H}_{SR}$ into the Born–Markov formula (\ref{eq:Born_master}) and evaluating the trace over the reservoir (using the cyclic property of the trace and the bosonic commutation relations) leads, after considerable algebra, to the Lindblad master equation: \begin{equation} \frac{\mathrm{d}}{\mathrm{d} t} \hat{\rho}_S = -\frac{i}{\hbar} [\hat{H}_S', \hat{\rho}_S] + \Gamma \left( \hat{\sigma}_- \hat{\rho}_S \hat{\sigma}_+ - \frac{1}{2} \{ \hat{\sigma}_+ \hat{\sigma}_-, \hat{\rho}_S \} \right). \label{eq:2level_Lindblad} \end{equation} Here, $\Gamma = \omega_0^3 |\mathbf{d}_{eg}|^2 / (3\pi\varepsilon_0\hbar c^3)$ is the spontaneous emission rate (the Einstein $A$ coefficient), $\hat{H}_S' = \hat{H}_S + \hbar\Delta\,\hat{\sigma}_+\hat{\sigma}_-$ includes the Lamb shift $\Delta$, and the single Lindblad operator is $\hat{L} = \sqrt{\Gamma}\,\hat{\sigma}_-$.
The Lindblad operator $\hat{\sigma}_- = | g \rangle \langle e |$ describes the quantum jump from the excited state to the ground state, accompanied by the emission of a photon. The term $\hat{\sigma}_- \hat{\rho}_S \hat{\sigma}_+$ transfers population from $| e \rangle$ to $| g \rangle$, while the anti-commutator term causes the decay of the off-diagonal elements (coherences).
3.2 Recovery of the Optical Bloch Equations
Writing $\hat{\rho}_S$ in the basis $\{ | e \rangle, | g \rangle \}$ as \begin{equation} \hat{\rho}_S = \begin{pmatrix} \rho_{ee} & \rho_{eg} \\ \rho_{ge} & \rho_{gg} \end{pmatrix}, \end{equation} and inserting into (\ref{eq:2level_Lindblad}) (with $\hat{H}_S' \approx \hat{H}_S$ for simplicity, neglecting the small Lamb shift), we obtain the equations of motion for the matrix elements: \begin{align} \dot{\rho}_{ee} &= -\Gamma \rho_{ee}, \label{eq:pop_e} \\ \dot{\rho}_{gg} &= +\Gamma \rho_{ee}, \label{eq:pop_g} \\ \dot{\rho}_{eg} &= -i\omega_0 \rho_{eg} - \frac{\Gamma}{2} \rho_{eg}, \label{eq:coh_eg} \\ \dot{\rho}_{ge} &= +i\omega_0 \rho_{ge} - \frac{\Gamma}{2} \rho_{ge}. \label{eq:coh_ge} \end{align} These are the optical Bloch equations for a two-level atom in the vacuum. The populations decay exponentially with rate $\Gamma$, while the coherences decay at half the rate, $\Gamma/2$. This factor-of-two relationship between the population decay rate ($T_1 = 1/\Gamma$) and the coherence decay rate ($T_2 = 2/\Gamma$ for pure radiative damping) is a hallmark of the Lindblad description.
[Figure 1: (a) A two-level system coupled to the electromagnetic vacuum reservoir. (b) The quantum jump process: the system makes a transition from $| e \rangle$ to $| g \rangle$ with rate $\Gamma$, emitting a photon. (c) Exponential decay of the excited state population $\rho_{ee}(t)$ and the coherence $|\rho_{eg}(t)|$.]
4. Three-Level Atomic System: Optical Bloch Equations
4.1 The Three-Level Model and Its Configurations
We now extend the formalism to a three-level atom. Depending on the relative energies of the three states and the optical selection rules, three-level systems are classified into three canonical configurations:
- Lambda ($\Lambda$) configuration: Two lower states $|1\rangle$ and $|2\rangle$ (which may be ground-state hyperfine or Zeeman sublevels) are optically coupled to a common excited state $|3\rangle$. Dipole transitions between $|1\rangle$ and $|2\rangle$ are forbidden (e.g., by parity or angular momentum selection rules).
- V configuration: A single lower state $|1\rangle$ is coupled to two upper states $|2\rangle$ and $|3\rangle$.
- Ladder (cascade) configuration: The states are arranged in energy order: $|1\rangle \to |2\rangle \to |3\rangle$, with allowed dipole transitions between adjacent levels.
We focus on the $\Lambda$-configuration, which is of paramount importance for EIT, CPT, and quantum memory applications. The three states have energies $E_1, E_2, E_3$, with $E_3 > E_1, E_2$. Two classical, monochromatic laser fields drive the transitions:
- A probe field of frequency $\omega_p$ couples $|1\rangle \leftrightarrow |3\rangle$ with Rabi frequency $\Omega_p$.
- A control (coupling) field of frequency $\omega_c$ couples $|2\rangle \leftrightarrow |3\rangle$ with Rabi frequency $\Omega_c$.
The transition $|1\rangle \leftrightarrow |2\rangle$ is dipole-forbidden. The atomic level scheme is illustrated in Figure 2.
[Figure 2: Three-level atom in the $\Lambda$-configuration. Two ground states $|1\rangle$ and $|2\rangle$ are coupled to an excited state $|3\rangle$ by a probe field $\Omega_p$ and a control field $\Omega_c$, respectively. The detunings are $\Delta_p = \omega_p - \omega_{31}$ and $\Delta_c = \omega_c - \omega_{32}$.]
4.2 Hamiltonian in the Rotating-Wave Approximation
The free atomic Hamiltonian (with the zero of energy set at state $|1\rangle$ for convenience) is \begin{equation} \hat{H}_0 = \hbar\omega_{21} |2\rangle\langle 2| + \hbar\omega_{31} |3\rangle\langle 3|, \end{equation} where $\hbar\omega_{ij} = E_i - E_j$. The classical electric fields are \begin{equation} \mathbf{E}_p(t) = \frac{1}{2} \mathbf{E}_{0p} (e^{-i\omega_p t} + e^{i\omega_p t}), \qquad \mathbf{E}_c(t) = \frac{1}{2} \mathbf{E}_{0c} (e^{-i\omega_c t} + e^{i\omega_c t}). \end{equation} In the electric dipole approximation, the interaction Hamiltonian is $\hat{H}_{\text{int}} = -\hat{\mathbf{d}} \cdot (\mathbf{E}_p + \mathbf{E}_c)$. The dipole operator has the matrix elements \begin{equation} \mathbf{d}_{31} = \langle 3 | \hat{\mathbf{d}} | 1 \rangle, \qquad \mathbf{d}_{32} = \langle 3 | \hat{\mathbf{d}} | 2 \rangle, \end{equation} with $\mathbf{d}_{12} = 0$. We define the (generally complex) Rabi frequencies \begin{equation} \Omega_p = -\frac{\mathbf{d}_{31} \cdot \mathbf{E}_{0p}}{\hbar}, \qquad \Omega_c = -\frac{\mathbf{d}_{32} \cdot \mathbf{E}_{0c}}{\hbar}. \label{eq:Rabi_def} \end{equation} Applying the rotating-wave approximation (dropping terms oscillating at $\pm 2\omega_p$, $\pm 2\omega_c$, and sum frequencies), the total Hamiltonian in the Schrödinger picture becomes \begin{align} \hat{H} &= \hbar\omega_{21} |2\rangle\langle 2| + \hbar\omega_{31} |3\rangle\langle 3| \\ &\quad + \frac{\hbar}{2} \left( \Omega_p e^{-i\omega_p t} |3\rangle\langle 1| + \Omega_p^* e^{i\omega_p t} |1\rangle\langle 3| \right) \nonumber \\ &\quad + \frac{\hbar}{2} \left( \Omega_c e^{-i\omega_c t} |3\rangle\langle 2| + \Omega_c^* e^{i\omega_c t} |2\rangle\langle 3| \right). \nonumber \label{eq:H_Lambda} \end{align}
4.3 Transformation to the Rotating Frame
To remove the explicit time dependence, we perform a unitary transformation to a rotating frame. We write the state vector as \begin{equation} |\tilde{\psi}(t)\rangle = \hat{U}(t) |\psi(t)\rangle, \end{equation} with the time-dependent unitary operator \begin{equation} \hat{U}(t) = e^{i(\omega_p t |3\rangle\langle 3| + (\omega_p - \omega_c)t |2\rangle\langle 2|)}. \end{equation} The Hamiltonian in the rotating frame is \begin{equation} \hat{\tilde{H}} = \hat{U} \hat{H} \hat{U}^\dagger - i\hbar \hat{U} \frac{\mathrm{d}}{\mathrm{d} t} \hat{U}^\dagger. \end{equation} A straightforward calculation yields the time-independent effective Hamiltonian \begin{equation} \boxed{\hat{\tilde{H}} = -\hbar \Delta_p |3\rangle\langle 3| - \hbar(\Delta_p - \Delta_c) |2\rangle\langle 2| + \frac{\hbar}{2} \left( \Omega_p |3\rangle\langle 1| + \Omega_c |3\rangle\langle 2| + \text{H.c.} \right)}, \label{eq:H_rot} \end{equation} where we have defined the one-photon detunings \begin{equation} \Delta_p = \omega_p - \omega_{31}, \qquad \Delta_c = \omega_c - \omega_{32}, \end{equation} and the two-photon detuning (Raman detuning) \begin{equation} \delta = \Delta_p - \Delta_c = \omega_p - \omega_c - \omega_{21}. \end{equation} The two-photon detuning $\delta$ quantifies how far the two-photon transition $|1\rangle \to |3\rangle \to |2\rangle$ is from Raman resonance. When $\delta = 0$, the two-photon resonance condition is satisfied.
It is conventional to re-define the zero of energy so that state $|1\rangle$ has zero energy in the rotating frame. The Hamiltonian then reads \begin{equation} \hat{\tilde{H}} = \hbar \begin{pmatrix} 0 & 0 & \frac{1}{2}\Omega_p^* \\ 0 & -\delta & \frac{1}{2}\Omega_c^* \\ \frac{1}{2}\Omega_p & \frac{1}{2}\Omega_c & -\Delta_p \end{pmatrix}, \label{eq:H_matrix} \end{equation} where the basis is ordered as $\{ |1\rangle, |2\rangle, |3\rangle \}$.
4.4 Spontaneous Emission and the Lindblad Operators
The excited state $|3\rangle$ can decay by spontaneous emission to both $|1\rangle$ and $|2\rangle$. Let $\Gamma_{31}$ and $\Gamma_{32}$ be the partial decay rates from $|3\rangle$ to $|1\rangle$ and $|2\rangle$, respectively. The total radiative decay rate of the excited state is \begin{equation} \Gamma_3 = \Gamma_{31} + \Gamma_{32}. \end{equation} In the Born–Markov approximation, each decay channel gives rise to a Lindblad jump operator: \begin{equation} \hat{L}_{31} = \sqrt{\Gamma_{31}} \, |1\rangle\langle 3|, \qquad \hat{L}_{32} = \sqrt{\Gamma_{32}} \, |2\rangle\langle 3|. \end{equation} The Lindblad superoperator for spontaneous emission is \begin{equation} \mathcal{L}_{\text{se}}[\hat{\rho}] = \sum_{j=1,2} \Gamma_{3j} \left( \hat{\sigma}_{1j} \hat{\rho} \hat{\sigma}_{1j}^\dagger - \frac{1}{2} \{ \hat{\sigma}_{1j}^\dagger \hat{\sigma}_{1j}, \hat{\rho} \} \right), \label{eq:Lindblad_3level} \end{equation} where $\hat{\sigma}_{1j} = |j\rangle\langle 3|$ (for $j=1,2$).
Additionally, there may be dephasing of the ground-state coherence between $|1\rangle$ and $|2\rangle$ due to collisions, magnetic field fluctuations, or laser linewidth. This is modeled by a pure dephasing Lindblad operator \begin{equation} \hat{L}_{\text{deph}} = \sqrt{\gamma_{12}} \, (|1\rangle\langle 1| - |2\rangle\langle 2|), \end{equation} where $\gamma_{12}$ is the ground-state dephasing rate. The total master equation in the rotating frame is \begin{equation} \boxed{\frac{\mathrm{d}}{\mathrm{d} t} \hat{\rho} = -\frac{i}{\hbar} [\hat{\tilde{H}}, \hat{\rho}] + \mathcal{L}_{\text{se}}[\hat{\rho}] + \mathcal{L}_{\text{deph}}[\hat{\rho}]}. \label{eq:master_3level} \end{equation}
4.5 Derivation of the Three-Level Optical Bloch Equations
The density operator in the basis $\{ |1\rangle, |2\rangle, |3\rangle \}$ is a $3 \times 3$ Hermitian matrix: \begin{equation} \hat{\rho} = \begin{pmatrix} \rho_{11} & \rho_{12} & \rho_{13} \\ \rho_{21} & \rho_{22} & \rho_{23} \\ \rho_{31} & \rho_{32} & \rho_{33} \end{pmatrix}, \qquad \rho_{ji} = \rho_{ij}^*, \quad \sum_{i=1}^3 \rho_{ii} = 1. \end{equation} Substituting this into the master equation (\ref{eq:master_3level}) and evaluating the commutators and Lindblad terms, we obtain the full set of optical Bloch equations for the three-level $\Lambda$-system:
4.5.1 Population Equations
For the excited state population: \begin{equation} \boxed{\dot{\rho}_{33} = -\Gamma_3 \rho_{33} - \frac{i}{2} \left( \Omega_p \rho_{13} - \Omega_p^* \rho_{31} \right) - \frac{i}{2} \left( \Omega_c \rho_{23} - \Omega_c^* \rho_{32} \right)}. \label{eq:OBE_33} \end{equation}
For the ground state $|1\rangle$: \begin{equation} \boxed{\dot{\rho}_{11} = +\Gamma_{31} \rho_{33} + \frac{i}{2} \left( \Omega_p \rho_{13} - \Omega_p^* \rho_{31} \right)}. \label{eq:OBE_11} \end{equation}
For the ground state $|2\rangle$: \begin{equation} \boxed{\dot{\rho}_{22} = +\Gamma_{32} \rho_{33} + \frac{i}{2} \left( \Omega_c \rho_{23} - \Omega_c^* \rho_{32} \right)}. \label{eq:OBE_22} \end{equation} Note that $\dot{\rho}_{11} + \dot{\rho}_{22} + \dot{\rho}_{33} = 0$, as required by probability conservation.
4.5.2 Optical Coherences (Coupled to $|3\rangle$)
For the coherence between $|1\rangle$ and $|3\rangle$ (probe transition): \begin{align} \dot{\rho}_{13} &= -\left( \frac{\Gamma_3}{2} + i\Delta_p \right) \rho_{13} + \frac{i}{2} \Omega_p^* (\rho_{33} - \rho_{11}) - \frac{i}{2} \Omega_c^* \rho_{12}. \label{eq:OBE_13} \end{align}
For the coherence between $|2\rangle$ and $|3\rangle$ (control transition): \begin{align} \dot{\rho}_{23} &= -\left( \frac{\Gamma_3}{2} + i\Delta_c \right) \rho_{23} + \frac{i}{2} \Omega_c^* (\rho_{33} - \rho_{22}) - \frac{i}{2} \Omega_p^* \rho_{21}. \label{eq:OBE_23} \end{align}
The Hermitian conjugates are $\dot{\rho}_{31} = (\dot{\rho}_{13})^*$ and $\dot{\rho}_{32} = (\dot{\rho}_{23})^*$.
4.5.3 Ground-State Coherence (Raman Coherence)
For the coherence between $|1\rangle$ and $|2\rangle$: \begin{align} \dot{\rho}_{12} &= -\left( \gamma_{12} + i\delta \right) \rho_{12} + \frac{i}{2} \Omega_p \rho_{32} - \frac{i}{2} \Omega_c^* \rho_{13}. \label{eq:OBE_12} \end{align} This equation is crucial. The two-photon detuning $\delta$ appears here, and the coherence $\rho_{12}$ is driven by a Raman-like process involving the product of the probe field and the control field. When $\delta = 0$ (Raman resonance) and the decay rates are negligible, the steady-state solution can yield $\rho_{12} \neq 0$—this is the basis of coherent population trapping and electromagnetically induced transparency.
Equations (\ref{eq:OBE_33})–(\ref{eq:OBE_12}), together with their complex conjugates, constitute the complete set of nine real-valued (or five complex, with trace condition) optical Bloch equations for the three-level $\Lambda$-system. These equations form the starting point for analyzing a vast array of phenomena in quantum optics and atomic physics.
[Figure 3: Flow diagram of the three-level Bloch equations. The populations $\rho_{11}, \rho_{22}, \rho_{33}$ are coupled to the optical coherences $\rho_{13}, \rho_{23}$, which in turn drive the ground-state coherence $\rho_{12}$. The decay rates $\Gamma_{31}, \Gamma_{32}, \gamma_{12}$ and the Rabi frequencies $\Omega_p, \Omega_c$ are shown on the corresponding arrows.]
5. Steady-State Solution and Key Phenomena
5.1 The Steady-State Equations
Setting the time derivatives to zero in the optical Bloch equations (\ref{eq:OBE_33})–(\ref{eq:OBE_12}) yields a system of linear algebraic equations for the steady-state density matrix elements. The solution can be found analytically, although the expressions are generally cumbersome. We summarize here the most important limits.
5.2 Weak Probe Limit and Electromagnetically Induced Transparency
Consider the weak probe limit: $|\Omega_p| \ll |\Omega_c|, \Gamma_3$. The control field is strong and treated to all orders, while the probe field is treated to first order. Initially, all population is in $|1\rangle$: $\rho_{11}^{(0)} = 1$, $\rho_{22}^{(0)} = \rho_{33}^{(0)} = 0$. To zeroth order in $\Omega_p$, the optical coherences $\rho_{13}$ and $\rho_{23}$ vanish. To first order, the probe-induced coherence $\rho_{13}^{(1)}$ determines the linear susceptibility $\chi_p \propto \rho_{13}^{(1)}/\Omega_p$.
Solving the steady-state equations to first order yields \begin{equation} \rho_{13}^{(1)} = \frac{i\Omega_p/2}{\Gamma_3/2 + i\Delta_p + \dfrac{|\Omega_c|^2/4}{\gamma_{12} + i\delta}}. \label{eq:rho13_weak} \end{equation} This is the fundamental result for EIT. The third term in the denominator, \begin{equation} \frac{|\Omega_c|^2/4}{\gamma_{12} + i\delta}, \end{equation} represents the control-field-induced modification of the probe response. On two-photon resonance ($\delta = 0$) and in the absence of ground-state dephasing ($\gamma_{12} = 0$), this term diverges, leading to $\rho_{13}^{(1)} \to 0$. Consequently, the probe absorption vanishes: the medium becomes transparent at the very line center where it would normally be most opaque. This is the phenomenon of electromagnetically induced transparency (EIT).
Physically, EIT arises from the destructive quantum interference between the two excitation pathways: the direct path $|1\rangle \xrightarrow{\Omega_p} |3\rangle$ and the indirect path $|1\rangle \xrightarrow{\Omega_p} |3\rangle \xrightarrow{\Omega_c^*} |2\rangle \xrightarrow{\Omega_c} |3\rangle$. The interference creates a "dark state"—a coherent superposition of $|1\rangle$ and $|2\rangle$ that does not couple to the excited state—and the probe field cannot be absorbed. This is accompanied by a steep positive dispersion, leading to slow light and enhanced nonlinear optical effects.
5.3 Coherent Population Trapping
When both fields are on exact resonance ($\Delta_p = \Delta_c = 0$) and $\delta = 0$, a remarkable solution to the full optical Bloch equations exists in which the excited state population is zero in the steady state: $\rho_{33} = 0$. The atom is trapped in a pure superposition of the two ground states. The density matrix is \begin{equation} \hat{\rho}_{\text{dark}} = \frac{1}{|\Omega_p|^2 + |\Omega_c|^2} \begin{pmatrix} |\Omega_c|^2 & -\Omega_p^*\Omega_c & 0 \\ -\Omega_p\Omega_c^* & |\Omega_p|^2 & 0 \\ 0 & 0 & 0 \end{pmatrix}. \label{eq:dark_state} \end{equation} The corresponding dark state vector is \begin{equation} |D\rangle = \frac{\Omega_c |1\rangle - \Omega_p |2\rangle}{\sqrt{|\Omega_p|^2 + |\Omega_c|^2}}. \end{equation} One verifies that $\hat{\tilde{H}} |D\rangle = 0$ (for $\Delta_p = \Delta_c = 0$): the dark state is a zero-energy eigenstate of the interacting Hamiltonian. This is the phenomenon of coherent population trapping (CPT), first observed by Alzetta et al. in 1976. CPT is the foundation of STIRAP, laser cooling below the recoil limit (velocity-selective CPT), and atomic clocks.
5.4 Stimulated Raman Adiabatic Passage (STIRAP)
STIRAP is a technique for transferring population between the two ground states $|1\rangle$ and $|2\rangle$ without ever populating the excited state $|3\rangle$. This is achieved by adiabatically varying the Rabi frequencies $\Omega_p(t)$ and $\Omega_c(t)$ in a "counter-intuitive" pulse sequence: the control field is turned on first, preparing the system in the dark state $|D(t)\rangle$, and the probe field is applied subsequently. By the adiabatic theorem, the system remains in the instantaneous dark state, which slowly rotates from $|1\rangle$ at early times to $|2\rangle$ at late times. Since the excited state is never populated, spontaneous emission losses are completely suppressed, enabling near-unit transfer efficiency. STIRAP has become a ubiquitous tool in atomic and molecular physics, quantum information, and chemical reaction dynamics.
[Figure 4: (a) The EIT transmission spectrum: probe absorption (Im$[\chi]$) vanishes on two-photon resonance, accompanied by steep normal dispersion (Re$[\chi]$). (b) The STIRAP pulse sequence: control ($\Omega_c$) precedes probe ($\Omega_p$). (c) The dressed-state picture showing the dark state $|D\rangle$ decoupled from the excited state.]
6. Numerical Considerations and Extensions
6.1 Solving the Bloch Equations Numerically
For arbitrary parameters, the three-level optical Bloch equations must be solved numerically. The standard approach is to vectorize the density matrix into a column vector $\vec{\rho}$ of length $N^2$ (for an $N$-level system) and write the master equation as a linear system: \begin{equation} \frac{\mathrm{d}}{\mathrm{d} t} \vec{\rho} = \mathcal{L} \vec{\rho}, \end{equation} where $\mathcal{L}$ is the $N^2 \times N^2$ Liouvillian superoperator. For the three-level system, $\mathcal{L}$ is a $9 \times 9$ matrix. The steady state is the eigenvector of $\mathcal{L}$ with eigenvalue zero. Standard numerical integration (e.g., Runge–Kutta) or matrix exponentiation ($e^{\mathcal{L}t}$) can be used for time-dependent problems.
6.2 Inclusion of Atomic Motion and Doppler Broadening
In atomic vapors, the thermal motion of atoms introduces Doppler shifts. For an atom moving with velocity $\mathbf{v}$, the laser frequencies are shifted to $\omega_p - \mathbf{k}_p \cdot \mathbf{v}$ and $\omega_c - \mathbf{k}_c \cdot \mathbf{v}$. The detunings become velocity-dependent, and the observed susceptibility is the convolution of the velocity-averaged Bloch equations over the Maxwell–Boltzmann distribution. Doppler broadening is often the dominant line-broadening mechanism at room temperature.
6.3 Multi-Level Extensions and the Full Density Matrix
The formalism developed here for three levels generalizes directly to systems with $N$ levels. The Hamiltonian and Lindblad operators are constructed from the atomic transition operators, and the resulting $N^2 - 1$ independent Bloch equations describe the complete open-system dynamics. This framework is the workhorse for modeling laser spectroscopy, quantum control, and quantum information processing in atomic and molecular systems.
References
- H.-P. Breuer and F. Petruccione, The Theory of Open Quantum Systems (Oxford University Press, Oxford, 2002). The definitive monograph on master equations, the Born–Markov approximation, and the Lindblad form. Essential reading for any serious student of open quantum systems.
- C. W. Gardiner and P. Zoller, Quantum Noise, 3rd ed. (Springer, Berlin, 2004). A comprehensive treatment of quantum stochastic processes and the input–output formalism, with detailed derivations of master equations for quantum optical systems.
- M. O. Scully and M. S. Zubairy, Quantum Optics (Cambridge University Press, Cambridge, 1997). Chapters 5 and 7 provide a pedagogical derivation of the master equation and its application to three-level systems, EIT, and correlated emission lasers.
- M. Fleischhauer, A. Imamoglu, and J. P. Marangos, "Electromagnetically induced transparency: Optics in coherent media," Reviews of Modern Physics 77, 633–673 (2005). An authoritative review of EIT and its applications.
- K. Bergmann, H. Theuer, and B. W. Shore, "Coherent population transfer among quantum states of atoms and molecules," Reviews of Modern Physics 70, 1003–1025 (1998). The classic review of STIRAP and adiabatic passage techniques.
- G. Lindblad, "On the generators of quantum dynamical semigroups," Communications in Mathematical Physics 48, 119–130 (1976). The original mathematical paper establishing the general form of completely positive Markovian master equations.
- V. Gorini, A. Kossakowski, and E. C. G. Sudarshan, "Completely positive dynamical semigroups of N-level systems," Journal of Mathematical Physics 17, 821–825 (1976). The independent derivation of the Lindblad form.
- C. Cohen-Tannoudji, J. Dupont-Roc, and G. Grynberg, Atom–Photon Interactions: Basic Processes and Applications (Wiley-VCH, Weinheim, 1998). Supplementary material on master equations and the optical Bloch equations for multi-level atoms.
- D. A. Steck, "Quantum and Atom Optics," available online at http://steck.us/teaching (revision 0.14.5, 2023). Comprehensive lecture notes with detailed derivations of the three-level Bloch equations and numerical solution methods.
- R. Loudon, The Quantum Theory of Light, 3rd ed. (Oxford University Press, Oxford, 2000). A classic text covering reservoir theory and the master equation approach to spontaneous emission.