Molecular electronics in junctions with energy disorder

Franz J Kaiser¹, Peter Hänggi and Sigmund Kohler

Theoretische Physik I, Institut für Physik der Universität Augsburg, Universitätsstraße 1, 86135 Augsburg, Germany

¹ Author to whom any correspondence should be addressed.

E-mail: franz.josef.kaiser@physik.uni-augsburg.de

Received 31 March 2008
Published 30 June 2008

Contents

• 1. Introduction
• 2. Wire-lead model and master equation
  • 2.1. Perturbation theory
  • 2.2. Floquet theory
  • 2.3. Master equation and current formula
• 3. Electron transport with slowly fluctuating energies
• 4. ac-driven disordered junctions
• 5. Conclusions
• Acknowledgments
• Appendix. Analytical solution for two levels
• References

1. Introduction

Chemical adsorption of sulfur atoms on a gold surface allows a stable bond between gold tips and thiol groups of molecules.This has been exploited for measuring the conductance and the current–voltage characteristics of gold–molecule–gold junctions [1]–[6]. Repeated measurements even with the same sample, however, reveal small but noticeable differences which possibly stem from environmental fluctuations that impact upon the effective molecule parameters. Moreover, the particular form of the gold tip can have a significant influence on the transport properties [7].

A present line of experimental research is the measurement of molecular conductance when the electrons are excited by electromagnetic waves. There one expects various phenomena ranging from photo-assisted transport to ratchet or non-adiabatic pump effects, i.e. the induction of dc currents by ac fields even in the absence of any voltage bias [8]–[10]. Moreover, it has been predicted that properly tailored laser pulses can give rise to short current pulses [11]–[14]. Since a dc current flows in one particular direction, a ratchet effect can occur only in `sufficiently asymmetric' systems [10]. In that respect, static disorder is sufficient to break the reflection symmetry of an individual realization and, thus, may support a ratchet effect.

The quantitative prediction of the current through a molecule is still a great challenge despite the significant progress achieved in recent years [15]–[18]. For a more qualitative understanding of the mechanisms involved in molecular transport, it is thus advantageous to employ for the molecule a rather generic tight-binding model [9, 10], [19]–[24]. Then a flexible method for the computation of transport properties is provided by master equations of the Bloch–Redfield type which allow one to include electron–electron and electron–phonon interactions, as well as time-dependent fields [10, 25]. Similar methods have also been used for describing incoherent transport [26, 27].

Here, we explore the role of slow fluctuations or static disorder for molecular conductance. Thereby we will assume that the relevant environmental fluctuations are so slow that they can be described as static disorder which defines an ensemble of wire Hamiltonians. Then a natural quantity of interest is the corresponding distribution of stationary currents. A setup for which this current distribution is also directly relevant is an array of molecular junctions that conduct in parallel. We employ a tight-binding model for the molecule and treat it with the Floquet master equation formalism derived in [25] which we review briefly in section 2. In section 3, we present results for a static model with a large voltage bias, while in section 4, we investigate pumping effects caused by an interplay of ac driving fields and disorder. The analytical derivation of the current distribution for a wire with two sites is deferred to the appendix.

2. Wire–lead model and master equation

The system of the driven molecular wire, the leads, and the coupling between the molecule and the leads, as sketched in figure 1, is described by the Hamiltonian

The wire is modelled by N tight-binding orbitals |n, n = 1,  ..., N, such that

with the tunnel matrix element Δ and the capacitive energy U. Each onsite energy E_n(t) + ξ_n contains a random contribution ξ_n that subsumes the influence of environmental fluctuations. We assume that these fluctuations are Gaussian distributed and so slow that we can treat them as static disorder. Thus, the probability that the onsite energy of orbital n lies in an interval of size dξ around E_n(t) + ξ_n reads

where the variance σ² is assumed to be position-independent. This implies that the energy fluctuations are spatially uncorrelated, such that ξ_nξ_n ' = σ²δ_nn '. The onsite energies E_n(t) = E_n + Ax_ncos(Ωt) are modulated by a harmonically time-dependent dipole force, where A denotes the electrical field amplitude multiplied by the electron charge and the distance between neighbouring sites, with the scaled position of site |n. Our goal will be to compute for many realizations of the wire Hamiltonian the resulting dc current which provides the current distribution P(I). The last term in equation (2) captures the electron–electron interaction within a capacitor model and the operator describes the number of excess electrons residing on the molecule. Below we shall assume that U is so large that only states with zero or one excess electron play a role. The first and the last sites of the molecule, |1 and |N, couple via the tunnelling Hamiltonian

to the respective lead. The operator ( ) creates an electron in the left (right) lead in state |Lq which is orthogonal to all wire states. The influence of the tunnelling Hamiltonian is fully characterized by the spectral density . If the lead states are dense and located at the centre of the conduction band, the spectral densities can be replaced by a constant, i.e. we assume Γ_ℓ() = Γ for both leads.

The leads are modelled as ideal Fermi gases

which are initially at thermal equilibrium with the chemical potential μ_L/R and, thus, are described by the density operator , where is the electron number operator for lead ℓ = L, R. Since a typical metal screens all electric fields with a frequency below the plasma frequency, we assume that the bulk properties of the leads are not affected by the laser irradiation.

2.1. Perturbation theory

The derivation of a master equation starts from the Liouville–von Neumann equation for the total density operator ρ(t) for which one obtains by standard techniques the approximate equation of motion [10], [28]–[31]

Here, the first term corresponds to the coherent dynamics of both the wire electrons and the lead electrons, while the second term describes resonant electron tunnelling between the leads and the wire. The tilde denotes operators in the interaction picture with respect to the molecule and the lead Hamiltonian without the molecule–lead coupling, , where U₀ is the propagator without the coupling. The net (incoming minus outgoing) electrical current through the left contact is given by minus the time-derivative of the electron number in the left lead multiplied by the electron charge –e. From equation (6) follows for the current in the wide-band limit the expression

where f_ℓ is the Fermi function of the respective lead and . Furthermore, ··· = tr_wireρ_wire··· denotes the expectation value with respect to the wire density operator. We emphasize that the expectation values in equation (7) depend directly on the dynamics of the isolated wire and are thus influenced by the driving.

2.2. Floquet theory

An important feature of the current formula (7) is its dependence on solely the wire operators while all lead operators have been eliminated. Therefore, it is sufficient to consider the reduced density operator ρ_wire = tr_leadsρ of the wire electrons. The effort necessary to calculate ρ_wire can be reduced significantly by exploiting the fact that the master equation (6) inherited from the total Hamiltonian a periodic time-dependence, which opens the way for a Floquet treatment.

One possibility for such a treatment is to use the Floquet states of the central system, i.e. the driven wire, as a basis. Thereby we also use the fact that in the wire Hamiltonian (2), the single-particle contribution commutes with the interaction term and, thus, these two Hamiltonians possess a complete set of common eigenstates. In analogy to the quasimomenta in Bloch theory for spatially periodic potentials, the quasienergies ε_α come in classes , of which all members represent the same physical solution of the Schrödinger equation. Thus we can restrict ourselves to states within one Brillouin zone such as, for example, .

For the computation of the current it is convenient to have an explicit expression for the interaction picture representation of the wire operators. It can be obtained from the (fermionic) Floquet creation and annihilation operators defined via the transformation , which reads in the interaction picture , with the important feature that the time difference t–t ' enters only via the exponential prefactor [10].

2.3. Master equation and current formula

In the following, we assume the interaction U to be the dominant energy scale in the system, therefore we can restrict the wire Hilbert space to the N + 1 dimensional subspace of states with zero or one electron, such that a basis for the decomposition of the reduced operator is {|0, c_α^† (t) |0}, where |0 denotes the wire state in the absence of excess electrons. Moreover, it can be shown [25] that at large times, the density operator of the wire becomes diagonal in the electron number . Therefore a proper ansatz reads

Note that we keep terms with α≠β, which means that we work beyond a rotating-wave approximation (RWA).

Following our evaluation of the master equation [25], we arrive at a set of N² coupled equations of motion for ρ_αβ(t) which in Fourier representation read

In an analogous manner we obtain for the dc current the expression

The results of this section allow us to numerically compute the dc current through a driven conductor as well as studying the undriven limit. The current distributions discussed below are obtained by computing the dc current for typically 10⁴ realizations of the wire Hamiltonian (2). Then these values are taken for a histogram with 150 bins which finally will be scaled such that we obtain a normalized probability density.

3. Electron transport with slowly fluctuating energies

We first address an undriven wire in the configuration sketched in figure 1 where the distribution of all wire levels is centred at energy E_n = 0. The transport voltage is so large that all eigenenergies lie well within the voltage window and, consequently, the transport is unidirectional. Then in the absence of onsite energy fluctuations (σ = 0), the current can be evaluated analytically within a RWA and reads , i.e. it decays with increasing wire length [32]. The index `max' refers to the fact that any modification of the onsite energies can only reduce the current—which is confirmed by our simulations. The physical reason for this is that for equal onsite energies, solely the kinetic energy determines the eigenstates which, consequently, are delocalized. Different onsite energies, by contrast, tend to `localize' the eigenstates. Thus in the limit of small disorder, the current distribution P(I) is expected to possess a clear peak at I = I_max and some minor contribution for lower values of I.

Figure 2 shows the simulated current distributions for two different variances. For a small variance (panel (a)), the distributions for short wires (N = 2,3) show the expected behaviour. With increasing wire length, the peak at I = I_max disappears and is eventually replaced by an apparently parabolic distribution. This length dependence can be understood in the following way: for a short wire, the probability that a level is out of resonance is rather low and, thus, most realizations of the wire Hamiltonian will allow resonant inter-site tunnelling. With an increasing number of levels, however, the probability of having at least one misaligned level increases and a current significantly smaller than I_max becomes more likely. The precise values will depend on the details and, consequently, we expect a broad distribution. This means that whenever a large number of levels plays a role, the transport through a molecule is extremely sensitive to even small disorder induced by environmental fluctuations.

With a larger variance, this scenario becomes even more pronounced, as can be seen in figure 2(b): then the peak at I_max is rather small and noticeable only for 2 and 3 sites. The most likely realization is a completely disordered wire with an accordingly low conductance. For N > 3, the distributions even possess a significant peak at I = 0 which corresponds to isolating behaviour. A closer inspection of the numerical data reveals that the crossover between conducting and isolating behaviour occurs when the effective disorder exceeds the tunnel matrix element Δ.

Interestingly enough, for N = 2, 3 the distribution turns out to even be non-monotonic, which means that one most likely finds either a current close to the theoretical maximum or a significantly smaller lower value. The non-monotonic behaviour for N = 2 can also be derived analytically. The derivation of the corresponding current distribution (A.3) can be found in the appendix. The excellent agreement of this analytical solution and the simulated distributions emphasize that the simulation with approximately 10⁴ realizations ensures good convergence.

4. ac-driven disordered junctions

In order to investigate the influence of an ac driving, we employ the same model as above, but now with an additional dipole driving modelled by time-dependent onsite energies E_n(t) = Ax_ncos(Ωt) as discussed in section 2. The driving frequency is chosen such that it matches the average splitting of the wire energies, while the amplitude A = Δ corresponds to intermediately strong driving. The solid line in figure 3 shows the current distribution in the absence of a voltage bias, V = 0. The reflection symmetry of the ensemble relates to the symmetric shape of the distribution, which implies that the current vanishes in the ensemble average. An individual realization of the wire, however, generally does not possess reflection symmetry because the random energy shifts are spatially uncorrelated. This asymmetry in combination with the non-adiabatic driving induces a coherent ratchet current, i.e. a dc current even in the absence of any net voltage bias. In the present case, the ratchet current is of the order of 10–20% of the current observed above in the large bias limit. This order of magnitude is typical when the driving frequency or a multiple of the driving frequency lies close to an internal resonance, while the intensity is moderate [25]. In addition to the broad distribution of ratchet currents, P(I) exhibits a peak at I = 0. This stems from realizations for which the driving is well out of resonance.

Figure 3. Current distribution for an ac-driven wire with N = 2 sites for various bias voltages. The fluctuations of the onsite energies are characterized by the standard deviation σ = 0.5Δ, the driving frequency and amplitude are A = Δ and , respectively. All the other parameters are as in figure 2.

For a bias voltage V > 0, the ensemble no longer possesses reflection symmetry and the current distribution is shifted towards positive values (see figure 3). For sufficiently small voltages V  Δ/e, non-adiabatic pumping against the voltage bias is still possible. Rather surprisingly, the peak at zero current remains. It now corresponds to realizations for which on the one hand, the driving is off-resonant while on the other hand, both levels lie outside the voltage window. With increasing bias voltage, the second condition is less frequently fulfilled, and eventually the distribution assumes a form similar to that found for a large voltage in the absence of driving. Already for V≈4Δ/e, the distribution is hardly distinguishable from the one shown in figure 2(a) for a wire with N = 2 sites in the absence of driving.

5. Conclusions

We have investigated the current through a molecular wire with disordered onsite energies. Such disorder can stem from the interaction with slow fluctuations of background charges in the substrate. In particular, we computed the resulting current distribution for two typical cases, namely an `open transport channel' and a driven molecular wire for which random energy shifts break reflection symmetry and, thus, the driving can induce a ratchet current.

The open transport channel is characterized by tight-binding levels with equal onsite energies, such that any misalignment stems from the disorder. Its main consequence is that as soon as the standard deviation of the onsite energies exceeds the tunnel matrix elements, the current distribution no longer peaks only at a finite value, but also at zero. For longer wires, only the peak at zero current remains. This isolating behaviour resembles Anderson localization which is found for electrons in a one-dimensional disordered lattice [33]. Note however, that we here considered short wires far from the scaling limit in which Anderson localization is usually studied.

Since the random energy shifts break reflection symmetry, driving the molecular wire with a laser field induces ratchet currents for which we found a relatively broad distribution. If the driving frequency is far from any molecular excitation energy, the ratchet current will be rather small, and we indeed found that this is the case for very many wire realizations. It has the consequence that the corresponding distribution possesses a spike at zero current. This means that non-adiabatic pumping of electrons against a voltage bias becomes generally impossible whenever the relevant wire energy levels lie well within the voltage window.

In conclusion, our results reveal that slow fluctuations or a static disorder can have a significant effect on molecular conduction. In various cases, the current distribution emerges to be rather flat, which means that even the magnitude of the current depends sensitively on environmental influences.

Acknowledgments

Financial support of the German Excellence Initiative via the `Nanosystems Initiative Munich (NIM)' and of the Elite Network of Bavaria via the International Doctorate Program `NanoBioTechnology' is gratefully acknowledged. This work has been supported by Deutsche Forschungsgemeinschaft through SFB 484 and SPP 1243.

Appendix. Analytical solution for two levels

A wire model with N = 2 sites represents an analytically solvable case for which one observes a non-monotonic current distribution and which can serve as test case for numerical implementations. Here, we consider a two-level system with on-site energies ξ_{1, 2}, i.e. with a bias 2η = ξ₁–ξ₂. Since the random energy shifts ξ_n are Gaussian distributed with variance σ², the bias 2η is also Gaussian distributed but with variance 2σ², i.e. its distribution function reads .

For the computation of the current, we restrict ourselves to the limit of a large transport voltage such that both eigenenergies of the two-level system lie within the voltage window. Then, the Fermi functions of the left and the right lead effectively become f_L = 1 and f_R = 0. In this case, transport can be described within the RWA which practically means that the reduced density operator of the wire is diagonal in energy representation [32]. Within the RWA follows from the master equation (9) the occupation probability ρ_αα = w_α¹/w_α² and, thus, . The coefficients w_αⁿ = |_α|n|² denote the overlap between the eigenstate |_α and the localized state |n (note that in the undriven case, all non-vanishing contributions have sideband index k = 0, such that here the sideband index k can be omitted). Inserting this solution into the current formula (10), we obtain .

The remaining task is now to diagonalize the single-particle Hamiltonian which provides the coefficients w_αⁿ. For bias 2η and tunnelling matrix element Δ, the Hamiltonian in pseudo-spin notation reads H = ησ_z + Δσ_x and possesses the eigenenergies ±δ = ±(η² + Δ²)^1/2. The corresponding eigenvectors _α are proportional to (δ + η, Δ) and (δ–η, Δ), respectively, such that w_α¹/w_α² = (δ±η)²/Δ². Then we obtain for the current the expression

which assumes its maximum in the unbiased limit η = 0.

The probability distribution for the current is related to w(η) via

where the summation considers all values of η that fulfil the condition I = I(η). After some straightforward algebra, we obtain by evaluating expression (A.2) the current distribution

which is defined and normalized on the interval [0, I_max].

References

[1]

Datta S, Tian W, Hong S, Reifenberger R, Henderson J I and Kubiak C P 1997 Phys. Rev. Lett. 79 2530 
CrossRef

[2]

Reed M A, Zhou C, Muller C J, Burgin T P and Tour J M 1997 Science 278 252 
CrossRef

[3]

Kergueris C, Bourgoin J P, Palacin S, Esteve D, Urbina C, Magoga M and Joachim C 1999 Phys. Rev. B 59 12505 
CrossRef

[4]

Cui X D, Primak A, Zarate X, Tomfohr J, Sankey O F, Moore A L, Moore T A, Gust D, Harris G and Lindsay S M 2001 Science 294 571 
CrossRefPubMed

[5]

Reichert J, Ochs R, Beckmann D, Weber H B, Mayor M and von Löhneysen H 2002 Phys. Rev. Lett. 88 176804 
CrossRefPubMed

[6]

Dadosh T, Gordin Y, Krahne R, Khivrich I, Mahalu D, Frydman V, Sperling J, Yacoby A and Bar-Joseph I 2005 Nature 436 677 
CrossRefPubMed

[7]

Yaliraki S N and Ratner M A 1998 J. Chem. Phys. 109 5036 
CrossRef

[8]

Lehmann J, Kohler S, Hänggi P and Nitzan A 2002 Phys. Rev. Lett. 88 228305 
CrossRefPubMed
Lehmann J, Kohler S, Hänggi P and Nitzan A 2003 J. Chem. Phys. 118 3283 
CrossRef

[9]

Platero G and Aguado R 2004 Phys. Rep. 395 1 
CrossRef

[10]

Kohler S, Lehmann J and Hänggi P 2005 Phys. Rep. 406 379 
CrossRef

[11]

Franco I, Shapiro M and Brumer P 2007 Phys. Rev. Lett. 99 126802 
CrossRefPubMed

[12]

Li G Q, Schreiber M and Kleinekathöfer U 2007 Europhys. Lett. 79 27006 
IOPscience

[13]

Kohler S and Hänggi P 2007 Nat. Nanotechnol. 2 675 
CrossRefPubMed

[14]

Fainberg B D, Jouravlev M and Nitzan A 2007 Phys. Rev. B 76 245329 
CrossRef

[15]

Di Ventra M, Pantelides S T and Lang N D 2000 Phys. Rev. Lett. 84 979 
CrossRefPubMed

[16]

Damle P, Ghosh A W and Datta S 2002 Chem. Phys. 281 171 
CrossRef

[17]

Heurich J, Cuevas J C, Wenzel W and Schön G 2002 Phys. Rev. Lett. 88 256803 
CrossRefPubMed

[18]

Evers F, Weigend F and Koentopp M 2004 Phys. Rev. B 69 235411 
CrossRef

[19]

Mujica V, Kemp M and Ratner M A 1994 J. Chem. Phys. 101 6849 
CrossRef

[20]

Segal D, Nitzan A, Davis W B, Wasielewski M R and Ratner M A 2000 J. Phys. Chem. B 104 3817 
CrossRef

[21]

Foa Torres L E F, Pastawski H M and Makler S S 2001 Phys. Rev. B 64 193304 
CrossRef

[22]

Čížek M, Thoss M and Domcke W 2004 Phys. Rev. B 70 125406 
CrossRef

[23]

del Valle M, Tejedor C and Cuniberti G 2005 Phys. Rev. B 71 125306 
CrossRef

[24]

del Valle M, Gutirrez R, Tejedor C and Cuniberti G 2007 Nat. Nanotechnol. 2 176 
CrossRefPubMed

[25]

Kaiser F J, Hänggi P and Kohler S 2006 Eur. Phys. J. B 54 201 
CrossRef

[26]

Petrov E G and Hänggi P 2001 Phys. Rev. Lett. 86 2862 
CrossRefPubMed

[27]

Lehmann J, Ingold G L and Hänggi P 2002 Chem. Phys. 281 199 
CrossRef

[28]

Bruder C and Schoeller H 1994 Phys. Rev. Lett. 72 1076 
CrossRefPubMed

[29]

Gurvitz S A and Prager Ya S 1996 Phys. Rev. B 53 15932 
CrossRef

[30]

Lehmann J, Kohler S, May V and Hänggi P 2004 J. Chem. Phys. 121 2278 
CrossRefPubMed

[31]

Brandes T 2005 Phys. Rep. 408 315 
CrossRef

[32]

Kaiser F J, Strass M, Kohler S and Hänggi P 2006 Chem. Phys. 322 193 
CrossRef

[33]

Anderson P W 1958 Phys. Rev. 109 1492 
CrossRef