Thermal Master Equations

A thermal master equation is a Markovian open-system equation for a system weakly coupled to an equilibrium heat bath. Its rates are constrained by bath temperature. A standard weak-coupling secular generator has the bare Gibbs state as a stationary state; convergence to that state additionally requires sufficient coupling between sectors and the absence of persistent modes.

The minimal thermal consistency test is:

$$\mathcal L(\rho_\beta)=0, \qquad \rho_\beta = \frac{e^{-\beta H_S}} {\operatorname{Tr}(e^{-\beta H_S})}.$$

For a standard weak-coupling secular derivation, this stationarity is enforced by Detailed Balance. A master equation can be in Lindblad form and still fail to be the correct thermal equation if its rates do not satisfy the thermal relations.

Helpful background. Detailed balance fixes the rate ratios, and the secular approximation separates resolved energy transitions. Steady states and relaxation distinguishes stationarity from attraction.

Thermal rates and Gibbs stationarity

A Markovian master equation is thermal when the environment is modeled as a stationary equilibrium reservoir and the reduced generator reflects that equilibrium.

The usual ingredients are:

• a system Hamiltonian $H_S$ with resolved energy gaps;
• a bath Hamiltonian $H_B$ and thermal state $\rho_B=e^{-\beta H_B}/Z_B$;
• weak system–bath coupling;
• bath correlations satisfying the KMS condition;
• Markov and secular approximations;
• rates satisfying detailed balance;
• Gibbs stationarity, with corrections only when calculated consistently beyond the approximation.

The word “thermal” should not be attached merely because a dissipator has positive rates. The rates must know the bath temperature.

Standard Weak-Coupling Form

Take a finite-dimensional system, or a justified infinite-dimensional generator domain, and Hermitian coupling operators. Start from an interaction

$$H_I = \sum_\alpha A_\alpha\otimes B_\alpha.$$

Decompose the system operators into Bohr-frequency components:

$$A_\alpha(\omega) = \sum_{\epsilon'-\epsilon=\hbar\omega} \Pi(\epsilon)A_\alpha\Pi(\epsilon').$$

With this convention, positive $\omega$ labels a system operator that lowers the system energy by $\hbar\omega$.

After Born, Markov, and secular approximations, a common thermal generator is

$$\frac{d\rho}{dt} = - \frac{i}{\hbar} [H_S+H_{\mathrm{LS}},\rho] + \sum_{\omega,\alpha,\beta} \gamma_{\alpha\beta}(\omega) \left( A_\beta(\omega)\rho A_\alpha^\dagger(\omega) - \frac12 \{A_\alpha^\dagger(\omega)A_\beta(\omega),\rho\} \right).$$

For dimensionless $A_\alpha$ and energy-valued $B_\alpha$, the rate matrix is $\gamma_{\alpha\beta}(\omega)=\hbar^{-2}\int_{-\infty}^\infty dt\,e^{i\omega t}\langle B_\alpha(t)B_\beta(0)\rangle$, with any small coupling included in $B_\alpha$. For an equilibrium bath, the KMS relation implies detailed-balance relations between the $\omega$ and $-\omega$ blocks. For one Hermitian bath operator, the relation is

$$\gamma(-\omega) = e^{-\beta\hbar\omega}\gamma(\omega), \qquad \omega>0.$$

This makes upward transitions thermally suppressed. Also $[H_{\mathrm{LS}},H_S]=0$ in the fully secular construction. Hence the Lamb shift does not change the bare Gibbs stationary state. Inserting shifted gaps into only the Gibbs weights, while retaining bare-gap rates, would be inconsistent.

Two-Level System

Let $E_e-E_g=\hbar\omega_0$ with $\omega_0>0$. A finite-temperature two-level thermal master equation is

$$\frac{d\rho}{dt} = - \frac{i}{\hbar}[H,\rho] + \Gamma_\downarrow\mathcal D[\sigma_-]\rho + \Gamma_\uparrow\mathcal D[\sigma_+]\rho,$$

Here $H$ is diagonal in this energy basis, with no transverse coherent drive, and the two rates have positive sum. The jump operators are

$$\sigma_-=\lvert g\rangle\langle e\rvert, \qquad \sigma_+=\lvert e\rangle\langle g\rvert.$$

Thermal detailed balance requires

$$\frac{\Gamma_\uparrow}{\Gamma_\downarrow} = e^{-\beta\hbar\omega_0}.$$

The excited-state population obeys

$$\dot p_e = -\Gamma_\downarrow p_e + \Gamma_\uparrow(1-p_e),$$

This is the two-state instance of a Pauli Rate Equation.

The steady state is therefore

$$p_e^{\mathrm{ss}} = \frac{\Gamma_\uparrow} {\Gamma_\uparrow+\Gamma_\downarrow}.$$

Using detailed balance,

$$p_e^{\mathrm{ss}} = \frac{1} {1+e^{\beta\hbar\omega_0}},$$

which is the Gibbs excited-state population for a two-level system.

Zero-temperature amplitude damping is the special case $\Gamma_\uparrow=0$. It is not a finite-temperature thermal model unless upward excitation is negligible. See Amplitude Damping Master Equation for the two-level relaxation generator and its zero-temperature limit.

Harmonic Oscillator

For the occupation factors and vacuum contribution, see Thermal and Vacuum Noise. For a harmonic oscillator coupled to a thermal bath with mean occupation

$$\bar n = \frac{1} {e^{\beta\hbar\omega}-1},$$

a standard thermal damping equation is

$$\frac{d\rho}{dt} = - \frac{i}{\hbar} [\hbar\omega a^\dagger a,\rho] + \kappa(\bar n+1)\mathcal D[a]\rho + \kappa\bar n\,\mathcal D[a^\dagger]\rho.$$

The $a$ term describes loss of one quantum; $a^\dagger$ describes absorption from the bath. With $\kappa>0$ and finite $\bar n$, the stationary distribution is geometric:

$$p_n=(1-r)r^n,\qquad r=\frac{\bar n}{\bar n+1}=e^{-\beta\hbar\omega}.$$

For $\bar n=1$, $p_n=2^{-n-1}$ and $\langle n\rangle=1$. This is not a two-level distribution: truncating to levels zero and one and renormalizing would give mean occupation $1/3$. A numerical oscillator cutoff must retain enough of the thermal tail. The moment equations below assume finite first occupation moment.

The mean occupation obeys

$$\frac{d}{dt}\langle n\rangle = -\kappa(\langle n\rangle-\bar n), \qquad n=a^\dagger a.$$

Thus the oscillator relaxes toward the thermal occupation:

$$\langle n(t)\rangle = \bar n + e^{-\kappa t} \left( \langle n(0)\rangle-\bar n \right).$$

For the reverse pair $n-1\leftrightarrow n$, both rates contain the same matrix-element factor $n$, so their ratio is

$$\frac{\kappa\bar n} {\kappa(\bar n+1)} = e^{-\beta\hbar\omega}.$$

Pure Dephasing Is Not Thermalization

A pure-dephasing term such as

$$\gamma_\phi\mathcal D[\sigma_z]\rho$$

may arise from thermal noise at zero frequency, but it does not by itself drive the system toward a Gibbs distribution. It damps coherences in the $\sigma_z$ basis while leaving populations unchanged.

Thermalization requires energy-exchange channels with rates that satisfy detailed balance. Dephasing can accompany thermal relaxation, but it is not a substitute for upward and downward transition terms.

Global Versus Local Thermal Equations

For a composite system with interacting parts, the safest weak-coupling thermal derivation usually uses the eigenbasis of the full system Hamiltonian $H_S$, not the bare Hamiltonians of the subsystems separately.

This is called a global master equation. It constructs jump operators from the Bohr frequencies of the complete $H_S$ and is the natural setting for detailed balance.

A local master equation uses dissipators such as

$$\mathcal D[a_1], \qquad \mathcal D[a_2],$$

on subsystems as if their couplings to the bath were independent of the interactions inside $S$. Local equations can be useful approximations, especially when internal couplings are weak or the modeling target is phenomenological, but they may violate detailed balance or predict heat currents in situations that should be equilibrium.

A fully secular global construction also needs resolved transition gaps. If internal splittings compete with the dissipative linewidth, neither an indiscriminately secular global equation nor an arbitrary local equation is automatically accurate; use a controlled partial-secular or coarse-grained derivation.

Two coupled qubits: which Gibbs state?

Let

$$H_S=\hbar\omega(n_1+n_2) +\hbar J(|10\rangle\langle01|+|01\rangle\langle10|), \qquad 0<J<\omega.$$

The single-excitation eigenstates $|\pm\rangle=(|10\rangle\pm|01\rangle)/\sqrt2$ have energies $\hbar(\omega\pm J)$. Their unequal Gibbs weights produce

$$\langle10|\rho_\beta|01\rangle =-\frac{e^{-\beta\hbar\omega}\sinh(\beta\hbar J)}{Z},$$

where $Z=1+2e^{-\beta\hbar\omega}\cosh(\beta\hbar J)+e^{-2\beta\hbar\omega}$.

In contrast, identical local thermal dissipators at the bare gap preserve the product of the two bare qubit Gibbs states. That product commutes with the exchange Hamiltonian but has zero $10,01$ coherence, so it differs from the global Gibbs state. For $\beta\hbar\omega=\ln3$ and $\beta\hbar J=\ln2$, the global coherence is $-9/70$ while the local product gives zero. This is a concrete thermal-consistency difference, not merely a choice of notation.

Strong-Coupling Caveat

At stronger system–bath coupling, the bare Gibbs state

$$\frac{e^{-\beta H_S}} {\operatorname{Tr}(e^{-\beta H_S})}$$

need not be the correct reduced equilibrium state. The equilibrium reduced state of the coupled total system is

$$\rho_S^{\mathrm{eq}} = \operatorname{Tr}_B \frac{ e^{-\beta(H_S+H_B+H_I)} }{ \operatorname{Tr} e^{-\beta(H_S+H_B+H_I)} }.$$

This can be represented using a Hamiltonian of mean force, not simply $H_S$. A weak-coupling thermal master equation may still be useful, but one should not demand exact relaxation to the bare Gibbs state outside its regime. See Strong Coupling for the broader warning.

Multiple Baths and Nonequilibrium Steady States

If several baths are present, the generator may be a sum

$$\mathcal L = \sum_r\mathcal L_r,$$

where each $\mathcal L_r$ satisfies detailed balance at its own inverse temperature $\beta_r$.

If the temperatures or chemical potentials differ, there is generally no single Gibbs steady state. The system may reach a nonequilibrium steady state with heat or particle currents.

This is not a failure of Lindblad form. It is a different physical situation from a single equilibrium thermal bath.

Slowly changing Hamiltonians

Quantum annealing applies thermal-rate models to changing energy levels. Equilibrium is then an instantaneous reference, not a guarantee that the state follows it. Compare the schedule timescale with bath memory, energy-gap resolution, and population relaxation before using a sequence of thermal generators.

Exercises

Two-level Gibbs population

Use detailed balance to show that the steady state of

$$\dot p_e = -\Gamma_\downarrow p_e + \Gamma_\uparrow(1-p_e)$$

has the Gibbs excited-state population.

Solution

The steady state is

$$p_e^{\mathrm{ss}} = \frac{\Gamma_\uparrow} {\Gamma_\uparrow+\Gamma_\downarrow}.$$

Detailed balance gives

$$\Gamma_\uparrow = e^{-\beta\hbar\omega_0}\Gamma_\downarrow.$$

Substitute:

$$p_e^{\mathrm{ss}} = \frac{e^{-\beta\hbar\omega_0}} {1+e^{-\beta\hbar\omega_0}} = \frac{1} {1+e^{\beta\hbar\omega_0}}.$$

Oscillator occupation relaxation

Given

$$\frac{d}{dt}\langle n\rangle = -\kappa(\langle n\rangle-\bar n),$$

solve for $\langle n(t)\rangle$.

Solution

This first-order equation has solution

$$\langle n(t)\rangle-\bar n = e^{-\kappa t} \left( \langle n(0)\rangle-\bar n \right).$$

Therefore

$$\langle n(t)\rangle = \bar n + e^{-\kappa t} \left( \langle n(0)\rangle-\bar n \right).$$

Zero-temperature limit

Take $\bar n\to0$ in the oscillator thermal damping equation. What remains?

Solution

When $\bar n=0$,

$$\kappa(\bar n+1)\mathcal D[a]\rho + \kappa\bar n\,\mathcal D[a^\dagger]\rho \longrightarrow \kappa\mathcal D[a]\rho.$$

Only loss remains. The oscillator relaxes toward the vacuum.

Why local can fail

Why can a local dissipator built from bare subsystem lowering operators fail detailed balance for an interacting composite system?

Solution

Detailed balance relates rates at the actual transition frequencies of the Hamiltonian whose Gibbs state is expected. If the subsystems interact, the eigenstates and energy gaps of the full $H_S$ need not match the bare subsystem transitions. A dissipator built from bare lowering operators can therefore assign rates to the wrong frequencies and fail to make the global Gibbs state stationary.

References

• R. Alicki and K. Lendi, Quantum Dynamical Semigroups and Applications, Springer (1987).
• H.-P. Breuer and F. Petruccione, The Theory of Open Quantum Systems, Oxford University Press (2002). Publisher.
• M. Cattaneo, G. L. Giorgi, S. Maniscalco, and R. Zambrini, “Local versus global master equation with common and separate baths: superiority of the global approach in partial secular approximation,” New Journal of Physics 21, 113045 (2019). Author manuscript.
• E. B. Davies, Quantum Theory of Open Systems, Academic Press (1976).
• C. W. Gardiner and P. Zoller, Quantum Noise, 3rd ed., Springer (2004). Publisher.
• Á. Rivas and S. F. Huelga, Open Quantum Systems: An Introduction, Springer (2012). Author manuscript.