Doublon lifetimes in dissipative environments

TL;DR

This paper investigates how different types of environmental noise affect the lifetime of doublons in a two-electron Hubbard model, revealing that charge noise can stabilize doublons while current noise destabilizes them, with implications for quantum dot systems.

Contribution

It provides a comparative analysis of charge and current noise effects on doublon decay, combining numerical and analytical methods for the first time.

Findings

Doublon lifetime increases with interaction under charge noise.

Doublons become less stable with larger interaction under current noise.

Predicted doublon lifetimes can reach up to 50 ns in quantum dot setups.

Abstract

We study the dissipative decay of states with a doubly occupied site in a two-electron Hubbard model, known as doublons. For the environment we consider charge and current noise which are modelled as a bosonic heat bath that couples to the onsite energies and the tunnel couplings, respectively. It turns out that the dissipative decay depends qualitatively on the type of environment as for charge noise, the life time grows with the electron-electron interaction. For current noise, by contrast, doublons become increasingly unstable with larger interaction. Numerical studies within a Bloch-Redfield approach are complemented by analytical estimates for the decay rates. For typical quantum dot parameters, we predict that the doublon life times up to 50 ns.

Equations

$$H_S$$

$$\equiv -JT+UD,$$

$$\frac{1}{\sqrt{2}}(c_{i\uparrow }^{†}{c}_{j\downarrow }^{†}-c_{i\downarrow }^{†}{c}_{j\uparrow }^{†})|0⟩,1\le i<j\le N,$$

$$c_{j\uparrow }^{†}{c}_{j\downarrow }^{†}|0⟩,j=1,\cdots ,N.$$

$$P_1=P_D-\frac{J}{U}(T^{+}+T^{-})+\mathcal{O}\left(\frac{J^2}{U^2}\right),$$

$$\dot{\rho }$$

$$\equiv -i[H_S,\rho ]+\mathcal{L}[\rho ].$$

$$Q_j$$

$$R_j$$

$$H_{SB}^Q=\sum _{j,\sigma }{n}_{j,\sigma }{\xi }_j,$$

$$\frac{⟨D{⟩}_{T_1}-⟨D{⟩}_{\infty }}{1-⟨D{⟩}_{\infty }}=\frac{1}{e}.$$

$$|\mathrm{tr}\left(P_1\rho \right)-\mathrm{tr}\left(D\rho \right)|\le$$

$$\simeq$$

$$⟨P_1⟩\simeq \Delta e^{-\Gamma t}+⟨P_1{⟩}_{\infty },$$

$$\Gamma \simeq -\frac{1}{\Delta }{\frac{d⟨P_1⟩}{dt}|}_{t=0}=-\frac{\mathrm{tr}\left(P_1\mathcal{L}[{\rho }_0]\right)}{⟨P_1{⟩}_0-⟨P_1{⟩}_{\infty }}.$$

$$\begin{gathered} \overline{\Gamma }=\frac{1}{N\Delta }\sum _j\mathrm{tr}\left(P_D[Q_j,[X_j,P_1]]\right) \\ -\mathrm{tr}\left(P_D\{R_j,[X_j,P_1]\}\right). \end{gathered}$$

$$Q_j$$

$$R_j$$

$${\overline{\Gamma }}_{\mathrm{HT}}\simeq \frac{4\pi \alpha J^2}{U^2\Delta }\left(2k_{B}T+U\right),$$

$${\overline{\Gamma }}_{\text{LT}}\simeq \frac{8\pi \alpha J^2}{U\Delta },$$

$$H_{SB}^I=\sum _{j,\sigma }(c_{j+1\sigma }^{†}{c}_{j\sigma }+c_{j\sigma }^{†}{c}_{j+1\sigma }){\xi }_j.$$

$$\Gamma \simeq -{\frac{d⟨P_1⟩}{dt}|}_{t=0}=-\mathrm{tr}\left(P_1\mathcal{L}[{\rho }_0]\right).$$

$${\overline{\Gamma }}_{\mathrm{HT}}$$

$${\overline{\Gamma }}_{\text{LT}}$$

$${\dot{\rho }}_{\alpha \beta }=-i({ϵ}_{\alpha }-{ϵ}_{\beta }){\rho }_{\alpha \beta }+\sum _{{\alpha }^{\prime }{\beta }^{\prime }}{\mathcal{L}}_{\alpha \beta ,{\alpha }^{\prime }{\beta }^{\prime }}{\rho }_{{\alpha }^{\prime }{\beta }^{\prime }}.$$

$$\begin{split}\mathcal{L}_{\alpha\beta,\alpha^{\prime}\beta^{\prime}}={}&\sum_{j}\Big{[}(\Gamma_{\beta^{\prime}\beta}+\Gamma_{\alpha^{\prime}\alpha})X^{(j)}_{\alpha\alpha^{\prime}}X^{(j)}_{\beta^{\prime}\beta}\\ &-\delta_{\beta\beta^{\prime}}\sum_{\beta^{\prime\prime}}\Gamma_{\alpha^{\prime}\beta^{\prime\prime}}X^{(j)}_{\alpha\beta^{\prime\prime}}X^{(j)}_{\beta^{\prime\prime}\alpha^{\prime}}\\ &-\delta_{\alpha\alpha^{\prime}}\sum_{\alpha^{\prime\prime}}\Gamma_{\beta^{\prime}\alpha^{\prime\prime}}X^{(j)}_{\beta^{\prime}\alpha^{\prime\prime}}X^{(j)}_{\alpha^{\prime\prime}\beta}\Big{]}\,,\end{split}$$

$$\Gamma (\omega )=\{\begin{array}{ll}J(\omega )(1+n_B(\omega ))& \omega >0\\ J(-\omega )n_B(-\omega )& \omega <0\end{array},$$

$$\begin{gathered} \overline{\Gamma }=\frac{1}{N\Delta }\sum _j\mathrm{tr}\left(P_1[Q_j,[X_j,P_1]]\right) \\ -\mathrm{tr}\left(P_1\{R_j,[X_j,P_1]\}\right). \end{gathered}$$

$$\overline{\Gamma}=-\frac{1}{\Delta}\frac{d\overline{\langle P_{1}\rangle}}{dt}\Big{|}_{t=0}=-\frac{1}{\Delta}\overline{\mathrm{tr}\left(P_{1}\mathcal{L}[\rho]\right)}\ ,$$

$$\overline{\Gamma }\Delta =-\frac{1}{N}\sum _{\alpha ,\beta \in I_1}{\mathcal{L}}_{\alpha \alpha ,\beta \beta }=\frac{1}{N}\sum _{\alpha \in I_0}\sum _{\beta \in I_1}{\mathcal{L}}_{\alpha \alpha ,\beta \beta },$$

$$P(c_1,\cdots ,c_N)=\frac{(N-1)!}{{\pi }^N}\delta \left(1-r^2\right),$$

$$\overline{c_n{c}_m^{\ast }}$$

$$\overline{c_n{c}_m^{\ast }{c}_{n^{\prime }}{c}_{m^{\prime }}^{\ast }}$$

$$\overline{\mathrm{tr}\left(\rho A\right)}$$

Full text

I Introduction

In recent years experiments with strongly-interacting cold atomic gases have attracted much attention.Bloch et al. (2008) A particular advantage of these systems is that their parameters can be controlled to a high degree either directly or via oscillating forces that lead to synthetic gauge fields.Creffield and Sols (2011); Dalibard et al. (2011) This allows a flexible engineering and simulation of many-body Hamiltonians. For a theoretical description, one frequently employs the Hubbard model. Despite its seeming simplicity, it captures a great variety of condensed-matter phenomena ranging from metallic behavior to insulators, magnetism, and superconductivity.

In the strongly interacting limit of the Hubbard model, particles occupying the same lattice site can bind together, even for repulsive interactions. This occurs when the onsite interaction is much larger than the tunneling such that energy conservation inhibits the decay into a state with two distant particles. In principle, both bosons Valiente and Petrosyan (2008); Compagno et al. (2017) and fermions Essler et al. (2005) can form such $N$-particle states. While the former allow any occupation number, for fermions with spin $s$, the occupation of one site is restricted to at most $2s+1$ particles. In particular, two spin-1/2 fermions may reside in a singlet spin configuration on one lattice site and, thus, form a doublon. Over the last years, they have been investigated both theoretically Hofmann and Potthoff (2012); Bello et al. (2016, 2017) and experimentally Winkler et al. (2006); Folling et al. (2007); Strohmaier et al. (2010); Preiss et al. (2015) with cold atoms in optical lattices.

In the context of solid-state based quantum information and quantum technologies, arrays of tunnel coupled quantum dots represent a recent platform for similar experiments with electrons.Puddy et al. (2015); Zajac et al. (2016) In comparison to optical lattices, however, these systems are way more sensitive to decoherence and dissipation stemming from the interaction with environmental degrees of freedom such as phonons or charge and current noise. Since environments may absorb energy, the separation of two electrons in a doublon state is no longer energetically forbidden. In this paper we cast some light on this issue by studying the life times of doublons in a one-dimensional lattice in the presence of charge and current noise, as is sketched in Fig. 1. For the environment we employ a Caldeira-Leggett model, Leggett et al. (1987); Hänggi et al. (1990) where depending on the type of noise, the bath couples locally to the onsite energies or to the tunnel matrix elements.

In Sec. II, we specify our model and sketch the derivation of a Bloch-Redfield master equation for the dissipative dynamics. Section III is devoted to the influence of charge noise, while the results for current noise are worked out in Sec. IV. Boundary effects and experimental consequences are discussed in Sec. V, while the appendix contains details of the master equation and the averaging of decay rates.

II Model and master equation

The Fermi-Hubbard model considers particles on a lattice with nearest neighbor tunneling and onsite interaction. For electrons, its Hamiltonian reads

[TABLE]

with the hopping matrix element $J$ and the interaction strength $U$. The fermionic operator $c^{\dagger}_{j\sigma}$ creates an electron with spin $\sigma$ on site $j$, while $n_{j\sigma}$ is the corresponding number operator. For convenience, we define the hopping operator between sites $j$ and $j+1$, as $T_{j}=\sum_{\sigma}c^{\dagger}_{j+1\sigma}c_{j\sigma}+\mathrm{H.c.}$ While the Hamiltonian (1) has open boundary conditions, we will also study the case of periodic boundary conditions (ring configuration) by adding the corresponding term for the hopping between the first and the last site.

Henceforth, we focus on the case of two fermions forming a spin singlet. Then we work in a Hilbert space that contains two types of states, single-occupancy states

[TABLE]

and the double-occupancy states, known as doublons,

[TABLE]

Both kinds of states are eigenstates of the operator $D$, which in the Hilbert space considered is equal to the projector onto the doublon states (3), in the following denoted as $P_{D}$.

While being different from the states in Eqs. (2) and (3), for sufficiently large values of $U$, the eigenstates of $H_{S}$ also discern into two groups, namely $N(N-1)/2$ states with energies $|\epsilon_{n}|\lesssim 4J$ and $N$ states, with energies $|\epsilon_{n}|\approx U$. We will refer to the two groups as the low-energy subspace $\mathcal{H}_{0}$, and the span of the latter as the high-energy subspace $\mathcal{H}_{1}$. In the strongly-interacting regime with $U\gg J$, treating the tunneling term as a perturbation, it is possible to express the projector onto the high-energy subspace $P_{1}$ as a power series in $J/U$, see Ref. MacDonald et al., 1988,

[TABLE]

where ${T^{+}=P_{D}T(\mathbb{I}-P_{D})}$ and ${T^{-}=(\mathbb{I}-P_{D})TP_{D}}$ comprise the hopping processes that increase and decrease the double occupancy respectively. $\mathbb{I}$ is the identity operator.

A key ingredient to our model is the coupling to environmental degrees of freedom described as $N$ independent baths of harmonic oscillators, Leggett et al. (1987); Hänggi et al. (1990) ${H_{B}=\sum_{j,n}\omega_{n}a^{\dagger}_{jn}a_{jn}}$. They couple to the Fermi-Hubbard chain via the Hamiltonian ${H_{SB}=\sum_{j}X_{j}\xi_{j}}$, where the $X_{j}$ are system operators that will be specified below. For ease of notation, we introduce the collective bath coordinates ${\xi_{j}=\sum_{n}g_{n}(a^{\dagger}_{jn}+a_{jn})}$. Moreover, we assume that all baths are equal and statistically independent, such that ${\langle\xi_{i}(t),\xi_{j}(t^{\prime})\rangle=2S(t-t^{\prime})\delta_{ij}}$.

Assuming weak coupling and Markovianity, the time evolution of the system’s density matrix $\rho$, can be suitably described by a master equation of the form Redfield (1957); Breuer and Petruccione (2007)

[TABLE]

with the operators

[TABLE]

The tilde denotes the interaction picture with respect to the system Hamiltonian, ${\tilde{X}_{j}(-\tau)=e^{-iH_{S}\tau}X_{j}e^{iH_{S}\tau}}$, while ${\mathcal{J}(\omega)=\pi\sum_{n}|g_{n}|^{2}\delta(\omega-\omega_{n})}$ is the spectral density of the baths and ${\mathcal{S}(\omega)=\mathcal{J}(\omega)\coth(\beta\omega/2)}$ is the Fourier transformed of the symmetrically-ordered equilibrium autocorrelation function $\langle\{\xi_{j}(\tau),\xi_{j}(0)\}\rangle/2$. $\mathcal{J}(\omega)$ and $\mathcal{S}(\omega)$ are independent of $j$ since all baths are identical. We will assume an ohmic spectral density ${\mathcal{J}(\omega)=\pi\alpha\omega/2}$, where the dimensionless parameter $\alpha$ characterizes the dissipation strength.

III Charge noise

Fluctuations of the background charges in the substrate essentially act upon the charge distribution of the chain. Therefore, we model it by coupling the occupation of each site to a heat bath, such that

[TABLE]

which means $X_{j}=n_{j}$. This fully specifies the master equation (5).

To get a qualitative impression of the decay dynamics of a doublon, let us start by discussing the time evolution of a doublon state in the strongly interacting regime shown in Fig. 2. For $\alpha=0$, i.e., in the absence of dissipation, the two electrons will essentially remain together throughout time evolution. This is due to energy conservation and the fact that kinetic energy in a lattice is bounded, it can be at most $2|J|$ per particle. Thus, particles forming a doublon cannot split, as they would not have enough kinetic energy on their own to compensate for the large $U$. However, since the doublon states are not eigenstates of the system Hamiltonian, we observe some slight oscillations of the double occupancy $\langle D\rangle$. Still the time average of this quantity stays close to unity, see Fig. 2(a).

On the contrary, if the system is coupled to a bath, doublons will be able to split releasing energy into the environment. Then the density operator eventually becomes the thermal state $\rho_{\infty}\propto e^{-\beta H_{S}}$. Depending on the temperature and the interaction strength, the corresponding asymptotic doublon occupancy $\langle D\rangle_{\infty}$ may still assume an appreciable value.

III.1 Numerical analysis

To gain quantitative insight, we decompose our master equation (5) into the system eigenbasis and obtain a form convenient for numerical treatment (for details, see Appendix A). A typical time evolution of the occupancy $\langle D\rangle$ is shown in Fig. 2(a). It exhibits an almost mono-exponential decay, such that the doublon life time $T_{1}$ can be defined as the $1/e$ time of the difference between initial and final value of $\langle D\rangle$,

[TABLE]

The corresponding decay rate $\Gamma=1/T_{1}$ is shown in Fig. 3 as a function of the temperature for different values of the dissipation strength $\alpha$. For small $\alpha$ and intermediate temperatures, $\Gamma$ increases with the temperature, reaching a maximum after which the tendency inverts. For sufficiently large temperatures, $\Gamma\propto(\alpha k_{B}T)^{-1}$.

III.2 Analytical estimate for the decay rate

An analytical estimate for the decay rates can often be gained from the behavior at the initial time $t=0$, i.e. from $\dot{\rho}(0)=-i[H_{S},\rho_{0}]+\mathcal{L}\rho_{0}$ with $\rho_{0}=\rho(0)$ being the pure initial state. In the present case, however, the calculation is hindered by the fast initial oscillations witnessed in Fig. 2(a). These oscillations stem from the mixing of the doublon states with the single-occupancy states. To circumvent this problem, we focus for the present purpose on the occupancy of the high-energy subspace, $\langle P_{1}\rangle$ shown in Fig. 2(b). It turns out that this quantity evolves more smoothly while it decays also on the time scale $T_{1}$. The reason for its lack of fast oscillations is that the projector $P_{1}$ commutes with the system Hamiltonian, so that it expectation value is determined solely by dissipation. Notice that the initial decay is temperature independent, while at a later stage, the decay is strongest for intermediate temperatures.

A formal way of understanding the similarity of the long time dynamics of $\langle D\rangle$ and $\langle P_{1}\rangle$ is provided by the estimate

[TABLE]

where the first lines follows from the Cauchy-Schwarz inequality for the inner product of operators, $(A,B)=\mathrm{tr}\left(A^{\dagger}B\right)$, while the second line stems from the perturbative expansion of $P_{1}$ given by Eq. (4). The result implies that when neglecting corrections of the order of $J/U$, we may determine $T_{1}$ and $\Gamma$ from either quantity. Nevertheless it is instructive to analytically evaluate $\Gamma$ for the decay of both $\langle D\rangle$ and $\langle P_{1}\rangle$.

Following our hypothesis of a mono-exponential decay, we expect

[TABLE]

therefore,

[TABLE]

This expression still depends slightly on the specific choice of the initial doublon state, in particular for open boundary conditions (see Sec. V.1, below). To obtain a more global picture, we consider an average over all doublon states, which can be performed analytically.Storcz et al. (2005) From Eq. (13), we find the average decay rate

[TABLE]

For details of the averaging procedure, see Appendix B.

For a further simplification, we have to evaluate the expressions (6) and (7) which is possible by approximating the interaction picture coupling operator as $\tilde{X}_{j}(-\tau)\simeq X_{j}-i\tau[H_{S},X_{j}]$. This is justified as long as the decay of the environmental excitations is much faster than the typical system evolution, i.e., in the high-temperature regime (HT). Inserting our approximation for $\tilde{X}_{j}$ and neglecting the imaginary part of the integrals, we arrive at

[TABLE]

With these expressions, Eq. (14) results in a temperature independent decay rate. Notice that any temperature dependence stems from the $Q_{j}$ in the first term of Eq. (14) which vanishes in the present case. While this observation agrees with the numerical findings in Fig. 2 for very short times, it does not reflect the temperature dependent decay of $\langle P_{1}\rangle$ at the more relevant intermediate stage.

This particular behavior hints at the mechanism of the bath-induced doublon decay. Let us notice that the coupling to charge noise, $X_{j}=n_{j}$, commutes with $D$. Therefore, the initial state is robust against the influence of the bath. Only after mixing with the single-occupancy states due to the coherent dynamics, the system is no longer in an eigenstate of the $n_{j}$, such that decoherence and dissipation become active. Thus, it is the combined action of the system’s unitary evolution and the effect of the environment which leads to the doublon decay.

An improved estimate of the decay rate, can be calculated by averaging the transition rate of states from the high-energy subspace to the low-energy subspace. Let us first focus on regime $k_{B}T\gtrsim U$ in which we can evaluate the operators $Q_{j}$ in the high-temperature limit. Then the average rate can be computed using expression (14) and replacing $P_{D}$ by $P_{1}$, see Appendix A. With the perturbative expansion of $P_{1}$ in Eq. (4) we obtain to leading order in $J/U$ the averaged rate

[TABLE]

valid for periodic boundary conditions. For open boundary conditions, the rate acquires an additional factor $(N-1)/N$. Notice that we have neglected back transitions via thermal excitations from singly occupied states to doublon states. We will see that this leads to some smaller deviations when the temperature becomes extremely large. Nevertheless, we refer to this case as the high-temperature limit.

In the opposite limit, for temperatures ${k_{B}T<U}$, the decay rate saturates at a constant value. To evaluate $\overline{\Gamma}$ in this limit, it would be necessary to find an expression for $\tilde{X}_{j}(-\tau)$ dealing properly with the $\tau$-dependence for evaluating the noise kernel, a formidable task that may lead to rather involved expressions. Nevertheless, one can make some progress by considering the transition of one initial doublon to one particular single-occupancy state. This corresponds to approximating our two-particle lattice model by the dissipative two-level system for which the decay rates in the Ohmic case can be taken from the literature,Weiss and Wollensak (1989); Makhlin et al. (2001) see Appendix C. Relating $J$ to the tunnel matrix element of the two-level system and $U$ to the detuning, we obtain from Eq. (41) the temperature-independent expression

[TABLE]

which formally corresponds to Eq. (17) with the temperature set to $k_{B}T=U/2$.

Figure 4 provides a comparison of these analytical findings with numerical results. The data in panel (a) reveal that the transition between the low-temperature regime and the high-temperature regime is rather sharp and occurs at $U\approx k_{B}T$. Panel (b) shows $\Gamma$ as a function of the temperature. For low temperatures, the numerical values saturate at $\overline{\Gamma}_{\text{LT}}$ obtained from the approximate mapping to a two-level system. For high temperatures, the analytical prediction $\overline{\Gamma}_{\text{HT}}$ seems slightly too large. The discrepancy stems from neglecting thermal excitations, as mentioned above.

IV Current noise

Fluctuating background currents mainly couple to the tunnel matrix elements of the system. Then the system-bath interaction is given by setting $X_{j}=T_{j}$ and reads

[TABLE]

Depending on the boundary conditions, the sum may include the term with $j=N$. The main qualitative difference of this choice is that in contrast to charge noise, $H_{SB}^{I}$ does not commute with the projector to the doublon subspace and, thus, generally ${\mathrm{tr}\left(D\mathcal{L}[\rho]\right)\neq 0}$. This enables a direct dissipative decay without the detour via an admixture of single-occupancy states to the doublon states. As a consequence, for the same value of the dimensionless dissipation parameter $\alpha$, the decay may be much faster. Also the temperature dependence of the decay changes significantly, as can bee seen in Fig. 5. While $\overline{\Gamma}$ is still proportional to $\alpha$, it now grows monotonically with the temperature.

As in the last section, we proceed by calculating analytical estimates for the decay rates. However, since the time evolution is no longer mono-exponential (not shown), we no longer start from the ansatz (14), but estimate the rate from the slope of the occupancy $\langle P_{1}\rangle$ at initial time,

[TABLE]

We again perform the average over all doublon states for $\rho_{0}$ in the limits of high and low temperatures. For periodic boundary conditions, we obtain to lowest order in $J/U$ the high and low temperature rates

[TABLE]

respectively, while open boundary conditions lead to the same expressions but with a correction factor $(N-1)/N$. In Fig. 6, we compare these results with the numerically evaluated ones as a function of the interaction [Fig. 6(a)] and the temperature [Fig. 6(b)]. Both show that the analytical approach correctly predicts the (almost) linear behavior at large values of $U$ and $k_{B}T$, as well as the saturation for small values. However, the approximation slightly overestimates the influence of the bath.

While the rates reflect the decay at short times, it is worthwhile to comment on the long time behavior under the influence of current noise. For open chains as well as for closed chains with an even number of sites, it is not ergodic as the long-time solution is not unique. The reason for this is the existence of a doublon state ${\lvert\Phi\rangle=\frac{1}{\sqrt{N}}\sum_{j=1}^{N}(-1)^{j}c^{\dagger}_{j\uparrow}c^{\dagger}_{j\downarrow}\lvert 0\rangle}$ which is an eigenstate of $H_{S}$ without any admixture of single-occupancy states. Since $T_{j}\lvert\Phi\rangle=0$ for all sites $j$, current noise may affect the phase of ${\lvert\Phi\rangle}$, but cannot induce its dissipative decay. For a closed chain with an odd number of sites, by contrast, the alternating phase of the coefficients of $\lvert\Phi\rangle$ is incompatible with periodic boundary conditions, unless a flux threads the ring. As a consequence, the chain eventually resides in the thermal state $\propto\exp(-\beta H_{S})$. The difference is manifest in the final value of the doublon occupancy at low temperatures. For closed chains with an odd number of sites, it will fully decay, while in the other cases, the population of $\lvert\Phi\rangle$ will survive.

V Discussion

V.1 Dimension and boundary effects

So far, we have considered decay rates as the averages of all possible initial doublon or high-energy states. While this is sufficient for a generic estimate of the life times, it ignores the fact that the behavior of individual states may differ significantly, in particular when the initial state is located at a boundary, which reduces the number of accessible decay channels. In Fig. 7 we present the decay rates for doublons as a function of the initial site. It reveals that in comparison to states at the center, an initial localization at the first or last site, may double the life time for charge noise and enhance by it by a factor three for current noise. The dashed lines in these plots marks the value for periodic boundary conditions, for which the value is practically the same as for a states in the center.

This knowledge about the role of boundaries and nearest neighbors provides some hint on the doublon life time in higher-dimensional lattices. Let us notice that Let the decay rates (13) and (20) contain one term for each single-occupancy state that is directly tunnel coupled to the initial site. Assuming that all terms are of the same order, we expect that $\overline{\Gamma}$ is by and large proportional to the coordination number of the lattice sites. Therefore the life time should decrease only moderately with the dimension, roughly as $T_{1}=\Gamma^{-1}\sim 2^{-D}$. From the data in Fig. 7(b), we can appreciate that for current noise, the difference between center and border is even larger. Thus, increasing dimensionality should have a slightly larger impact on the doublon life times.

V.2 Experimental implications

A current experimental trend is the fabrication of larger arrays of quantum dots,Puddy et al. (2015); Zajac et al. (2016) which triggered our question on the feasibility of doublon experiments in solid-state systems. While the size of these arrays would be sufficient for this purpose, their dissipative parameters are not yet fully known. For an estimate we therefore consider the values for GaAs/InGaAs quantum dots which have been determined recently via Landau-Zener interference.Forster et al. (2014); ona Notice, that for the strength of the current noise, only an upper bound has been reported. We nevertheless use this value, but keep in mind that it leads to a conservative estimate. In contrast to the former sections, we now compute the decay for the simultaneous action of charge noise and current noise.

Figure 8(a) shows the dissipative time evolution for a doublon initially localized at the center of a chain with 5 sites. The dynamics exhibits a few coherent oscillations in which the doublon evolves into a superposition of the kind $|2,0,0\rangle+|0,0,2\rangle$, which represents an example of a NOON state.Lee et al. (2002) Each component propagates to one end of the chain, where is it reflected such that subsequently the initial states revives. In Fig. 8(b), we depict the evolution of the corresponding doublon occupancy and the purity. Both quantities decay rather smoothly. This agrees to the finding found in Sec. IV for pure current noise which obviously dominates. It is also consistent with the values for the respective analytical decay rates in the low-temperature limit. Figure 9 shows the $T_{1}$ times for two different interaction strengths. It reveals that for low temperatures $T\lesssim J/k_{B}T$, the life time is essentially constant, while for larger temperatures, it decreases moderately until $k_{B}T$ comes close to the interaction $U$. For higher temperatures, $\Gamma$ starts to grow linearly. On a quantitative level, we expect life times of the order $T_{1}\sim 5$ ns already for a moderately low temperatures $T\lesssim 100$ mK. Since we employed the value of the upper bound for the current noise, the life time might be even larger.

Considering the analytical estimates for the decay rates at low temperatures, Eqs. (18) and (22), separately, lets us conclude that for smaller values of $U$, current noise becomes less important, while the impact of charge noise grows. Therefore, a strategy for reaching larger $T_{1}$ times is to design quantum dots arrays with smaller onsite interaction, such that the ratio $U/J$ becomes more favorable. The largest $T_{1}$ is expected in the case in which both low-temperature decay rates are equal, $\overline{\Gamma}_{\text{LT,charge}}=\overline{\Gamma}_{\text{LT,current}}$, which for the present experimental parameters is found at $U\sim 10J$ (while our data is for $U\sim 100J$). This implies that in an optimized device, the doublon life times could be larger by one order of magnitude to reach values of $T_{1}\sim 50$ ns, which is corroborated by the data in the inset of Fig. 9.

VI Conclusions

We have investigated the life times of double-occupancy states or doublons in a one-dimensional Hubbard model under the influence of dissipating environments. While in optical lattices, the resulting dissipative decay may be of minor influence, for quantum dot arrays, it will be a limiting factor.

We have considered two different couplings between the system and its environment, which physically correspond to the impact of charge noise and current noise, respectively. Within a Bloch-Redfield formalism, this model can be treated with a master equation, which allows one to numerically determine the life times from the time evolution of the reduced density operator. Moreover, it provides analytical estimates for the initial decay rates. It turned out that the striking difference between the two couplings is that the impact of charge noise decreases with the interaction, while current noise becomes increasingly relevant.

For present quantum dots, the doublon life time is expected to be of the order 5 ns, which would limit the coherent dynamics to only a few periods. However, our analytical estimates suggest that for quantum dot arrays with smaller onsite interaction, an extension by one order of magnitude should be feasible. Thus, the recent trend towards arrays with ever more coherently coupled quantum dots will allow the experimental realization of effects that so far have been measured only in optical lattices.

Appendix A Master equation in the system eigenbaasis

To bring the master equation (5) into a form that is suitable for a numerical implementation, we have to evaluate the $\tau$-integrals in Eqs. (6) and (7). This is possible after a decomposition into the system eigenbasis $\{\lvert\phi_{\alpha}\rangle\}$ with $H_{S}\lvert\phi_{\alpha}\rangle=\epsilon_{\alpha}\lvert\phi_{\alpha}\rangle$. Then the transformation to the interaction picture provides phase factors yielding a Dirac delta function and a principal value integral. Neglecting the latter, as it usually consist in a renormalization of the free system parameters, and using the notation ${\rho_{\alpha\beta}\equiv\langle\phi_{\alpha}\rvert\rho\lvert\phi_{\beta}\rangle}$ and ${X^{(j)}_{\alpha\beta}\equiv\langle\phi_{\alpha}\rvert X_{j}\lvert\phi_{\beta}\rangle}$, the master equation becomes

[TABLE]

The generalized golden-rule rates

[TABLE]

are determined by the transition matrix elements of the system operator that couples to the bath and the factors ${\Gamma_{\alpha\beta}\equiv\Gamma(\epsilon_{\alpha}-\epsilon_{\beta})}$ with

[TABLE]

and the thermal bosonic occupation number $n_{B}(\omega)=(e^{\beta\omega}-1)^{-1}$.

The Bloch-Redfield equation allows the direct computation of decay rates averaged over all possible initial states, which in our case are doublon states or high-energy states. To this end, we distinguish, those from a set $I_{1}$ labeling the high-energy states and $I_{0}$ for the low-energy states. With the formulas for the averages derived in the Appendix B and the projector to the high-energy subspace $P_{1}$, we arrive at

[TABLE]

Notice that the factor $\Delta$ accounts for the finite final value of the decay in Eq. (12). Therefore,

[TABLE]

where the bar denotes the average over all pure states belonging to the high-energy subspace, instead of the doublon subspace, see Appendix B. An alternative form for this quantity is

[TABLE]

where the last equality follows from the trace preserving property of the master equation, $\sum_{\alpha}\mathcal{L}_{\alpha\alpha,\beta\beta}=0$.

Appendix B Average over pure initial states

As an ensemble of pure states, we consider any linear combination $\lvert\psi\rangle=\sum_{n=1}^{N}c_{n}\lvert n\rangle$ of orthonormal basis states $|n\rangle$, ${n=1,\ldots,N}$. As a minimal requirement, we postulate that $\lvert\psi\rangle$ is normalized and invariant under unitary transformations. Then its probability density reads

[TABLE]

where $r^{2}=\sum_{n=1}^{N}|c_{n}|^{2}$. Then averages of the kind $\overline{c_{n}c^{*}_{m}}$ or $\overline{c_{n}c^{*}_{m}c_{n^{\prime}}c^{*}_{m^{\prime}}}$ can be expressed as integrals of polynomials over the $(2N-1)$-dimensional unit sphere. Following Ref. Folland, 2001, we find

[TABLE]

which implies

[TABLE]

To compute average rates for the transitions between two groups of states, cf. Eq. (28), the initial linear combination $\lvert\psi\rangle$ is restricted to the doublon subspace which has dimension $N_{D}$. Therefore we have to replace the prefactor $N$ by $N_{D}$ and the operators $A$ and $B$ by their projections to the subspace, $P_{D}AP_{D}$ and $P_{D}BP_{D}$.

Appendix C Two-level system

For completeness, we summarize the Bloch-Redfield result for the decay rates of the two-level system coupled to an Ohmic bath.Weiss and Wollensak (1989); Makhlin et al. (2001) For the notation used in the main text, it is defined by the Hamiltonian

[TABLE]

with the tunnel matrix element $\Delta$ and the detuning $\epsilon$. The bath coupling is specified by (i) $X=\sigma_{z}$ for charge noise and (ii) $X=\sigma_{x}$ for current noise, respectively. To establish a relation to our Hubbard chain, we identify the detuning by the interaction, $\epsilon\simeq U$, and $\Delta=\sqrt{8}J$. Note that replacing charge noise by current noise corresponds to interchanging $\epsilon$ and $\Delta$. Therefore, we can restrict the derivation of the decay rate to case (i).

It is straightforward to transform the Hamiltonian into the eigenbasis of the two-level system, where it reads

[TABLE]

with $E=\sqrt{\epsilon^{2}+\Delta^{2}}$, while the system-bath coupling becomes

[TABLE]

In the interaction picture, it is

[TABLE]

Again ignoring the imaginary part of the integral in (7), the noise kernel can be written as

[TABLE]

The projector to the initial state is $P_{1}=(\sigma_{0}+\sigma_{z})/2$, so that the decay rate can be found as

[TABLE]

where for an Ohmic spectral density

[TABLE]

Accordingly, we find for case (ii) the rate

[TABLE]

which provides the analytical high-temperature result (18) for charge noise.