Dissipation assisted Thouless pumping in the Rice-Mele model

TL;DR

This paper shows that dissipation from a thermal environment can enhance the robustness of topological charge pumping in the Rice-Mele model, making it more quantized and experimentally feasible.

Contribution

It demonstrates that dissipation can improve topological pumping quantization in the Rice-Mele model, using a Floquet framework and quantum master equations.

Findings

Dissipation increases Floquet band population with topological winding.

Dissipative dynamics lead to more accurate quantized charge pumping.

Quantum master equation matches exact MPS results at low temperature.

Abstract

We investigate the effect of dissipation from a thermal environment on topological pumping in the periodically-driven Rice-Mele model. We report that dissipation can improve the robustness of pumping quantisation in a regime of finite driving frequencies. Specifically, in this regime, a low-temperature dissipative dynamics can lead to a pumped charge that is much closer to the Thouless quantised value, compared to a coherent evolution. We understand this effect in the Floquet framework: dissipation increases the population of a Floquet band which shows a topological winding, where pumping is essentially quantised. This finding is a step towards understanding a potentially very useful resource to exploit in experiments, where dissipation effects are unavoidable. We consider small couplings with the environment and we use a Bloch-Redfield quantum master equation approach for our numerics: Comparing these results with an exact MPS numerical treatment we find that the quantum master equation works very well also at low temperature, a quite remarkable fact.

Equations

$${\hat{H}}_{\mathrm{tot}}(t)={\hat{H}}_S(t)+{\hat{H}}_{SB}+{\hat{H}}_B$$

$$\widehat{H}_{S}(t)=-\sum_{j=1}^{N}\Big{(}J_{1}(t){\hat{c}^{\dagger}}_{j,\mathrm{\scriptscriptstyle B}}{\hat{c}^{\phantom{\dagger}}}_{j,\mathrm{\scriptscriptstyle A}}+J_{2}(t){\hat{c}^{\dagger}}_{j+1,\mathrm{\scriptscriptstyle A}}{\hat{c}^{\phantom{\dagger}}}_{j,\mathrm{\scriptscriptstyle B}}+\mathrm{H.c.}\Big{)}+\Delta(t)\,\sum_{j=1}^{N}\Big{(}{\hat{c}^{\dagger}}_{j,\mathrm{\scriptscriptstyle A}}{\hat{c}^{\phantom{\dagger}}}_{j,\mathrm{\scriptscriptstyle A}}-{\hat{c}^{\dagger}}_{j,\mathrm{\scriptscriptstyle B}}{\hat{c}^{\phantom{\dagger}}}_{j,\mathrm{\scriptscriptstyle B}}\Big{)}$$

$${\hat{H}}_S(t)=\sum _k^{\mathrm{BZ}}\left[\begin{array}{cc}{{\hat{c}}^{†}}_{k,\mathrm{A}}& {{\hat{c}}^{†}}_{k,\mathrm{B}}\end{array}\right]{\hat{\mathcal{H}}}_S^k(t)\left[\begin{array}{c}{{\hat{c}}^{}}_{k,\mathrm{A}}\\ {{\hat{c}}^{}}_{k,\mathrm{B}}\end{array}\right]=\sum _k^{\mathrm{BZ}}\left[\begin{array}{cc}{{\hat{c}}^{†}}_{k,\mathrm{A}}& {{\hat{c}}^{†}}_{k,\mathrm{B}}\end{array}\right]\mathbf{R}(k,t)\cdot \hat{\mathbf{\sigma }}\left[\begin{array}{c}{{\hat{c}}^{}}_{k,\mathrm{A}}\\ {{\hat{c}}^{}}_{k,\mathrm{B}}\end{array}\right]$$

$$\mathbf{R}(k,t)=(-J_1(t)-J_2(t)\cos (ka),-J_2(t)\sin (ka),\Delta (t){)}^{\mathrm{T}}.$$

$${\hat{H}}_{\mathrm{tot}}(t)$$

$${\hat{\mathcal{H}}}_k(t)$$

$$\hat{\mathrm{J}}(t)=\frac{1}{L}\sum _k^{\mathrm{BZ}}\left[\begin{array}{cc}{{\hat{c}}^{†}}_{k,\mathrm{A}}& {{\hat{c}}^{†}}_{k,\mathrm{B}}\end{array}\right]{\hat{\mathcal{J}}}_k(t)\left[\begin{array}{c}{{\hat{c}}^{}}_{k,\mathrm{A}}\\ {{\hat{c}}^{}}_{k,\mathrm{B}}\end{array}\right]$$

$$\hat{\mathcal{J}}_{k}(t)=\frac{ea}{2\hbar}\bigg{(}J_{2}(t)\sin(ka)\hat{\sigma}^{x}+\Big{(}J_{1}(t)-J_{2}(t)\cos(ka)\Big{)}\hat{\sigma}^{y}\bigg{)}\;,$$

$$Q_{m}=\int_{-\frac{\pi}{a}}^{\frac{\pi}{a}}\frac{\mathrm{d}k}{2\pi}\int_{(m-1)\tau}^{m\tau}\mathrm{d}t\;{\mathrm{Tr}}\Big{(}\hat{\mathcal{J}}_{k}(t)\,\hat{\rho}_{\scriptscriptstyle\mathrm{S}}^{k}(t)\Big{)}\;.$$

$$Q_m\stackrel{m\to \infty }{⟶}\overline{Q}\equiv \underset{M\to \infty }{\lim }\frac{1}{M}\sum _{m=1}^M{Q}_m.$$

$$|{\psi }_n(t)⟩=e^{-i{ϵ}_{n}t/\hslash }|u_n(t)⟩$$

$${\hat{\rho }}_{\mathrm{S}}^k(t)=\sum _{n,n^{\prime }}{\rho }_{n{n}^{\prime }}^k(t)|{\psi }_n^k(t)⟩⟨{\psi }_{n^{\prime }}^k(t)|$$

$$\overline{Q}=\underset{M\to \infty }{\lim }\frac{1}{M}\sum _{n,n^{\prime }}{\int }_{-\frac{\pi }{a}}^{\frac{\pi }{a}}\frac{\mathrm{d}k}{2\pi }{\int }_0^{M\tau }\mathrm{d}t{e}^{-\frac{i}{\hslash }({ϵ}_n^k-{ϵ}_{n^{\prime }}^k)t}{F}_{n{n}^{\prime }}^k(t)$$

$$F_{n{n}^{\prime }}^k(t)={\rho }_{n{n}^{\prime }}^k(t){\mathbb{J}}_{n^{\prime }n}^k(t)$$

$${\mathbb{J}}_{n^{\prime }n}^k(t)=⟨u_{n^{\prime }}^k(t)|{\hat{\mathcal{J}}}_k(t)|u_n^k(t)⟩$$

$${\rho }_{n{n}^{\prime }}^{k,\mathrm{coh}}=⟨u_n^k(0)|{\psi }_k(0)⟩⟨{\psi }_k(0)|u_{n^{\prime }}^k(0)⟩.$$

$$\overline{Q}=Q_{\mathrm{d}}^{\mathrm{coh}}={\int }_{-\frac{\pi }{a}}^{\frac{\pi }{a}}\frac{\mathrm{d}k}{2\pi }\sum _n{\rho }_{nn}^{k,\mathrm{coh}}{\int }_0^{\tau }\mathrm{d}t{\mathbb{J}}_{nn}^k(t).$$

$$Q_m\stackrel{m\to \infty }{⟶}\overline{Q}=Q_{\mathrm{d}}^{\mathrm{diss}}={\int }_{-\frac{\pi }{a}}^{\frac{\pi }{a}}\frac{\mathrm{d}k}{2\pi }\sum _n{\int }_0^{\tau }\mathrm{d}t{\rho }_{nn}^k(t){\mathbb{J}}_{nn}^k(t).$$

$$Q_{\mathrm{d}}^{\mathrm{diss}}=\frac{1}{\hslash \omega }{\int }_{-\frac{\pi }{a}}^{\frac{\pi }{a}}\mathrm{d}k\sum _n{\overline{\rho }}_{nn}^k{\partial }_k{ϵ}_n^k,$$

$$\dot{\mathbf{r}}=\frac{2}{\hslash }\mathbf{R}\times \mathbf{r}-{\mathbb{A}}_{\mathrm{diss}}\cdot \mathbf{r}-\mathbf{b}.$$

$${\gamma }_{\phi }$$

$${\gamma }_{\mathrm{R}}$$

$$\gamma (\omega )={\int }_{-\infty }^{+\infty }\mathrm{d}t{\text{e}}^{i\omega t}⟨\hat{X}(t)\hat{X}(0){⟩}_0.$$

$${\mathbb{A}}_{\mathrm{diss}}$$

$$\mathbf{b}$$

$${\gamma }_{\mathrm{D},i}$$

$${\gamma }_{ij}$$

$${\hat{\mathcal{H}}}_k(t)=\mathbf{R}(k,t)\cdot \hat{\mathbf{\sigma }}+(\mathbf{n}\cdot \hat{\mathbf{\sigma }})\otimes {\int }_0^1\mathrm{d}x\lambda (x)({{\hat{b}}^{†}}_k(x)+{{\hat{b}}^{}}_k(x))+{\int }_0^1\mathrm{d}x\hslash \omega (x){{\hat{b}}^{†}}_k(x){{\hat{b}}^{}}_k(x)$$

$${\mathcal{J}}_{\beta }(\omega )=\frac{\mathrm{sgn}(\omega )\mathcal{J}(|\omega |)}{2}\left(1+\coth (\beta \omega /2)\right),$$

$${{\hat{a}}^{†}}_{k,j}={\int }_{-1}^1\mathrm{d}x{U}_j(x){{\hat{b}}^{†}}_k(x)\text{with}{U}_j(x)={\lambda }_{\beta }(x)p_j(x)$$

$$⟨p_j,p_{j^{\prime }}⟩={\int }_{-1}^1\mathrm{d}x{\lambda }_{\beta }^2(x)p_j(x)p_{j^{\prime }}(x)={\delta }_{j,j^{\prime }}.$$

$${\pi }_{j+1}(x)=(x-{\alpha }_j){\pi }_j(x)-{\beta }_j{\pi }_{j-1}(x)j=0,1,2,...$$

Full text

1 Introduction

Quantised adiabatic transport in insulators, discovered by Thouless [1] in 1983, is a topic of current interest in the field of topological insulators, recently fostered by the experimental realization of a Thouless pump with ultra-cold atoms [2, 3].

The strict quantisation of the pumped charge requires the quantum dynamics to be, in principle, adiabatic and unitary. However, in concrete experimental realizations these two requirements cannot be perfectly fulfilled. The study of non-adiabatic effects on a specific example of Thouless pump, the so-called Rice-Mele model [4], showed the emergence of quadratic corrections in the driving frequency whenever the system is initially prepared in the initial Hamiltonian ground state [5]. Similar non-adiabatic effects have been discussed in the context of the quantisation of the Hall conductivity in the Harper-Hofstadter model [6]. The effect of disorder on the topological pumping of the Rice-Mele model has been recently discussed in Ref. [7]. An important point is also the effect of interactions, which make the system non integrable and lead the system to heat up $T=\infty$ in the long-time limit. In this case, adiabatic pumping is asymptotically washed out and can occur only as a transient condition [8].

Dissipative effects may also disrupt pumping quantisation. However, the impact and role of dissipation in the performance of a Thouless pump remains an open question. Some studies have taken into account thermal effects by using a thermal initial state, instead of the Hamiltonian ground state [9, 5], followed by a unitary dynamics. Within such a framework, it was found [9] that charge quantisation is robust against non-zero temperatures in the initial state when a single pumping cycle is considered. Moreover, in the limit of an infinite number of pumping cycles [5], thermal corrections were found to be exponentially small for low enough temperatures but increasingly relevant when the temperature approached the insulating gap. Concerning papers where a genuine dissipative dynamics is considered, we mention a study [10] of the Qi-Wu-Zhang model through a Lindblad Markovian quantum master equation, where it is found that the pumped charge, starting from the quantisation value, decreases monotonically to zero with increasing noise.

In this work, we aim at understanding how topological quantum pumping is affected by the interaction with a bosonic thermal bath. To do so, we study the time evolution of the Rice-Mele Hamiltonian [4] and analyse how charge pumping is affected both by dissipation and by non-adiabatic effects. The remarkable finding of our study is that — in appropriate conditions, i.e. if the temperature of the bath is low enough, and in the limit of infinite cycles — a dissipative dynamics may be beneficial to a Thouless pump, fighting against non-adiabatic effects, and leading to a pumped charge which is closer to the quantised value.

We rationalise our findings by analysing the dissipative results in terms of the Floquet states of the unitary dynamics. Interestingly, we find that the charge pumped at stationarity is still expressed, exactly as in the unitary case, only in terms of the populations of the Floquet bands, with no role for the quantum coherences. Quantised pumping is essentially related to the perfect occupation of what one might call the “lowest-energy Floquet band”, i.e. the Floquet band constructed by choosing, for each momentum, the Floquet mode with (period-averaged) lowest-energy expectation: this is the quasi-energy band “closest” to the instantaneous ground state, that shows the correct winding [11] across the Brillouin zone (BZ) and leads to a quantised charge transport, up to non-adiabatic corrections which are exponentially small in the driving period, see Refs. [12, 5]. The dissipation-induced improved pumping is found to be strictly related to the fact that dissipation brings to a higher population for such a lowest-energy Floquet band. We might regard it as a dissipative preparation of a topological Thouless pump, somewhat similar in spirit to the results of Ref. [13], which however deals with a two-dimensional non-driven system with a specifically engineered dissipation. It is also worth to mention Ref. [14], dealing with non-topological pumping in superconducting nanocircuits. The authors show that a zero-temperature bath can make the pumped charge at finite frequency closer to its coherent adiabatic counterpart. This phenomenon occurs also in our case.

The paper is organised as follows: in Sec. 2 we introduce our dissipative version of the Rice-Mele model. In Sec. 3 we define the pumped charge, i.e. the observable displaying quantisation in the non-dissipative adiabatic limit, and we describe it within the Floquet theory. In Sec. 4 we illustrate how we compute the dissipative time evolution of the system. In particular, we benchmark the reliability of the Bloch-Redfield Quantum Master Equation against a non-perturbative approach, which exploits a unitary chain-mapping transformation in order to carry out time evolution with Matrix Product States. We find a very good agreement, especially at T=0 where the quantum master equation was not expected to work [15].

In Sec. 5 we present our numerical results for the pumped charge and discuss the pumping enhancement due to a low-temperature bath. We deepen the understanding of this effect in the Floquet framework, which provides insights in terms of populations of the Floquet bands. Finally, in Sec. 6 we summarize the main results and draw our conclusions.

2 The model

We study here a dissipative version of the Rice-Mele model, where the non-unitary dynamics comes from coupling the system to a bosonic bath at thermal equilibrium. The Hamiltonian is

[TABLE]

where the three terms on the r.h.s. are, respectively, the system, system-bath interaction and bath Hamiltonians.

The system Hamiltonian consists of a bipartite lattice on which spinless fermions hop on nearest-neighbor sites, according to the Rice-Mele Hamiltonian [4]:

[TABLE]

where $N$ is the number of diatomic cells, ${\hat{c}^{\dagger}}_{j,\mathrm{\scriptscriptstyle A}(\mathrm{\scriptscriptstyle B})}$ creates a fermion on site $A(B)$ of the $j^{th}$ cell, $J_{1}(t)$ and $J_{2}(t)$ are respectively intra-cell and inter-cell hopping terms and $\Delta(t)$ modulates the on-site energies. We assume periodic boundary conditions (PBC), so that translational invariance allows us to Fourier transform the fermionic operators for the $A$ and $B$ sites separately, ${\hat{c}^{\phantom{\dagger}}}_{j,\mathrm{\scriptscriptstyle A}(\mathrm{\scriptscriptstyle B})}=\sum_{k}e^{ikaj}{\hat{c}^{\phantom{\dagger}}}_{k,\mathrm{\scriptscriptstyle A}(\mathrm{\scriptscriptstyle B})}/\sqrt{N}$, where the sum over the discrete wave-vectors $k=2\pi n/(Na)$, with $n=0,1,\dots,N-1$ and $a$ the cell length, runs inside the first BZ. By applying this transformation to Eq. (2), we block-diagonalise the Hamiltonian in sectors of different $k$:

[TABLE]

where $\hat{\mathcal{H}}_{S}^{k}(t)={\bf R}(k,t)\cdot\hat{\boldsymbol{\sigma}}$ is a 2-dimensional operator, conveniently parameterised, using the Bloch-vector notation, with the vector of Pauli matrices $\hat{\boldsymbol{\sigma}}$ and an effective magnetic field ${\bf R}(k,t)$

[TABLE]

The energy bands $E_{\pm}(k,t)=\pm|{\bf R}(k,t)|$ of $\hat{\mathcal{H}}_{S}^{k}(t)$ never “touch” except at a single “metallic” point with $\Delta=0$ and $J_{1}=J_{2}$, where ${\bf R}=0$ (for $k=\pi/a$). We completely fill the lowest band with $N$ particles — the half-filling condition — and realise a quantised adiabatic pumping by driving the band insulator around the metallic point with a schedule $J_{1}(t)=J_{0}+\delta_{0}\cos\omega t$, $J_{2}(t)=J_{0}-\delta_{0}\cos\omega t$ and $\Delta(t)=\Delta_{0}\sin\omega t$. Here $\omega=2\pi/\tau$ is the driving frequency, associated to a period $\tau$, so that $\widehat{H}_{S}(t+\tau)=\widehat{H}_{S}(t)$.

To account for dissipation in the simplest and most practical way, we choose to include identical but independent harmonic baths for each $k$-subsector, coupled in the usual Caldeira-Leggett spin-boson [16] fashion

[TABLE]

where ${\bf n}$ specifies a unit vector for the bath coupling, $\hat{X}_{k}=\sum_{l}\lambda_{l}({\hat{b}^{\dagger}}_{k,l}+{\hat{b}^{\phantom{\dagger}}}_{k,l})$, $\lambda_{l}$ are coupling constants, and $[{\hat{b}^{\phantom{\dagger}}}_{k,l},{\hat{b}^{\dagger}}_{k^{\prime},l^{\prime}}]=\delta_{k,k^{\prime}}\delta_{l,l^{\prime}}$. In terms of the original fermions, a $\hat{\sigma}^{z}$-bath coupling, given by ${\bf n}=(001)$, would correspond to a term $({\hat{c}^{\dagger}}_{k,\mathrm{\scriptscriptstyle A}}{\hat{c}^{\phantom{\dagger}}}_{k,\mathrm{\scriptscriptstyle A}}-{\hat{c}^{\dagger}}_{k,\mathrm{\scriptscriptstyle B}}{\hat{c}^{\phantom{\dagger}}}_{k,\mathrm{\scriptscriptstyle B}})\otimes\hat{X}_{k}$. This will be our standard choice unless otherwise specified. Our dissipative Rice-Mele model can then be effectively regarded as a collection of $N$ independent dissipative 2-level systems, one for each momentum $k$ in the BZ, with a system-bath coupling effectively acting on $\hat{\sigma}^{z}$. The interaction between system and environment is encoded in the bath spectral function ${\mathcal{J}}(\omega)=\sum_{l}\lambda_{l}^{2}\delta(\omega-\omega_{l})$. We will consider a standard Ohmic dissipation [16], modelled in the frequency continuum limit as ${\mathcal{J}}(\omega)=2\alpha\hbar^{2}\omega F_{\scriptscriptstyle\mathrm{cutoff}}(\omega,\omega_{c})$, where $\alpha$ is the (dimensionless) coupling strength, and $F_{\scriptscriptstyle\mathrm{cutoff}}(\omega,\omega_{c})$ implements a suitable frequency cutoff $\omega_{c}$. In the following we will adopt either the exponential form $F_{\scriptscriptstyle\mathrm{cutoff}}=\exp(-\omega/\omega_{c})$ or the hard-cutoff $F_{\scriptscriptstyle\mathrm{cutoff}}=\Theta(\omega_{c}-\omega)$.

3 The pumped charge

The current density operator $\hat{\mathrm{J}}(t)$ is obtained as a derivative of the system Hamiltonian with respect to an external flux $\Phi$ passing through the hole of the PBC ring, $\hat{\mathrm{J}}=\frac{\partial\widehat{H}_{S}(\Phi)}{\partial\Phi}\vline_{\Phi=0}$. In the present case, the current density operator can be expressed as

[TABLE]

where $L=Na$ is the system length,

[TABLE]

and $e$ is the charge of the electron. Given the density matrix $\hat{\rho}_{\scriptscriptstyle\mathrm{S}}^{k}(t)$, the pumped charge during the $m^{th}$ driving period, $Q_{m}$, is thus given by integrating the current density over the period. In the thermodynamic limit $L\to\infty$ we can write

[TABLE]

After a sufficiently large number of cycles, the average pumped charge is expected to converge to an asymptotic value

[TABLE]

In the next section we are going to use Floquet theory to study in detail how this convergence occurs.

3.1 Floquet analysis for the pumped charge

We introduce here some important notions of Floquet theory applied to charge pumping, extending the coherent evolution treatment presented in Ref. [5]. Since the driving is periodic, i.e. $\widehat{H}_{S}(t)=\widehat{H}_{S}(t+\tau)$, from Floquet theory [17, 18, 19] we know that the solutions to the time-dependent Schrödinger equation for the closed system have the following form

[TABLE]

where ${n}$ labels the possible solutions, $|\psi_{n}(t)\rangle$ are called Floquet states, $|u_{n}(t)\rangle$ are called Floquet modes and are $\tau-$periodic, and $\epsilon_{n}$ are the quasi-energies. In the present case, using $k$ as a quantum number and $n=\pm$ for two Floquet states $|\psi_{n}^{k}(t)\rangle$ at each $k$, we can rewrite the density matrix $\hat{\rho}_{\scriptscriptstyle\mathrm{S}}^{k}(t)$ in the coherent Floquet basis

[TABLE]

where $\rho_{n{n^{\prime}}}^{k}(t)=\langle\psi_{n}^{k}(t)|\hat{\rho}_{\scriptscriptstyle\mathrm{S}}^{k}(t)|\psi_{{n^{\prime}}}^{k}(t)\rangle$.

In this framework, the infinite-time average pumped charge, Eq. (9), can be written as

[TABLE]

where

[TABLE]

and

[TABLE]

is the matrix element of the current operator between Floquet modes, hence an explicitly periodic quantity. In the coherent evolution case [5], the density matrix $\rho_{{n}{n^{\prime}}}^{k}(t)$ turns out to be time-independent and related to the initial state $|\psi_{k}(0)\rangle$ as

[TABLE]

In turn, if the quasi-energies are non degenerate, when ${n}\neq{n^{\prime}}$, the $k$-integral will vanish in the limit $t\to\infty$, since the oscillating phase factors $\textrm{e}^{-i(\epsilon_{{n}}^{k}-\epsilon_{{n^{\prime}}}^{k})t/\hbar}$ will lead to destructive interference cancellations. More formally, this is a consequence of the Riemann-Lebesgue lemma applied to the $k$-integration, as explained in detail in Refs. [20, 21].

Using this result in Eq. (12), and exploiting the infinite-time integration, it follows that only the populations $\rho_{{n}{n}}^{k,\mathrm{coh}}=|\langle u_{{n}}^{k}(0)|\psi_{k}(0)\rangle|^{2}$ of the Floquet bands come into play, and one arrives at the so-called Floquet diagonal ensemble [20]

[TABLE]

With a very similar application of the Riemann-Lebesgue lemma it is also possible to see that, in the thermodynamic limit, $Q_{m}$ defined in Eq. (8) tends to $Q_{\rm d}^{\mathrm{coh}}$ when $m\to\infty$.

In the dissipative case, $\rho_{{n}{n^{\prime}}}^{k}(t)$ is generally time dependent. What we find — and explicitly discuss in Sec. 5 and Fig. 4 — is that after a certain transient, because of dissipative effects, $\rho_{{n}{n^{\prime}}}^{k}(t)$ becomes $\tau$-periodic. Hence, we can again apply the Riemann-Lebesgue lemma as done in Ref. [21] and show that only the diagonal terms contribute, arriving at the dissipative version of the Floquet diagonal ensemble formula for the average pumped charge

[TABLE]

As in Ref. [5], the asymptotic pumped charge can be related to the properties of the Floquet quasienergies. In order to do that, a result coming from our numerics is crucial: if we approximate $\rho_{{n}{n}}^{k}(t)$ at stationarity with its average value on one period $\overline{\rho}_{{n}{n}}^{k}$, we get corrections to the pumped charge of the order $\sim 10^{-6}$. So, a very good approximation [22] consists in replacing $\rho_{{n}{n}}^{k}(t)$ with its time-average $\overline{\rho}_{{n}{n}}^{k}$. With this approximation, using arguments similar to those of Ref. [5], we find that the asymptotic pumped charge can be written in the form

[TABLE]

which looks very similar to Eq. (4) in Ref. [5], obtained for a coherent dynamics. Because the derivative with respect to $k$ can be recast as the derivative with respect to an external flux [5], this equation is also strictly analogous to Eq. (19) of Ref. [23], where pumping in a dissipative superconducting nanocircuit was considered.

4 Numerical methods

In this section we introduce and compare two different methods for computing the dissipative dynamics of our model. The first method discussed in Sec. 4.1 is the standard Bloch-Redfield Quantum Master Equation (QME) [24, 25] for a two-level system, which is fast, but valid only within the weak-coupling Born-Markov approximation. The second method, described in Sec. 4.2, is based on a unitary mapping of the harmonic bath into a bosonic chain interacting with the system  [26, 27, 28, 29]. The full system-plus-bath time evolution is then obtained by using Matrix Product States (MPS) and the 2-site Time-Dependent Variational Principle [30, 31, 32]. To the best of our knowledge, we present here the very first results obtained by this method in a non-equilibrium scenario, i.e. in a system with explicit time-dependence. This second method is reliable at any coupling strength and accounts for non-Markovian effects, at the expense of a higher computational complexity. We therefore used it here just to benchmark the QME results in some specific regimes. We remark that the reliability of a QME was thoroughly tested in Ref. [33] against the numerically-exact QUasi-Adiabatic Path Integral (QUAPI) method [34, 35] for the dissipative Landau-Zener model. We perform here a similar check with the present driving scheme, see Sec. 4.3.

All the results presented in Sec. 5 have been obtained using the QME, since this allows us to check the thermodynamic limit ($N\to\infty$), in a reasonable computational time.

4.1 Bloch-Redfield equation

According to the Hamiltonian in Eq. (5) and under the assumptions of weak system-bath coupling and Born-Markov approximation [24, 25, 36, 37, 38, 33], we can write a QME to describe the reduced density matrix $\hat{\rho}_{\scriptscriptstyle\mathrm{S}}^{k}(t)$ for the system for each $k$-vector. There are many slightly different ways of writing down the relevant QME, depending, for instance, on whether or not one adopts a rotating wave approximation (RWA). The approach used below makes use of a RWA and is essentially equivalent to the “double-sided adiabatic QME” in Lindblad form explained in Ref. [39]. Observe that this treats the dissipative dynamics assuming an adiabatic driving, but we will also apply it in regimes where this condition is slightly violated.

Writing the system density matrix in the Bloch-vector notation $\hat{\rho}_{\scriptscriptstyle\mathrm{S}}^{k}(t)=(\mathbb{1}+{\bf r}_{k}(t)\cdot\hat{\boldsymbol{\sigma}})/2$, the equations to solve for the dynamics in Schrödinger representation read

[TABLE]

where we dropped the $k$ and $t$ labels from ${\bf r}$ and related quantities. From the Lindblad form derived using the rotating-wave-approximation (RWA) [39], we can express the relevant ingredients appearing in ${\mathbb{A}}_{\mathrm{diss}}$ and ${\bf b}$ in terms of pure-dephasing $\gamma_{\varphi}$ and relaxation $\gamma_{\scriptscriptstyle\mathrm{R}}$ rates [40]

[TABLE]

where $S_{X}(\omega)=\gamma(\omega)+\gamma(-\omega)$ is the symmetrised Fourier transform of the free thermal bath correlation function

[TABLE]

For the case of Ohmic dissipation with exponential cutoff, ${\mathcal{J}}(\omega)=2\alpha\hbar^{2}\omega\textrm{e}^{-\omega/\omega_{c}}$, we find that $S_{X}(0)=8\pi\hbar\alpha/\beta$ and $S_{X}(2E/\hbar)=2\pi{\mathcal{J}}(2E/\hbar)\coth(\beta E)$, with $E=|{\bf R}|$, $\beta=(k_{B}T)^{-1}$ and $k_{B}$ the Boltzmann constant. In terms of these quantities, the dissipation matrix ${\mathbb{A}}_{\mathrm{diss}}$ and the vector ${\bf b}$ have the form

[TABLE]

where

[TABLE]

These equations agree with those discussed in Ref. [40], obtained, with similar approximations, from a perturbative diagrammatic approach. As one can easily verify, if ${\bf R}$ is time independent, the final steady state value of ${\bf r}(t\to\infty)$ correctly describes the thermal density matrix for each momentum $k$. To compute the system’s dynamics, we solved Eqs. (19) through a standard fourth-order Runge-Kutta method.

4.2 Chain mapping based time evolution using Matrix Product States

In order to benchmark the QME we employ a numerically exact technique using MPS, recently introduced by Tamascelli et al. [29]. The method is essentially based on two key ideas. First, a more efficient representation of the finite temperature bath is obtained by modifying the spectral density. In particular, the thermal state will be represented by the vacuum state of a transformed bosonic environment, which allows us to write the state of system-plus-bath as a pure state rather than a density matrix. Second, the spin-boson (star-like) Hamiltonian is transformed to a chain geometry, using orthogonal polynomials. This chain Hamiltonian can be simulated efficiently using MPS. For our benchmarks we use the Ohmic spectral density with hard cutoff: ${\mathcal{J}}(\omega)=2\alpha\hbar^{2}\omega\Theta(\omega_{c}-\omega)$.

Let us now discuss the method in more detail. We start from the continuum version of the Hamiltonian in Eq. (5) for a given wave vector $k$. Parameterizing the modes by the dimensionless variable $x$ we have [27, 28]

[TABLE]

where now the bosonic operators have a (frequency) continuum normalisation in terms of Dirac delta: $[{\hat{b}^{\phantom{\dagger}}}_{k}(x),{\hat{b}^{\dagger}}_{k^{\prime}}(x^{\prime})]=\delta_{k,k^{\prime}}\delta(x-x^{\prime})$. Here, the functions $\lambda(x)$ and $\omega(x)$ define the interaction and dispersion, respectively, and are connected to the spectral density through the relation $\lambda^{2}(x)=\omega^{\prime}(x){\mathcal{J}}(\omega(x))$ [26, 28]. The freedom in choosing those functions is exploited by setting $\omega(x)=\omega_{c}x$, which is useful in the algebra of orthogonal polynomials performed below. At this point the bath is described by the density matrix corresponding to the thermal state, requiring large local dimensions at high temperatures. Now the new insight of Ref. [29] — see also Ref. [27] for an alternative route through a thermal Bogoliubov transformation — is that the thermal state can be represented by the bosonic vacuum, if the spectral density ${\mathcal{J}}(\omega)$ is replaced by the temperature-dependent spectral density

[TABLE]

extending its support to negative frequencies. As a consequence, the integrals in Eq. (24) have an extended domain $[-1,1]$, and $\lambda(x)$ has to be replaced by $\lambda_{\beta}(x)$, which is defined through the thermal spectral density: $\lambda_{\beta}^{2}(x)=\omega_{c}{\mathcal{J}}_{\beta}(\omega(x))$.

Next, following Refs. [41, 26, 28], we perform a mapping from the star-like spin-boson geometry to a chain geometry. We do that by defining the unitary transformation

[TABLE]

and inverse transformation ${\hat{b}^{\dagger}}_{k}(x)=\sum_{j=0}^{\infty}U_{j}(x){\hat{a}^{\dagger}}_{k,j}$. The (real) polynomials $p_{j}(x)$ of degree $j$ are orthonormal with respect to the inner product

[TABLE]

Moreover, the corresponding monic polynomials $\pi_{j}(x)$ — obtained when dividing each polynomial $p_{j}$ by the coefficient of the leading power — satisfy the three term recurrence relation [26]

[TABLE]

with initial polynomial $\pi_{-1}\equiv 0$. While those polynomials are not needed explicitly, we need to find the recurrence coefficients $\alpha_{j}$ and $\beta_{j}$. For zero temperature they are known for some special spectral densities [26], but in most cases we need to calculate them numerically [42, 43]. Exploiting this recurrence relation, the transformation in Eq. (26) yields the chain geometry Hamiltonian

[TABLE]

with system-bath coupling $\kappa_{0}=\sqrt{\int_{-\omega_{c}}^{\omega_{c}}\mathrm{d}\omega\ {\mathcal{J}}_{\beta}(\omega)}$, hopping amplitudes $\kappa_{j}=\hbar\omega_{c}\sqrt{\beta_{j}}$ and on-site energies $\Omega_{j}=\omega_{c}\alpha_{j}$. This infinite chain can be truncated at some finite length $L_{\mathrm{chain}}$, such that the excitations reflected at the border of the truncated chain do not reach the system [28]. We then carry out the time evolution using the 2-site Time Dependent Variational Principle [30, 31, 32], which has proven to combine high accuracy and efficiency [45]. Let us now briefly discuss the sources of errors in this method. First, we need to find the recurrence coefficients for the chain mapping, which is done numerically using the package in Refs.[42, 43]. The error, however, is under full control and any target accuracy can be reached. Once we have our chain Hamiltonian, we carry out time-evolution. Here, additional errors are due to the splitting of the differential equation in the Time-Dependent Variational Principle and due to limited representation capabilities of the MPS. While the splitting error can be controlled through the time step, the representation capabilities of the MPS are enhanced by increasing the bond dimension and the local bosonic dimension, which has to be truncated to some finite value. In summary, all errors are well under control and can be made arbitrarily small.

4.3 Benchmarking the QME

We benchmark here the QME data against the numerically-exact MPS approach. Since the latter technique is computationally heavy, we consider here a chain of $N=16$ diatomic cells. Moreover, the most interesting results, in the following, will be related to a low-temperature bath, which is also the regime were the QME is believed to be less reliable [15]. For this reason, we fix here $T=0$. The bath spectral function is again Ohmic with hard cutoff, as described in Sec. 4.2, with $\alpha=0.001$ and $\hbar\omega_{c}=10J_{0}$. The pumping is performed over $M=12$ cycles and for two different driving frequencies $\hbar\omega/J_{0}=0.2$ and $\hbar\omega/J_{0}=0.4$. In Fig. 1 we plot the total current density as a function of time, focusing on the first and the last period, where stationarity is already attained. We observe an excellent agreement between the QME data (orange lines) and the MPS results (blue circles) in all regimes. The corresponding pumped charges agree within an error of the order of the precision we use to compute the $k$-integrals. Indeed, for $\hbar\omega/J_{0}=0.2$ we obtain $\overline{Q}=0.99999$ from the QME and $\overline{Q}=0.99998$ from the MPS, while for $\hbar\omega/J_{0}=0.4$ we get $\overline{Q}=1.00066$ from the QME and $\overline{Q}=1.00068$ from the MPS.

5 Dissipative pumping results

We present here how dissipation affects the pumped charge at different driving frequencies. In all the following numerical analysis we will approximate the integral over $k$ with a discrete sum in the first BZ. All the calculations are performed with sizes $N$ which we have verified to be large enough to be representative of the thermodynamic limit: in practice, $N\sim 100$ is enough in the presence of dissipation. We choose the cutoff frequency $\omega_{c}$ in the spectral function to be much bigger than the widest spectral gap, $\omega_{c}=1000J_{0}/\hbar$ (we comment upon different choices of $\omega_{c}$ in A). We start considering the bath coupling strength to be $\alpha=0.001$; as we will see in Sec. 5.2, increasing it will have quite interesting effects on the pumping dynamics, stabilizing the quantised pumping against the driving for $\hbar\omega>0.4J_{0}$. On the opposite, for smaller frequencies, results are insensitive on $\alpha$.

Let us start considering the behaviour of the pumped charge after a single cycle. In Fig. 2(a) we plot the charge pumped after a single cycle, $Q_{1}$, versus the driving frequency $\omega$. On the one hand, at larger values of the frequency, the bath has almost no effect, and the behaviour at all temperatures remains almost identical to the coherent one, which coincides with that reported in Ref. [5]. On the other hand, at smaller values of the frequency, the system has enough time to “feel” thermal effects and in general moves away from the ideal quantised pumping, here corresponding to $Q=1$. The charge converges to a finite value which depends on the bath temperature $T$, see the inset of Fig. 2(a). We will further comment on this point later.

These results change remarkably when pumping over a larger number of cycles. We find, and this is one of the main results of the paper, that the charge pumped over the $m$-th cycle $Q_{m}$ can overcome the corresponding coherent result in presence of a thermal bath of sufficiently low temperature. This is shown in Fig. 2(b). Observe that dissipation makes the convergence to the infinite-time average much faster than the coherent case. Notice also that the infinite-time average results are precisely described by the Floquet diagonal ensemble formulas, Eqs. (16) and (17), shown by horizontal dashed lines.

Figure 2(c) shows the infinite-time pumped charge $\overline{Q}$ versus the driving frequency $\omega$, for both coherent and dissipative evolutions at different $T$. The dissipative results are obtained either from Eq. (9) (with $M=100$) or equivalently from Eq. (17); the coherent results are obtained from Eq. (16). Observe that at $T=0$ the dissipative results are always well above the coherent ones. We remark that dissipation at $T=0$ restores a nearly quantised pumped charge $\overline{Q}=1$ away from the strict adiabatic limit $\omega\to 0$. $\overline{Q}$ starts to deviate from one only for frequencies around $\hbar\omega\simeq 0.5J_{0}$, as discussed in Sec. 5.2. This dissipative improvement of the pumped charge persists also at finite $T$, for large enough $\omega$. We define this phenomenon dissipation assisted Thouless pumping. This finding is supported by the benchmark done in Sec. 4.3, where the QME and the MPS results are found to be in excellent agreement. Our finding is also qualitatively independent of whether one assumes in the QME some form of RWA or not, as shown in B. Notice, however, that some of the quantitative aspects of our results do depend on the details of the QME chosen.

Interestingly, in the small frequency regime we observe that $\overline{Q}(\omega\rightarrow 0)\equiv Q_{m}(\omega\rightarrow 0)$ for any $m\geq 1$, i.e. the $\omega\to 0$ limit is independent of the number of driving periods. This is because in the $\omega\to 0$ limit the dissipative transient induced by the bath occurs within a single driving period.

Notice that previous discussions of thermal effects in the Rice-Mele model [9, 5, 46] have considered a coherent Schrödinger dynamics starting from an initial thermal state. We compared the results obtained from such an approach with the genuinely dissipative evolution studied in the present paper. In general, we observed completely different results, both in the short and in the intermediate frequency ranges of study, as illustrated in Fig. 3. This is definitely not surprising, but worth mentioning. Fig. 3 is important also in another respect: we see that in the dissipative case, if we take very different initial conditions, we get the same asymptotic regime. This is not at all surprising in a dissipative system and marks the difference with the asymptotic regime of the coherent case [20].

Focusing again on Fig. 2, observe the flat $\omega$-behaviour of $\overline{Q}$ for the dissipative evolution at $T=0$. It would be interesting to pin down if the corrections to the strict adiabatic limit $\overline{Q}(\omega\to 0)=1$ change from power-law [5] to exponentially small in the presence of zero-temperature dissipation. Unfortunately, this question is extremely difficult to answer numerically. Indeed, this aspect of the analysis is highly sensitive to the type of weak-coupling approximation performed in the QME. Although the QME results obtained with different flavours of RWA and without RWA are in good qualitative agreement, they are quantitatively different in that respect (see B). A similar question might be asked concerning the behaviour of $\overline{Q}(\omega\to 0,T)$ as a function of the bath coupling temperature $T$, a question that is once again numerically elusive and rather sensitive to the details of the QME used. A detailed MPS analysis might help in clarifying these points, but we suspect that the high numerical accuracy necessary to answer some of these questions would make these issues extremely delicate.

5.1 Floquet analysis of dissipative results

It is insightful to understand the dissipative improvement shown in Fig. 2(c) within a Floquet framework. Let us therefore go back to the discussion of Sec. 3.1 and study the dynamics of the system density matrix in the Floquet basis.

A crucial point will be to numerically show that the the populations and the coherences in the Floquet basis $\rho_{nn^{\prime}}^{k}(m\tau)$ converge to an asymptotic $\tau$-periodic behaviour, in order to justify the application of the Riemann-Lebesgue lemma and the validity of Eq. (17).

Another point will be to understand the dissipation assisted Thouless pumping in terms of these quantities. In order to do that, we consider the lowest-energy Floquet state. This state is the Floquet state with maximum overlap with the instantaneous ground state. The instantaneous ground state gives rise to perfectly quantised pumping in the full adiabatic limit $\omega\to 0$. For finite $\omega$ the lowest-energy Floquet state does the same (up to corrections exponentially small in $\omega$) [12, 5]. We can construct the lowest-energy Floquet state by choosing, for each $k$, the Floquet state with (period-averaged) lowest-energy expectation $|\psi_{-}^{k}(t)\rangle$. On the opposite, choosing for each $k$ the Floquet state with period-averaged highest energy expectation we construct the highest energy Floquet state $|\psi_{+}^{k}(t)\rangle$.

It is interesting to understand if the dissipative improvement of the pumped charge is induced by the environment making this special state more populated. In order to answer to this question we will focus on the population of the highest-energy Floquet state $\rho_{++}^{k}(m\tau)$ and show that it is reduced by the environment in comparison to the unitary case. This means that the population of the lowest-energy state $\rho_{--}^{k}(m\tau)=1-\rho_{++}^{k}(m\tau)$ is increased.

We analyze the two points raised above in Fig. 4. We start focusing on the point about convergence. We show the stroboscopic dynamics of the coherences absolute values $|\rho_{-+}^{k}(m\tau)|$ between the two Floquet bands in Fig. 4(a) and of the excited-band populations $\rho_{++}^{k}(m\tau)$ in Fig. 4(c). We see that both quantities converge to stationary values for $m\to\infty$. The phases of the coherences (not shown) also converge to fixed values. After stroboscopic stationarity is reached, the intra-period behaviour of these quantities is illustrated in panels (b) and (d), respectively: observe a $\tau/2$-periodicity in both cases. This justifies the application of the Riemann-Lebesgue lemma in Sec. 3.1 and the validity of Eq. (17).

Now let us focus on the point about the population of the highest-energy Floquet state. In Fig. 4 (d) we compare the populations of the highest-energy Floquet band in the coherent and dissipative cases, showing that $\rho_{++}^{k}$ is generally reduced by several orders of magnitude in presence of dissipation at $T=0$. Hence $\rho_{--}^{k}$, the population of the lowest-energy Floquet state, is increased and then the topological pumping is improved at finite frequencies. In the main figure we fix $k$ and look at the dependence on time, while in the inset we plot the period-averaged $k^{th}$-population $\overline{\rho}_{++}^{k}$ versus the momentum $k$. In the inset we note, incidentally, the presence of a value of $k$ where the coherent value shows an irregularity and the dissipative value shows a peak. That peak corresponds to a Floquet quasi-resonance which gives rise to a non-adiabaticity and increases the asymptotic dissipative population of the highest-energy Floquet band [47].

We conclude that dissipation moves the system towards the lowest-energy Floquet state, and that is the way in which the bath improves the topological pumping. This asymptotic condition where the Floquet state with maximum overlap with the instantaneous ground state is populated is very different from what happens in the high-frequency limit where the asymptotic condition is given by the Floquet-Gibbs states [48, 49].

5.2 Effect of the bath coupling strength

We focus here on how the stationary pumped charge $\overline{Q}$ changes as the interaction $\alpha$ is changed over different orders of magnitude. We have used the QME to study this behaviour and we show in Fig. 5 our results for $T=0$ and for the three frequencies $\hbar\omega/J_{0}=0.4,0.5,0.6$. Notice that we show results up to $\alpha=0.1$ for completeness, but we consider our QME reliable up to $\alpha\simeq 10^{-2}$, see Ref. [33]. Data for $\hbar\omega\leq 0.4J_{0}$ are almost constant in the range $\alpha\in[10^{-5},10^{-2}]$ (we just show $\hbar\omega=0.4J_{0}$). However, at larger frequencies things get more interesting: The curves are non monotonic and there is a maximum of $\overline{Q}$ at around $\alpha\simeq 6\times 10^{-3}$. There is therefore an interval of $\alpha$ where a stronger dissipation stabilizes the pumping against the driving.

Focusing on the data in Fig. 5 for $\alpha=0.001$, we can furthermore notice that the trend $\overline{Q}=1$ observed previously in Fig. 2(c) starts to deviate from the quantised value at frequencies $\hbar\omega\simeq 0.5J_{0}$.

6 Conclusions

We analysed the role of dissipation from a somewhat idealised thermal environment — coupling independent baths to each fermionic $k$-mode — on Thouless pumping in the Rice-Mele model. We found that a low temperature bath can assist against undesired (inevitable) non-adiabatic effects. Indeed, at fixed finite driving frequency, the pumped charge obtained from dissipative evolution can be closer to the quantised value with respect to the one obtained from purely coherent dynamics [5]. Dissipation induces this improvement because it increases the population of the lowest-energy Floquet band. Indeed, the pumped charge would be essentially quantised — up to exponentially small terms — when this band is completely filled. This is somewhat in line, in a non-topological context, with the finding of improved pumping in a three-site fermionic chain [50] subject to dissipation.

Our findings are qualitatively independent of the system-bath coupling chosen as long as we stay in a weak coupling regime, as discussed in C. They are furthermore validated by benchmarking against a numerically-exact method that does not rely upon any weak-coupling nor Born-Markov approximation.

The fact that thermal effects can be beneficial is remarkable and interesting for future experimental realizations. We stress that the effect is not related to a bath engineering, exploited in the literature for other topological models[13].

A further step towards a deeper understanding would be to study more realistic couplings to the environment, e.g. via operators acting on sites in real space, which break the entanglement in physical space. However, this analysis requires more sophisticated approaches [51, 52, 53], and is left to future studies.

Another interesting direction would be to study topological measures, such as the Uhlmann phase [54] and the Ensemble Geometric Phase (EGP) [46]. In particular, it would be interesting to inquire if the Uhlmann phase of the asymptotic time-periodic effective density matrix has a relation with the pumped charge, in analogy with the Berry or the Aharonov-Anandan geometric phase in the coherent cyclic case [23, 1].

We acknowledge fruitful discussions with R. Fazio, S. Barbarino, L. Privitera and M. Wauters. Research was partly supported by EU FP7 under ERC-ULTRADISS, Grant Agreement No. 834402. GES acknowledges that his research has been conducted within the framework of the Trieste Institute for Theoretical Quantum Technologies (TQT).

Appendix A Dependence on the cutoff

Let us focus on the issue of choosing the cutoff $\omega_{c}$ in the bath spectral function. Generally $\omega_{c}$ is taken to be the largest energy scale of the system, so that the dynamics becomes insensitive to the detail of this parameter. In the present case, since the system energy gap is always of the order of $J_{0}$ and we consider temperatures $T\leq J_{0}$, we require $\hbar\omega_{c}\gg J_{0}$. The behavior of $\overline{Q}$ vs the cutoff $\omega_{c}$, see Fig. 6, shows the range of cutoff frequencies for which we observe a convergence of $\overline{Q}$. We therefore selected $\hbar\omega_{c}=1000J_{0}$. Fig. 6 is also useful to illustrate the effect of some basic dissipation mechanisms. If $\omega_{c}$ is much smaller than the minimum system energy gap, the probability of having jumps between energy levels is negligible and the result tends to become again insensitive to the cutoff value. Then, the only relevant dissipation mechanism comes from pure dephasing, given by $\gamma_{\varphi}$ in Eq. (20a). Notice however that $\gamma_{\varphi}\sim T$ and hence it vanishes at $T=0$. This is consistent with what we observe in Fig. 6: for $\hbar\omega_{c}\ll J_{0}$, $\overline{Q}$ is insensitive to the cutoff; moreover, for $T=0$, we recover precisely the coherent result, pinpointed by the horizontal dashed line.

Appendix B Different approximations in the Bloch-Redfield QME

In this appendix, starting again from the Bloch-Redfield QME, we derive and employ two sets of equations alternative to Eq. (19) for the study of the steady state pumped charge.

The first one simply consists in the Bloch-Redfield QME without any additional RWA. The equations have again the general form written in Eq. (19), with the coherent part unchanged, while the incoherent part, for ${\bf n}=(001)$ system-bath coupling, reads

[TABLE]

where

[TABLE]

and the vector ${\bf b}$ is

[TABLE]

The second approach makes use of the Bloch-Redfield QME expanded in the coherent system Floquet basis, $\{|\psi_{n}^{k}(t)\rangle\}_{{n}=\pm}$, defined by Eq. (10). We will consider here a single dissipative two-level system at fixed momentum and we will omit the $k$ label in all the operators for clarity. Following Refs. [55, 19, 23, 47], it is possible to perform a RWA according to the quasi-energies (analogous to the standard one done on the system’s energies). Eventually, one would see that the equations for the coherences decouple from the ones for the populations. It can be shown that the coherences go to zero after a finite time, so that the steady state is diagonal in the Floquet basis. The populations $\rho_{{n}{n}}(t)=\langle\psi_{{n}}(t)|\hat{\rho}_{\scriptscriptstyle\mathrm{S}}(t)|\psi_{{n}}(t)\rangle$ can be determined by the rate equation [55, 19, 23, 47]

[TABLE]

where $\rho_{++}=1-\rho_{--}$ and the rates are given by

[TABLE]

where we defined $A_{{n^{\prime}}{n},l}$ as the $l^{th}$ Fourier coefficient of the $\tau$-periodic function $\langle u_{n^{\prime}}(t)|({\bf n}\cdot\boldsymbol{\sigma})|u_{n}(t)\rangle$, while $\gamma(\omega)$ is the Fourier transform of the bath correlation function and $\Delta_{{n^{\prime}}{n},l}=(\epsilon_{n^{\prime}}-\epsilon_{n})/\hbar-2\pi l/\tau$. The steady state is then very easily determined by setting Eq. (33) to zero, which leads to

[TABLE]

We employed both approaches to compute the pumped charge at stationarity $\overline{Q}$ vs the driving frequency. These results are compared to the ones shown in the main text in Fig. 2(c) for the dynamics induced by Eq. (19). In Fig. 7, we provide the outcomes from the three approaches: In the plot, “$\textrm{RWA}_{\textrm{energy}}$” refers to the RWA according to the instantaneous system’s energies, “no RWA” points to the data obtained from the Bloch-Redfield QME without RWA, while “$\textrm{RWA}_{\textrm{Floquet}}$” refers to the RWA performed on the QME written in the system’s Floquet basis. We observe that the improvement over the coherent curve is obtained in all the three cases, giving further support to our claims. Furthermore, the results obtained from the two versions of the RWA seem to match quite well, especially at smaller frequencies. Nevertheless, the various approached give results which are quantitatively different. For example, at $T=0.75J_{0}$ and in the frequency range studied, one might get or not an improvement over the coherent case depending on the approach used. It is hard to say which version of the QME is best approximating the true dissipative time-evolution, until a careful benchmarking against the numerically exact MPS approach is performed for an extensive set of parameters. This is left to a future study.

Appendix C Dependence on the coupling operator

To test the generality of our findings, we study also system-bath coupling operators different from $\hat{\sigma}_{z}$. We focus here on the case in which each two-level system is coupled to the reservoir via the $\hat{\sigma}_{x}$ operator, which would correspond to choosing ${\bf n}=(100)$ in Eqs. (20a) and (20b). Fig. 8 shows the result for $\overline{Q}$ vs the driving frequency. We see that there is no qualitative difference with the result shown in Fig. 2(c), corresponding to a coupling via $\hat{\sigma}_{z}$. We tried also other coupling operators (for generic ${\bf n}$) and we obtained qualitatively similar results.