A Landauer Formula for Bioelectronic Applications
Eszter Papp, Dávid P. Jelenfi, Máté T. Veszeli, Gábor Vattay
AA L
ANDAUER F ORMULA FOR B IOELECTRONIC A PPLICATIONS
A P
REPRINT
Eszter Papp, Dávid P. Jelenfi, Máté T. Veszeli and Gábor Vattay
Department of Physics of Complex Systems,Eötvös Loránd UniversityH-1117 Budapest, Pázmány Péter sétány 1/A, Hungary [email protected]
September 4, 2019 A BSTRACT
Recent electronic transport experiments using metallic contacts attached to proteins identified some’stylized facts’ which contradict conventional wisdom that increasing either the spatial distancebetween the electrodes or the temperature suppresses conductance exponentially. These include nearlytemperature independent conductance over the protein in the 30-300K range, distance independentconductance within a single-protein in the 1-10 nm range and an anomalously large conductancein the 0.1-10 nS range. In this paper we develop a generalization of the low temperature Landauerformula which can account for the joint effects of tunneling and decoherence and can explain thesenew experimental findings. We use novel approximations which greatly simplify the mathematicaltreatment and allow us to calculate the conductance in terms of a handful macroscopic parametersinstead of the myriads of microscopic parameters describing the details of an atomic level quantumchemical computation. The new approach makes it possible to get predictions for the outcomes of newexperiments without relying solely on high performance computing and can distinguish importantand unimportant details of the protein structures from the point of view of transport properties. K eywords Landauer fromula · conductance of biomolecules · metallic contacts Electron transport measurements via metallic contacts attached to proteins show anomalous properties relative toelectron transfer in homologous structures[1, 2]. Borrowing the concept of stylized facts from economics[3], we canintroduce here three simplified presentations of empirical findings: • Conductance measured between metallic electrodes attached well to large protein structure is unexpectedlyhigh. It falls into the nano Siemens scale even over distances of several nano meters[4, 5, 6]. • The conductance does not show significant decay by increasing the distance of the electrodes[7, 8, 9]. • The conductance remains nearly constant when temperature is changed from tens of Kelvins to ambienttemperatures[10].Bioelectronic measurements with metallic contacts chemically bound to molecules can be regarded as molecularjunctions and the Landauer-Büttiker (LB) formula is one of the best theoretical tools to describe quantum conductance atzero temperature in such systems[11]. It expresses the conductance in terms of the scattering matrix elements betweenmetallic leads. In the simplest case, only a single scattering channel is open in a narrow lead and a single transmission T ( E F ) at the Fermi energy E F determines the conductance G = 2 e h T ( E F ) , (1) a r X i v : . [ c ond - m a t . d i s - nn ] S e p PREPRINT - S
EPTEMBER
4, 2019where the unit of the quantum conductance e /h ≈ nS . At high temperatures this formula is not applicable andelectron transfer is usually treated in semiclassical Marcus theory[12] (MT) k et = 2 π (cid:126) | H AB | √ πλk B T exp (cid:18) − ( λ + ∆ G ◦ ) λk B T (cid:19) , (2)where the electron transfer rate k et is expressed in terms of the electronic coupling between the initial and final states | H AB | , the reorganization energy λ , and the total Gibbs free energy change ∆ G ◦ . Then Nitzan et al. showed[13, 14, 15]that the conductance of molecular junctions is proportional with the electron transfer rate within the same approximation.Since electronic states in biomolecules are highly localized, overlap between distant electronic states decay fast andboth LB and MT yields exponentially decaying conductance G ∼ exp( − βl ) , where l is the distance of the electrodesand /β is about Å. Temperature dependence is also exponential due to the Arrhenius factor in (2). Both LB and MTare limiting cases only and the electron-vibrational (electron-phonon) interactions should be treated more carefully inthe intermediate regime. Recently, in Ref.[16] this derivation has been carried out for a molecular junction modelledas a single electronic level coupled with a collection of normalized vibrational modes. Using a generalized quantummaster equation, it has been shown that LB and MT can be viewed as two limiting cases of this more general expression.The current impasse in interpreting experimental results is coming from the fact that charge transport through molecularjunctions is described either as a purely coherent or a purely classical phenomenon. In recent years it became clear thatdecoherence plays an important role in biological energy transfer processes [17, 18] and these effects are not coveredby the semiclassical approximation[19]. In this paper we show that decoherence due to strong coupling to vibrationalmodes plays an important role in electron transport processes as well. We generalize the LB formula for conditionsrelevant in bioelectronic systems operating at strong decoherence. We capture new physics, which is absent in bothlimiting cases but plays a significant role when a metallic electrode is attached to the molecule and the chemical bondingis strong between the metal and the nearest localized electronic state of the molecule while direct tunneling between thetwo distant localized electronic states is exponentially suppressed.Our starting point is the derivation of low temperature LB formula for molecules by Datta et al. in Refs.[20, 21] whichwe summarize here briefly. The molecule is coupled to a left and a right electrode. The discrete levels of the molecule ε n are non-resonantly coupled to to the left and right electrodes with coupling strengths Γ Ln and Γ Rn respectively. Thepresence of contacts broadens the levels and can be described with a Lorentzian density of states d n ( E ) = 12 π Γ n ( E − ε n ) + Γ n / , (3)where Γ n = Γ Ln + Γ Rn is the broadening due to the contacts. If the level ε n were in equilibrium with the left contactthen the number of electrons N Ln , occupying the level would be given by N Ln = 2 (cid:90) + ∞−∞ d n ( E ) f ( E, µ L ) dE, (4)where µ L is the chemical potential in the left lead, f e ( E, µ ) = (1 + e ( E − µ ) /kT ) − is the Fermi distribution and thefactor 2 stands for spin degeneracy. A similar expression is valid for N Rn when the molecule is in equilibrium with theright lead. Under non-equilibrium conditions the number of electrons N n will be somewhere in between N Ln and N Rn and we can write the net current at the left junction as I Ln = e Γ Ln (cid:126) ( N Ln − N n ) , (5)where Γ Ln / (cid:126) is the escape rate from the level to the left lead. Similarly, for the right junction I Rn = e Γ Rn (cid:126) ( N n − N Rn ) . (6)Steady state requires I Ln = I Rn yielding N n = (Γ Ln N Ln + Γ Rn N Rn ) / (Γ Ln + Γ Rn ) . The current through the level is then I n = I Ln = I Rn = e (cid:126) Γ Ln Γ Rn Γ Ln + Γ Rn ( N Ln − N Rn ) . (7)Using the density of states the current trough the molecule can be written as I = (cid:88) n I n = (cid:88) n eh (cid:90) + ∞−∞ Γ Ln Γ Rn ( E − ε n ) + (Γ Ln + Γ Rn ) / f ( E, µ L ) − f ( E, µ R )) dE. (8)2 PREPRINT - S
EPTEMBER
4, 2019In the linear regime and at low temperatures kT → the chemical potential in the left and right electrodes is µ L/R = E F ± eU/ , the difference f ( E, µ L ) − f ( E, µ R ) ≈ δ ( E − E F ) eU and then I = 2 e h T ( E F ) · U, (9)where the transmission is given by the Breit-Wigner formula[11] T ( E F ) = (cid:88) n Γ Ln Γ Rn ( E F − ε n ) + (Γ Ln + Γ Rn ) / . (10)The formula is valid when the Fermi energy E F is close to an eigenenergy of the isolated molecule and the level spacingof the isolated molecule is larger than Γ Ln + Γ Rn . When energy of the isolated molecule ε n is above the Fermi energy E F the expression Γ Ln Γ Rn ( E F − ε n ) + (Γ Ln + Γ Rn ) / . (11)describes electron transmission. An electron tunneling trough an unoccupied orbital of the molecule from the left leadto the right one when electric field is switched on in that direction. When ε n is below the Fermi energy E F the orbitalis occupied and E F − ε n is positive. This case describes hole transport. A positively charged hole tunneling troughthe molecule from the right to the left electrode with negative tunneling energy ε n − E F . This way, both processescontribute to the net current with the same sign. In the next section we generalize Datta’s result for finite temperatures. The derivation of the LB formula in the Breit-Wigner approximation is an especially suitable starting point forgeneralization to include vibronic effects. When such effects are present, the electron (or hole) which tunnels intoan orbital of the molecule is able to transit to another orbital of the molecule since the energy difference betweenthe orbitals can be taken away (or supplied) by the interaction with the vibrational modes. Note that even at zerotemperature the electron can hop to lower energy so vibronic effects modify the LB formula even in that case. Thesteady state condition I Ln = I Rn should be modified to account for hoping in and out of an orbital. Electrons can hopbetween the (nearly unoccupied) electronic states above the Fermi energy while holes can hop between the (nearlyoccupied) states below it. Accordingly, we should treat electrons and holes separately. For brevity, we derive the resultsfor electrons in detail and then give the analogous expressions for holes.Quantum master equations are the most convenient way to describe the transition between electronic states. Theyare in general non-Markovian, but for practical purposes can be approximated with Markovian equations such as theRedfield equation[22]. The reduced density matrix elements of the molecule in the energy basis (cid:37) nm then satisfy alinear equation ∂ t (cid:37) nm = i (cid:126) ( ε m − ε n + i Γ n / i Γ m / (cid:37) nm + (cid:88) kl R nmkl (cid:37) kl + J n δ nm , (12)where J n is the external current, ε n + i Γ n / is the broadened level and R nmkl is the Bloch-Redfield tensor describingtransitions due to the couplings to the phononic vibrations. The tetradic matrix R mnkl is the transfer rate from (cid:37) kl to (cid:37) mn and can be expressed as R mnkl = Γ lmnk + Γ ∗ knml − δ ml (cid:88) p Γ nppk − δ nk (cid:88) q Γ ∗ mqql , (13)where Γ mnkl = 1 (cid:126) (cid:90) ∞ dτ e − i ( ε k − ε l ) τ/ (cid:126) (cid:104) V nm ( τ ) V kl (0) (cid:105) b , are Fourier-Laplace transforms of correlation functions of matrix elements V ij of the system-bath coupling operatorbetween system eigenstates i and j and the brackets represent a trace over the thermalized bath. We note that in thislevel of description the electron and hole states do not mix. Electrons can hop on states above, while holes below theFermi energy. Consequently, there are two separate Redfield equations, one for the electrons and one for the holes. Thefour indices of R nmkl and the two indices of (cid:37) nm should be either all electron or hole states.We normalize this equation such a way that the diagonal elements of the density matrix can correspond to the occupations (cid:37) nn = N n introduced in the previous section. In this case the external (material) current becomes J n = Γ Ln (cid:126) (cid:37) Ln + Γ Rn (cid:126) (cid:37) Rn , (14)3 PREPRINT - S
EPTEMBER
4, 2019where (cid:37)
L/Rn = 2 (cid:90) + ∞−∞ d n ( E ) f ( E, µ
L/R ) dE, (15)is the occupation of the levels when the molecule is in equilibrium with the left/right lead.In absence of an electric field ( U = 0 ) the system is in equilibrium and µ L/R = E F . The Fermi energy is betweenthe HOMO and the LUMO energies ε N/ = ε HOMO < E F < ε N/ = ε LUMO . The number of electrons N in themolecule is given by the sum of occupancies N = (cid:80) ∞ n =1 f e ( ε n , E F ) . This can be written also as N/ (cid:88) n =1 (1 − f ( ε n , E F )) = ∞ (cid:88) n = N/ f e ( ε n , E F ) , (16)and with the Fermi distribution for holes f h ( E, µ ) = 1 − f e ( E, µ ) = (1 + e ( µ − E ) /kT ) − it can be simplified to N/ (cid:88) n =1 f h ( ε n , E F ) = ∞ (cid:88) n = N/ f e ( ε n , E F ) , (17)meaning that the number of holes below the Fermi energy is the same as the electrons above it and the molecule ischarge neutral. An important aspect of bioelectronic systems is that they have highly localized electronic states anda large HOMO-LUMO gap ( ∼ eV ) in accordance with standard Density Functional Theory (DFT) calculations[23].Based on this we can assume that e ( E F − ε n ) /kT (cid:28) for electronic states and also e ( ε n − E F ) /kT (cid:28) for hole statesand we can replace the Fermi distribution with the Boltzmann distribution in both cases f e ( E, µ ) = e − ( E − µ ) /kT and f h ( E, µ ) = e − ( µ − E ) /kT . Accordingly, the equilibrium occupancy for electrons is given by the Boltzmann distribution (cid:37) L/Rn = (cid:26) e − ( ε n − E F ) /kT for electrons , e − ( E F − ε n ) /kT for holes . (18)Using (17) in the Boltzmann approximation we can introduce the partition function Z ( T ) = N/ (cid:88) n =1 e − ( E F − ε n ) /kT = ∞ (cid:88) n = N/ e − ( ε n − E F ) /kT , which is the same for electrons and holes.In absence of electric field the steady state solution of the Redfield equation (12) is also the Boltzmann distribution (cid:37) nn = (cid:26) e − ( ε n − E F ) /kT for electrons , e − ( E F − ε n ) /kT for holes . (19)Then in equilibrium each term vanishes separately in I = e (cid:88) n Γ Ln (cid:126) ( (cid:37) Ln − (cid:37) nn ) = 0 , (20)and there is no electric current.When the electric field is switched on ( U (cid:54) = 0 ) the system is out of equilibrium. In the linear regime we can expand thedeviation from the equilibrium (cid:37) L/Rn = 2 (cid:90) + ∞−∞ dEd n ( E )[ f ( E, E F ) ∓ f (cid:48) ( E, E F ) eU/ ... ] ≈ (cid:37) L/Rn ± D n ( E F , T ) eU, (21)where D n ( E F , T ) = − (cid:90) + ∞−∞ f (cid:48) ( E, E F ) d n ( E ) dE. (22)We can introduce the deviation of the density matrix elements from their equilibrium value (cid:37) (cid:48) nm = (cid:37) nm − (cid:37) nm andusing (12) and (14) we can write the steady state equation − Γ Ln − Γ Rn (cid:126) D n ( E F , T ) eU δ nm = i (cid:126) ( ε m − ε n + i Γ n / i Γ m / (cid:37) (cid:48) nm + (cid:88) kl R nmkl (cid:37) (cid:48) kl , (23)4 PREPRINT - S
EPTEMBER
4, 2019where we grouped the external current to the left hand side. The current through the molecule is I = e (cid:88) n Γ Ln (cid:126) ( D n ( E F , T ) eU − (cid:37) (cid:48) nn ) . (24)To get the general solution of (23) we can introduce the tetradic matrix L nmkl = ( i/ (cid:126) )( ε m − ε n + i Γ n / i Γ m / δ nk δ ml + R nmkl and can write − Γ Ln − Γ Rn (cid:126) D n ( E F , T ) eU δ nm = (cid:88) kl L nmkl (cid:37) (cid:48) kl . (25)The solution of this equation can be given in terms of the inverse matrix (cid:37) (cid:48) nm = − (cid:88) k L − nmkk Γ Lk − Γ Rk (cid:126) D k ( E F , T ) eU, (26)where the inverse satisfies the relation (cid:80) pq L nmpq L − pqkl = δ nk δ ml . Substituting this solution into (24) we get thegeneralized Landauer-Büttiker formula I = e U (cid:126) (cid:88) n D n ( E F , T ) (cid:34) Γ Ln + 1 (cid:126) (cid:88) k Γ Lk L − kknn (Γ Ln − Γ Rn ) (cid:35) . (27)This formula can be brought (see Appendix A) to a form which reflects the left-right symmetry G = e (cid:126) (cid:88) n D n ( E F , T ) (cid:34) Γ Ln + Γ Rn + 1 (cid:126) (cid:88) k (Γ Lk − Γ Rk ) L − kknn (Γ Ln − Γ Rn ) (cid:35) , (28)which is our main result. We note that L − kknn = 0 unless both k and n are electron or hole states, consequently theconductance can be split to an electron and a hole part G = G e + G h , where G e = e (cid:126) ∞ (cid:88) n = N/ D n ( E F , T ) Γ Ln + Γ Rn + 1 (cid:126) ∞ (cid:88) k = N/ (Γ Lk − Γ Rk ) L − kknn (Γ Ln − Γ Rn ) , (29)and G h = e (cid:126) N/ (cid:88) n =1 D n ( E F , T ) Γ Ln + Γ Rn + 1 (cid:126) N/ (cid:88) k =1 (Γ Lk − Γ Rk ) L − kknn (Γ Ln − Γ Rn ) , (30)just like in the zero temperature LB formula discussed before. The present formalism allows us to calculate the electron transfer along the same lines. For specificity, the left electrodeplays the role of donor and the right electrode the acceptor site. The electron charge on the donor and acceptor sitesfollows the Fermi distribution, which can be approximated by the Boltzmann distribution due to the large HOMO-LUMO gap. The electron-hole picture is useful here as well. Electrons traversing the molecule via almost unoccupiedorbitals above the Fermi energy contribute to the electron part, while transfer via almost fully occupied orbitals belowthe Fermi level can be regarded as hole transport. Introducing (cid:37) D for the total density on the left electrode and (cid:37) A forthe right electrode, the left and right densities become (cid:37) L/Rn = (cid:37) D/A p Bn , where p Bn = (cid:26) e − ( ε n − E F ) /kT /Z for electrons , e − ( E F − ε n ) /kT /Z for holes . (31)Then, given the external current J n = Γ Ln (cid:126) (cid:37) Ln + Γ Rn (cid:126) (cid:37) Rn , (32)we have to solve the Redfield equation − J n δ nm = (cid:88) kl L nmkl (cid:37) kl , (33)5 PREPRINT - S
EPTEMBER
4, 2019to get the total material current J = (cid:88) n Γ Ln (cid:126) ( (cid:37) Ln − (cid:37) nn ) , (34)which leads to J = 1 (cid:126) (cid:37) D (cid:88) n p Bn Γ Ln (cid:34) (cid:126) (cid:88) k Γ Lk L − kknn (cid:35) + 1 (cid:126) (cid:37) A (cid:88) n p Bn Γ Ln (cid:34) (cid:126) (cid:88) k Γ Rk L − kknn (cid:35) . (35)Electrons can escape from the acceptor site with escape rate κ = J/(cid:37) A and the electron transfer rate is the ratio of thematerial current and the density at the donor site k ET = J/(cid:37) D . We can express the transfer rate with the escape rate k ET = (cid:80) n p Bn Γ Ln (cid:2) (cid:126) (cid:80) k Γ Lk L − kknn (cid:3) / (cid:126) /κ ) (cid:80) n p Bn Γ Ln (cid:2) (cid:126) (cid:80) k Γ Lk L − kknn (cid:3) / (cid:126) , (36)where we used the result of Appendix A to show that (cid:88) n p Bn Γ Ln (cid:34) (cid:126) (cid:88) k Γ Rk L − kknn (cid:35) = − (cid:88) n p Bn Γ Ln (cid:34) (cid:126) (cid:88) k Γ Lk L − kknn (cid:35) . (37) The strength of the contacts relative to the thermal energy plays a crucial role in the conductance properties of thesesystems. When the contacts are weak the electrons and holes can enter the molecule from the lead via thermal excitation.In this case the conductance is intimately related to electron transfer as it has been shown in Refs.[14, 24, 15]. Here, wederive an exact formula between electron transfer and electron conductance.When the contacts are weak compared to the thermal energy Γ n (cid:28) kT the density of states consists of delta peaks d n ( E ) ≈ δ ( E − ε n ) , the Fermi distribution is well approximated with the Boltzmann and we get D n ( E F , T ) ≈ − (cid:90) + ∞−∞ f (cid:48) ( E, E F ) δ ( E − ε n ) dE = − f (cid:48) ( ε n , E F ) ≈ (cid:26) e − ( ε n − E F ) /kT /kT for electrons ,e − ( E F − ε n ) /kT /kT for holes . (38)This can be written in the more compact form using the normalized Boltzmann distribution D n ( E F , T ) =( Z ( T ) / kT ) p Bn . Substituting this into (27) and using the sum rule derived in Appendix D we can eliminate Γ n and get the form G = e Z ( T ) kT (cid:126) (cid:88) n p Bn Γ Ln (cid:34) (cid:126) (cid:88) k Γ Lk L − kknn (cid:35) . (39)The sum in this expression appears in (36) as well so that the electron transfer rate can be expressed with the conductancedirectly k ET = ( kT /e Z ( T )) G kT /e Z ( T )) G/κ . (40)When the escape from the acceptor is strong, we can neglect /κ and the conductance is proportional with the electrontransfer rate G = e kT Z ( T ) · k ET . (41)In biomolecules where the HOMO-LUMO gap ∆ HL = E LUMO − E HOMO is large compared to the thermal energy kT the partition sum is dominated by the gap Z ( T ) ≈ e − ( E LUMO − E F ) /kT ≈ e − ∆ HL / kT and the conductance is G = e kT e − ∆ HL / kT k ET , (42)where e − ∆ HL / kT ∼ − making these systems practically insulators. In case of bridged molecular systemsconsidered by Nitzan in Ref.[24] the partition function is Z ≈ e − ∆ E/kT , where ∆ E = E B − E F is the differencebetween the Fermi energy of the metallic leads and the average energy of the bridged system and we recover Nitzan’sformula G = e kT e − ∆ E/kT k ET . (43)6 PREPRINT - S
EPTEMBER
4, 2019Finally, in case of a coherent one dimensional bridged molecule system where the gap is zero ∆ E = 0 considered inRef.[13] we can use the classical partition function of a free particle in a one dimensional box: Z = (cid:90) L M (cid:90) + ∞−∞ dxdph e − p / mkT = L M √ πmkTh , (44)where L M is the length of the molecule and m is the effective mass of the electron in the molecule. In this case we getthe formula G = 2 e h L M (cid:114) πm kT k ET , (45)which is slightly different from the heuristically derived result in Ref.[13] G = e h L M (cid:113) m E F k ET . The differencecomes from the fact that in Ref.[13] it is assumed that the electron traverses the molecule with effective velocity v F = (cid:112) E F /m and energy E F , while in reality the electron transfer happens at the thermal energy scale kT and witheffective velocity v T = (cid:112) kT /m . The two expressions differ in a factor of v F /v T ≈ only, which is hard to verifyexperimentally. In the opposite case, when the contacts are strong, electrons and holes enter the molecule via tunneling. This is anew regime not covered by the previous studies and we show that the relation between electron transfer rate and theconductance breaks down.When the thermal energy is small compared to the strengths of the contacts kT (cid:28) Γ n the density of states is smoothand the thermal distribution is approximately − f (cid:48) ( E, E F ) ≈ δ ( E − E F ) so that D n ( E F , T ) = (cid:90) + ∞−∞ d n ( E ) δ ( E − E F ) dE = 12 π Γ n ( E F − ε n ) + Γ n / . (46)Substituting this into (29) and (30) we get the electron G e = e h ∞ (cid:88) n = N/ Γ n ( E F − ε n ) + Γ n / Γ Ln + Γ Rn + 1 (cid:126) ∞ (cid:88) k = N/ (Γ Lk − Γ Rk ) L − kknn (Γ Ln − Γ Rn ) , (47)and the hole conductance G h = e h N/ (cid:88) n =1 Γ n ( E F − ε n ) + Γ n / Γ Ln + Γ Rn + 1 (cid:126) N/ (cid:88) k =1 (Γ Lk − Γ Rk ) L − kknn (Γ Ln − Γ Rn ) . (48)In this case, tunneling populates the levels and the relation with the Boltzmann distribution breaks down. It is no longerpossible to connect electron transfer and conductance with a simple formula.The other equally crucial factor in the conductance of these systems is the strength of the coupling to the environmentthrough the vibrational degrees of freedom, which is encoded in the matrix L − kknn . In general it is complicated tocalculate this quantity since the microscopic details of the couplings between the electron and vibration degrees offreedom can play an important role. To get an insight we consider the two limiting cases, when the coupling tothe vibrations is negligible (coherent case) and the opposite case, when coupling to the environment dominates (fulldecoherence). Surprisingly, in the latter case the details of the coupling drop out and the conductance depends on thecontact strengths and the energy spectrum of the molecule as we show next.In the coherent case, when the Redfield tensor elements describing the coupling to the heath bath are small ( | R | (cid:28) Γ )we can neglect them and the inverse operator matrix elements become L − nnmm = − (cid:126) δ nm Γ Ln + Γ Rn . (49)Substituting this into (29) and (30) we get the electron G e = 2 e h ∞ (cid:88) n = N/ D n ( E F , T ) Γ Ln Γ Rn Γ Ln + Γ Rn , (50)7 PREPRINT - S
EPTEMBER
4, 2019and the hole conductance G h = 2 e h N/ (cid:88) n =1 D n ( E F , T ) Γ Ln Γ Rn Γ Ln + Γ Rn . (51)In case of strong contacts and ( | R | , kT (cid:28) Γ ) we can substitute this and recover the Landauer-Büttiker formula in theBreit-Wigner approximation G = G e + G h = 2 e h (cid:88) n Γ Ln Γ Rn ( E F − ε n ) + (Γ Ln + Γ Rn ) / . (52) It is an important experimental situation[5, 6], when one of the contacts is strong and the other one is weak. In somecases, the strongly coupled electrode forms a covalent bond with a specific atom of the molecule, while the otherelectrode is coupled non-specifically with weak coupling. Interestingly, to meet the condition of the strong contact case kT (cid:28) Γ n = Γ Ln + Γ Rn it is sufficient if only one of the contacts is strong. For example, if the left contact is strong kT (cid:28) Γ Ln and the right contact is weak we can neglect Γ Rn in (47) and (48) and arrive at the expression for the electron G e = e h ∞ (cid:88) n = N/ Γ Ln ( E F − ε n ) + (Γ Ln ) / Γ Ln + 1 (cid:126) ∞ (cid:88) k = N/ Γ Lk L − kknn Γ Ln , (53)and for the hole conductance G h = e h N/ (cid:88) n =1 Γ Ln ( E F − ε n ) + (Γ Ln ) / Γ Ln + 1 (cid:126) N/ (cid:88) k =1 Γ Lk L − kknn Γ Ln . (54)This means that high conductance can arise not just between two strong contacts but also in the single strong contactcase. We discuss this possibility further in the next sections. When the Bloch-Redfield terms are small compared to the couplings to the leads we can expect just some moderatedeviations from the LB formula. Interesting new physics arises in the opposite case ( Γ (cid:28) | R | ), when the electrons orholes arriving from the leads are strongly mixed in the molecule. In absence of the coupling to the leads the matrix L nmkl = ( i/ (cid:126) )( ε m − ε n ) δ nk δ ml + R nmkl , (55)describes the isolated molecule. The steady state solution of the density matrix of this system is the Boltzmanndistribution. The normalized equilibrium density matrix is (cid:37) nn = e − ( ε n − E F ) /kT /Z ( T ) for electrons and (cid:37) nn = e − ( E F − ε n ) /kT /Z ( T ) for holes. The inverse operator can be calculated perturbatively in the limit of small coupling Γ n → . For the details of the calculation see Appendix B and the diagonal elements of (56) become L − nnmm ≈ − (cid:126) (cid:37) nn (cid:80) p Γ p (cid:37) pp , (56)where all the indices should be either electrons or holes. Substituting this into (28) yields the hole conductance G h = e h N/ (cid:88) n =1 Γ n ( E F − ε n ) + Γ n / Rn (cid:104) Γ Lh ( T ) (cid:105) + Γ Ln (cid:104) Γ Rh ( T ) (cid:105)(cid:104) Γ Lh ( T ) (cid:105) + (cid:104) Γ Rh ( T ) (cid:105) , (57)and the electron conductance G e = e h ∞ (cid:88) N/ Γ n ( E F − ε n ) + Γ n / Rn (cid:104) Γ Le ( T ) (cid:105) + Γ Ln (cid:104) Γ Re ( T ) (cid:105)(cid:104) Γ Le ( T ) (cid:105) + (cid:104) Γ Re ( T ) (cid:105) , (58)such that the total conductance is the sum of the two G = G h + G e . (59)8 PREPRINT - S
EPTEMBER
4, 2019We introduced the weighted sums of the left and the right coupling strengths for electron (cid:104) Γ L/Re ( T ) (cid:105) = ∞ (cid:88) n = N/ Γ L/Rn e − ( ε n − E F ) /kT , (60)and for the hole states (cid:104) Γ L/Rh ( T ) (cid:105) = N/ (cid:88) n =1 Γ L/Rn e − ( E F − ε n ) /kT . (61)We can carry out the summation and get the conductance in the vibration dominated regime G = e h (cid:20) T Rh ( E F ) (cid:104) Γ Lh ( T ) (cid:105) + T Lh ( E F ) (cid:104) Γ Rh ( T ) (cid:105)(cid:104) Γ Lh ( T ) (cid:105) + (cid:104) Γ Rh ( T ) (cid:105) + T Re ( E F ) (cid:104) Γ Le ( T ) (cid:105) + T Le ( E F ) (cid:104) Γ Re ( T ) (cid:105)(cid:104) Γ Le ( T ) (cid:105) + (cid:104) Γ Re ( T ) (cid:105) (cid:21) , (62)where we introduced the sums for holes T L/Rh ( E F ) = N/ (cid:88) n =1 Γ n Γ L/Rn ( E F − ε n ) + Γ n / , (63)and for electrons T L/Re ( E F ) = ∞ (cid:88) n = N/ Γ n Γ L/Rn ( E F − ε n ) + Γ n / . (64)A remarkable property of this new conductance formula is that it is independent of the details of the vibrational processand relies solely on the equilibrium distribution and the couplings to the leads. The ratio P Lh ( T ) = (cid:104) Γ Lh ( T ) (cid:105)(cid:104) Γ Lh ( T ) (cid:105) + (cid:104) Γ Rh ( T ) (cid:105) , (65)in the expression of the conductance can be interpreted as the probability that a hole entering anywhere into the moleculeleaves it towards the left lead. The sum T Rh ( E F ) is the probability that a hole tunnels into the molecule from the rightlead. The product T Rh ( E F ) (cid:104) Γ Lh ( β ) (cid:105) is the probability that a hole entering from the right lead leaves the molecule troughthe left lead. The four terms in the formula G = e h (cid:2) T Rh ( E F ) P Lh ( T ) + T Lh ( E F ) P Rh ( T ) + T Re ( E F ) P Le ( T ) + T Le ( E F ) P Re ( T ) (cid:3) , (66)represent the four scenarios in which electrons and holes can generate current. It has a very modest dependence ontemperature and on the distance between the contacts as we show in the next sections.Finally, here we can also discuss the sub-case when the left electrode forms a strong specific bond with the moleculewhile the other electrode is coupled weakly or non-specifically. We can then neglect T Re ( E F ) and T Rh ( E F ) in (66) andget the simplified expression G = e h (cid:2) T Lh ( E F ) P Rh ( T ) + T Le ( E F ) P Re ( T ) (cid:3) . (67)This means that electrons and holes can tunnel into the molecule via the strong left contact and some of them can leavetrough the weak contact with the equilibrium probabilities P Re ( T ) and P Rh ( T ) . The temperature dependence of the conductance is coming from the probabilities. In Appendix C we show that theyare temperature independent if kT (cid:28) ε HOMO − ε HOMO − and kT (cid:28) ε LUMO +1 − ε LUMO , which is usually holdsin proteins, where the level spacings are typically in the order of . − . eV and the experimental temperatures arein the kT = 0 . − . eV range. In certain cases, due to the fluctuation of level spacing it can happen that ε HOMO − ε HOMO − or ε LUMO +1 − ε LUMO is somewhat lower accidentally, therefore we keep the temperaturedependent terms in leading order to account for these effects. Using the probabilities derived in Appendix C we get thefollowing expression for the temperature independent part of the conductance G = e h (cid:20) T Rh ( E F )Γ LHOMO + T Lh ( E F )Γ RHOMO Γ LHOMO + Γ
RHOMO + T Re ( E F )Γ LLUMO + T Le ( E F )Γ RLUMO Γ LLUMO + Γ
RLUMO (cid:21) , (68)9 PREPRINT - S
EPTEMBER
4, 2019and for the temperature dependence in leading order G T = e h ( T Rh ( E F ) − T Lh ( E F )) Γ LHOMO Γ RHOMO (Γ LHOMO + Γ
RHOMO ) (cid:34) Γ LHOMO − Γ LHOMO − Γ RHOMO − Γ RHOMO (cid:35) e − ( ε HOMO − ε HOMO − ) /kT + e h ( T Re ( E F ) − T Le ( E F )) Γ LLUMO Γ RLUMO (Γ LLUMO + Γ
RLUMO ) (cid:34) Γ LLUMO +1 Γ LLUMO − Γ RLUMO +1 Γ RLUMO (cid:35) e − ( ε LUMO +1 − ε LUMO ) /kT , and G = G + G T . The sign of the temperature dependent part is determined by the combined effect of the sign of Γ LHOMO − / Γ LHOMO − Γ RHOMO − / Γ RHOMO , Γ LLUMO +1 / Γ LLUMO − Γ RLUMO +1 / Γ RLUMO , T Rh ( E F ) − T Lh ( E F ) and T Re ( E F ) − T Le ( E F ) , which depends not just on whether the left or the right side is coupled stronger, but also from thedetails of the couplings to the HOMO vs. HOMO-1 and LUMO vs. LUMO-1. Looking at the formula (66) we can realize that neither T nor P has a systematic dependence on the distance of theelectrodes. The orbitals of large molecules are localized. Assuming that the electrodes are far from each other, they arecoupled the most strongly to some localized orbital of the molecule near the electrode. These orbitals do not overlapand direct tunneling is negligible. The electrons are transported due to the strong vibrational effect. Due to the strongmixing inside of the molecule, the electron moves ergodically inside and looses information about its point of arrival.The exit direction (left or right) is determined solely by the escape rates from the molecule. The probability that wefind the electron on an orbital is given by the Boltzmann distribution e − ( E F − ε n ) /kT /Z ( T ) and the rate of exit fromthis state to the left electrode is the escape rate multiplied with the probability (Γ Ln / (cid:126) ) e − ( E F − ε n ) /kT /Z ( T ) . The totalrate of escape to the left electrode is (cid:80) n (Γ Ln / (cid:126) ) e − ( E F − ε n ) /kT /Z ( T ) and the probability that the electron leaves themolecule trough the left exit is the ratio of the total rate of exit to the left divided by the total rate of exit (left+right) (cid:80) n (Γ n / (cid:126) ) e − ( E F − ε n ) /kT /Z ( T ) . The result is independent of the distance of the electrodes. Looking at one of thesums when the electrodes are far away from each other T Lh ( E F ) = N/ (cid:88) n =1 Γ n Γ Ln ( E F − ε n ) + Γ n / , (69)we can realize that the strongest contributions come from large Γ Ln . However, in this case the corresponding Γ Rn is verysmall, since the left and the right electrodes can’t couple strongly to the same localized state at the same time. So, whenwe calculate this sum, we can drop the terms related to the right electrode and get T Lh ( E F ) ≈ N/ (cid:88) n =1 (Γ Ln ) ( E F − ε n ) + (Γ Ln ) / , (70)which depends only on the couplings of the left electrode and is independent of the relative position of the two electrodes.Similar arguments are true for the right electrode, and it is distance independent as well. In summary, none of the termsin formula (66) show any systematic distance dependence. We can make a rough estimate of the conductance in a typical arrangement. The probabilities P are in order ofunity and the typical value of T determines the order of magnitude. The magnitude of | E F − ε n | is bigger than halfof the HOMO-LUMO gap and can be − eV typically. The couplings Γ n are in the . eV range. The ratios (Γ n / | E F − ε n | ) are then typically − − − and e /h ≈ , nS . The resulting conductances then should betypically in the . − . nS range. It is instructive to calculate the temperature dependence of the electron transfer rate in the strong decoherence case.Using the same approximation yields k ET = Z ( T ) (cid:126) (cid:20) (cid:104) Γ Le ( T ) (cid:105)(cid:104) Γ Re ( T ) (cid:105)(cid:104) Γ Le ( T ) (cid:105) + (cid:104) Γ Re ( T ) (cid:105) + (cid:104) Γ Lh ( T ) (cid:105)(cid:104) Γ Rh ( T ) (cid:105)(cid:104) Γ Lh ( T ) (cid:105) + (cid:104) Γ Rh ( T ) (cid:105) (cid:21) , (71)10 PREPRINT - S
EPTEMBER
4, 2019where we took the large escape rate κ → ∞ limit in (36). Since the HOMO-LUMO gap is much larger than the thermalenergy, we can again keep only the leading terms and get k ET ≈ e − ∆ HL /kT (cid:126) (cid:20) Γ LLUMO Γ RLUMO Γ LLUMO + Γ
RLUMO + Γ
LHOMO Γ RHOMO Γ LHOMO + Γ
RHOMO (cid:21) , (72)where ∆ HL is the HOMO-LUMO gap. This expression shows that the electron transfer rate has an Arrehnius typetemperature dependence. It contains the products Γ LLUMO Γ RLUMO and Γ LHOMO Γ RHOMO which describe tunnelingtrough the entire molecule and decay exponentially with the size of the molecule. It is obvious that the electron transferrate is not proportional with the conductance, which - unlike the transfer rate - shows no exponential dependence ondistance or temperature.
In this section we show how the present theory explains some of the key experimental findings of protein conductancein the presence of strongly coupled electrodes. C u rr e n t d e n s i t y [ A / c m ] -1 ] Native m-Mbapo MbReconstituted m-Mb Figure 1: Current density curves of experiment Ref[25] reconstructed using (73). Start and end points of curves visuallyextracted from original figure. For activation energy the gap between HOMO and HOMO-1 energies have been used inaccordance with (68).In Ref.[25] it has been found that the current through the system increases with the coupling strength if metallicelectrodes are attached to various Myoglobin structures. The change of the current with the temperature decreaseswith increasing strength and the current becomes temperature independent for small temperature. In all cases theyfound only mild temperature dependence inconsistent with the large HOMO-LUMO gap of Myoglobin. In Fig 1 wereproduced Figure 3. of the original article with the formula I = I + I T exp( − ∆ E/kT ) , (73)where the parameter values are shown in Table 1. In Table 2 we show the calculated values of the energies of molecularorbitals of Myoglobin. The observed temperature dependence is not consistent with the large HOMO-LUMO gap andactivation energy ∆ E = ( ε LUMO − ε HOMO ) / . eV , therefore we can exclude all traditional explanations,which rely on the thermal excitation trough the gap. On the other hand (73) is fully consistent with our formula (68) and(69). According to our formula the conductance becomes temperature independent for low temperatures and the weaktemperature dependence is governed by the smallest of the HOMO and HOMO-1 energy difference or LOMO andLOMO+1 energy difference. These are ε LUMO +1 − ε LUMO = 0 . eV and ε HOMO − ε HOMO − = 0 . eV PREPRINT - S
EPTEMBER
4, 2019for Myoglobin. The HOMO and HOMO-1 difference dominates and ∆ E = ε HOMO − ε HOMO − reproduces theexperimental results correctly. We note that in Myoglobin hole transport dominates the temperature dependent part andthat hole transport in general is often disregarded in intuition driven theoretical studies of electron transfer.Table 1: Parameter values reproducing the Myoglobin measurement results of Ref.[25]. I in A/cm I T in A/cm Native m-Mb . · − . · − apo Mb . · − . · − Reconstituted m-Mb . · − . · − Table 2: Myoglobin energies near the HOMO-LUMO gap. The energies have been calculated with the semiempiricalextended Hückel method implemented in the YaEHMOP package http://yaehmop.sourceforge.net . Myoglobinstructure taken from RCSB PDB . energy in eV LUMO+1 -8.9052LUMO -9.5185HOMO -9.9637HOMO-1 -10.0282Our second example is the measurement of Cytochrome C in Ref[26]. We reproduced two measurement curves shown inFigure 1.c. of the original paper and we show them in Fig 2. We used again the fit (73) and the parameters are in Table 3. -16-15-14-13-12 5 10 15 20 25 l n ( J - . V ) -1 ]Electrostatic binding (WT)Covalent binding (E104C) Figure 2: Current density curves of experiment Ref[26] reconstructed using (73). Start and end points of curves visuallyextracted from original figure. In case of covalent bonding for activation energy the gap between LUMO and LUMO+1energies have been used in accordance with (68), while for the electrostatic bonding case the experimentally foundvalue has been used.For the electrostatic case we used the activation energy ∆ E = 0 . eV found experimentally[26]. In Table 4 we showthe calculated orbital energies of Cytochrome C. This experimentally found activation energy is in reasonable agreementwith the value calculated from the numerical HOMO-LUMO gap ∆ E = ( ε LUMO − ε HOMO ) / . eV . Accordingto our formula the conductance becomes temperature independent for low temperatures and the weak temperaturedependence is governed by the smallest of the HOMO and HOMO-1 energy difference or LOMO and LOMO+1 energydifference, which are ε LUMO +1 − ε LUMO = 0 . eV and ε HOMO − ε HOMO − = 0 . eV for Cytochrome C. In12 PREPRINT - S
EPTEMBER
4, 2019this case the LUMO and LUMO+1 difference dominates and ∆ E = ε LUMO +1 − ε LUMO reproduces the experimentalresults correctly. We note that in Cytochrome C electron transport dominates the temperature dependent part.Table 3: Parameter values reproducing the Cytochrome C measurement results of Ref.[26] . I in A/cm I T in A/cm ∆ E in eV Covalent binding (E104C) . · − . · − . Electrostatic binding (WT) . · − . · − . Table 4: Cytochrome C energies near the HOMO-LUMO gap. The energies have been calculated with the semiempiricalextended Hückel method implemented in the YaEHMOP package http://yaehmop.sourceforge.net . CytochromeC structure taken from RCSB PDB . energy in eV LUMO+1 -9.9457LUMO -9.9252HOMO -9.6021HOMO-1 -9.0685
In Ref.[27] it has been found that in certain peptide nuclear acid structures both electron transfer rates and conductancedecays exponentially with the length of the structure, but the two exponents differ considerably. While electron transferrates decay as e − βl , where β ET ≈ . Å − like in most biological structures, conductance decay is about two thirdsslower β G ≈ . β ET . The authors attributed this to a possible new power law scaling relation between the two G ∼ k . ET in place of the linear relation G ∼ k ET . As we have seen, there is no such a simple relation between thesetwo quantities in general, but it is possible that the two quantities decay with different exponents.We note, that in the experiment the contact connecting the PNA structure to the Au electrode is weak. This hasbeen concluded in Ref.[28] suggesting that charge transfer between PNA and the Au substrate may be difficult dueto the required amounts of energy to overcome the injection barriers. The metal electrode used in the conductancemeasurement has been strong indicated by the high values of measured conductance as well. Electron transfer rate inthis situation is well described by (72) and conductance by (67) which can be simplified further by keeping the leadingterm in case of a large HOMO-LUMO gap G = e h (cid:20) T Lh ( E F ) Γ RHOMO Γ LHOMO + T Le ( E F ) Γ RLUMO Γ LLUMO (cid:21) . (74)The products of couplings show exponential dependence Γ LLUMO Γ RLUMO ∼ Γ LHOMO Γ RHOMO ∼ e − β ET l , where l isthe length of the molecule since they describe tunneling across the molecule. The HOMO and LUMO orbitals aresomewhere midway in the molecule and they are located very close to each other. It is then reasonable to assumethat Γ LHOMO ∼ Γ LLUMO ∼ e − β ET xl and Γ RHOMO ∼ Γ RLUMO ∼ e − β ET (1 − x ) l , where xl and (1 − x ) l is the distanceof the left and right electrodes from the location of the HOMO-LUMO orbitals respectively. Inserting this into theconductance and assuming that T Lh ( E F ) and T Le ( E F ) don’t change with the distance as we discussed before we cansee that G ∼ e − β ET (1 − x ) l . (75)Assuming x ≈ . can explain the relation of exponents observed in the experiment. This suggests that the HOMO andLUMO are closer to the metal contact, which is due to the Ferrocene redox center attached to the end of the molecule incontact with the metal. This example clarifies that in case of a strong and a weak contact distance dependence can stillbe observed. Only two strong contacts guarantee distance independence. In Ref.[9] the distribution of conductance between metallic contacts attached to various proteins has been investigated.Here we reconstruct the experiment where conductance have been measured between thiolated Lysines of a Streptavidinmolecule. Streptavidin is the smallest protein used in Ref.[9] and it is computationally feasible to calculate its electronic13
PREPRINT - S
EPTEMBER
4, 2019structure using semiempirical methods. For the calculations the structure has been used. Biotine molecules have been removed from the structure. The energies and orbitals have beencalculated with the semiempirical extended Hückel method implemented in the YaEHMOP package http://yaehmop.sourceforge.net . Here we attempt to reproduce the experimental result shown in Figure 2A (T-T), in Figure 3A(lower part) and in Figure 4A of the original paper. Original data has been kindly provided by the Authors. We assumethat the thiolated sites have been randomly coupled to the substrate which and to the STM tip. We select pairs ofthiolated sites using geometric information so that they lie on opposite sites of the molecule. We took all possible pairswith equal weight. For the coupling strengths the standard formula
Γ = V | Φ Mk | has been used, where V is thestrength of the coupling and | Φ Mk | is the quantum mechanical probability to find the electron in orbital M on atom k . We assumed that the strengths are random and Gaussian distributed. We also assumed that the STM tip is stronglycoupled with an average coupling strength V = 0 . eV with variance . eV , while we assumed a weaker coupling tothe substrate with average coupling strength V = 0 . eV and variance . eV . These concrete values are the resultsof optimization carried out by visually comparing the result of the calculation with the experimentally obtained data.In Fig3 we show the result of this calculation. The number of parameters in simulating the experimental situation is P r o b a b ili t y d e n s i t y log G (nano Siemens) Theorythiolated-Streptavidin experiment Figure 3: Comparison of measurement data from Figure 4A of Ref.[9] (yellow points) and our simulation (red line).Horizontal axis is logarithmic and conductance is in units of nS .enormous and the details of the structure of the protein and the quantum mechanical calculation contain a wide range ofapproximations and errors. Never the less, it is obvious that (66) gives conductances in the correct order of magnitudewith physically realistic coupling strengths. Both the average and the variance and the general shape of the resultingdistribution is compatible with the experimental finding. Note, that the inputs of (66) are the energies and the wavefunctions of the molecule. The theory presented here explains other aspects of this experiment such as the tip-substratedistance independence of the conductance distribution. The generalization of the Landauer-Büttiker formula revealed that decoherence plays an important role in the electrontransport properties of proteins and other biological molecules. Strong coupling to vibration modes and the resultingdecoherence explains the high conductance inside of these structures. When accessed via weak electrostatic links thesestructures look like insulators since electron transfer rate decaying with distance and temperature governs transportproperties. When accessed via strong covalent or nearly covalent bonds, the same structures show good conductanceproperties over long distances and at high temperatures. Beyond explaining novel experiments, what can be thebiological role of these effects? In 1941 Nobel Prize winner biochemist Albert Szent-Györgyi put forward[29] manyexamples when electrons travel over large distances very fast within a biomolecule or across the entire cell. Most ofhis problems are still open as Marcus theory strongly suppresses long distance transport over energy barriers. Strong14
PREPRINT - S
EPTEMBER
4, 2019decoherence via coupling of electronic and vibrational degrees of freedom can open up new possibilities to understandthe special electronic wiring of biological systems.
Electronic properties of medium size proteins (Myoglobin, Cytochrome C and Streptavidin) have been calculated withthe semiempirical extended Hückel method implemented in the YaEHMOP package http://yaehmop.sourceforge.net . Structures have been taken from RCSB PDB . Computer codes for conductancecalculations in Matlab and Phyton are provided as Supplementary material.
G.V. thanks Stuart Lindsay introduction and guidance in the subject and for sharing the data of his group’s experiments.This research was supported by the National Research Development and Innovation Office of Hungary (Project No.2017-1.2.1-NKP-2017-00001). A We start with the definition of the inverse operator (cid:88) pq L nmpq L − pqkl = δ nk δ ml , (76)and separate the operator to a coupling strength dependent part plus the operator introduced in (55) (cid:88) pq (( − / (cid:126) )(Γ n / m / δ np δ mq + L nmpq ) L − pqkl = δ nk δ ml . (77)We take the diagonal elements ( n = m ) and get − Γ n (cid:126) L − nnkl + (cid:88) pq L nnpq L − pqkl = δ nk δ nl . (78)We can summ up the left and right hand sides for n and get − (cid:88) n Γ n (cid:126) L − nnkl + (cid:88) pq (cid:88) n L nnpq L − pqkl = δ kl . (79)Then (cid:80) n L nnkl = 0 due to the probability conservation law (cid:80) n R nnkl = 0 , which can be be verified using (13). Thensetting k = l we get the sum rule (cid:88) n Γ n (cid:126) L − nnkk = − . (80)The the expression (27) of the current can be written also as I = e U (cid:126) (cid:88) n D n ( E F , T ) (cid:34) Γ Ln + 1 (cid:126) (cid:88) k (Γ k − Γ Rk ) L − kknn (Γ Ln − Γ Rn ) (cid:35) , (81)and using the sum rule we can carry out the summation for k for the first part and get I = e U (cid:126) (cid:88) n D n ( E F , T ) (cid:34) Γ Rn − (cid:126) (cid:88) k Γ Rk L − kknn (Γ Ln − Γ Rn ) (cid:35) . (82)This expression is the same as (27) just the L and R indices are interchanged. To get a formula manifestly symmetric inthese indices we have to add up (27) and (82) and divide by two to get (28).15 PREPRINT - S
EPTEMBER
4, 2019 B We can start from the operator L nmkl = ( − / (cid:126) )(Γ n / m / δ nk δ ml + L nmkl . (83)The first term containing Γ can be regarded as a perturbation to the second term. The unperturbed system has a zeroeigenvalue, since it corresponds to a closed system in equilibrium (cid:88) kl L nmkl (cid:37) kl = 0 , (84)where the right eigenvector is the equilibrium density matrix. Since electrons and holes are not mixed, we have two sepa-rate right eigenvectors (cid:37) nn = e − ( ε n − E F ) /kT /Z ( T ) for electrons (and zeros for holes) and (cid:37) nn = e − ( E F − ε n ) /kT /Z ( T ) for holes (and zeros for electrons), so there is an eigenvalue λ = 0 both in the electron and in the hole sector. The lefteigenvector corresponding to the zero eigenvalue is δ nm (assuming n and m are both electron or hole states), which wecan verified by direct substitution (cid:88) nm δ nm L nmkl = (cid:88) n L nnkl = 0 , (85)where we used the probability conservation law (cid:80) n R nnkl = 0 , which can be be verified using (13). In the first orderof perturbation theory the new eigenvalue replacing the zero eigenvalue λ = 0 becomes the perturbation operatorsandwiched between the left and the right unperturbed eigenvectors λ = (cid:88) nmkl δ nm ( − / (cid:126) )(Γ n / m / δ nk δ ml (cid:37) kl = − (cid:88) n Γ n (cid:126) (cid:37) nn . (86)The eigenvalues of the inverse operator L − nmkl are the reciprocals or the eigenvalues of the operator. Since the perturbedeigenvalue is very small, its reciprocal will dominate the inverse operator in leading order and we can write it in termsof the corresponding left and right eigenvectors L − nmkl ≈ − δ nm (cid:37) nn δ kl (cid:80) p Γ p (cid:126) (cid:37) pp , (87)where all four indices and the summation should be done for electron and hole states separately. C The left and right probabilities are dominated by the orbital close to the Fermi energy. The first two most importantcontributions in case of holes come from the HOMO orbital and from HOMO-1 (the orbital below the HOMO) P Lh = Γ LHOMO e − ( E F − ε HOMO ) /kT + Γ LHOMO − e − ( E F − ε HOMO − ) /kT + ... (Γ LHOMO + Γ
RHOMO ) e − ( E F − ε HOMO ) /kT + (Γ LHOMO − + Γ RHOMO − ) e − ( E F − ε HOMO − ) /kT + ... . (88)This can be written as P Lh = Γ LHOMO Γ LHOMO + Γ
RHOMO
LHOMO − / Γ LHOMO ) e − ( ε HOMO − ε HOMO − ) /kT + ... LHOMO − + Γ RHOMO − ) / (Γ LHOMO + Γ
RHOMO )) e − ( ε HOMO − ε HOMO − ) /kT + ... , (89)and then we can expand in the small parameter e − ( ε HOMO − ε HOMO − ) /kT and get in leading order and group the terms P Lh = Γ LHOMO Γ LHOMO + Γ
RHOMO + Γ
LHOMO Γ RHOMO (Γ LHOMO + Γ
RHOMO ) (cid:34) Γ LHOMO − Γ LHOMO − Γ RHOMO − Γ RHOMO (cid:35) e − ( ε HOMO − ε HOMO − ) /kT + ... . (90)The sign of the temperature dependent part is determined by the ratios Γ LHOMO − / Γ LHOMO and Γ RHOMO − / Γ RHOMO which are the relative strengths of couplings of the left and right electrodes to the HOMO and HOMO-1 orbitals. Byexchanging left and right we can get the probability P Rh = Γ RHOMO Γ LHOMO + Γ
RHOMO − Γ LHOMO Γ RHOMO (Γ LHOMO + Γ
RHOMO ) (cid:34) Γ LHOMO − Γ LHOMO − Γ RHOMO − Γ RHOMO (cid:35) e − ( ε HOMO − ε HOMO − ) /kT + ... . (91)16 PREPRINT - S
EPTEMBER
4, 2019In case of electrons the LUMO and LUMO+1 orbitals play the same role and we can get the analogous expression P Le = Γ LLUMO Γ LLUMO + Γ
RLUMO + Γ
LLUMO Γ RLUMO (Γ LLUMO + Γ
RLUMO ) (cid:34) Γ LLUMO +1 Γ LLUMO − Γ RLUMO +1 Γ RLUMO (cid:35) e − ( ε LUMO +1 − ε LUMO ) /kT + ..., (92)and P Re = Γ RLUMO Γ LLUMO + Γ
RLUMO − Γ LLUMO Γ RLUMO (Γ LLUMO + Γ
RLUMO ) (cid:34) Γ LLUMO +1 Γ LLUMO − Γ RLUMO +1 Γ RLUMO (cid:35) e − ( ε LUMO +1 − ε LUMO ) /kT + ... . (93) D We start with the definition of the inverse operator (cid:88) pq L − nmpq L pqkl = δ nk δ ml , (94)and separate the operator to a coupling strength dependent part plus the operator introduced in (55) (cid:88) pq L − nmpq (( − / (cid:126) )(Γ k / l / δ pk δ ql + L pqkl ) = δ nk δ ml . (95)We take the diagonal elements ( k = l ) and get − L − nmkk Γ k (cid:126) + (cid:88) pq L − nmpq L pqkk = δ nk δ mk . (96)We can multiply both sides with the Boltzmann distribution p Bk , sum up for k and get − (cid:88) k L − nmkk Γ k (cid:126) p Bk + (cid:88) pq L − nmpq (cid:88) k L pqkk p Bk = δ nm p Bn (97)Then (cid:80) k L pqkk p Bk = 0 since the Boltzmann distribution is the steady state solution. Then setting m = n we get thesum rule (cid:88) k L − nmkk Γ k (cid:126) p Bk = − p Bn . (98) References [1] Christopher D Bostick, Sabyasachi Mukhopadhyay, Israel Pecht, Mordechai Sheves, David Cahen, and DavidLederman. Protein bioelectronics: a review of what we do and do not know.
Reports on Progress in Physics ,81(2):026601, 2018.[2] Nadav Amdursky, Debora Marchak, Lior Sepunaru, Israel Pecht, Mordechai Sheves, and David Cahen. Electronictransport via proteins.
Advanced Materials , 26(42):7142–7161, 2014.[3] Nicholas Kaldor. Capital accumulation and economic growth. In
The theory of capital , pages 177–222. Springer,1961.[4] Juan M Artés, Ismael Díez-Pérez, and Pau Gorostiza. Transistor-like behavior of single metalloprotein junctions.
Nano letters , 12(6):2679–2684, 2011.[5] Bintian Zhang, Weisi Song, Pei Pang, Yanan Zhao, Peiming Zhang, István Csabai, Gábor Vattay, and StuartLindsay. Observation of giant conductance fluctuations in a protein.
Nano futures , 1(3):035002, 2017.[6] Bintian Zhang, Weisi Song, Pei Pang, Huafang Lai, Qiang Chen, Peiming Zhang, and Stuart Lindsay. Role ofcontacts in long-range protein conductance.
Proceedings of the National Academy of Sciences , 116(13):5886–5891,2019.[7] Si-Youl Yoo, Taek Lee, Yong-Ho Chung, Junhong Min, and Jeong-Woo Choi. Fabrication of biofilm in nanoscaleconsisting of cytochrome f/2-maa bilayer on au surface for bioelectronic devices by self-assembly technique.
Journal of nanoscience and nanotechnology , 11(8):7069–7072, 2011.17
PREPRINT - S
EPTEMBER
4, 2019[8] Katalin V Korpany, Pinky Langat, Dong Myeong Kim, Neil Edelman, Daniel R Cooper, Jay Nadeau, andAmy Szuchmacher Blum. Conductance switching in the photoswitchable protein dronpa.
Journal of the AmericanChemical Society , 134(39):16119–16122, 2012.[9] Bintian Zhang and Stuart Lindsay. Electronic decay length in a protein molecule.
Nano letters , 2019.[10] Lior Sepunaru, Israel Pecht, Mordechai Sheves, and David Cahen. Solid-state electron transport across azurin:From a temperature-independent to a temperature-activated mechanism.
Journal of the American ChemicalSociety , 133(8):2421–2423, 2011.[11] CJ Lambert. Basic concepts of quantum interference and electron transport in single-molecule electronics.
Chemical Society Reviews , 44(4):875–888, 2015.[12] Rudolph A Marcus. On the theory of oxidation-reduction reactions involving electron transfer. i.
The Journal ofChemical Physics , 24(5):966–978, 1956.[13] Dvira Segal, Abraham Nitzan, William B Davis, Michael R Wasielewski, and Mark A Ratner. Electron transferrates in bridged molecular systems 2. a steady-state analysis of coherent tunneling and thermal transitions.
TheJournal of Physical Chemistry B , 104(16):3817–3829, 2000.[14] Abraham Nitzan. Electron transmission through molecules and molecular interfaces.
Annual review of physicalchemistry , 52(1):681–750, 2001.[15] Abraham Nitzan. A relationship between electron-transfer rates and molecular conduction.
The Journal ofPhysical Chemistry A , 105(12):2677–2679, 2001.[16] Jakub K Sowa, Jan A Mol, G Andrew D Briggs, and Erik M Gauger. Beyond marcus theory and the landauer-büttiker approach in molecular junctions: A unified framework.
The Journal of chemical physics , 149(15):154112,2018.[17] Mohan Sarovar, Akihito Ishizaki, Graham R Fleming, and K Birgitta Whaley. Quantum entanglement inphotosynthetic light-harvesting complexes.
Nature Physics , 6(6):462, 2010.[18] Gregory S Engel, Tessa R Calhoun, Elizabeth L Read, Tae-Kyu Ahn, Tomáš Manˇcal, Yuan-Chung Cheng, Robert EBlankenship, and Graham R Fleming. Evidence for wavelike energy transfer through quantum coherence inphotosynthetic systems.
Nature , 446(7137):782, 2007.[19] Masoud Mohseni, Patrick Rebentrost, Seth Lloyd, and Alan Aspuru-Guzik. Environment-assisted quantum walksin photosynthetic energy transfer.
The Journal of chemical physics , 129(17):11B603, 2008.[20] Ferdows Zahid, Magnus Paulsson, and Supriyo Datta. Electrical conduction through molecules. In
AdvancedSemiconductor and Organic Nano-Techniques , pages 1–41. Elsevier, 2003.[21] Supriyo Datta.
Electronic transport in mesoscopic systems . Cambridge university press, 1997.[22] Heinz-Peter Breuer, Francesco Petruccione, et al.
The theory of open quantum systems . Oxford University Presson Demand, 2002.[23] Greg Lever, Daniel J Cole, Nicholas DM Hine, Peter D Haynes, and Mike C Payne. Electrostatic considera-tions affecting the calculated homo–lumo gap in protein molecules.
Journal of Physics: Condensed Matter ,25(15):152101, 2013.[24] Abraham Nitzan. The relationship between electron transfer rate and molecular conduction 2. the sequentialhopping case.
Israel journal of chemistry , 42(2-3):163–166, 2002.[25] Sara Raichlin, Israel Pecht, Mordechai Sheves, and David Cahen. Protein electronic conductors: hemin–substratebonding dictates transport mechanism and efficiency across myoglobin.
Angewandte Chemie , 127(42):12556–12560, 2015.[26] Nadav Amdursky, Doron Ferber, Carlo Augusto Bortolotti, Dmitry A Dolgikh, Rita V Chertkova, Israel Pecht,Mordechai Sheves, and David Cahen. Solid-state electron transport via cytochrome c depends on electroniccoupling to electrodes and across the protein.
Proceedings of the National Academy of Sciences , 111(15):5556–5561, 2014.[27] Emil Wierzbinski, Ravindra Venkatramani, Kathryn L Davis, Silvia Bezer, Jing Kong, Yangjun Xing, Eric Borguet,Catalina Achim, David N Beratan, and David H Waldeck. The single-molecule conductance and electrochemicalelectron-transfer rate are related by a power law.
Acs Nano , 7(6):5391–5401, 2013.[28] Matthäus A Wolak, Alexander Balaeff, Sebastian Gutmann, Harry J Helmrich, Ruan Vosloo, Martin M Beerbom,Emil Wierzbinski, David H Waldeck, Silvia Bezer, Catalina Achim, et al. Electronic structure of self-assembledpeptide nucleic acid thin films.
The Journal of Physical Chemistry C , 115(34):17123–17135, 2011.[29] Albert SZENT-GYORGYI. The study of energy–levels in biochemistry.