View source for Open Quantum Systems/Projection Operator Techniques

{{lecture}}

==== Decay into two bands ====

We consider a two-level systems coupled to an environment consisting of two distinct energy bands. The system is described by the Hamiltonian
:<math>H = H_S + H_E + H_I.</math>
The system Hamiltonian <math>H_S</math> is an ordinary two-level system with an energy difference of <math>\Delta E</math>, according to
:<math>H_S = \Delta E |e\rangle\langle e| = \Delta E\sigma_+\sigma_-,</math>
where we have set the energy of the ground state <math>|g\rangle</math> to zero. For the energy bands of the environment, we assume a linear spectrum, i.e,
:<math>H_E = \sum\limits_{n_1=1}^{N_1} \frac{\delta\varepsilon}{N_1} n_1 |n_1\rangle \langle n_1| + \sum\limits_{n_2=1}^{N_2} \left(\Delta E +\frac{\delta\varepsilon}{N_2} n_2\right) |n_2\rangle \langle n_2|,</math>
where <math>N_1</math> and <math>N_2</math> are the total number of states in each band. The interaction Hamiltonian is given by
:<math>H_I = \sum\limits_{n_1=1}^{N_1}\sum\limits_{n_2=1}^{N_2} \lambda_{n_1,n_2} \big(\sigma_+ \otimes 1_E \big) (1_S \otimes |n_1\rangle\langle n_2|) + \text{H.c.},</math>
i.e., whenever the two-level system is de-excited, an excitation is created in the environment, and vice versa. The matrix elements <math>\lambda_{n1,n2}</math> are random Gaussianly distributed variables with a variance <math>\lambda^2</math>. In the following, we assume a separation of energy scales according to <math>\Delta E \gg \delta E \gg \lambda</math>.

This model has some very interesting properties. First of all, the size of the environment is clearly finite, meaning that excitations in the environment may not decay fast enough and the Markov approximation may not be applicable. This is even further strengthened when taking into account that the interaction will lead to very strong correlations between system and environment. At the same time, the weak interaction naturally gives rise to a peturbative expansion.

We proceed by defining a correlated projection operator for the relevant part of the dynamics <ref name="Breuer2007">{{cite book|last1=Breuer |first1=H.-P. |last2 = Petruccione | first2= F. |title=The Theory of Open Quantum Systems |publisher=Oxford University Press |year=2002}}</ref>, according to
:<math>\mathcal{P}\rho = \text{Tr}\left\{\Pi_1\rho\right\}\frac{1}{N_1}\Pi_1 + \text{Tr}\left\{\Pi_2\rho\right\}\frac{1}{N_2}\Pi_2 = P_e \frac{1}{N_1}\Pi_1 + P_g \frac{1}{N_2}\Pi_2,</math>
where we have introducted the projectors
:<math>\begin{align}\Pi_1 &= \sum\limits_{n_1=1}^{N_1} 1_S \otimes |n_1 \rangle\langle n_1|\\
\Pi_2 &= \sum\limits_{n_2=1}^{N_2} 1_S \otimes |n_2 \rangle\langle n_2|\end{align}</math>
and <math>P_e</math> and <math>P_g</math> denote the probability to find the two-level system in its excited or ground state, respectively. Because of conservation of probability, we have <math>P_g = 1-P_e</math>, meaning that the relevant part of the dynamics consists of a single variable!

In the following, we will derive its equation of motion. The second order TCL master equation is given by
:<math>\mathcal{P}\frac{d}{dt}\rho = \int\limits_0^t ds \mathcal{P}\mathcal{L}(t)\mathcal{L}(s)\mathcal{P}\rho,</math>
which in our case can be written as
:<math>\dot{P}_e \frac{1}{N_1}\Pi_1 + (1-\dot{P}_e) \frac{1}{N_2}\Pi_2 = \int\limits_0^t ds \text{Tr}\left\{\Pi_1\mathcal{L}(t)\mathcal{L}(s)\mathcal{P}\rho\right\}\frac{1}{N_1}\Pi_1 + \text{Tr}\left\{\Pi_2\mathcal{L}(t)\mathcal{L}(s)\mathcal{P}\rho\right\}\frac{1}{N_2}\Pi_2.</math>
since the projection and the time derivative commute. Multiplying with <math>\Pi_1</math> and taking the trace yields
:<math>\begin{align}\dot{P}_e &= \int\limits_0^t ds\text{Tr}\left\{\Pi_1\mathcal{L}(t)\mathcal{L}(s)\mathcal{P}\rho\right\}\\
&=-\int\limits_0^t ds\text{Tr}\left\{\Pi_1\left[H_I(t),\left[H_I(s),P_e \frac{1}{N_1}\Pi_1 + (1-P_e) \frac{1}{N_2}\Pi_2\right]\right]\right\}.\end{align}</math>
Expanding the commutator and making use of the relation <math>\Pi_i\Pi_j=\delta_{ij}\Pi_i</math> results in
:<math>\dot{P}_e = -\int\limits_0^t ds\text{Tr}\left\{P_e \frac{1}{N_1} \Pi_1 H_I(s)H_I(t) - (1-P_e)\frac{1}{N_2}\Pi_1 H_I(t)H_I(s) + \text{H.c.}\right\}.</math>
We can now introduce two relaxation rates <math>\gamma_{1,2}</math> describing the transition rates between the two bands. Specifically, we have
:<math>\gamma_i = \frac{1}{N_i}\int\limits_0^t ds\sum\limits_{n1}\langle n_1| H_I(s)H_I(t)+\text{H.c} |n_1\rangle.</math>
As the interaction Hamiltonian couples does not contain terms that leave the band index unchanged, we may write the matrix elements as
:<math>\langle n_1| H_I(s)H_I(t)+\text{H.c} |n_1\rangle = \sum\limits_{n_2} \langle n_1| H_I(s)|n_2\rangle\langle n_2|H_I(t)|n_1\rangle + \text{c.c.}.</math>
As we work in the interaction picture, we have
:<math>\langle n_1| H_I(s)|n_2\rangle\langle n_2|H_I(t)|n_1\rangle = \exp(i\omega_{12}s)\exp(-i\omega_{12}t) |\lambda_{n1,n2}|^2,</math>
where we have introduced the transitions frequencies <math>\omega_{12} = \delta\varepsilon(n_1/N_1-n_2/N_2)</math>.
Using the statistical properties of the interaction Hamiltonian and performing the integration over <math>s</math> finally yields <ref name="Weimer2008"></ref>
:<math>\gamma_i = \frac{2\lambda^2}{N_i}\sum\limits_{n1,n2} \frac{\sin(\omega_{12}t)}{\omega_{12}}.</math>
For the next step, we note that the sinc function appearing in the relaxation rates is a representation of the Dirac <math>\delta</math> distribution, i.e.,
:<math>\pi\delta_t(\omega_{12}) = \lim_{t\to\infty}\frac{\sin(\omega_{12}t)}{\omega_{12}}.</math>
Then, at large enough times and small enough couplings <ref name="Schwabl2010"></ref>, we may approximate the rates by
:<math>\begin{align}\gamma_1 &\approx \frac{2\pi\lambda^2}{N_1}\sum_{k,l}\delta(E_{n_1}-E_{n_2}) = 2\pi\lambda^2 \frac{N_2}{\delta\varepsilon}\\\gamma_2 &\approx 2\pi\lambda^2 \frac{N_1}{\delta\varepsilon}.\end{align}</math>
The equation of motion for the excited state probability then reads
:<math>\dot{P}_e = -(\gamma_2+\gamma_1) P_e + \gamma_2.</math>
Assuming the two-level system is initially in its excited state, the stationary state given by <math>P_e^s = \gamma_2/(\gamma_1+\gamma_2)</math> is approached in an exponential decay according to
:<math> P_e(t) = P_e^s + \frac{\gamma_1}{\gamma_1+\gamma_2}\exp[-2(\gamma_1+\gamma_2)t].</math>
This solution is in excellent agreement with numerical simulations of the full Schrödinger equation <ref name="Breuer2006"></ref>.

It is important to stress that this efficient description of such a complex system is only possible because our projection operator correctly captures the essential features of the dynamics. If one choses a different projection operator, say, one without correlations between system and environments, according to
:<math>\mathcal{P}\rho = \text{Tr}_E\{\rho\}\otimes\rho_E,</math>
then the second order TCL master equation no longer correctly describes the dynamics, as higher orders are divergent <ref name="Breuer2006" />.