Quantum mechanics often deals with the description of single state vectors, but when dealing with ensembles of quantum systems it becomes more convenient to employ the density operator formalism described in this section.
The density operator
The conventional approach of using a state-vector description of quantum mechanics is slightly different from the use of a density operator or density matrix, \(\require{newcommand} \newcommand{\vc}[1]{\mathbf{#1}} \newcommand{\desc}[1]{\begin{quote}\emph{#1}\end{quote}} \newcommand{\p}{\partial} \newcommand{\bra}[1]{\langle #1 |\,} \newcommand{\ket}[1]{\,| #1 \rangle} \newcommand{\braket}[2]{\langle #1 | #2 \rangle} \newcommand{\matel}[3]{\langle #1 | #2 | #3 \rangle} \newcommand{\mexp}[1]{\langle #1 \rangle} \newcommand{\kb}[2]{| #1 \rangle\langle #2 |} \DeclareMathOperator{\Tr}{Tr} \rho = \ket{\psi}\bra{\psi}\), as \(\rho\) describes an ensemble of spin systems rather than a single instance of the system. For a system in thermal equilibrium, the density operator is diagonal and can be written in the form:
\begin{equation} \rho = \sum_i p_i \ket{i}\bra{i}, \label{eq_density_operator} \end{equation}where \(p_i\) is the probability that a spin system from the ensemble occupies the state \(\ket{i}\). An ensemble with \(p_j = 1\) and \(p_i = 0\) for all \(i \neq j\) is called a pure ensemble; otherwise the ensemble is mixed. Note that entangled quantum states (superpositions) is an entirely different and unrelated concept to that of a mixed ensemble - but entangled states are in fact considered non-equilibrium states, and in order to handle such states a more general definition of the density matrix is required:
\begin{equation} \rho = \sum_{i,j} p_{ij} \ket{i}\bra{j} \, . \label{eq_density_operator_noneq} \end{equation}This is the most general definition of \(\rho\), describing an ensemble that is out of equilibrium. The off-diagonal matrix elements, \(p_{ij}\), are often referred to as coherences as they describe the coherence or entanglement of the various states in the ensemble. When the quantum system is coupled to a thermal reservoir, all the off-diagonal elements of \(\rho\) will vanish over time in a process called decoherence, which describes the relaxation of the system towards thermal equilibrium and the form in Eq. (\ref{eq_density_operator}). Note that the density operator is hermitian and a basis that diagonalizes the operator always exists - the concept of coherences therefore depends on the chosen basis; usually coherences are discussed relative to the basis that diagonalizes the density operator at thermal equilibrium (i.e. enforced by the environment), and an ensemble with coherences with respect to that basis, therefore, is often referred to as a non-equilibrium ensemble.
Since the density operator always may be written in the form of Eq. (\ref{eq_density_operator}) it is evident that its eigenvalues are probabilities. This also means that the trace of the density operator must be 1, since the probabilities should sum up to unity. When measuring the state of a spin system instance from the ensemble, the eigenvalues provides the probabilities for measuring the various states. Note that performing a measurement does not turn the ensemble into a pure ensemble, as a mixed ensemble is a statistical mixture of pure quantum states. If the measurement operator is not diagonal in the basis that diagonalizes \(\rho\), the measurement will make the state of a spin system instance from the ensemble collapse to an eigenstate of the measurement operator. While a measurement affects the composition of the ensemble it does not, in general, make the ensemble pure. Whether an ensemble is pure can be checked by evaluating \(\Tr(\rho^2)\):
\begin{equation} \Tr(\rho^2) = \sum_i p_i^2 \leq 1, \end{equation}where equality holds if and only if the ensemble is pure, according to the definition of a pure ensemble.
The expectation value of an observable, \(\vc{A}\), is given by \(\matel{\psi}{\vc{A}}{\psi}\) for a system in the quantum state \(\ket{\psi}\). If an ensemble of systems is considered instead, an ensemble average has to be made:
\begin{equation} \overline{\mexp{\vc{A}}} = \sum_i p_i \matel{i}{\vc{A}}{i} = \sum_i \sum_n p_i \matel{i}{\vc{A}}{n}\braket{n}{i} = \sum_n \matel{n}{\vc{A}\rho}{n} = \Tr(\vc{A}\rho) \, . \end{equation}Here \(\overline{\mexp{\vc{A}}}\) denotes the ensemble averaged expectation value.
The trace of the density operator is only ensured to be equal to 1 for closed quantum systems, but it is often necessary to deal with open quantum systems; assume for example that we consider a radical pair embedded in a molecular system, say a protein. In this case the spins of the radical pair will be coupled to external degrees of freedom, and the density matrix describing the physical system can be defined as \(\rho_{AB} = \rho_A \otimes \rho_B\), where \(\rho_A\) describes the spin system itself - the part of the system we are interested in, in terms of spin dynamics - and \(\rho_B\) describes all the external degrees of freedom. The density operator \(\rho_{AB}\) retains a trace of 1, but the trace of \(\rho_A\) or \(\rho_B\) may be different. In the following sections the partial trace over the environmental degrees of freedom, \(\rho \equiv \Tr_B(\rho_{AB}) = \rho_A \cdot \Tr(\rho_B)\), will be referred to as the density operator for the system, since we do not need to concern ourselves with the environment, \(\rho_B\), or the combined (total) density operator, \(\rho_{AB}\), unless we discuss decoherence processes, i.e. spin relaxation. Thus the trace of the density operator is allowed to be less than 1 in the remainder of the article series, as the spin system may participate in processes involving the environment.
Dynamics: The Liouville-von Neumann equation
The time-evolution of the density operator follows directly from the Schrödinger equation:
\begin{equation} \frac{d\rho}{dt} = \sum_j\left(-\frac{ip_j}{\hbar}\left[\vc{H},\ket{j}\bra{j}\right] + \frac{dp_j}{dt} \ket{j}\bra{j}\right), \label{eq_preLvN} \end{equation}where it is usually assumed that \(dp_i = 0\). This leads to the equation known as the Liouville-von Neumann equation:
\begin{equation} \frac{d\rho}{dt} = -\frac{i}{\hbar}[\vc{H},\rho] \, . \label{eq_LiouvilleVonNeumann} \end{equation}The solution to the Liouville-von Neumann equation also follow from the solution to the Schrödinger equation, and if the Hamiltonian is time-independent the time-evolution of the density operator is, therefore, given by:
\begin{equation} \rho(t) = \exp\left(-\frac{i}{\hbar}\vc{H}t\right) \rho(0) \exp\left(\frac{i}{\hbar}\vc{H}t\right) \, . \label{eq_LvN_solution} \end{equation}The situation is more complicated when the spin Hamiltonian is time-dependent, and time-dependent interactions are treated in another article. Also note that additional reaction terms may be added to Eq. (\ref{eq_LiouvilleVonNeumann}) as described in section article about reaction operators.
Note that the last term in Eq. (\ref{eq_preLvN}) is the only way to change composition of the ensemble mixture, and while it is assumed to vanish in order to obtain the Liouville-von Neumann equation, Eq. (\ref{eq_LiouvilleVonNeumann}), it is often reintroduced in the form of reaction operators or spin relaxation operators.
State vector versus density operator formalism
The state vector and density operator formalisms describe the same physics, but the density operator formalism is superior when a mixed ensemble of systems is considered. As an example, consider a radical pair with a single spin-\(\frac{1}{2}\) nucleus on one of the radicals, and assume that the radical pair is created in the singlet state of the unpaired electronic spins. Since nothing is known about the nuclear spin, it must be assumed that a mixed ensemble of radical pairs is obtained with 50% radical pairs created in the \(\ket{S\alpha}\) state and 50% in the \(\ket{S\beta}\) state. The \(\alpha\) and \(\beta\) refers to the spin up and the spin down state of the nuclear spin, respectively. Note that the mixture would not be exactly 50\%/50\% at thermal equilibrium in the presence of an external magnetic field; a difference in the population of the states is crucial for NMR spectroscopy. If a property such as the expectation value of the electronic \(\vc{S}_z\) operator has to be measured for the initial states, it would be done as:
\begin{equation} \mexp{\vc{S}_z(t=0)} = \underbrace{0.5\matel{S\alpha}{\vc{S}_z}{S\alpha} + 0.5\matel{S\beta}{\vc{S}_z}{S\beta}}_\text{state vector formalism} = \underbrace{\Tr(\vc{S}_z\rho(t=0))}_\text{density operator formalism} \, . \end{equation}The density operator formalism already looks simpler, but then consider calculating the expectation value at a later time instance, \(\mexp{\vc{S}_z(t=T)}\); this requires the application of the Schrödinger equation twice in order to propagate both \(\ket{S\alpha}\) and \(\ket{S\beta}\) in time, while the density operator formalism only needs a single application of the Liouville-von Neumann equation, as the Liouville-von Neumann equation propagates the entire ensemble in one go. This is a big deal if you have a large system with, for example, 8 spin-\(\frac{1}{2}\) nuclei in the radical pair; in this case the ensemble would be a mixture of not just 2 but \(2^8 = 256\) different singlet states, which means that the state vector approach would require solving the Schrödinger equation 256 times while the density operator approach still only needs a single use of the Liouville-von Neumann equation.
