When simulating the dynamics of open quantum systems, one often assumes a pure initial state of the system part. The validity of the assumption, of course, depends on the experimental realization. In most of the excitation energy transfer experiments, a pure initial state is a faithful representation of the system because the coherent light that excites molecules makes the wave function of the system collapse to a certain state, hence factorizing the system-environment density matrix. The factorized density matrix is $\rho(t=0)=\rho_\mathrm{sys}(0)\otimes\rho_\mathrm{env}(0)$ and the density matrix of the system at time zero (initial state)$\rho_\mathrm{sys}(t=0)$ is a pure-state in this case.

The lights in nature, however, are not always coherent. For example, the sun’s light, the energy source of light-harvesting activities in plants, is not coherent. Incoherent lights don’t always generate a pure state [1]. Therefore It is relevant to study the relaxation dynamics of the system starting from a mixed initial state. A very useful formalism called thermo-field dynamics can handle mixed states by mapping them to pure states. Once a mixed state is mapped to a pure state, all the conventional numerical methods designed for the Schroedinger equation can be applied to simulate the relaxation dynamics. In the following, thermo-field formalism is introduced.

Consider an arbitrary time-dependent density matrix $\rho(t)$. The basis of the density matrix is eigenstates $\ket{n}$ of the system Hamiltonian $H$. The expectation value of an observable $A$ can be obtained from the trace formula $\mathrm{Tr}(\rho A)$. There exists a composite pure state $$\ket{\psi_\rho(t)}=\rho^{1/2}\otimes \mathrm{I}\cdot \sum_{n,n'}\ket{n}\otimes \ket{n'}$$ that can produce the same expectation value by the normal average formula $\braket{\psi_\rho(t)|A|\psi_\rho(t)}$. Here $\hat{\mathrm{I}}$ is the identity operator in the space spanned by $\{\ket{n'}\}$ and $\ket{n'}$ is a copy of state $\ket{n}$.

What’s left is to study the time evolution of the state $\ket{\psi_\rho(t)}=\rho^{1/2}\otimes \mathrm{I}\cdot \sum_{n=n'}\ket{n}\otimes \ket{n'}$. We first note the dynamics of the density matrix are governed by the von Neumann equation $$i\hbar \frac{\partial \rho}{\partial t} = H\rho - \rho H.$$ Then we directly differentiate the state $\ket{\psi_\rho(t)}$ with respect to the time $t$ and put the von Neumann equation in the result. We have [2]

$$ i\hbar \frac{\partial \ket{\psi_\rho(t)}}{\partial t} = (H\rho^{1/2} - \rho^{1/2} H)\sum_{n=n'}\ket{n}\otimes \ket{n'} .$$

Expanding the equation above, the first term $H\rho^{1/2}\sum_{n=n'}\ket{n}\otimes \ket{n'}$ is equal to $H \ket{\psi_\rho(t)}$ by definition. The second term is $-\rho^{1/2} H\sum_{n,n'}\ket{n}\otimes \ket{n'}=- \rho^{1/2} \sum_{n,n'}E_n\ket{n}\otimes \ket{n'}$. It can be rewritten in new way. Noting the fact the $\ket{n'}$ is a copy of $\ket{n}$, one can replace $E_n$ by $E_{n'}$ without affecting the result. Thus,
\begin{align} -\rho^{1/2} H\sum_{n=n'}\ket{n}\otimes \ket{n'}=&- \sum_{n=n'}E_n(\rho^{1/2}\ket{n})\otimes \ket{n'}\\ =&- \sum_{n=n'}E_{n'}(\rho^{1/2} \ket{n})\otimes \ket{n'}\\ =&- \sum_{n=n'}(\rho^{1/2} \ket{n})\otimes (E_{n'}\ket{n'})\\ =&- \sum_{n=n'}(\rho^{1/2} \ket{n})\otimes (\tilde{H}\ket{n'})\\ =&- \tilde{H}\rho^{1/2}\sum_{n=n'} \ket{n}\otimes \ket{n'}\\ =&- \tilde{H}\ket{\psi_\rho(t)} \end{align}
where $\tilde{H}$ is a copy of $H$ and it operates on the duplicate Hilbert spanned by $\{\ket{n'}\}$.

Summarizing all the results above, we can conclude that the pure state $\ket{\psi_\rho(t)}$, which is a pure counterpart of the midex state $\rho(t)$, evolves according to the new Schroedinger equation,
\begin{align} i\hbar \frac{\partial \ket{\psi_\rho(t)}}{\partial t} = (H-\tilde{H}) \ket{\psi_\rho(t)} \end{align}
in the composite Hilbert space. Solving this equation can provide us with the necessary information to obtain the expectation values.

Although the density matrix $\rho$ is mixed, the formalism we introduced above can map it to a pure state in a doubled Hilbert space and avoid working with density matrices.

References

[1] Brumer, P. (2018). Shedding (incoherent) light on quantum effects in light-induced biological processes. The journal of physical chemistry letters, 9(11), 2946-2955.
[2] Suzuki, M. (1991). Emerging broken symmetry in space and time. In Thermal Field Theories (pp. 5-15). Elsevier.