Non-renewal statistics in quantum transport through the eyes of first-passage and waiting time distributions

TL;DR

This paper compares first-passage and waiting time distributions in quantum transport, demonstrating that first-passage times better capture non-renewal behavior and backtunneling effects in nanoscale systems.

Contribution

It introduces the use of first-passage time distributions to accurately analyze non-renewal behavior in quantum transport, especially when waiting times are insufficient.

Findings

First-passage times correctly identify non-renewal behavior.

Waiting times fail to incorporate backtunneling effects.

First-passage analysis applies across multiple transport regimes.

Abstract

The waiting time distribution has, in recent years, proven to be a useful statistical tool for characterising transport in nanoscale quantum transport. In particular, as opposed to moments of the distribution of transferred charge, which have historically been calculated in the long-time limit, waiting times are able to detect non-renewal behaviour in mesoscopic systems. They have failed, however, to correctly incorporate backtunneling events. Recently, a method has been developed that can describe unidirectional and bidirectional transport on an equal footing: the distribution of first-passage times. Rather than the time between successive electron tunnelings, the first-passage refers to the first time the number of extra electrons in the drain reaches $+1$. Here, we demonstrate the differences between first-passage time statistics and waiting time statistics in transport scenarios where the waiting time either cannot correctly reproduce the higher order current cumulants or cannot be calculated at all. To this end, we examine electron transport through a molecule coupled to two macroscopic metal electrodes. We model the molecule with strong electron-electron and electron-phonon interactions in three regimes: (i) sequential tunneling and cotunneling for a finite bias voltage through the Anderson model, (ii) sequential tunneling with no temperature gradient and a bias voltage through the Holstein model, and (iii) sequential tunneling at zero bias voltage and a temperature gradient through the Holstein model. We show that, for each transport scenario, backtunneling events play a significant role; consequently, the waiting time statistics do not correctly predict the renewal and non-renewal behaviour, whereas the first-passage time distribution does.

Equations

$$\dot{\mathbf{P}}(t)$$

$$\mathbf{P}(t)=e^{\mathbf{L}t}\mathbf{P}(0).$$

$$\dot{\mathbf{P}}(n,t)$$

$$\mathbf{P}(\chi ,t)$$

$$\mathbf{P}(n,t)$$

$$\dot{\mathbf{P}}(\chi ,t)$$

$$\mathbf{L}(\chi )$$

$$\mathbf{P}(\chi ,t)$$

$$M(\chi ,t)$$

$$\left(\mathbf{I},\mathbf{P}(n,t)\right)$$

$$M(\chi ,t)$$

$$F$$

$$=\frac{\mathcal{S}(0)}{2⟨I⟩}.$$

$$w_F(\tau )$$

$$\Pi (\tau )$$

$$\tilde{w}(z)$$

$$⟨⟨{\tau }^k⟩⟩$$

$$⟨I⟩=\frac{1}{⟨\tau {⟩}_F}-\frac{1}{⟨\tau {⟩}_B},$$

$$⟨I{⟩}_F$$

$$⟨{\tau }_1+{\tau }_2+\dots +{\tau }_k⟩$$

$$=\sum _{p=1}^k\int d{\tau }_p{\tau }_{p}w({\tau }_p)$$

$$=k⟨\tau ⟩,$$

$$⟨I{⟩}_F$$

$$R$$

$$p$$

$$⟨{\tau }_1{\tau }_2⟩$$

$$=\frac{\left(\mathbf{I},\mathbf{J}{\mathbf{L}}_0^{-2}\mathbf{J}{\mathbf{L}}_0^{-2}\mathbf{J}\bar{\mathbf{P}}\right)}{\left(\mathbf{I},\mathbf{J}\bar{\mathbf{P}}\right)}.$$

$$\mathbf{P}(n,t)$$

$$=\mathbf{T}(n|t)\mathbf{P}(0).$$

$$\tilde{\mathbf{T}}(n|z)$$

Full text

I Introduction

Molecular electronic devices generally consist of some quantum system, such as a molecule or quantum dot, sandwiched between two macroscopic metal electrodes Cuevas and Scheer (2010); Xiang et al. (2016). One of the main objectives of mesoscopic physics is to use the particle transport through such devices to understand the nature of these quantum systems. Directly measuring transitions internal to the quantum system destroys any quantum coherence it may exhibit; instead, one must reconstruct the system behaviour from transport properties available from measurements of bath variables, such as the stationary particle current $\langle I\rangle$. Molecular devices are small enough, however, that thermal effects and stochastic single electron transport cause fluctuations in $\langle I\rangle$. Such fluctuations hold a dual interest in mesoscopic electron transport. On the one hand they impede the performance of electronic components that require reliable currents; but on the other hand they are largely determined by the molecule’s underlying dynamics and can thus be used to probe quantum behaviour, which is difficult to observe. Consequently, multiple theoretical tools describing fluctutations in quantum electron transport have been developed Bagrets and Nazarov (2003); Brandes (2008); Saito and Dhar (2016); Budini et al. (2014); Esposito et al. (2010).

The full counting statistics (FCS) calculates cumulants $\langle\langle I^{k}\rangle\rangle$ of the distribution of total transferred charge $n$ in a fixed time interval. Levitov et al. (1996); Bagrets and Nazarov (2003); Hassler et al. (2008); Nazarov (2003) The second current cumulant, the zero-frequency noise, can reveal the internal dynamics of the molecule Blanter and Büttiker (2000); Mahé et al. (2010); Bednorz and Belzig (2010). For example, super-Poissonian shot noise has many origins: it can arise from telegraphic switching due to spin-polarised leads Ptaszyński (2017, 2017) or inelastic cotunneling Thielmann et al. (2005a); Kaasbjerg and Belzig (2015); negative differential resistance due to asymmetric coupling Thielmann et al. (2005b); the dynamical channel blockade Belzig (2005); or avalanching electrons due to interactions with a vibrational mode Koch and Von Oppen (2005). Higher order cumulants are necessary to fully characterise the transport when the current distribution is non-Gaussian Bruderer et al. (2014); Reulet et al. (2003); Bomze et al. (2005). The major limitation of fixed time tools like the FCS is that they have historically been calculated in the infinite time, or zero-frequency, limit Budini et al. (2014). Although finite frequency, and thus finite time, FCS have been developed Emary et al. (2007); Ubbelohde et al. (2012); Galaktionov et al. (2003); Nagaev et al. (2004); Pilgram et al. (2004); Galaktionov and Zaikin (2011); Marcos et al. (2010); Salo et al. (2006), it is numerically and analytically easier to work in the zero-frequency regime. In recent years, however, statistical tools that operate in fluctuating time but fixed $n$ have arisen, motivated by the long time limitations of FCS.

The most prominent flucutating time tool is the waiting time distribution (WTD) $w(\tau)$: the conditional probability density that, given an electron has tunneled to the drain electrode, the next electron will tunnel to the drain after a time interval $\tau$ Brandes (2008). Since the inception of waiting times in mesoscopic electron transport, a common question has been whether the WTD contain information distinct from the FCS. It has been shown that the large deviation rate functions describing zero-frequency FCS are equivalent to those for waiting times Budini et al. (2014); the WTD cumulants $\langle\langle\tau\rangle\rangle$ can be combined to exactly reproduce the zero-frequency current cumulants. Budini (2011); Esposito et al. (2010); Albert et al. (2011) These relations only hold, however, when the renewal assumption is satisfied: $w(\tau_{1},\tau_{2})=w(\tau_{1})w(\tau_{2})$, where $w(\tau_{1},\tau_{2})$ is the joint probability for successive waiting times. When subsequent waiting times are correlated the renewal assumption is violated and these relationships no longer hold; the WTD contains information separate to the FCS. Interestingly, even in Markovian systems non-renewal behaviour can arise from telegraphic switching Ptaszyński (2017, 2017), superconductivity Albert et al. (2016), quantum coherence Dasenbrook et al. (2015), and the presence of strong interactions Rudge and Kosov (2018); Kosov (2017a). WTD measurements are experimentally accessible via sophisticated real-time single-electron detection techniques Gustavsson et al. (2006, 2009, 2008); Zby (2015); Lu et al. (2003); or alternatively the WTD may be reconstructed from low-order current correlation functions Haack et al. (2015).

The WTD in quantum electron transport can be calculated via scattering theory Albert et al. (2011, 2012); Chevallier et al. (2016); Albert and Devillard (2014); Albert et al. (2016); Haack et al. (2014); Dasenbrook et al. (2015) and non-equilibrium Green’s functions Tang et al. (2014); Seoane Souto et al. (2015). In order to connect it with other fluctuating time statistics, however, we focus on WTDs calculated from quantum master equations, which rely on quantum jump operators defined from the Liouvillian Brandes (2008); Rudge and Kosov (2018, 2016a, 2016b); Kosov (2017b); Thomas and Flindt (2013); Sothmann (2014); Kosov (2017a); Ptaszyński (2017); Rajabi et al. (2013); Walldorf et al. (2018); Potanina and Flindt (2017); Welack et al. (2009); Ptaszyński (2017); Kosov (2018). Unfortunately, WTDs calculated in this manner are so far only able to reproduce $\langle\langle I^{k}\rangle\rangle$ ($k>1$), as according to renewal theory, for unidirectional transport. For bidirectional transport, one may define WTDs between tunnelings to and from the drain; the total average current, regardless of whether the transport is renewal, is then $\langle I\rangle=\frac{1}{\langle\tau\rangle_{F}}-\frac{1}{\langle\tau\rangle_{B}}$, where $\langle\tau\rangle_{F}$ is the average waiting time between forward tunnelings, conditioned upon the initial probability of a tunneling to the drain, and $\langle\tau\rangle_{B}$ is the corresponding backward tunneling average waiting time. It is not clear, however, how to obtain the higher order total current cumulants from these two WTDs. In some cases, furthermore, the processes that contribute to the current do not appear in the Liouvillian; for example, elastic cotunneling transfers electrons from the source to the drain in the same quantum process without altering the state of the intermediate system Braggio et al. (2006); Kaasbjerg and Belzig (2015); Utsumi et al. (2006); Emary (2009). To fully describe the transport in these cases, one must use the $n$-resolved master equation Gurvitz (1998); but waiting times calculated from this method Thomas and Flindt (2013); Potanina and Flindt (2017); Seoane Souto et al. (2015) are only defined for unidirectional transport Rudge and Kosov (2018). Bidirectional transitions play an important role in mesoscopic electron transport Fujisawa et al. (2006), so it is desirable to have a fluctuating time statistic applicable to this regime.

Consider instead the first-passage time distribution (FPTD): the distribution of times until the net number of transferred electrons in the drain first reaches some number $n$. First-passage times have recently gained traction in non-equilibrium thermodynamics; first-passage time theory has been developed for fluctuations in entropic variables describing stationary Markovian processes Bauer and Cornu (2014); Roldán et al. (2015); Neri et al. (2017); Garrahan (2017); Saito and Dhar (2016); Ptaszyński (2018). This work, and subsequent experimental and theoretical research from Singh et al. Singh et al. (2018), has predicted and verified fluctutation relations for first-passage times. Saito and Dhar Saito and Dhar (2016) have developed first-passage time notation for Markovian quantum systems, which is suitable for mesoscopic electron transport described by $n$-resolved master equations. Ptaszynski Ptaszyński (2018) has also recently brought first-passage times into line with FCS and WTDs, by developing a method to calculate correlations between subsequent first-passage times; deriving exact relationships between the FPTD cumulants and the FCS cumulants; and demonstrating that these relationships disappear when the renewal assumption is violated. The expectation now is that first-passage times are able to analyse the non-renewal behaviour of molecular systems similarly to how WTDs have been used, except that first-passage times are applicable to bidirectional transport. Considering these developments, it is interesting to compare the WTD and the FPTD in systems and regimes where one expects bidirectional transitions.

This paper’s purpose is to demonstrate the usefulness of FPTDs by comparing them to the WTD in regimes where the WTD either cannot be calculated or does not correctly recreate the FCS. We demonstrate this in three transport scenarios. First, cotunneling through an Anderson impurity is considered. Several authors have theoretically Thielmann et al. (2005b); Kaasbjerg and Belzig (2015) and experimentally Okazaki et al. (2013); Zhang et al. (2007); Onac et al. (2006) shown that for a particular voltage range super-Poissonian noise arises from telegraphic switching induced by inelastic cotunneling processes. In this voltage range backtunneling processes can occur; so until now it has remained impenetrable to waiting times. One would expect that telegraphic switching would be accompanied by positive correlations in the subsequent waiting times, as the transport switches between channels with different characteristic $\langle\tau\rangle$; however, we demonstrate that the transport in this regime is actually characterised by small negative correlations in the first-passage times. Second, we examine sequential tunneling through the Holstein model due to a voltage bias. Waiting times have already been theoretically calculated for this scenario Kosov (2017a, b); we demonstrate that those calculations do not correctly reproduce the FCS and that the correlations between the first-passage times are slightly smaller than those between the waiting times. Finally, we consider sequential tunneling through the Holstein model with zero bias voltage but with a temperature gradient: a scenario in which backtunneling events majorly contribute to the transport De and Muralidharan (2018).

The paper is organised as follows. Section II outlines the three relevant statistical tools: FCS, the WTD, and the FPTD. Section III describes the dynamics of the three transport scenarios, written in the form of Markovian master equations. Section IV contains results from the statistical tools calculated for the three transport scenarios, while Section V summarises the work and the Appendix details cotunneling calculations and vibration jump operators.

Throughout the manuscript we use natural units: $\hbar=k_{B}=e=1$.

II Statistical tools

II.1 $n$-resolved master equation

We consider a general quantum system weakly coupled to a number of thermal reservoirs, which is described by a Markovian master equation for the reduced system density matrix:

[TABLE]

where $\{\mathbf{P}\}_{l}(t)$ is the probability that the system is in discrete state $l$ and $\mathbf{L}$ is the Liouvillian containing the system dynamics. Throughout the paper, we refer to this as the standard master equation, to distinguish it from the $n$-resolved master equation. We ignore coherence between system states, so that Eq.(1) reduces to a rate equation written in the system basis. The system moves from state $l$ to state $k$ via tunneling of particles or quasi-particles from or to the thermal reservoirs, and the transition rate for this process is $\Gamma_{kl}$. The Liouvillian off-diagonals are $[\mathbf{L}]_{kl}=\Gamma_{kl}$ and the diagonals are $[\mathbf{L}]_{kk}=-\sum\limits_{l\neq k}\Gamma_{lk}$.

The system starts in some state $\mathbf{P}(0)$, so that the general solution of Eq.(1) is

[TABLE]

We assume that the system has a unique stationary state $\bar{\mathbf{P}}$, which is the null vector satisfying $\mathbf{L}\bar{\mathbf{P}}=0$.

In order to calculate the FCS, FPTD, and WTD, it is necessary to resolve the master equation on the net number of particles $n$ transferred to the drain electrode:

[TABLE]

where $\{\mathbf{P}(n,t)\}_{k}$ now represents the probability that, at time $t$, the net number of extra electrons in the drain is $n$ and the system is in state $k$. It is important to note that despite the term “extra electrons” $n$ may be zero, a negative integer, or a positive integer. From here, following Ptaszynski’s notation Ptaszyński (2018), we will refer to $n$ as the jump number. Sequential and cotunneling processes add or remove at most one electron to or from the drain, so $n-n^{\prime}=0,\pm 1$. The Liouvillian can then be split into two jump operators containing all the transitions that move electrons forward to the drain $\mathbf{J}_{F}$ and all the transitions that move electrons back from the drain $\mathbf{J}_{B}$, and a part containing the remaining dynamics $\mathbf{L}_{0}=\mathbf{L}-\mathbf{J}_{F}-\mathbf{J}_{B}$. The $n$-resolved master equation is then

[TABLE]

The Fourier transform of $\mathbf{P}(n,t)$

[TABLE]

conveniently introduces a counting field $\chi$ that transforms Eq.(3) from an infinite set of coupled differential equations to a solvable problem Bagrets and Nazarov (2003); Nazarov (1999):

[TABLE]

Now, $\mathbf{P}(\chi,t)$ time-evolves according to the solution of Eq.(7), with the initial condition $\mathbf{P}(\chi,0)=\bar{\mathbf{P}}$ since it is assumed that counting begins at $t=0$ when the system has already reached the steady state:

[TABLE]

II.2 Full counting statistics

The generating function for the cumulants of transferred charge is

[TABLE]

where $P(n,t)=\operatorname{Tr}[\mathbf{P}(n,t)]$ is the probability that $n$ extra electrons are transferred to the drain in the time interval $[0,t]$. Throughout the paper we use the superoperator notation $\operatorname{Tr}[\mathbf{A}]=\left(\mathbf{I},\mathbf{A}\right)$, where $\mathbf{A}$ is a column vector of length $m$ and $\mathbf{I}=(1,1,1,...)$ is a row vector of length $m$. Consequently, Eq.(10) transforms to

[TABLE]

with the inverse Fourier transform

[TABLE]

Examining the trace of Eq.(6), one sees that the part corresponding to the cumulant generating function is

[TABLE]

Theoretically, one can now calculate all current cumulants in any time interval by direct differentiation of $M(\chi,t)$ in the limit $\chi\rightarrow 0$: \langle\langle I^{k}\rangle\rangle=\frac{1}{t}(-i)^{k}\frac{\partial^{k}}{\partial\chi^{k}}M(\chi,t)\Big{|}_{\chi=0}, where $t$ is the measurement time. The second cumulant, for example, is $\langle I^{2}\rangle-\langle I\rangle^{2}$. It is well known that, except for the simplest systems, these calculations remain too difficult at finite $t$. Bagrets and Nazarov have shown that in the long-time limit, the cumulant generating function reduces to

[TABLE]

where $\Lambda_{max}$ is the eigenvalue of $\mathbf{L}(\chi)$ with the largest real part Bagrets and Nazarov (2003); Flindt et al. (2008).

We focus on the Fano factor, the second cumulant scaled by the average current, which measures the zero-frequency noise power $\mathcal{S}(0)$ compared to that of a Poissonian process:

[TABLE]

For Poissonian processes $F=1$; sub-Poissonian and super-Poissonian processes are characterised by $F<1$ and $F>1$, respectively.

II.3 Waiting time distribution

If instead one wishes to examine the distribution of waiting times between subsequent tunnelings to the drain, then one may use the formula developed by Brandes Brandes (2008):

[TABLE]

The waiting time distribution between tunnelings from the drain, $w_{B}(\tau)$, is similarly defined.

There are tunneling events that do not appear in the standard master equation, yet change the drain number $n$. To include these events in the WTD, one must use the $n$-resolved master equation and the definition of the WTD from the idle time probability Thomas and Flindt (2013); Rudge and Kosov (2018); Albert et al. (2012); van Kampen (1992):

[TABLE]

where $p$ is the probability for observing a tunneling to the drain and $\Pi(\tau)$ is the idle time probability: the probability that no tunneling events to the drain are observed by time $\tau$.

In the case of unidirectional transport, when $\mathbf{J}_{B}=0$ and all backtunneling transitions are safely neglected from the Liouvillian, $p=\langle I\rangle$ and the idle time probability is obtained from $n$-resolved master equation:

[TABLE]

Now, $\mathbf{L}(\chi)=\mathbf{L}_{0}+\mathbf{J}_{F}e^{i\chi}$ and no backtunneling transition rates are included in $\mathbf{L}_{0}$. Evidently, backtunneling events will be accompanied by an $e^{-i\chi}$ factor that will diverge in the limit $i\chi\rightarrow-\infty$. Of course, if the transport is unidirectional in the backward direction as opposed to the forward direction, then the same definition may be used with the limit $i\chi\rightarrow\infty$; either way only one direction of tunneling events may be included. Using Eq.(19) in the waiting time definition from Eq.(18) yields the same form of the WTD as developed by Brandes, except now the exponent $\mathbf{L}(0)-\mathbf{J}_{F}$ is obtained from the $n$-resolved master equation; when the transport is unidirectional the two methods coincide. The only difference between them is that when the jump operators can be defined in terms of a standard master equation Brandes’ method can calculate both $w_{F}(\tau)\text{ and }w_{B}(\tau)$, whereas Eq.(18) cannot.

The Laplace transform of the WTD is

[TABLE]

from which the cumulants are calculated:

[TABLE]

We have included the notation $\mathbf{L}_{0}=\mathbf{L}-\mathbf{J}$ as shorthand for the exponent in Eq.(17).

The average current can be reconstructed from the average waiting times:

[TABLE]

where $\langle\tau\rangle_{F}\text{ and }\langle\tau\rangle_{B}$ are the average forward and backward waiting times, respectively. Waiting times calculated from an $n$-resolved only satisfy Eq.(22) for unidirectional transportRudge and Kosov (2018). The equality in Eq.(22) holds regardless of whether the transport is renewal, which follows from examining the forward and backward currents separately. The forward current, reconstructed from the WTD, is

[TABLE]

where $P(k)$ is the probability that the forward current is $\frac{k}{\langle\tau_{1}+\tau_{2}+\ldots+\tau_{k}\rangle_{F}}$. The average waiting times are defined via the joint probability distribution:

[TABLE]

where we have made no assumption on the renewal nature of the transport. Eq.(23) then reduces to

[TABLE]

assuming that $\{P(k)\}$ spans the entire current space. The same reasoning follows for the backward current.

The higher order current cumulants, in comparison, are only reconstructable from the WTD for unidirectional transport; for example, the randomness parameter

[TABLE]

is equal to the Fano factor defined in Eq.(15) for unidirectional transort satisfying the renewal assumption. The Pearson correlation coefficient measures the strength of correlations between successive waiting times:

[TABLE]

When $p$ is nonzero, the renewal assumption is violated and the equality between Eq.(29) and Eq.(15) no longer applies.

The first moment of the joint distribution $\langle\tau_{1}\tau_{2}\rangle$ is calculated similarly to the cumulants:

[TABLE]

II.4 First-passage time distribution

Unlike the WTD, which calculates the time interval between successive electron tunnelings to the drain, the FPTD calculates the time interval between a tunneling to the drain and when the jump number next reaches a specified integer $n$. Necessarily, the FPTD includes backward tunneling events, and it is defined from the $n$-resolved master equation. To derive the FPTD we follow the method neatly encapsulated by Saito and Dhar in the framework of Markovian master equations Saito and Dhar (2016), and explicitly explained by Ptaszynski Ptaszyński (2018).

We return to the $\chi$-dependent Liouvillian defined in Eq.(8) and start the time-evolution in some arbitrary state $\mathbf{P}(0)$. The $n$-resolved probability vector is written using Eq.(7) and Eq.(9):

[TABLE]

Now, $\mathbf{T}(n|t)=\frac{1}{2\pi}\int_{0}^{2\pi}d\chi e^{-in\chi}e^{\mathbf{L}(\chi)t}$ is the evolution operator for $\mathbf{P}(n,t)$, which in Laplace space is

[TABLE]

The matrix element $\left[\mathbf{T}(n,t)\right]_{kl}$ is the conditional probability that, given the system is initially in state $l$ at time $t=0$, it will be in state $k$ at time $t$ and the jump number will have reached $n$ for the first time. Saito and Dhar Saito and Dhar (2016), PtaszynskiPtaszyński (2018), and Singh et al.Singh et al. (2018) have all found exact expressions for $\tilde{\mathbf{T}}(n|z)$ for simple single reset systems, in which a tunneling to the drain leaves the molecule empty. In general, though, evaluating Eq.(35) requires a numerical contour integration, which is the approach we take.

Alternatively, $\mathbf{P}(n,t)$ is obtainable by considering the column vector of first-passage time probabilities $\mathbf{F}(n|\tau)$. Here, $\left[\mathbf{F}(n|\tau)\right]_{k}$ is the probability density that at time $\tau$ the system is in state $k$ and $n$ extra electrons have been recorded for the first time in the drain. Premultiplication by the evolution operator $\mathbf{T}(0,t-\tau)\mathbf{F}(n,\tau)$ transitions the system in the interval $[\tau,t]$ to a final probability vector, during which no total change to $n$ is made. This contributes to the probability vector $\mathbf{P}(n,t)$; the final step integrates over all possible first-passage times:

[TABLE]

From here it is easier to work in Laplace space. The RHS of Eq.(36) is a convolution, and using Eq.(34) one obtains

[TABLE]

In Laplace space, the total first-passage time probability is

[TABLE]

As Ptaszynski Ptaszyński (2018) notes, in order to connect the FPTD to the WTD, the initial state $\mathbf{P}(0)$ needs to be the normalised system state just after a forward jump to the drain has occurred:

[TABLE]

Combining this choice of $\mathbf{P}(0)$ with Eq.(38), we obtain

[TABLE]

where $F(n|\tau)$, defined as the inverse Laplace transform of Eq.(40) $\mathcal{L}^{-1}\left\{\tilde{F}(n|z)\right\}$, is the probability that, given an initial tunneling to the drain occurred at some arbitrary time in the stationary state, the jump number first reached $n$ after a time interval $\tau$. For unidirectional transport, $F(1|\tau)$ exactly reproduces the WTD defined in Eq.(17).

The FPTD in time-space is obtained from the numerical inverse Laplace transform of Eq.(40):

[TABLE]

We focus, however, on the cumulants of the FPTD, which are again calculated from $\tilde{F}(n|z)$:

[TABLE]

where the limit $z\rightarrow 0^{+}$ is required as $\mathbf{L}(\chi)$ is singular for $\chi=\{0,2\pi\}$ Ptaszyński (2018). We have introduced the ∗ notation for FPTD cumulants to distinguish them from WTD cumulants, and used Ptaszynski’s notation of $\langle\langle\tau^{k}_{n}\rangle\rangle^{*}$ referring to the $\text{k}^{th}$ cumulant of $F(n|\tau)$.

Ptaszynski has recently shown that when the transport is renewal the FPTD cumulants are linearly related:

[TABLE]

as well as deriving, similarly to the WTD, renewal relations between the FCS and the first-passage time cumulants Ptaszyński (2018). The average current is then

[TABLE]

Unlike the WTD, where the equality in Eq.(22) is satisfied regardless of renewal behaviour, the equality in Eq.(45) is only true for renewal processes; the logic applied in Eq.(26) is not available to first-passage times since, for non-renewal processes, $F(2|\tau)\neq F(1|\tau;1|\tau^{\prime})$.

One can also define a randomness parameter $R^{*}_{n}$ for each distribution $F(n|\tau)$:

[TABLE]

For renewal systems they all reduce to $R^{*}_{1}$, which is then equal to the Fano factor. Unlike $R$, we expect that the renewal relationships will hold in bidirectional transport. The Pearson correlation coefficient between subsequent first-passage times can be calculated using $F(2|\tau)$ Ptaszyński (2018):

[TABLE]

Throughout the majority of the paper we will not use the subscript $n$ for the cumulants, with the implication being that we are referring to cumulants of $F(1|\tau)$ unless explicitly stated.

III Model systems

In this section we introduce the three transport scenarios under consideration and define the rate equations describing their dynamics. Quantum rates are calculated using the T-matrix method, which is a suitable replacement for more rigorous approaches when the transport is Markovian; that is, when $T_{\alpha},V_{SD}\gg\gamma$, where $\gamma$ is the broadening of the molecular energy levels, $T_{\alpha}$ is the temperature of electrode $\alpha$, and $V_{SD}$ is the source-drain voltage bias Thielmann et al. (2005b). In order to use rate equations as opposed to comprehensive master equations Gurvitz et al. (1996); Gurvitz and Prager (1996); Gurvitz (1998); Li et al. (2005); Harbola et al. (2006); Nazarov (1993); Davies et al. (1993) we also operate under the secular approximation; coherence between molecular basis states are ignored. All scenarios describe the transport of electrons through a molecule weakly coupled to two macroscopic metal electrodes. The total Hamiltonian is

[TABLE]

where $H_{M}$ is the Hamiltonian of the nanoscale quantum system; $H_{S}\text{ and }H_{D}$ are the source and drain Hamiltonians

[TABLE]

and $H_{T}$ is the system-electrode interaction

[TABLE]

The operator $a_{\alpha,k,\sigma}^{\dagger}$ creates an electron in electrode $\alpha$ in the single-particle state $k$ with spin $\sigma=\left\{\uparrow,\downarrow\right\}$ and free energy $\varepsilon_{\alpha,k,\sigma}$, while $a_{\alpha,k,\sigma}$ annihilates an electron with spin $\sigma$ from the same state in electrode $\alpha$, and $t_{\alpha,k,\sigma}$ are the tunneling amplitudes between the quantum system and the electrodes.

The molecule is described by the Anderson-Holstein model; an energy level, which can be spin split, is coupled to a vibrational mode such that electrons occupying the molecule experience electron-electron and electron-phonon interactions.

[TABLE]

The operators $a^{\dagger}_{\sigma}\text{ and }a_{\sigma}$ create and annihilate an electron in the single-particle state $\sigma$ with free energy $\varepsilon_{\sigma}$, while $b^{\dagger}\text{ and }b$ create and annihilate phonons with frequency $\omega$. The particle number operators are $n_{\sigma}=a^{\dagger}_{\sigma}a_{\sigma}$ and $n=\sum\limits_{\sigma}a^{\dagger}_{\sigma}a_{\sigma}$. The $\lambda(b^{\dagger}+b)n$ term describes the electron-phonon interaction, where $\lambda$ is the coupling strength. Likewise, $U$ is the Coulomb repulsion strength.

The transition rate from the electrode-system many-body state $i$ to many-body state $f$ is

[TABLE]

where the T-matrix is

[TABLE]

The leading order approximation is Fermi’s golden rule, and describes sequential tunneling:

[TABLE]

Fermi’s generalised golden rule Bruus and Flensberg (2002); Timm (2008) details the higher order terms; for example, cotunneling events are described by going to fourth order in $H_{T}$.

III.1 Anderson impurity

For a rigid molecule with no electron-phonon coupling $\lambda=0$, Eq.(51) reduces to the Anderson model:

[TABLE]

The molecular orbital is then in one of four system states; it can be empty $|0\rangle$, occupied by a single electron of either spin $|\sigma\rangle$, or occupied by two electrons of opposite spin $|2\rangle$, with probability vector $\mathbf{P}(\chi,t)=\left[P_{0}(\chi,t),P_{\uparrow}(\chi,t),P_{\downarrow}(\chi,t),P_{2}(\chi,t)\right]^{T}$. Including cotunneling and sequential tunneling events, the Liouvillian minus the forward and backward jump operators is

[TABLE]

with the forward and backward jump operators defined as

[TABLE]

The total sequential rates are

[TABLE]

the total elastic cotunneling rates are

[TABLE]

and the total inelastic cotunneling rates are

[TABLE]

The sequential rates are

[TABLE]

where \gamma^{\alpha}=2\pi\big{|}t_{\alpha}\big{|}^{2}\rho(\varepsilon_{\alpha}); $\rho(\varepsilon_{\alpha})$ is the density of states for the source and drain electrodes, which is assumed to be constant; and $n_{F}(\varepsilon-\mu_{\alpha})$ is the Fermi-Dirac distribution:

[TABLE]

During cotunneling, the electron undergoes two tunneling events in the same quantum process through an intermediate “virtual” energy level in the molecule. Elastic cotunneling leaves the molecule in the same energetic state as before the process, while inelastic cotunneling leaves the molecular orbital occupied by the same number of electrons but in an excited or de-excited state. The notation $\Gamma_{kl}^{\alpha\beta}$ implies that an electron cotunneled from electrode $\alpha$ to electrode $\beta$, changing the system from state $l$ to state then $k$; if the process is elastic $l=k$. Elastic cotunneling processes do not appear in Eq.(1) as they do not affect the molecular probabilites $\mathbf{P}(t)$, while elastic cotunneling events involving the same electrode $\alpha=\beta$ are not included even in the $n$-resolved master equation as they neither change the jump number or the molecular state.

The cotunneling rates are calculated by going to fourth order in $H_{T}$ in Eq.(52) and Eq.(53). To avoid divergences in the T-matrix calculated cotunneling rates, and overcounting of sequential processes, one must use the regularisation procedure noted first by Averin Averin (1994) and detailed separately by Koch et al. Koch et al. (2006, 2004) and Turek and Matveev Turek and Matveev (2002). The details and specific form of the lengthy cotunneling rates for an Anderson impurity are in Appendix A.

III.2 Holstein model

In the limit $U\rightarrow\infty$, the molecule cannot be doubly occupied and behaves as a single level coupled to a vibrational mode: the Holstein model, with the Hamiltonian

[TABLE]

After the Lang-Firsov transformation Lang and Firsov (1963) the molecular Hamiltonian reduces to

[TABLE]

with $\varepsilon=\varepsilon_{0}-\frac{\lambda^{2}}{\omega}$ and new fermionic $\tilde{a}^{\dagger}\tilde{a}$ and bosonic $\tilde{b}^{\dagger}\tilde{b}$ operators obtained from the original operators via the canonical transformation. The molecular eigenstate is now described by two quantum numbers $|nq\rangle$, where $n=\left\{0,1\right\}$ is the electron occupation number and $q=\left\{0,1,...,+\infty\right\}$ is the vibrational quanta occupation number. The associated eigenenergy is $E_{nq}=\varepsilon n+\omega q$. The Fourier transformed vector is then

[TABLE]

which is described by the master equation Mitra et al. (2004)

[TABLE]

Evidently, there can be any number of interacting vibrations and $\mathbf{P}(\chi,t)$ is a vector of infinite length. In practice, however, we only include a finite number $N$ of vibrations in the transport calculations, so that $\mathbf{P}(\chi,t)$ and $\mathbf{I}$ are vectors of length $2(N+1)$. Likewise, the Liouvillian and two jump operators are $2(N+1)\times 2(N+1)$ matrices. $\mathbf{J}_{F}$ and $\mathbf{J}_{B}$ are explicitly written in Appendix B.

The rates are

[TABLE]

where the Franck-Condon matrix elements are

[TABLE]

We consider two distinct vibrational scenarios. So far we have described transport when the phonons are unequilibrated, but when the phonons are in equilibrium with an external bath at temperature $T_{V}$, the rate equation reduces to that of a single resonant level, with Liouvillian

[TABLE]

and jump operators

[TABLE]

where the probabilities are defined with the ansatz

[TABLE]

the total transition rates are $T_{kl}=\sum\limits_{\alpha}T_{kl}^{\alpha}$; and the individual transition rates are

[TABLE]

IV Results

In this section we compare the full counting statistics, waiting time distribution, and first-passage time distribution for the three transport scenarios.

IV.1 Sequential tunneling and cotunneling through an Anderson impurity

In a previous paper, we examined, via waiting times, the effect cotunneling has on strong electron-electron interactions in an Anderson impurity Rudge and Kosov (2018). In order to include elastic cotunneling events, it is necessary to use the $n$-resolved master equation to calculate the waiting times. Unfortunately, when bidirectional transitions are included the definition of the idle-time probability in Eq.(19) is no longer satisfied Thomas and Flindt (2013). One may of course treat the transport as unidirectional and calculate the waiting times in either direction; Fig.(1a) shows that doing this, in a regime where backtunneling processes are present, does not reproduce the correct cotunneling current from $\frac{1}{\langle\tau\rangle_{F,co}}-\frac{1}{\langle\tau\rangle_{B,co}}$. The first-passage time prediction $\frac{1}{\langle\tau\rangle_{co}^{*}}$, on the other hand, reproduces the cotunneling current; the black and red dashed lines in Fig.(1a) coincide. Over the voltage range $2-3\text{meV}$ there is a small numeric difference between $\langle I\rangle_{co}$ and $\frac{1}{\langle\tau\rangle_{co}^{*}}$, since the transport is non-renewal, but it is not visible to the naked eye. In comparison, the sequential current is exactly reproduced by $\frac{1}{\langle\tau\rangle_{F,seq}}-\frac{1}{\langle\tau\rangle_{B,seq}}$, since sequential tunneling is described completely by a standard master equation, and $\frac{1}{\langle\tau\rangle_{seq}^{*}}$, since the sequential transport is always renewal. One can see this in Fig.(1a), as all solid lines coincide.

These difficulties previously kept us to the large bias voltage regime. It is well known in the literature, however, that in a certain voltage range inelastic cotunneling processes cause telegraphic switching in an Anderson impurity, which is now accessible to fluctuating time statistics via the FPTD. The Fano factor in Fig.(1b) is a reproduction of Fig.(3) from Ref.[Thielmann et al., 2005a], with the first-passage and waiting time randomness parameters included as well. The molecule is in the Coulomb blockade regime $U\gg V_{SD}$, so that transport through one $\sigma$ level blocks transport through the other $\bar{\sigma}$ level. At small voltages, for sequential tunneling, Belzig Belzig (2005) has shown that the molecule goes through cycles of transferring $n-1$ electrons quickly through the $\uparrow$ level until it gets “stuck” in the $\downarrow$ level. The cycle completes when the $\downarrow$ empties via thermal broadening of the drain occupation function. The noise is thus characteristic of a sum of Poissonian processes, each transferring $n$ electrons, and weighted by $(\frac{1}{2})^{n}$. In contrast, when cotunneling processes are included, super-Poissonian noise arises from a different mechanism. Note that from here we have not included cumulants from the forward cotunneling WTD, such as $R$ and $p$, since in this voltage regime the transport is bidirectional.

At low bias $V_{SD}/2<\varepsilon_{\uparrow}$ the only transport mechanism is elastic cotunneling through the empty $\uparrow$ level or the occupied $\downarrow$ level; the noise is Poissonian and $F=1$. When $\varepsilon_{\downarrow}>V_{SD}/2\geq\varepsilon_{\uparrow}$ the $\uparrow$ level may be occupied via sequential tunneling, and there are now two main transport channels: elastic cotunneling through the $\downarrow$ level and sequential tunneling through $\uparrow$ level. Telegraphic switching between these two processes occurs when inelastic cotunneling induces a spin flip $\sigma\rightarrow\bar{\sigma}$ or the molecule relaxes from $\uparrow\rightarrow\downarrow$ via sequential tunneling. Super-Poissonian noise is a symptom of the switching between two Poissonian processes, which is the main observation of Fig.(3) from Ref.[Thielmann et al., 2005a] and is evident in Fig.(1b). In this voltage regime the transport is evidently non-renewal as $R^{*}$ does not match the Fano factor; there should thus be correlations between subsequent first-passage times.

Positive correlations in the waiting times are known to occur for systems undergoing telegraphic switchin: for example, in systems with spin-polarised electrodes Ptaszyński (2017). The current switches between two transport channels, each with a characteristic waiting time. When in the faster channel, waiting times between successive tunnelings are more likely to be shorter, and vice versa in the slower channel. Consequently, waiting times between successive tunnelings to the drain are positively correlated; one would expect that there would be similar positive correlations between the first-passage times. In Fig.(2), however, we report intriguing negative correlations between successive first-passage times. That is, if the jump number reaches $+1$ after a short time-interval, it is more likely that the time until the jump number next reaches $+1$ will be long. The values are small enough, however, that the first-passage times are still not considered correlated. A possible explanation is that the time it takes for the transport to switch between the two mechanisms is comparable to the average first-passage times for the two channels, so that the system does not spend long enough in either channel for subsequent first-passage times to become positively correlated.

IV.2 Vibrations with a finite bias voltage

One of the authors (Kosov) has previously reported positive correlations in successive waiting times for transport through the Holstein model Kosov (2017a), which we now compare to first-passage times.

We consider transport through the Holstein model with fixed temperature and over a range of $V_{SD}$. In such a scenario, Koch and von Oppen have previously explained the different behaviour of equilibrated and unequilibrated current at low voltages Koch and Von Oppen (2005), which originates from $X_{qq^{\prime}}$ being suppressed on diagonal and close-to-diagonal elements under a strong electron-phonon coupling $\lambda\gg 1$. At small $T$ and eV, when the electrons are not energetic enough to reach highly excited phonon states, the current is thus suppressed, which is evident more for equilibrated phonons than out-of-equilibrium phonons, as multiple unequilibrated transitions can produce the required excited phonon states. At higher voltages, the current also displays the characteristic steplike structure, with steps occurring when the voltage is a multiple of $2$, since $\omega=1$.

The total electronic current is correctly predicted from the WTD in Eq.(22) for both equilibrated and unequilibrated phonons; the blue and black dashed and solid lines in Fig.(3a) exactly overlap. In comparison, $\frac{1}{\langle\tau\rangle^{*}}$ correctly predicts the current over all voltages for equilibrated but not unequilibrated phonons, implying that at the voltages where they differ the transport is non-renewal. One can see that the red dashed line in Fig.(3a) coincides with both the blue and black dashed lines, but that there is a difference between the red solid line and the blue and black solid lines. Unlike the current plots for the other transport scenarios we consider, the discrepancy between the unequilibrated $\frac{1}{\langle\tau\rangle^{*}}\text{ and }\langle I\rangle$ due to renewal processes is large enough to be visually recognisable. As Fig.(4) shows, however, the largest correlations occur at voltages where the current, and thus the difference between $\frac{1}{\langle\tau\rangle^{*}}\text{ and }\langle I\rangle$, is negligible; consequently, we turn to the randomness parameters to analyse the non-renewal behaviour.

As found in Ref.[Kosov, 2017a], when phonons are unequilibrated, there exist strong correlations between successive waiting times, due to an elastic “shortcut” channel that opens in a small voltage range. In Fig.(4) we report similar correlations in the first-passage times over the same voltage range, albeit with a smaller peak. This is expected; due to backtunneling events the first-passage time can never be smaller than the waiting time $\tau^{*}\geq\tau$, and so correlations between first-passage times must be smaller. The Fano factor and randomness parameters in Fig.(3b) support this non-renewal behaviour. At low voltages $R$ fails to capture the divergence of the Fano factor and $R^{*}$. Between the $1\rightarrow 2\text{eV}$ voltage range $F$ and $R^{*}$ differ, correctly predicting the non-renewal behaviour in Fig.(4). In comparison, $R$ actually crosses the $F$ at $\sim 1.9\text{eV}$, implying renewal behaviour, while at the same voltage the waiting time Pearson coefficient is still a non-zero $p\approx 0.05$. The correlations disappear when the phonons are forced to relax to equilibrium immediately, as the excited phonon state providing the elastic shortcut channel is no longer occupied.

IV.3 Vibrations with a finite temperature gradient

At finite temperatures, and when the electrode Fermi energies are approximately in resonance with $\varepsilon+\omega q$, the bidirectional transport is dominated by thermal noise. A recent experimental paper has shown, however, that electronic noise distinct from the thermal noise can be generated from a temperature gradient $\Delta T$ in the absence of a bias voltage, which they label delta-T noise Lumbroso et al. (2018). We finish our investigation by exploring the effect a temperature gradient has on the first-passage times for sequential tunneling through the Holstein model.

Since the transport in this scenario is described by a standard master equation, the exact current matches the prediction from the WTD in Eq.(22), which is shown in Fig.(5a); all black and blue lines coincide. The FPTD current prediction from Eq.(45), on the other hand, differs from $\langle I\rangle$ at low $\Delta T$, implying non-renewal behaviour in this regime; otherwise, the red FPTD lines coincide with all blue and black lines. Interestingly, the equilibrium current is larger than the non-equilibrium current, which is understood by examining $X_{qq^{\prime}}$ and the phonon occupation probabilities. For $\lambda=1$ and $\omega=1$ the Franck-Condon factor is largest for elastic $q=0$ transitions. These transitions are more accessible to phonons in equilibrium, as the $q=0$ state has a larger relative occupation probability for equilibrium phonons than non-equilibrium phonons:

[TABLE]

We note that for large $\Delta T$ the total current is the same order as that induced by a voltage bias in Fig.(3a).

In Fig.(5b) $F$ diverges as $\Delta T\rightarrow 0$, which is expected since thermal processes dominate in this regime. At larger $\Delta T$ the noise no longer diverges but is still super-Poissonian, and distinct from the shot noise generated by a voltage bias. We attribute this to the relatively large temperatures; even the smallest drain temperature we consider is large enough to cause backtunnelings from the drain, hence the super-Poissonian partition noise normally associated with a voltage bias.

It is evident that $R^{*}\neq F$ at small $\Delta T$, implying that the transport is non-renewal in this temperature range. Despite this non-renewal behaviour, we have chosen not to display the Pearson correlations: they are non-zero but small. At larger $\Delta T$ when the renewal assumption is satisfied, $R^{*}$ successfully predicts $F$, while $R$ clearly does not as it only includes unidirectional transitions.

V Conclusion

In this paper we have compared the waiting time statistics to the first-passage time statistics for a molecule with strong electron-electron or electron-phonon interactions. All methods rely on a form of the master equation, in which the on-site interactions are treated exactly and the coupling between the molecule and electrode is treated perturbatively under the Born-Markov approximation. In bidirectional transport scenarios we have demonstrated the effectiveness of using first-passage times as opposed to waiting times. In particular, we compared the total current to the current computed from the WTD and FPTD, the Fano factor to the WTD and FPTD randomness parameters, and the correlation between subsequent waiting and first-passage times. We first examine sequential tunneling and cotunneling through an Anderson impurity. When a full transport description requires use of the $n$-resolved master equation, the WTD can be calculated in the forward tunneling regime only. Such difficulties arise when including elastic cotunneling, but do not apply to first-passage times; we demonstrate that there are unexpected, minor negative correlations between the first-passage times even though the transport displays telegraphic switching. Next, we analysed the Holstein model at finite bias voltage, in which correlations between successive first-passage times are smaller than correlations between successive waiting times when backtunneling events are present. Finally, we considered transport through the Holstein model at zero bias voltage and a finite temperature gradient. The temperature gradient induces a current, but also significant backtunneling events. Hence, the forward waiting time $\tau_{F}$ was incapable of accurately predicting the Fano factor, as opposed to the first-passage time $\tau^{*}$.

Appendix A Cotunneling rates

Cotunneling rates calculated from the T-matrix approach are regularised according to the procedure outlined by Koch et al Koch et al. (2006, 2004). Calculating cotunneling rates for transport through an Anderson impurity, in particular, has been thoroughly covered by in Ref.[Kaasbjerg and Belzig, 2015; Rudge and Kosov, 2018], to which we direct the reader for a comprehensive overview.

The elastic cotunneling rates are

[TABLE]

where $E_{1}\text{ and }E_{2}$ are derived from the cotunneling pathways involved in the rate, and the $\pm$ is only positive for elastic tunneling through an initially empty or initially doubly occupied system. Similarly, inelastic cotunneling rates have the general form

[TABLE]

where the broadening of the intermediate virtual state $i\eta$ in Eq.(116) and Eq.(117)was first noted by Averin.

The rates are then split into a sum via:

[TABLE]

All that remains is to compute the two types of integrals:

[TABLE]

Here, $E^{\alpha}_{\pm,k}=1/2\>\pm\>i\frac{\mu_{\alpha}-E_{k}}{2\pi T}$; $\psi(x)$ and $\psi^{(1)}(x)$ are the digamma and trigamma functions, respectively; and the Bose-Einstein distribution is

[TABLE]

Appendix B Jump operators for the Holstein model

The jump operators defined from the master equation in Eq.(100) and Eq.(101) are

[TABLE]