18 Aspects of time evolution and density matrices
Now that we know the basics of how to work with density matrices, we’ll explore some interesting phenomena that occur once we study how they evolve in time. In particular, considering how measurement interacts with time evolution will lead us to the quantum Zeno effect, while a brief detour into open quantum systems will let us briefly consider the phenomenon of decoherence.
Before we see these more interesting topics, we first need to understand how time evolution acts on density matrices without any of these complications. If we think about the time evolution of our state vectors, then the density matrix will inherit time dependence as well. Working in the Schrödinger picture, the time-dependent density matrix will be \hat{\rho}(t) = \sum_i p_i \ket{\psi_i(t)} \bra{\psi_i(t)} \\ = \sum_i p_i \hat{U}(t) \ket{\psi_i(0)} \bra{\psi_i(0)} \hat{U}^\dagger(t) \\ = \hat{U}(t) \hat{\rho}(0) \hat{U}^\dagger(t). Taking a time derivative gives us the equation of motion \frac{\partial \hat{\rho}}{\partial t} = \frac{1}{i\hbar} [\hat{H}, \hat{\rho}]. This equation is known as the quantum Liouville-von Neumann equation. Note that this is almost the Heisenberg equation for the density matrix, except that it has the wrong sign (in the Heisenberg equation for \hat{O}, the commutator on the right is [\hat{O}, \hat{H}].) This is because we’re working in Schrödinger picture, and \hat{\rho}(t) is defined directly from our state vectors \ket{\psi(t)}.
18.1 Measurement, revisited
Let’s return one more time to our discussion of measurements to add a couple of more observations. When we’re dealing with mixed states, we have two sources of uncertainty: classical uncertainty in the state, and then the usual quantum uncertainty in the outcome of a given measurement. For the latter, we know that measuring eigenvalue a of observable \hat{A} causes the state to collapse via the corresponding projection operator \hat{P}_a. In fact, this is relying on the fact that we can write the operator \hat{A} itself in terms of its own eigenvalues and projectors as \hat{A} = \sum_a a \hat{P}_a. This is known as the spectral decomposition of \hat{A}.
It turns out that we can use the same projection operator to conveniently handle both sources of uncertainty at once. This isn’t really an update of the measurement postulate, it just follows from what we’ve defined so far:
Given an observable \hat{A}, the probability of measuring outcome a in a system described by density matrix \hat{\rho} is p(a) = {\rm Tr} (\hat{\rho} \hat{P}_a), where \hat{P}_a is the projection operator corresponding to outcome a. After measurement, the system collapses to \hat{\rho} \rightarrow \frac{1}{p(a)} \hat{P}_a \hat{\rho} \hat{P}_a.
As with wavefunctions, the collapse of our density matrix through projection is a violent and non-unitary way for our system to evolve in response to measurements. However, this actually opens up some interesting physical effects for us to look at. For one thing, since most measurements have multiple outcomes, saying what happens to our system after a measurement is complicated by keeping track of the multiple possible states after collapse. But with density matrices, we can deal with the situation where we measure \hat{A} but keep all of the outcomes together in our system; the result is the mixed state
\hat{\rho} \rightarrow \sum_a p(a) \left( \frac{\hat{P}_a \hat{\rho} \hat{P}_a}{p(a)} \right) = \sum_a \hat{P}_a \hat{\rho} \hat{P}_a.
18.2 The Quantum Zeno effect
This brings us to the idea of considering the impacts of measurement itself on our system, not just thinking about the outcomes. If we consider a system with both unitary time evolution and applied measurements, the full time evolution of the system will be a result of both effects.
Suppose that we consider a specific measurement \hat{B} that collapses to an non-degenerate eigenstate \ket{b}, which is not an eigenstate of the Hamiltonian. If we measure once and find \ket{b}, and then let the system evolve in time, it will evolve in some unitary way away from \ket{b} according to its energy-eigenstate decomposition. But now suppose we measure repeatedly and find the result b over and over. After every measurement, the state collapses back to \ket{b} and the time evolution has to start over.
The basic idea of slowing down the time evolution using measurements to “reset” the state is known as the quantum Zeno effect, after Zeno’s paradox (in which a turtle takes successively smaller and smaller steps and is thus never able to reach its destination.) However, realization of the Zeno effect requires the measurement outcome to stay the same every time we measure, so we return to the same \ket{b}. But due to the unitary time evolution not preserving this state, there will be a non-zero chance that the system just jumps to another eigenstate when we measure - not a desirable effect if we’re looking to use the Zeno effect to keep the system in a certain state.
We can get a better understanding of the quantum Zeno effect, one which is closer to how it is actually realized in experiment, by using density matrices instead. Density matrices allow us to work with the system as a whole, no matter the measurement outcomes. We’ll do a simple two-state system example to illustrate what happens. Suppose our two-state system has Hamiltonian \hat{H} = \frac{\hbar \omega}{2} \left( \begin{array}{cc} 1 & 0 \\ 0 & -1 \end{array} \right), and we take the initial state to be spin-up in the x direction, \ket{\psi(0)} = \ket{\rightarrow} = \frac{1}{\sqrt{2}} (\ket{\uparrow} + \ket{\downarrow}). If we let this state evolve in time, we know the answer is \ket{\psi(t)} = \frac{1}{\sqrt{2}} (e^{-i\omega t/2} \ket{\uparrow} + e^{+i\omega t/2} \ket{\downarrow}), so that \left\langle \psi(0) | \psi(t) \right\rangle = \cos \frac{\omega t}{2}. We see that as the system evolves in time, the overlap of \ket{\psi(t)} with the \ket{\rightarrow} state decreases (until things turn around and return to the original state for t > \pi/\omega, at least in an ideal situation.) Can we slow this down with the quantum Zeno effect by making some measurements?
This question will be much easier to answer using density matrices, so let’s switch over to that. The original density matrix is \hat{\rho}(0) = \ket{\psi(0)} \bra{\psi(0)} = \frac{1}{2} \left( \begin{array}{cc} 1 & 1 \\ 1 & 1 \end{array} \right) \\ = \frac{1}{2} \hat{1} + \frac{1}{2} \hat{\sigma}_x, writing a decomposition in terms of Pauli matrices which will be useful to work with: in general, we can write \hat{\rho} = \rho_1 \hat{1} + \rho_x \hat{\sigma}_x + \rho_y \hat{\sigma}_y + \rho_z \hat{\sigma}_z and then \rho_1(0) = \rho_x(0) = 1/2. For the time evolution, we have the Liouville equation \frac{d\hat{\rho}}{dt} = \frac{1}{i\hbar} [\hat{H}, \hat{\rho}] \\ = -\frac{i\omega}{4} [\hat{\sigma}_z, \hat{\rho}] \\ = \frac{\omega}{2} (\rho_x \hat{\sigma}_y - \rho_y \hat{\sigma}_x), from which we read off that \dot{\rho}_1 = \dot{\rho}_z = 0, and the coupled system \dot{\rho}_x = -\tfrac{\omega}{2} \rho_y, \\ \dot{\rho}_y = +\tfrac{\omega}{2} \rho_x. Assuming only unitary evolution and no measurement yet, this gives us a time-dependent density matrix of the form \hat{\rho}_U(t) = \frac{1}{2}\left( \begin{array}{cc} 1 & e^{-i\omega t/2} \\ e^{+i\omega t/2} & 1 \end{array}\right). Instead of looking at the overlap with the initial state, let’s just consider the expectation value of x-direction spin, \left\langle \hat{\sigma}_x \right\rangle. This contains the same information since it is equal to +1 in the initial state and decreases as the system evolves in time. As a check on our density matrix, we have \left\langle \hat{\sigma}_x \right\rangle = {\rm Tr} (\hat{\sigma}_x \hat{\rho}(t)) \\ = {\rm Tr} \left( \begin{array}{cc} 1/2 e^{-i\omega t/2} & 1/2 \\ 1/2 & 1/2 e^{+i\omega t/2} \end{array} \right) = \cos \frac{\omega t}{2}.
Now we add in our measurement. Suppose that we measure \hat{\sigma}_x, but don’t separate things out based on the outcome, keeping our system unchanged. Then from above, the effect of projection will be \hat{\rho} \rightarrow \hat{P}_{\rightarrow} \hat{\rho} \hat{P}_{\rightarrow} + \hat{P}_{\leftarrow} \hat{\rho} \hat{P}_{\leftarrow}, where the projection matrices are \hat{P}_\rightarrow = \ket{\rightarrow} \bra{\rightarrow} = \frac{1}{2} \left(\begin{array}{cc} 1 & 1 \\ 1 & 1 \end{array} \right) \\ = \frac{1}{2} (\hat{1} + \hat{\sigma}_x), \\ \hat{P}_\leftarrow = \ket{\leftarrow} \bra{\leftarrow} = \frac{1}{2} \left(\begin{array}{cc} 1 & -1 \\ -1 & 1 \end{array} \right) \\ = \frac{1}{2} (\hat{1} - \hat{\sigma}_x). If you work this out using some Pauli matrix algebra, you will find that in the parameterization given above, the effect of this projection is to set \rho_y, \rho_z = 0 while leaving the other two components of the density matrix untouched.
Do the algebra and demonstrate that projection due to the unfiltered \hat{\sigma}_x measurement does as described above.
Answer:
In terms of Pauli matrices, we have \hat{\rho} \rightarrow \frac{1}{4} (\hat{1} + \hat{\sigma}_x) \hat{\rho} (\hat{1} + \hat{\sigma}_x) + \frac{1}{4} (\hat{1} - \hat{\sigma}_x) \hat{\rho} (\hat{1} -\hat{\sigma}_x) \\ = \frac{1}{2} \hat{\rho} + \frac{1}{2} \hat{\sigma}_x \hat{\rho} \hat{\sigma}_x, cancelling off the terms with one \hat{\sigma}_x. We can rewrite this further using a commutator as \hat{\rho} \rightarrow \hat{\rho} + \frac{1}{2} \hat{\sigma}_x [\hat{\rho}, \hat{\sigma}_x]. Now we need to find that commutator. Let’s expand out \hat{\rho} in terms of our matrix basis, and just keep the non-vanishing commutators: [\hat{\rho}, \hat{\sigma}_x] = [\rho_y \hat{\sigma}_y + \rho_z \hat{\sigma}_z, \hat{\sigma}_x] \\ = 2i\rho_y \epsilon_{yxz} \hat{\sigma}_z + 2i \rho_z \epsilon_{zxy} \hat{\sigma}_y \\ = 2i (\rho_z \hat{\sigma}_y - \rho_y \hat{\sigma}_z). Finally, we multiply by \hat{\sigma}_x on the left, \hat{\sigma}_x [\hat{\rho}, \hat{\sigma}_x] = 2i (\rho_z \hat{\sigma}_x \hat{\sigma}_y - \rho_y \hat{\sigma}_x \hat{\sigma}_z) \\ = 2i (\rho_z i\epsilon_{xyz} \hat{\sigma}_z - \rho_y i \epsilon_{xzy} \hat{\sigma}_y) \\ = -2 \rho_y \hat{\sigma}_y - 2\rho_z \hat{\sigma}_z. Plugging back in above, this cancels the y and z components off exactly, so we have simply \hat{\rho} \rightarrow \rho_1 \hat{1} + \rho_x \hat{\sigma}_x.
Note that this doesn’t change the expectation value \left\langle \hat{\sigma}_x \right\rangle at all - as we should expect, since without any time evolution repeated measurement of this operator should always give the same result. But with unitary time evolution added, if we track \left\langle \hat{\sigma}_x(t) \right\rangle, measurement won’t cause any sudden jumps, but it will alter the time evolution by zeroing out the \rho_y component.
At this point we need to resort to numerical solution, since the measurements aren’t described by a nice analytic function, they are impulses applied at specific times. But we have all the analytic pieces above needed to implement a numerical solution in Mathematica. Here is the result:
We can clearly see the Zeno effect; as we apply 3, 6, and 12 measurements within the timespan of a single half-period t = \pi / (\omega/2), the evolution becomes slower and slower, so that the system remains closer to the initial state \ket{\rightarrow}.
18.3 Open quantum systems
We have seen a fair bit of formalism for how density matrices can be used to describe bipartite quantum systems \mathcal{H} = \mathcal{H}_A \otimes \mathcal{H}_B, testing for entanglement and so on. Now I’d like to consider a specific bipartite system, which is an open quantum system, \mathcal{H} = \mathcal{H}_S \otimes \mathcal{H}_E. Here \mathcal{H}_S is a quantum system that we’re interested in, and \mathcal{H}_E represents the environment. Any realistic quantum experiment that we want to do always exists in the presence of an environment, so this is a useful setup for thinking about it. Typically we would have d(\mathcal{H}_E) \gg d(\mathcal{H}_S), that is, the environment is an enormous Hilbert space, like a thermal reservoir if we were doing statistical mechanics.
Suppose that we initialize our system \mathcal{H}_S into a pure state \ket{\psi}. If we only had the system itself (i.e. if there was no contact with the environment), then it would evolve ideally according to the Schrödinger equation. However, now we immerse the system in an environment, which means that the initial state is larger, \ket{\Psi(0)} = \ket{\psi}_S \otimes \ket{\phi}_E (we’ll assume that the whole system begins in a pure state, to begin with, and with no entanglement between the system and the environment.) Again, if the Hamiltonian itself is purely bipartite, \hat{H} = \hat{H}_S \otimes \hat{H}_E, then the time evolution of \ket{\psi}_S just follows the Schrödinger equation according to the system Hamiltonian. But this is only true if the environment is perfectly shielded; in general, there will be some off-diagonal parts of the Hamiltonian connecting \mathcal{H}_S and \mathcal{H}_E, which will complicate the time evolution.
Since we are just interested in what happens to the system and want to deal with the environment as a “nuisance”, it makes sense to “integrate out” the environment by partial tracing and just study the reduced density matrix of the system, \hat{\rho}_S(t) = {\rm Tr} _E \ket{\Psi(t)}\bra{\Psi(t)}. If we compute the full density matrix, then it obeys the Liouville-von Neumann equation. But the system density matrix evolves a bit differently. Let’s put in the time evolution of the state: \hat{\rho}_S(t) = {\rm Tr} _E \left( \hat{U}(t) \ket{\Psi(0)} \bra{\Psi(0)} \hat{U}^\dagger(t) \right) \\ = \sum_{i,j} \sum_\mu \ket{s_i} (\bra{s_i} \otimes \bra{e_\mu}) \hat{U}(t) \ket{\Psi(0)} \bra{\Psi(0)} \hat{U}^\dagger(t) (\ket{s_j} \otimes \ket{e_\mu}) \bra{s_j} This is (by construction) always a density matrix over the smaller system space, since it takes the form \sum_{i,j} \ket{s_i}\bra{s_j} (...). But the time evolution operator \hat{U}(t) operates on the full bipartite space. By relying on the assumption that the initial state is bipartite, we can sort of decompose what the time evolution is doing into components that correspond to each environmental basis vector \ket{e_\mu}. To do so, we define the following object: \hat{M}_\mu(t) \equiv \sum_{i,k} \ket{s_i} (\bra{s_i} \otimes \bra{e_\mu}) \hat{U}(t) (\ket{s_k} \otimes \ket{\phi}_E) \bra{s_k}, where \ket{\phi} is the initial state of the environment. The object \hat{M}_\mu(t) is known as a Kraus operator. It takes a bit of algebra, but you can verify from what we have so far that the time-evolved density matrix can be written in terms of the Kraus operators as \hat{\rho}_S(t) = \sum_\mu \hat{M}_\mu(t) \hat{\rho}_S(0) \hat{M}_{\mu}^\dagger(t). As we have already seen, mixtures of pure density operators are generally not pure, so this is an evolution that can take our system from a pure state to a mixed state, in principle. (This is not surprising; the evolution is pure with regards to the bipartite system, which we’ve already seen generally leads to a mixed reduced density matrix for the subsystems.)
The mapping from one density matrix \hat{\rho}_S(0) to another \hat{\rho}_S(t) is known more generally as a quantum map. Quantum maps that satisfy certain nice properties are further known as quantum channels; time evolution happens to satisfy all of these nice properties, so time evolution through the Kraus operators as we’ve derived here is in fact a quantum channel. Quantum channels are of general interest in the context of quantum information, but I won’t say any more about them here. Instead, I want to use the basic framework to highlight a fascinating and very important type of time evolution that appears in a system interacting with its environment.
18.4 Decoherence
Now we’ll do an example, which shows off an effect known as phase-damping; the specific example we consider appears in the lecture notes of David Tong and John Preskill on this topic. The setup is a system consisting of a single qubit, coupled to an environment with three basis states \{\ket{e_0}, \ket{e_1}, \ket{e_2}\}. We assume the time-evolution operator acts in the following way: \hat{U} (\ket{\uparrow}_S \otimes \ket{e_0}) = \ket{\uparrow}_S \otimes (\sqrt{1-p} \ket{e_0} + \sqrt{p} \ket{e_1}), \\ \hat{U} (\ket{\downarrow}_S \otimes \ket{e_0}) = \ket{\downarrow}_S \otimes (\sqrt{1-p} \ket{e_0} + \sqrt{p} \ket{e_2}), \\ In other words, the spin state is not affected at all by the time evolution. Physically, we can imagine that there is a large energy splitting between the system states \ket{\uparrow}_S and \ket{\downarrow}_S, while the interactions with the environment involve very small energy splittings that are much easier to excite. Both references note that one can think of this as a heavy particle interacting with a very light background gas, for example. Note that in this simple example there is no explicit time dependence; that will come with repeated application of the operator.
To proceed, we need the Kraus operators, which are easy to compute: \hat{M}_0(dt) = \sum_{i,k} \ket{s_i} (\bra{s_i} \otimes \bra{e_0}) \hat{U} (\ket{s_i} \otimes \ket{e_0}) \bra{s_k} I’m writing this as an infinitesmal interval dt since we’re imaginging this interaction is only happening once. The operator \hat{U} doesn’t change the spin, so only the diagonal elements of this object are non-zero: \hat{M}_0(dt) = \ket{\uparrow} (\bra{\uparrow} \otimes \bra{e_0}) \ket{\uparrow} \otimes (\sqrt{1-p} \ket{e_0} + \sqrt{p} \ket{e_1}) \bra{\uparrow} \\ + \ket{\downarrow} (\bra{\downarrow} \otimes \bra{e_0}) \ket{\downarrow} \otimes (\sqrt{1-p} \ket{e_0} + \sqrt{p} \ket{e_2}) \bra{\downarrow} \\ = \ket{\uparrow}\bra{\uparrow} \sqrt{1-p} + \ket{\downarrow}\bra{\downarrow} \sqrt{1-p} \\ = \sqrt{1-p} \hat{1}. It’s straightforward to see similarly that for the other two operators, only one of the four terms in the sum survives, and we have \hat{M}_1(dt) = \sqrt{p} \ket{\uparrow}\bra{\uparrow} = \left(\begin{array}{cc} \sqrt{p} & 0 \\ 0 & 0 \end{array} \right), \\ \hat{M}_2(dt) = \sqrt{p} \ket{\downarrow}\bra{\downarrow} = \left(\begin{array}{cc} 0 & 0 \\ 0 & \sqrt{p} \end{array} \right). It’s helpful to rewrite these in terms of Pauli matrices as \hat{M}_1(dt) = \frac{\sqrt{p}}{2} (\hat{1} + \hat{\sigma}_z), \\ \hat{M}_2(dt) = \frac{\sqrt{p}}{2} (\hat{1} - \hat{\sigma}_z). Continuing to the time-evolved density matrix, we thus have \hat{\rho}_S(dt) = \sum_\mu \hat{M}_\mu(t) \hat{\rho}_S(0) \hat{M}_\mu^\dagger(t) \\ = (1-p) \hat{\rho}_S(0) + \frac{p}{2} (\hat{\rho}_S(0) + \hat{\sigma}_z \hat{\rho}_S(0) \hat{\sigma}_z) \\ = (1 - \frac{p}{2}) \hat{\rho}_S(0) + \frac{p}{2} \hat{\sigma}_z \hat{\rho}_S(0) \hat{\sigma}_z Writing out the second term, letting the original density matrix be completely arbitrary, \hat{\sigma}_z \hat{\rho}_S(0) \hat{\sigma}_z = \left( \begin{array}{cc} 1 & 0 \\ 0 & -1 \end{array} \right) \left( \begin{array}{cc} \rho_{++} & \rho_{+-} \\ \rho_{-+} & \rho_{--} \end{array} \right) \left( \begin{array}{cc} 1 & 0 \\ 0 & -1 \end{array} \right) \\ = \left( \begin{array}{cc} \rho_{++} & -\rho_{+-} \\ -\rho_{-+} & \rho_{--} \end{array} \right) So the p dependence cancels completely on the diagonal, and we find \hat{\rho}_S(dt) = \left( \begin{array}{cc} \rho_{++} & (1-p) \rho_{+-} \\ (1-p) \rho_{-+} & \rho_{--} \end{array} \right)
So the effect of the interaction with the environment on the system is: diagonal elements of the density matrix are unchanged, but off-diagonal matrix elements are suppressed. To really see the time-dependence appearing, we have to apply the map repeatedly, representing continuous interaction with the environment. (I haven’t gone deep enough into the formalism to argue why this makes sense, but it does - see the referenced lecture notes above for more detail.) We define the probability of scattering per unit time to be \Gamma, so that if we wait for time dt, then p = \Gamma dt. Then we wait for total time t = N dt. With N interactions in total, the density matrix becomes \hat{\rho}_S(t) = \left( \begin{array}{cc} \rho_{++} & (1-p)^N \rho_{+-} \\ (1-p)^N \rho_{-+} & \rho_{--} \end{array} \right) = \left( \begin{array}{cc} \rho_{++} & e^{-\Gamma t} \rho_{+-} \\ e^{-\Gamma t} \rho_{-+} & \rho_{--} \end{array} \right) using (1-p)^N = (1-\frac{\Gamma t}{N})^N \approx e^{-\Gamma t}. So the off-diagonal elements decay exponentially as we wait. To give a better picture for what this means, take the initial state in the system to be \ket{\psi} = \alpha \ket{\uparrow}_S + \beta \ket{\downarrow}_S. Then, \hat{\rho}_S(t) = \left( \begin{array}{cc} |\alpha|^2 & e^{-\Gamma t} \alpha^\star \beta \\ e^{-\Gamma t} \beta^\star \alpha & |\beta|^2 \end{array} \right)
So the off-diagonal components that are being suppressed are precisely the components that allow our state to exist in a quantum superposition of spins. As time goes on, we’re left with either a pure state which is only \ket{\uparrow} or \ket{\downarrow}, or a classical mixture of the two.
The phenomenon that we’re seeing a glimpse of here is known as decoherence. It provides an explanation as to why quantum features such as superposition don’t survive up to the classical, macroscopic world; the larger the object is, the more it interacts with the environment around it, which causes any superposition present to rapidly decay away. This gives an explanation as to why Schrödinger’s cat will always be a thought experiment and not a real experiment; a cat (or any other macroscopic object) will simply decohere too rapidly to exist in a superposition state for any length of time, even if we were able to prepare one.