Disorder-Induced Rectification in a Molecular System

Abstract

In this work, We study on the electron transport of an ensemble of coupled sites that simulates an array of quantum dots or a molecular system. By using the Green’s function technique, we calculate current and shot noise for linear and disordered site arrays. While in the linear case the characteristic I–V curve reveals no current rectification, in the disordered configurations a robust rectification is found, thus indicating an operational regime typical of molecular diodes. Additionally, a negative differential resistance is observed due to the drop of the bias voltage along the structure, which yields to an energy mismatch of neighboring sites. Finally, the Fano factor reveals a stronger transport correlation for positive than for negative bias voltages in the disordered site configuration.

1 Introduction

One of the major goals in nanoelectronics is the development of diodes and transistors on the same scale of nanocontacts size [1, 2]. In order to achieve this, a variety of manufacturing techniques has been developed in the past few years. For instance, advanced lithographic techniques have allowed the design of quantum dots in semiconductor-based structures, which can operate as a single electron transistor [3]. Artificial molecules based on self-assembled semiconductor quantum dots have also been designed [4]. Additionally, in the context of organic junctions, a variety of manufacturing techniques such as break junction [5], scanning tunneling microscope (STM) [6–9], and Langmuir–Blodgett [10] have improved a great deal in the ability to manipulate and measure single-molecule systems.

Since the proposal of Aviram and Ratner [11] for a single-molecule diode, the search for molecules with rectifying characteristics has motivated several studies towards the molecular electronics. Due to its many-body interacting nature, the molecular systems do not admit exact solutions, even for a simple parabolic confining potential [12, 13]. However, methods based on first principles and density functional theory (DFT) calculations have allowed obtaining approximately the electric current in the few-electron interacting devices [14]. Bednarek et al., using an effective interaction show that it is possible to find an exact solution for a two-electron artificial molecule in coupled quantum dots [15].

In this work, we use a phenomenological approach to account for the quantum transport through an artificial molecule conceptually built by coupling single-level quantum dots, as illustrated in Fig. 1a–d.¹ An external bias (source–drain) voltage is applied to drive the system out of equilibrium, thus generating a tunneling current. By considering the same coupling to the right and left leads, we exploit the dot arrangement, which mimics a molecular geometry and affects the current. This assumption contrasts to those of previous works in the literature that found current rectification due to asymmetric bindings to the leads [5, 16]. It is well known that asymmetric couplings of the central molecule to the electrodes result in an uneven distribution of the voltage drop between the contacts, which induces rectification. Here, we show that even for a symmetric coupling, rectification can occur due to the internal degrees of freedom of the molecular site distribution. The rectification we found comes from the interplay between the potential drop along the structure and the site disorder. In order to gain further insights on the underlining quantum physics of the electronic transport, we also exploit the shot-noise signal. Shot noise, an unavoidable fluctuation of the current due to the charge quantization, has been extensively used to access additional information not contained in the average current [17, 18]. A few examples in this regard encompass the fractional charge observations in the quantum effect of Hall [19], correlations by Pauli and Coulomb [20], suppression of the Kondo effect due to ferromagnetic leads [21], and spin-flip effects [22]. Based on Green function technique [23], we find a strong current and noise rectification arising from a disordered geometry.² The Fano factor also reveals a robust Pauli correlation in the charge flow for these configurations. Additionally, a negative differential resistance (NDR) is found for both geometries, which was identified as due to misalignment of the site energies when a source–drain bias voltage is applied.

Fig. 1

(Colored online) Sketch of the systems considered in this study. A molecular structure is attached to both left and right reservoirs. Four distinct site configurations are illustrated: a linear and b–d three distinct disordered distributions. From b to d, the site concentration increases to the left side of the system. This results in an amplification of the rectification characteristics. The parameter η _i denotes the bias voltage drop along the structure for the linear drop profile. For the abrupt drop at the barriers, we have η _i = 1/2 throughout the system. The molecular coupling Γ₀ is the same in both leads

2 Formulation

2.1 Tight-Binding Model

The system illustrated in Fig. 1 will be described separately for each spin component by the hamiltonian H _σ = H _Lσ + H _Rσ + H _Mσ + H _{T σ} , where H _Lσ and H _Rσ give the free electron energy to the left (L) and right (R) leads. Explicitly, we have

$$ H_{\eta\sigma}=\sum\limits_{\mathbf{k}} \varepsilon_{\mathbf{k}\sigma\eta} c_{\mathbf{k}\sigma\eta}^{\dagger} c_{\mathbf{k}\sigma\eta} $$

(1)

where k is the free electron wave vector, σ is the electron spin, and η is the lead index, such that η = L (η = R) for left (right) lead. The operator c _{k σ η} ($c_{\mathbf {k}\sigma \eta }^{\dagger }$) annihilates (creates) one electron with wave vector k and spin σ in the reservoir η. The molecular structure is described by the term H _Mσ which reads

$$ H_{\rm {M}\sigma} = \sum\limits_{i=1}^{N} \varepsilon_{i\sigma} d_{i\sigma}^{\dagger} d_{i\sigma} + \frac{1}{2} \displaystyle\sum\limits_{i,j=1}^{N}\left(t_{ij} d_{i\sigma}^{\dagger} d_{j\sigma} + t_{ji} d_{j\sigma}^{\dagger} d_{i\sigma}\right) $$

(2)

where d _{i σ} $(d_{i\sigma }^{\dagger })$ annihilates (creates) one electron with spin σ at site i. The matrix elements t _{i j} give the hopping between sites i and j. For simplicity, we assume all nonzero t _{i j} equal to t. The nonvanishing t _{i j} elements correspond to the neighboring sites linked by dashed lines in Fig. 1b–d. In particular to the linear configuration, the Hamiltonian H _Mσ assumes a tridiagonal matrix form, with the sites energy ε _{i σ} in the main diagonal and constant parameters t _{i j} = t in the first diagonals below and above the main diagonal. Several disordering configurations were also investigated, as illustrated in Fig. 1b–d. Finally, H _Tσ gives the tunneling coupling between the central system and the reservoirs:

$$ H_{\rm {T}\sigma} = \sum\limits_{\mathbf{k}\eta i} \left( t_{\mathbf{k}\sigma i}^{\eta} c_{\mathbf{k}\sigma\eta}^{\dagger} d_{i\sigma} + t_{\mathbf{k}\sigma i}^{\eta^ *} d_{i\sigma}^{\dagger} c_{\mathbf{k}\sigma\eta} \right) $$

(3)

with $t_{\mathbf {k}\sigma i}^{\eta }$ being the matrix element that couples the site orbital i with lead η. A source–drain bias voltage is applied across the system in order to drive it away the equilibrium, thus generating a current. The voltage drop along the structure is modeled by ε _{i σ} = ε ₀ − η _i e V,³ where ε ₀ is the localized level energy without bias, e V is the source–drain bias voltage, and η _i is a dimensionless parameter that gives the drop of the voltage inside the chain at site i. The left and right chemical potentials are related to each other accordingly to μ _L − μ _R = e V. In particular, we take μ _L = 0 as our reference of energy. In particular, for the linear chain configuration, we assume η _i = i / (N + 1), where N is the total number of sites. For the more disordered structures sketched in Fig. 1b–d, the parameters η _i will depend mainly on the distance between site i and the electrodes. For a linear drop of the bias voltage along the disordered arrays, we assume η ₁ = 1/6, η ₂ = η ₃ = 2/6, η ₄ = η ₅ = η ₆ = 3/6, η ₇ = η ₈ = 4/6, and η ₉ = η ₁₀ = 5/6, which means that sites with approximately the same distance from the electrodes experience the same potential drop. For completeness, we also consider abrupt drop of the bias voltages at the contacts. For this case, we have simply η _i = 1/2 for all i values.

2.2 Current and Shot noise

The electrical current is calculated by using the standard Landauer–Büttiker equation [23]:

$$ I_{\sigma} = \frac{e}{\hbar}\int \frac{d\varepsilon}{2\pi} T_{\sigma}(\varepsilon) [f_{\rm L}(\varepsilon)-f_{\rm R}(\varepsilon)], $$

(4)

where T _σ (ε) = Tr{T}. The coefficient transmission matrix can be written as $\mathbf {T}=\{{\Gamma }_{\rm {L}\sigma } G_{\rm {M}\sigma }^{r} {\Gamma }_{\rm {R}\sigma } {G}_{\rm {M}\sigma }^{a}\}$.⁴ The function f _L(R)(ε) is the Fermi distribution function of the left (right) reservoir. The shot noise can also be computed from the following expression [17, 22]:

$$\begin{array}{rll} S_{\sigma} &=&\frac{2 e^{2}}{\hbar}\int \frac{d\varepsilon}{2\pi} \text{Tr}\{ [f_{\rm L}(1-f_{\rm L})+f_{\rm R}(1-f_{\rm R})]\mathbf{T} \\ && \qquad\quad +(f_{\rm L}-f_{\rm R})^{2} \mathbf{T}[1-\mathbf{T}]\}. \end{array} $$

(5)

Due to the fact that only the first and the last sites of the chain couple to their respective electrodes, a further simplification to the noise expression can be done, thus resulting in a factorized form as follows:

$$\begin{array}{rll} S_{\sigma} &=&\frac{2 e^{2}}{\hbar}\int \frac{d\varepsilon}{2\pi} \{ [f_{\rm L}(1-f_{\rm L})+f_{\rm R}(1-f_{\rm R})] T_{\sigma}(\varepsilon) \\ && \qquad\quad +(f_{\rm L}-f_{\rm R})^{2} T_{\sigma}(\varepsilon)[1-T_{\sigma}(\varepsilon)]\}, \end{array} $$

(6)

where T _σ (ε) = Tr{T}. A detailed derivation of this result is presented in Appendix 1.

The retarded Green function appearing in the transmission coefficient is given by

$$ \mathbf{G}^{r}_{\rm {M}\sigma}(\varepsilon) =\left[\varepsilon \mathbf{I} -\mathbf{H}_{\rm {M}\sigma} -\mathbf{\Sigma}_{\sigma}^{r}\right]^{-1}, $$

(7)

where $\mathbf {\Sigma }_{\sigma }^{r} = \mathbf {\Sigma }^{r}_{\rm {L}\sigma } + \mathbf {\Sigma }^{r}_{\rm {R}\sigma }$ is the self-energy due to the leads. More explicitly, we have

$$ \mathbf{\Sigma}^{r}_{\eta\sigma} = \mathbf{\Lambda}_{\eta\sigma} -\frac{i}{2} \mathbf{\Gamma}_{\eta\sigma}, $$

(8)

with Λ _{η σ} and Γ _{η σ} being the real and imaginary parts of the self-energy. We assume the wideband limit [23],⁵ where the real part can be neglected, and the matrix elements of the imaginary part are given by

$$ [\mathbf{\Gamma}_{\eta\sigma}]_{ij} = 2 \pi \sum\limits_{\mathbf{k}} t_{\mathbf{k}\sigma i}^{\eta^{\star}} t_{\mathbf{k}\sigma j}^{\eta} \delta(\epsilon - \epsilon_{\mathbf{k} \sigma \eta}) $$

(9)

In our particular case, we attach only the first site to the left lead and the last site to the right lead. So, the only nonzero elements of Γ _{η σ} are [Γ_Lσ ]₁₁ = Γ_Lσ and [Γ_Rσ ] _{N N} = Γ_Rσ . The advanced Green function $\mathbf {G}^{a}_{\rm {M}\sigma }(\varepsilon )$ is given by the Hermitian conjugate of $\mathbf {G}^{r}_{\rm {M}\sigma }(\varepsilon )$. Finally, we define the Fano factor as γ = S / (2e I), which provides information regarding transport correlations [22], where S = S _↑ + S _↓ and I = I _↑ + I _↓. Since the present model is spin-independent, we obtain simply I _↑ = I _↓ and S _↑ = S _↓.

In what follows, we assume Γ_Lσ = Γ_Rσ = Γ₀ and take Γ₀ as our energy scale. In order to keep the generality of our results, we present the current and the noise in units of I ₀ and S ₀, respectively, which are quantified by a general energy parameter Γ₀. More specifically, we have $I_0=e\Gamma _0/\hbar$ and $S_0=2e^{2}\Gamma _0/\hbar$. From the experimental point of view, we typically find Γ₀ ∼ 1. 0 meV in molecular systems [26], which results in I ₀ ∼ 0. 24 μA and S ₀ ∼ 0. 77 × 10⁻²⁵ A ²/Hz. Moreover, in semiconductor quantum dot systems, we observe Γ₀ ∼ 1–100 μeV [27] which turns in I ₀ ∼ 0. 24–24 nA and S ₀ ∼ 0. 77 × 10⁻²⁸–77 × 10⁻²⁸ A ²/Hz. The additional parameters adopted in the present study are the site energy Γε ₀, the tunneling rate Γ₀, and the hopping matrix elements t _{i j} = t. In order to optimize the diode effect, the set of parameters were chosen as in Table 1.

Table 1 Phenomenological parameters adopted in the numerical calculation

3 Results

Figure 2 shows the transmission coefficient against energy for all configurations illustrated in Fig. 1. The left panels correspond to positive bias voltage (e V > 0), while the right panels, to negative bias voltage (e V < 0). The vertical dotted lines indicate the chemical potentials μ _L and μ _R, with μ _L = 0 set as our energy reference. For positive bias voltage, the conduction window (CW) corresponds to the energy interval [μ _R, μ _L]. In contrast, for negative voltages, the conduction window is the interval [μ _L, μ _R].

Fig. 2

(Colored online) Transmission coefficient against energy for fixed bias. a–b Linear configuration. c–h Disordered cases. Both direct (eV > 0) and reverse (eV < 0) voltages are presented. In the linear configuration, the transmission peaks are equally distributed inside the conduction window (CW ) for both direct and reverse bias. In contrast, asymmetric peak structure emerges inside the CW for direct and reverse bias voltages in the disordered system. This turns in rectification characteristics in the transport. Parameters: εε ₀ = 5Γ₀, t= 2Γ₀, eV = ±7Γ₀, k _B T = 0.01Γ₀

Starting with the linear chain configuration, one can clearly see in Fig. 2a, b that the transmission coefficient profiles inside the conduction window (CW) are exactly the same for both direct and reverse bias voltages. In both cases, we find a small peak inside CW range and a larger peak at μ _L or μ _R, for e V > 0 and e V < 0, respectively. The additional peaks of the transmission coefficient lie outside the CW, so they do not contribute to the electronic transport at this specific voltage. These peaks can be moved into the CW by just changing ε ₀ via a gate potential. Due to this symmetric behavior inside the CW region for both e V < 0 and e V > 0, we have no rectification.

In contrast, the disordered configurations introduce asymmetries in the transmission coefficient. Figure 2c, d show the transmission coefficient for the geometry sketched in Fig. 1b. In the same way, Fig. 2e–h presents the transmission coefficient for the geometries in panels c and d of Fig. 1, respectively. This yields a completely distinct profile of T(ε) in the range [μ _R, μ _L] when the sign of the bias voltage is switched. By comparing Fig. 2c, d, we note that for negative bias voltage (e V = − 7Γ₀), only a small peak of T(ε) lies inside the CW. In contrast, for positive bias voltage (e V = 7Γ₀), a couple of peaks arise in the conduction window. In particular, a much higher peak emerges in the CW for e V > 0 than for e V < 0. This contrasting behavior is responsible for the rectification observed in the characteristic I–V curve, as will be described later. Similar asymmetries, responsible for the rectification, can be observed for the others disordered geometries (Fig. 2e–h).

In order to have a more clear picture of the conduction channels in the present systems, we show in Fig. 3a–d a two-dimensional map of the transmission coefficient as a function of bias voltage and energy. The relation between Figs. 2 and 3 is as follows: Figs. 2a, b and 3a correspond to the linear configuration (Fig. 1a); Figs. 2c, d and 3b to the disordered configuration illustrated in Fig. 1b; Figs. 2e, f and 3c to configuration in Fig. 1c; and finally Figs. 2g, h, and 3d to configuration in Fig. 1d.

Fig. 3

(Colored online) Color map of the transmission coefficient as a function of bias voltage and energy for both a linear and b–d disordered configurations. The horizontal dashed line indicates μ _L, while the diagonal dashed line, μ _R. The conduction window in the energy–voltage map is labeled by CW . In the linear case, the transmission presents a more symmetric profile, while the disordered configurations give rise to asymmetric patterns. This asymmetry results in an enhancement of the number of channels that cross the chemical potential μ _L for positive bias. Parameters: ε ₀ = 5Γ₀, t = 2Γ₀, k _B T = 0.01Γ₀

In this plot, the horizontal dashed line denotes the left chemical potential at zero (energy reference), while diagonal dashed line corresponds to μ _R = μ _L − e V. The conduction windows are labeled by CW.

To the linear chain (Fig. 3a), the transmission coefficient reveals a symmetric peak distribution. Close to e V = 0, we find exactly ten peaks due to the split of the localized levels. The pattern is centered at energy ε ₀ = 5Γ₀. Our particular choice of the sites energy at ε ₀ = 5Γ₀ is to have all the conduction channels above the Fermi level of the leads at small bias, as seen in Fig. 3a. When the bias is turned on, these channels start gradually to get into the conduction window and contribute to the current. Larger values of ε ₀ would result at same transport features, however, at higher bias voltages. Additionally, increasing the hopping parameter t, we find a larger split of the channels. Whenever a channel matches μ _L or μ _R, a step in the current arises. The peaks lose amplitude as |e V| is enhanced. This is related to the mismatch of the neighboring energy levels as |e V| enlarges, which suppresses the resonant hopping and consequently the transmission coefficient.

In the disordered distributions (Fig. 3b–d), the transmission coefficient loses symmetry, and the overall pattern is deformed. This is so because the disorder lifts the degeneracy of the site levels, thus introducing anti-crossings in the transmission channels. The loss of symmetry in T(ε) to yields to a higher number of transmission channels matching the conduction window for positive than negative voltages, which favors resonant tunneling for positive biases. We should point out that by changing ε ₀, one can rigidly move all the pattern to higher or lower energies, thus changing the number of levels in both left and right conduction windows.

Figure 4a–h shows the current and shot noise for all the sites configurations. To the linear array (Fig. 4a), the current presents steps whenever a conduction channel crosses μ _L. Similar features are seen for both e V < 0 and e V > 0 with no current rectification, I(V) = − I (− V). This result is in agreement with Figs. 2a, b and 3a. After each step, the current is slightly suppressed, thus presenting a NDR. This is so because for increasing biases, each site level tends to lose resonance with its neighbors, thus decreasing the probability of intersite hopping. In Fig. 4b, we also present the current for linear arrays, but for an abrupt bias, voltage drops at the contacts. No NDR is seen for this particular case because all the site levels ε _i remain aligned to each other as |e V| increases, contrasting to the linear drop regime. Each step in the current arises whenever a transmission channel crosses the Fermi level of the source lead. In the disordered case, we find a much richer behavior. All the left curves present a diode-like effect with I(V) < I(−V) on average. Interestingly, the current attains higher values for e V > 0 when more sites became closer to the left electrode (compare Fig. 4c, e, g). In order to make our geometric analysis more quantitative, we define a linear density of sites λ to the three disordered cases (see insets). To quantify λ, we count the number of sites per unit length along the molecule. The unity length is indicated by the parallel dashed lines in Fig. 1b–d. In Fig. 4c, we have a broader λ. As we move to Fig. 4e, g, λ becomes more localized close to the left lead.

Fig. 4

(Colored online) a–h Current and shot noise against bias voltage for both linear and disordered configurations. The left panels correspond to the linear drop of the bias voltage along the structures, while the right panels, to the abrupt drop of the bias voltage at the contacts. The current and the shot noise present steps whenever a transmission channel crosses the reservoir chemical potential for both potentials profiles. However, after each step in the linear drop case, the current reveals a negative differential resistance. Additionally, the disordered configuration presents a transport rectification, with current and shot noise being much larger for positive than negative bias voltage. Similar rectification is not seen for the abrupt drop of the bias voltage. In the insets in panels c, e, and g, we show the linear density of sites along the junction. Parameters: ε ₀ = 5Γ₀, t = 2Γ₀, k _B T = 0.01Γ₀

The disordering induced rectification can be understood by looking at Fig. 3. Due to the distortion of the transmission channels pattern imposed by the disordering, more channels go into CW area for e V > 0 than e V < 0. This is the main mechanics responsible for the rectification. In the disordered case with an abrupt drop of the bias voltage at the contacts, no rectification is observed. The current simply presents the step-like behavior already found to the linear chain configuration. In addition to the current, in Fig. 4a–h, we present the shot noise, which reveals similar behavior to those observed in the current. In the linear site configuration, we find symmetric shot noise curves, S(V) = S(−V). In contrast, the disordered case reveals an asymmetry with S(V) ≠ S(−V). However, to the abrupt bias voltage drop, no distinction between S(V) and S(−V) is found.

Finally, in Fig. 5a–h, we plot the Fano factor (γ) for all system configurations and the two bias voltage profiles (linear and abrupt drops). In the linear chain, we observe a dip in γ whenever a resonant transport takes place (steps in the current). When a conduction channel aligns with the Fermi level of the source, it increases the transmission coefficient. As a result, it fills a scattering state and consequently suppresses the Fano Factor.⁶ The role of Pauli correlations on transport becomes more relevant as the number of channels in the conduction window enlarges for increasing bias voltage; consequently, the dips at the Fano factor attain lower values. In contrast, in the disordered configurations, the Fano factor is strongly asymmetric. It can attain lower values for e V > 0 than e V < 0 in the linear bias voltage drop regime. In contrast, the Fano Factor is symmetric with respect to e V for any geometry if the bias presents an abrupt drop at the contacts (Fig. 5d, f, h).

Fig. 5

(Colored online) a–h The Fano factor against bias voltage for both linear and disordered cases. The left panels correspond to the linear bias voltage applied along the structures, and the right panels, to the abrupt drop of the bias voltage at the contacts. While the Fano factor in the linear configuration (panel a) exhibits symmetric dips for positive and negative voltages, it shows asymmetric structures for the disordered systems (panels c, e, g). In particular, in the disordered cases, the Fano factor is strongly suppressed for direct bias. In the case of an abrupt bias voltage drop at the contacts (panels b, d, f, h), no asymmetries are observed in the Fano factor. Parameters: ε ₀ = 5Γ₀, t = 2Γ₀, k _B T = 0.01Γ₀

4 Conclusion

In summary, we report on a theoretical analysis of current and shot noise in an array of sites in four distinct configurations: (a) linear and (b) three disordered distributions. We analyze two regimes of the bias voltage: linear drop of the bias voltage along the structures and abrupt drop of the bias voltage at the contacts. In the linear geometry, for both regimes of the bias voltage drop, the current, the shot noise and the Fano factor present similar behavior for both e V > 0 and e V < 0 bias voltages. In the disordered case with an abrupt drop of the bias voltage at the contacts, no rectification is observed. However, in the disordered configuration with a linear drop of the bias voltage, we observe a strong rectification in the transport properties. The current and the shot noise attain much higher values for positive than negative voltages, thus behaving as a diode. Additionally, as can be more clearly seen in Fig. 5c, g, the Fano factor tends to be more suppressed for e V > 0 than e V < 0, which indicates that the correlation in the charge flow is asymmetric with respect to the bias voltage. The present diode is directly related to the internal molecular geometry and not to asymmetries in the left and right coupling to the leads. This opens the potentiality to future development of geometry-based diodes.

Appendices

Appendix 1

One may note that the transmission coefficient element can be written as follows (in the demonstration below, we omit the spin label):

$$ T_{il}= \sum\limits_{mnk} \Gamma^{L}_{im} G^{r}_{mn} \Gamma^{R}_{nk} G^{a}_{kl} $$

(10)

where $\Gamma ^{\rm L}_{im}=\Gamma _{0} \delta _{im} \delta _{i1}$ gives the tunneling rate between site i = 1 and the left (L) lead and $\Gamma ^{\rm R}_{nk}=\Gamma _{0} \delta _{nk} \delta _{nN}$ corresponds to the tunneling between the last site (N = 10) and the right (R) lead. $G^{r}_{mn}$ and $G^{a}_{kl}$ are the matrix elements of the retarded and advanced Green functions. Using explicitly $\Gamma ^{\rm L}_{im}$ and $\Gamma ^{\rm R}_{nk}$ in Eq. (6), we find

$$ T_{il}= \sum\limits_{mnk} \Gamma_{0} \delta_{im} \delta_{i1} G^{r}_{mn} \Gamma_{0} \delta_{nk} \delta_{nN} G^{a}_{kl} $$

(11)

Making use of the Kronecker delta function, we obtain

$$ T_{il}= \Gamma_{0} \delta_{i1} G^{r}_{iN} \Gamma_{0} G^{a}_{Nl} . $$

(12)

Using Eq. (12) into Tr{T T}, we have

$$\begin{array}{rll} \text{Tr} \{ \mathbf{T}\mathbf{T}\} &=& \sum\limits_{i}\sum\limits_{l} \Gamma_{0} \delta_{i1} G^{r}_{iN} \Gamma_{0} G^{a}_{Nl} \Gamma_{0} \delta_{l1} G^{r}_{lN} \Gamma_{0} G^{a}_{Ni} \\ &=& \Gamma_{0} G^{r}_{1N} \Gamma_{0} G^{a}_{N1} \Gamma_{0} G^{r}_{1N} \Gamma_{0} G^{a}_{N1} \\ &=& \left(\Gamma_{0} G^{r}_{1N} \Gamma_{0} G^{a}_{N1} \right)^2. \end{array} $$

(13)

On the other hand, one can easily show that

$$ \text{Tr} \{ \mathbf{T}\} = \Sigma_{i} T_{ii} = \Gamma_{0} G^{r}_{1N} \Gamma_{0} G^{a}_{N1}. $$

(14)

So comparing Eqs. (13) and (14), we conclude that

$$ \text{Tr} \{ \mathbf{T} \mathbf{T}\} = [ \text{Tr} \{ \mathbf{T}\} ]^{2} = T^2. $$

(15)

Appendix 2

It is still possible to find an analytical expression for the transmission coefficient. To do that, we start writing the component $G_{1,N}^{r}$ in terms of the Dyson equation:

$$ G_{1,N}^{r}=g_{1,N}^{r} + \displaystyle\sum\limits_{l,m}^{N} g_{1,l}^{r} \Sigma_{l,m}^{r} G_{m,N}^{r} $$

(16)

where $g_{1,l}^{r}$ is electron retarded green function for the molecule only (no leads), and $\Sigma _{l,m}^{r}=\frac {-i}{2} ( \Gamma ^{L} \delta _{lm} \delta _{l1} + \Gamma ^{R} \delta _{lm} \delta _{lN})$ is the element (l, m) for the tunneling self-energy. With this, we have

$$\begin{array}{rll} G_{1,N}^{r} &=& g_{1,N}^{r} + \displaystyle\sum\limits_{l,m}^{N} g_{1,l}^{r} \left[\delta_{lm} \delta_{l1} \Sigma_{l,m}^{r} + \delta_{lm} \delta_{lN} \Sigma_{l,m}^{r}\right] G_{m,N}^{r} \\ &=& \frac{ g_{1,N}^{r}\left[ 1 + \Sigma_{N,N}^{r} G_{N,N}^{r}\right]}{1 - g_{1,1}^{r} \Sigma_{1,1}^{r}}. \end{array} $$

(17)

To solve this equation, we need to find the Green’s function $G_{N,N}^{r}$. Using the same reasoning developed above, we obtain

$$ G_{N,N}^{r}= \frac{ g_{N,N}^{r} + g_{N,1}^{r} \Sigma_{1,1}^{r} G_{1,N}^{r}}{1 - g_{N,N}^{r} \Sigma_{N,N}^{r}}, $$

(18)

which results in

$$\left[1 - g_{1,1}^{r} \Sigma_{1,1}^{r}\right]G_{1,N}^{r}=g_{1,N}^{r}\left[ 1 + \frac{\Sigma_{N,N}^{r} g_{N,N}^{r} + \Sigma_{N,N}^{r} g_{N,1}^{r} \Sigma_{1,1}^{r} G_{1,N}^{r}}{1 - g_{N,N}^{r} \Sigma_{N,N}^{r}} \right],$$

(19)

or

$$ G_{1,N}^{r}=\frac{g_{1,N}^{r}}{\left[1- g_{1,1}^{r} \Sigma_{1,1}^{r}\right]\left[1- g_{N,N}^{r} \Sigma_{N,N}^{r}\right]- g_{1,N}^{r} \Sigma_{N,N}^{r} g_{N,1}^{r} \Sigma_{1,1}^{r}} $$

(20)

Noting that

$$ T=\text{Tr} \{ \mathbf{T}\} = \Gamma_{0} G^{r}_{1N} \Gamma_{0} G^{a}_{N1}=\Gamma_{0}^{2} \arrowvert G^{r}_{1N} \arrowvert ^{2}, $$

(21)

we finally arrive at

$$ T= \Gamma_{0}^{2} \arrowvert \frac{g_{1,N}^{r}}{\left[1- g_{1,1}^{r} \Sigma_{1,1}^{r}\right]\left[1- g_{N,N}^{r} \Sigma_{N,N}^{r}\right]- g_{1,N}^{r} \Sigma_{N,N}^{r} g_{N,1}^{r} \Sigma_{1,1}^{r}} \arrowvert ^2. $$

(22)