Lindblad equation (also called the GKSL equation, after Gorini, Kossakowski, Sudarshan, and Lindblad) is the most general Markovian master equation for the density matrix $\rho$ of an open quantum system. It extends the von Neumann equation by adding a dissipator term that models incoherent processes — energy loss, dephasing, and other forms of coupling to an environment.

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = -\frac{i}{\hbar}[H,\,\rho] + \sum_k\!\left(L_k\rho L_k^\dagger - \frac{1}{2}\{L_k^\dagger L_k,\,\rho\}\right)$$

The first term $-\frac{i}{\hbar}[H, \rho]$ is the von Neumann term. It describes the unitary, reversible part of the evolution driven by the system Hamiltonian $H$ — the same dynamics that the Schrödinger equation captures for pure states. When the dissipator is zero (no environment coupling), the Lindblad equation reduces exactly to the von Neumann equation.

The sum $\sum_k\!\left(L_k\rho L_k^\dagger - \frac{1}{2}\{L_k^\dagger L_k, \rho\}\right)$ is the dissipator. Each term corresponds to one incoherent channel indexed by $k$. The $L_k\rho L_k^\dagger$ part is the quantum jump: the environment has “observed” the system through channel $k$, collapsing it toward a new state. The $-\frac{1}{2}\{L_k^\dagger L_k, \rho\}$ part — with $\{A, B\} = AB + BA$ the anticommutator — is the no-jump correction that keeps $\text{tr}(\rho) = 1$ by accounting for the probability that no jump has occurred. Together, the two parts form a completely positive trace-preserving (CPTP) map on $\rho$.

The operators $L_k$ are called jump operators (or Lindblad operators). Their choice encodes the physics of the environment — different environments couple to the system in different ways, and each coupling mechanism gets its own $L_k$. The six examples below cover the most common cases for qubits and bosonic modes.

Relaxation is spontaneous decay from the excited state $\lvert 1\rangle$ down to the ground state $\lvert 0\rangle$. Physically this is what happens when a two-level atom emits a photon, or when a superconducting qubit dissipates energy into its microwave environment. The jump operator is the lowering operator $\sigma^-$, which maps $\lvert 1\rangle \to \lvert 0\rangle$ and annihilates $\lvert 0\rangle$:

$$L = \sqrt{\gamma}\,\sigma^-, \qquad \sigma^- = \lvert 0\rangle\langle 1\rvert = \begin{pmatrix}0 & 1 \\ 0 & 0\end{pmatrix}$$

Its Hermitian adjoint is the raising operator $\sigma^+ = (\sigma^-)^\dagger = \lvert 1\rangle\langle 0\rvert$:

$$\sigma^+ = \begin{pmatrix}0 & 0 \\ 1 & 0\end{pmatrix}$$

Note that $\sigma^+\sigma^- = \lvert 1\rangle\langle 1\rvert$ is the projector onto the excited state. With $L^\dagger L = \gamma\,\sigma^+\sigma^-$, the special-case Lindblad equation (single channel, $H = 0$) is:

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \gamma\!\left(\sigma^-\rho\,\sigma^+ - \frac{1}{2}\{\sigma^+\sigma^-,\,\rho\}\right)$$

In components this drives $\dot\rho_{11} = -\gamma\rho_{11}$ (excited population decays at rate $\gamma = 1/T_1$) and $\dot\rho_{01} = -\frac{\gamma}{2}\rho_{01}$ (coherences also decay, at half the rate). The ground-state population grows accordingly: $\dot\rho_{00} = +\gamma\rho_{11}$.

Pure dephasing destroys off-diagonal coherences with no energy exchange — the populations $\rho_{00}$ and $\rho_{11}$ stay fixed while the off-diagonals $\rho_{01}$ and $\rho_{10}$ decay. Physically this comes from random fluctuations in the qubit's transition frequency (e.g. charge noise in superconducting qubits, or collisions in atomic systems) that scramble the relative phase between $\lvert 0\rangle$ and $\lvert 1\rangle$.

$$L = \sqrt{\frac{\gamma_\phi}{2}}\,\sigma_z, \qquad \sigma_z = \begin{pmatrix}1 & 0 \\ 0 & -1\end{pmatrix}$$

Because $\sigma_z^2 = I$, we have $L^\dagger L = (\gamma_\phi/2)\,I$, so the anticommutator term simplifies: $\frac{1}{2}\{L^\dagger L, \rho\} = \frac{\gamma_\phi}{2}\rho$. The special-case Lindblad equation becomes:

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \frac{\gamma_\phi}{2}\!\left(\sigma_z\rho\,\sigma_z - \rho\right)$$

In components: $\dot\rho_{00} = 0$, $\dot\rho_{11} = 0$, $\dot\rho_{01} = -\gamma_\phi\,\rho_{01}$, $\dot\rho_{10} = -\gamma_\phi\,\rho_{10}$. The coherences decay exponentially at rate $\gamma_\phi = 1/T_\phi$. Combined with relaxation, the total dephasing rate is $1/T_2 = 1/(2T_1) + 1/T_\phi$.

At nonzero temperature the environment contains thermal photons, which can excite the qubit from $\lvert 0\rangle$ to $\lvert 1\rangle$. This is the reverse process from relaxation. The jump operator is the raising operator $\sigma^+$:

$$L = \sqrt{\gamma_\uparrow}\,\sigma^+, \qquad \sigma^+ = \lvert 1\rangle\langle 0\rvert = \begin{pmatrix}0 & 0 \\ 1 & 0\end{pmatrix}$$

The thermal excitation rate is $\gamma_\uparrow = \gamma\,\bar{n}$, where $\bar{n} = (e^{\hbar\omega/k_BT} - 1)^{-1}$ is the mean photon number of the bath at the qubit frequency $\omega$ and temperature $T$. Note that $\sigma^-\sigma^+ = \lvert 0\rangle\langle 0\rvert$ is the projector onto the ground state. The special-case Lindblad equation is:

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \gamma_\uparrow\!\left(\sigma^+\rho\,\sigma^- - \frac{1}{2}\{\sigma^-\sigma^+,\,\rho\}\right)$$

In components: $\dot\rho_{00} = -\gamma_\uparrow\rho_{00}$ (ground population depletes) and $\dot\rho_{11} = +\gamma_\uparrow\rho_{00}$. Together with the relaxation channel ($L = \sqrt{\gamma}\,\sigma^-$, rate $\gamma = \gamma_\downarrow$), the steady state is the thermal equilibrium $\rho_\infty = \text{diag}(1/(1+\bar{n}),\; \bar{n}/(1+\bar{n}))$, i.e. the qubit thermalises to the bath temperature.

Depolarizing noise applies random Pauli errors in all three directions with equal probability, driving $\rho$ toward the maximally mixed state $I/2$. It requires three separate jump operators, one per Pauli:

$$L_x = \sqrt{\frac{p}{3}}\,\sigma_x, \qquad L_y = \sqrt{\frac{p}{3}}\,\sigma_y, \qquad L_z = \sqrt{\frac{p}{3}}\,\sigma_z$$

The three Pauli matrices in full:

$$\sigma_x = \begin{pmatrix}0 & 1 \\ 1 & 0\end{pmatrix}, \qquad \sigma_y = \begin{pmatrix}0 & -i \\ i & 0\end{pmatrix}, \qquad \sigma_z = \begin{pmatrix}1 & 0 \\ 0 & -1\end{pmatrix}$$

All three satisfy $\sigma_j^2 = I$, so $\frac{1}{2}\{L_j^\dagger L_j, \rho\} = \frac{p}{3}\rho$ for each $j$. The special-case Lindblad equation (H = 0, three channels):

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \frac{p}{3}\sum_{j\in\{x,y,z\}}\!\left(\sigma_j\rho\,\sigma_j - \rho\right)$$

Using the Pauli identity $\sum_j \sigma_j\rho\,\sigma_j = 2\,\text{tr}(\rho)\,I - \rho = 2I - \rho$ (which holds for any $2\times 2$ matrix $\rho$ with $\text{tr}(\rho) = 1$), this simplifies to:

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \frac{p}{3}(2I - \rho - 3\rho) = \frac{4p}{3}\!\left(\frac{I}{2} - \rho\right)$$

The fixed point is $\rho = I/2$: any initial state is driven toward the completely mixed state at rate $4p/3$.

For a bosonic mode (cavity field, mechanical oscillator, LC circuit), photon loss into the environment is described by the bosonic annihilation operator $a$. In the Fock (photon-number) basis $\{\lvert 0\rangle, \lvert 1\rangle, \lvert 2\rangle, \ldots\}$, $a$ acts as $a\lvert n\rangle = \sqrt{n}\,\lvert n-1\rangle$ — it removes one photon and lowers the mode by one excitation level.

$$L = \sqrt{\kappa}\,a, \qquad a = \begin{pmatrix}0 & \sqrt{1} & 0 & 0 & \cdots \\ 0 & 0 & \sqrt{2} & 0 & \cdots \\ 0 & 0 & 0 & \sqrt{3} & \cdots \\ \vdots & & & & \ddots\end{pmatrix}$$

The rate $\kappa = 1/T_{\text{cavity}}$ is the photon loss rate (inverse cavity lifetime). The Hermitian adjoint $a^\dagger$ is the creation operator; $a^\dagger a = \hat{n}$ is the photon number operator. The special-case Lindblad equation (H = 0, single channel):

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \kappa\!\left(a\rho\, a^\dagger - \frac{1}{2}\{a^\dagger a,\,\rho\}\right)$$

In the Fock basis, the matrix element $\rho_{mn} = \langle m\rvert\rho\lvert n\rangle$ evolves as:

$$\dot\rho_{mn} = \kappa\!\left(\sqrt{(m+1)(n+1)}\;\rho_{m+1,\,n+1} - \frac{m+n}{2}\,\rho_{mn}\right)$$

The mean photon number decays exponentially: $d\langle\hat{n}\rangle/dt = -\kappa\langle\hat{n}\rangle$. The steady state is the vacuum $\lvert 0\rangle\langle 0\rvert$.

Bosonic dephasing preserves photon number populations $\rho_{nn}$ but destroys coherences between Fock states with different photon numbers. The jump operator is the number operator $\hat{n} = a^\dagger a$:

$$L = \sqrt{\kappa_\phi}\,\hat{n} = \sqrt{\kappa_\phi}\,a^\dagger a, \qquad a^\dagger a = \begin{pmatrix}0 & 0 & 0 & \cdots \\ 0 & 1 & 0 & \cdots \\ 0 & 0 & 2 & \cdots \\ \vdots & & & \ddots\end{pmatrix}$$

Physically this arises when the mode frequency fluctuates randomly (e.g. due to a dispersive coupling to a fluctuating qubit, or dielectric loss that is proportional to photon number). The special-case Lindblad equation (H = 0, single channel):

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \kappa_\phi\!\left(\hat{n}\,\rho\,\hat{n} - \frac{1}{2}\{\hat{n}^2,\,\rho\}\right)$$

In the Fock basis the result is particularly clean. Since $\hat{n}\lvert n\rangle = n\lvert n\rangle$ is diagonal, every matrix element picks up a simple factor:

$$\dot\rho_{mn} = -\frac{\kappa_\phi}{2}(m - n)^2\,\rho_{mn}$$

Diagonal elements $(m = n)$ are exactly preserved. Off-diagonal elements decay at a rate proportional to the squared difference in photon number, so states separated by many quanta dephase much faster than adjacent Fock states.

A qubit is the simplest open quantum system: a two-level system with states $\lvert 0\rangle$ (ground) and $\lvert 1\rangle$ (excited). Its density matrix is a $2\times 2$ matrix:

$$\rho = \begin{pmatrix}\rho_{00} & \rho_{01} \\ \rho_{10} & \rho_{11}\end{pmatrix}$$

The diagonal entries $\rho_{00}$ and $\rho_{11}$ are real, non-negative populations summing to one ($\rho_{00} + \rho_{11} = 1$). The off-diagonal entries $\rho_{01} = \rho_{10}^*$ are complex coherences. A pure state $\lvert\psi\rangle = \alpha\lvert 0\rangle + \beta\lvert 1\rangle$ corresponds to $\rho_{00} = |\alpha|^2$, $\rho_{11} = |\beta|^2$, $\rho_{01} = \alpha\beta^*$.

The standard qubit Hamiltonian describes a two-level system with energy splitting $\hbar\omega$ between ground and excited state:

$$H = \frac{\hbar\omega}{2}\sigma_z = \frac{\hbar\omega}{2}\begin{pmatrix}1&0\\0&-1\end{pmatrix}$$

To write the full Lindblad equation for a realistic qubit, use two jump operators — one for T1 relaxation and one for pure dephasing:

$$L_1 = \sqrt{\gamma}\,\sigma^- = \sqrt{\gamma}\begin{pmatrix}0&1\\0&0\end{pmatrix}, \qquad L_2 = \sqrt{\frac{\gamma_\phi}{2}}\,\sigma_z = \sqrt{\frac{\gamma_\phi}{2}}\begin{pmatrix}1&0\\0&-1\end{pmatrix}$$

Now work out each term of $\dot\rho = -\frac{i}{\hbar}[H,\rho] + \sum_k\!\left(L_k\rho L_k^\dagger - \frac{1}{2}\{L_k^\dagger L_k,\rho\}\right)$ in matrix form.

Coherent term. The commutator $[H, \rho]$:

$$[\sigma_z, \rho] = \begin{pmatrix}1&0\\0&-1\end{pmatrix}\begin{pmatrix}\rho_{00}&\rho_{01}\\\rho_{10}&\rho_{11}\end{pmatrix} - \begin{pmatrix}\rho_{00}&\rho_{01}\\\rho_{10}&\rho_{11}\end{pmatrix}\begin{pmatrix}1&0\\0&-1\end{pmatrix} = \begin{pmatrix}0&2\rho_{01}\\-2\rho_{10}&0\end{pmatrix}$$

$$-\frac{i}{\hbar}[H,\rho] = -\frac{i\omega}{2}\begin{pmatrix}0&2\rho_{01}\\-2\rho_{10}&0\end{pmatrix} = \begin{pmatrix}0&-i\omega\rho_{01}\\i\omega\rho_{10}&0\end{pmatrix}$$

The coherent term drives oscillation of the coherences at frequency $\omega$ (Larmor precession). It leaves populations untouched.

T1 dissipator. With $L_1^\dagger L_1 = \gamma\,\sigma^+\sigma^- = \gamma\begin{pmatrix}0&0\\0&1\end{pmatrix}$:

$$\sigma^-\rho\,\sigma^+ = \begin{pmatrix}0&1\\0&0\end{pmatrix}\begin{pmatrix}\rho_{00}&\rho_{01}\\\rho_{10}&\rho_{11}\end{pmatrix}\begin{pmatrix}0&0\\1&0\end{pmatrix} = \begin{pmatrix}\rho_{11}&0\\0&0\end{pmatrix}$$

$$\frac{1}{2}\{L_1^\dagger L_1,\rho\} = \frac{\gamma}{2}\left\{\begin{pmatrix}0&0\\0&1\end{pmatrix},\rho\right\} = \frac{\gamma}{2}\begin{pmatrix}0&\rho_{01}\\\rho_{10}&2\rho_{11}\end{pmatrix}$$

$$L_1\rho L_1^\dagger - \frac{1}{2}\{L_1^\dagger L_1,\rho\} = \gamma\begin{pmatrix}\rho_{11}&0\\0&0\end{pmatrix} - \frac{\gamma}{2}\begin{pmatrix}0&\rho_{01}\\\rho_{10}&2\rho_{11}\end{pmatrix} = \gamma\begin{pmatrix}\rho_{11}&-\tfrac{1}{2}\rho_{01}\\-\tfrac{1}{2}\rho_{10}&-\rho_{11}\end{pmatrix}$$

Dephasing dissipator. With $L_2^\dagger L_2 = (\gamma_\phi/2)\,\sigma_z^2 = (\gamma_\phi/2)\,I$ (since $\sigma_z^2 = I$):

$$\sigma_z\rho\,\sigma_z = \begin{pmatrix}\rho_{00}&-\rho_{01}\\-\rho_{10}&\rho_{11}\end{pmatrix}, \qquad \frac{1}{2}\{L_2^\dagger L_2,\rho\} = \frac{\gamma_\phi}{2}\rho$$

$$L_2\rho L_2^\dagger - \frac{1}{2}\{L_2^\dagger L_2,\rho\} = \frac{\gamma_\phi}{2}(\sigma_z\rho\,\sigma_z - \rho) = \begin{pmatrix}0&-\gamma_\phi\rho_{01}\\-\gamma_\phi\rho_{10}&0\end{pmatrix}$$

Full equation. Sum all three contributions:

$$\dot\rho = \begin{pmatrix}0&-i\omega\rho_{01}\\i\omega\rho_{10}&0\end{pmatrix} + \gamma\begin{pmatrix}\rho_{11}&-\tfrac{1}{2}\rho_{01}\\-\tfrac{1}{2}\rho_{10}&-\rho_{11}\end{pmatrix} + \begin{pmatrix}0&-\gamma_\phi\rho_{01}\\-\gamma_\phi\rho_{10}&0\end{pmatrix}$$

$$\dot\rho = \begin{pmatrix}\gamma\rho_{11} & \left(-i\omega - \tfrac{\gamma}{2} - \gamma_\phi\right)\rho_{01} \\ \left(i\omega - \tfrac{\gamma}{2} - \gamma_\phi\right)\rho_{10} & -\gamma\rho_{11}\end{pmatrix}$$

Reading off the four equations for the matrix entries:

$$\dot\rho_{00} = +\gamma\,\rho_{11}$$

$$\dot\rho_{11} = -\gamma\,\rho_{11}$$

$$\dot\rho_{01} = \left(-i\omega - \frac{1}{T_2}\right)\rho_{01}, \qquad \frac{1}{T_2} = \frac{\gamma}{2} + \gamma_\phi$$

$$\dot\rho_{10} = \left(+i\omega - \frac{1}{T_2}\right)\rho_{10}$$

These are the optical Bloch equations. The solutions are direct:

$$\rho_{11}(t) = \rho_{11}(0)\,e^{-\gamma t} = \rho_{11}(0)\,e^{-t/T_1}$$

$$\rho_{01}(t) = \rho_{01}(0)\,e^{(-i\omega - 1/T_2)\,t}$$

The excited population decays exponentially at rate $\gamma = 1/T_1$. The coherences oscillate at the qubit frequency $\omega$ while their magnitude decays at rate $1/T_2 = \gamma/2 + \gamma_\phi$. Since $T_2 \leq 2T_1$, dephasing is always at least as fast as the T1-limited floor.

The Lindblad equation is valid under the Markov approximation: the environment is assumed to be memoryless, forgetting the interaction faster than the system evolves. This is analogous to an interrupt-driven peripheral that completes its transaction well within one main-loop cycle — the system never needs to track the peripheral's internal state. When the bath correlation time $\tau_{\text{bath}}$ is comparable to the system timescale, memory effects matter and the full non-Markovian dynamics must be used instead.

The Lindblad equation can be derived from the quantum channel picture. Over a short time step $dt$, a quantum channel is described by a set of Kraus operators $\{K_k\}$ satisfying the completeness relation $\sum_k K_k^\dagger K_k = I$. The density matrix evolves as:

$$\rho(t + dt) = \sum_k K_k\,\rho(t)\,K_k^\dagger$$

This is the most general evolution that is linear, completely positive, and trace-preserving. For a single dissipation channel described by a jump operator $L$, use two Kraus operators:

$$K_0 = I - \frac{i}{\hbar}H\,dt - \frac{1}{2}L^\dagger L\,dt, \qquad K_1 = \sqrt{dt}\;L$$

$K_1$ is the jump Kraus operator: it fires with probability proportional to $dt$, representing one quantum jump. $K_0$ is the no-jump Kraus operator: it captures the coherent Hamiltonian evolution plus the accumulated “penalty” from not jumping.

First, verify that these satisfy the completeness relation to first order in $dt$ (anything proportional to $dt^2$ is dropped). Take the dagger of $K_0$:

$$K_0^\dagger = I + \frac{i}{\hbar}H\,dt - \frac{1}{2}L^\dagger L\,dt$$

Then:

$$K_0^\dagger K_0 = \left(I + \frac{i}{\hbar}H\,dt - \frac{1}{2}L^\dagger L\,dt\right)\!\left(I - \frac{i}{\hbar}H\,dt - \frac{1}{2}L^\dagger L\,dt\right)$$

Expand and drop $dt^2$ terms (e.g. $H^2 dt^2$ and $H L^\dagger L\,dt^2$ all vanish). The surviving terms are:

$$K_0^\dagger K_0 = I - L^\dagger L\,dt$$

The jump term gives:

$$K_1^\dagger K_1 = \left(\sqrt{dt}\;L\right)^\dagger\!\left(\sqrt{dt}\;L\right) = dt\;L^\dagger L$$

Sum them:

$$K_0^\dagger K_0 + K_1^\dagger K_1 = I - L^\dagger L\,dt + L^\dagger L\,dt = I \quad\checkmark$$

The completeness relation holds to first order. Now expand the full evolution $\rho(t+dt) = K_0\rho K_0^\dagger + K_1\rho K_1^\dagger$.

Start with the no-jump term $K_0\rho K_0^\dagger$. Multiply out, keeping only first-order terms in $dt$:

$$K_0\rho K_0^\dagger = \left(I - \frac{i}{\hbar}H\,dt - \frac{1}{2}L^\dagger L\,dt\right)\rho\left(I + \frac{i}{\hbar}H\,dt - \frac{1}{2}L^\dagger L\,dt\right)$$

$$= \rho + \frac{i}{\hbar}\rho H\,dt - \frac{1}{2}\rho L^\dagger L\,dt - \frac{i}{\hbar}H\rho\,dt - \frac{1}{2}L^\dagger L\rho\,dt + \mathcal{O}(dt^2)$$

Recognise the commutator $[H, \rho] = H\rho - \rho H$ and the anticommutator $\{L^\dagger L, \rho\} = L^\dagger L\rho + \rho L^\dagger L$:

$$K_0\rho K_0^\dagger = \rho - \frac{i}{\hbar}[H,\,\rho]\,dt - \frac{1}{2}\{L^\dagger L,\,\rho\}\,dt$$

Now the jump term $K_1\rho K_1^\dagger$:

$$K_1\rho K_1^\dagger = \sqrt{dt}\;L\;\rho\;L^\dagger\sqrt{dt} = L\rho L^\dagger\,dt$$

Add both terms together:

$$\rho(t+dt) = \rho - \frac{i}{\hbar}[H,\,\rho]\,dt - \frac{1}{2}\{L^\dagger L,\,\rho\}\,dt + L\rho L^\dagger\,dt$$

Subtract $\rho(t)$ from both sides and divide by $dt$:

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = -\frac{i}{\hbar}[H,\,\rho] + L\rho L^\dagger - \frac{1}{2}\{L^\dagger L,\,\rho\}$$

That is the Lindblad equation for a single jump channel. For multiple independent channels, each with its own Kraus pair $(K_0^{(k)}, K_1^{(k)})$ and jump operator $L_k$, the no-jump Kraus operator picks up one $-\frac{1}{2}L_k^\dagger L_k\,dt$ term per channel:

$$K_0 = I - \frac{i}{\hbar}H\,dt - \frac{1}{2}\sum_k L_k^\dagger L_k\,dt$$

Repeating the same expansion yields the general form:

$$\frac{\mathrm{d}\rho}{\mathrm{d}t} = -\frac{i}{\hbar}[H,\,\rho] + \sum_k\!\left(L_k\rho L_k^\dagger - \frac{1}{2}\{L_k^\dagger L_k,\,\rho\}\right)$$

• Lindblad equation (QuTiP)