A Self Energy Model of Dephasing in Molecular Junctions - The

Subscriber access provided by United Arab Emirates University | Libraries Deanship

Article

A Self Energy Model of Dephasing in Molecular Junctions Gabriele Penazzi, Alessandro Pecchia, Verena Gupta, and Thomas Frauenheim J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.6b04185 • Publication Date (Web): 05 Jul 2016 Downloaded from http://pubs.acs.org on July 9, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry C is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

A Self Energy Model of Dephasing in Molecular Junctions Gabriele Penazzi,∗,† Alessandro Pecchia,‡ Verena Gupta,¶ and Thomas Frauenheim¶ †Bremen Center for Computational Materials Science,Universit¨ at Bremen, Am Fallturm 1, 28359 Bremen, Germany ‡CNR-ISMN, Via Salaria km 29.300, 00017 Monterotondo, Roma, Italy ¶Bremen Center for Computational Materials Science, Universit¨ at Bremen, Am Fallturm 1, 28359 Bremen, Germany E-mail: [email protected]

1

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Quantum decoherence plays an important role in the charge transport characteristics of molecular wires at room temperature. In this paper we propose a generalization of an electron-phonon dephasing model to non orthogonal LCAO basis. We implemented the model in combination with a Density Functional based Tight Binding (DFTB) theory framework and utilized it to model charge transport characteristics of an anthraquinone (AQ) based molecular wire. We demonstrate a modulation of Quantum Interference (QI) effects compatible with experiments and confirm the robustness of QI signatures with respect to dephasing. An analysis of the spatial localization of the dephasing process reveals that both the QI and the dephasing process are localized in the AQ region, hence justifying the general robustness of the transmission temperature dependence in different AQ based systems.

Introduction The enormous progress in the experimental investigation of single molecule mechanically controlled break junctions (MCBJ) 1 and ultra-thin self-assembled monolayers 2 shed new light on the physics of charge transport in molecular junctions. Quantum transport plays a central role in these systems: in single molecule junctions the transport occurs mostly via coherent tunneling assisted by the closest molecular frontier orbital and it can be often described in terms of Simmons 3 or resonant tunneling models. However in the last years we witnessed an increasing amount of experimental evidence breaking this simple paradigm and stressing the importance of temperature and dynamic fluctuations. Charge transport regime transition and intermediate transport regimes between coherent tunneling and thermally activated hopping have been observed in oligomers 4–6 and short DNA strands; 7 several independent experiments demonstrated that Quantum Interference (QI) effects can influence dramatically the conductance of molecular junctions. 8–12 Such effects have been suggested as a possible way to engineer efficient thermoelectric devices. 13,14 In both tunneling/hopping 2


Page 2 of 32

Page 3 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


transition and QI modulation the partial or total loss of quantum coherence plays a pivotal role: the investigation and understanding of partially coherent transport both theoretically and experimentally is indeed recognized as one of the pressing challenges in future molecular electronics research. 15 Coherent transport is commonly modeled by means of Landauer-B¨ uttiker theory or Non Equilibrium Green’s Function (NEGF) theory. 16–18 These approaches are combined with model Hamiltonian, empirical methods or at a more sophisticated level with Density Functional Theory (DFT) or Density Functional Tight Binding (DFTB) theory. DFT offers the advantage to deliver a solid parameter free ab-initio framework which is particularly suitable in the study of chemically heterogeneous systems and it can be therefore considered the workhorse of molecular electronic modeling. DFTB is a fast approximate method based on DFT, 19 which offers a trade off between empirical and ab-initio methods: it provides results reasonably close to DFT in a very wide range of biological and organic systems and it can be combined with NEGF in a way methodologically very close to DFT. 20 Even though molecular systems are typically rather small, fast approximate methods can also play a primary role if we consider that the statistical nature of those experiments or the presence of static and dynamic disorder calls for statistical analysis of hundreds of configurations 15 or averaging along molecular dynamics trajectories. 21,22 The problem of decoherence in molecular junctions has attracted so far the interest of several authors. 23–27 At the level of Landauer theory decoherence can be introduced by adding additional virtual electrodes, called B¨ uttiker probes. 28 The B¨ uttiker model is conceptually simple and offers a physically intuitive picture and for this reason it is commonly used to describe model Hamiltonian systems. It has been applied to the study of molecular wires and DNA chains either in combination with Green’s function or transfer matrix formalism. 7,24,29,30 The model was originally formulated to be used in combination with model systems and in linear response in the applied bias, making its generalization to LCAO basis and self consistent NEGF not straightforward. The virtual probes can be generalized be-

3



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

yond the linear response, the generalization is not unique and can be formulated in a way which may or may not introduce dissipation, as recently elaborated by other authors. 31,32 Alternatively, electron-phonon interactions are typically responsible for decoherence can be rigorously introduced at the level of NEGF within a many body perturbation theory approach: this has been implemented in the context of DFT 33,34 and DFTB, 18,35,36 among the others. The method has the advantage of a strong physical foundation and it is successfully used to reproduce Inelastic Electron Tunneling Spectra at low temperatures, however such rigorousness come at the price of a much larger computational expense with respect to coherent NEGF. In this work we will first introduce a dephasing model formulated in the language of NEGF Keldysh theory, extending the work of previous authors 37–39 to the non orthogonal LCAO basis case. The model describes dephasing in a computationally efficient way through a single phenomenological coupling parameter, similar in spirit to B¨ uttiker probes, nevertheless offering an advantage to rely on well defined approximations to the full inelastic electron-phonon model. Therefore it can be more easily generalized to be combined with arbitrary basis or to recover the more complex physics of the full treatment by relaxing some approximations. In the methodological sections, we will discuss the formulation and define different implementations based on the locality of the dephasing model. The numerical implementation has been developed in the open source DFTB+ code, 40 but most of the arguments are applicable at a general level of LCAO-DFT. In the second part of the paper we apply the model to a characteristic QI system: an anthraquinone based molecular wire. This system is well characterized experimentally and has been the object of intensive DFT studies, 10,12,25,41 allowing for a direct comparison to experimental data and different level of theory.

4


Page 4 of 32

Page 5 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Method We start by writing the non interacting electronic Hamiltonian for a molecular junction according to the partitioning commonly used in NEGF methods 18 as

H = HL + HR + HD + VLD + VRD

(1)

where HL , HR represent the projection of the total Hamiltonian on a left and right electrode respectively, HD the projection on a device region comprising the molecule and possibly electrode buffers and VLD , VRD the interaction between the device region and the left and right electrode respectively. We assume then that the device region is interacting with a phonon bath

HD = Hel + Helph

(2)

The electronic component on the molecular region is written as a non interacting single particle Hamiltonian

Hel =

X

εiµ a†iµ aiµ +

X

εiµ,jν a†iµ ajν

(3)

iµ6=jν

iµ

where we assume for generality an atom based local orbital (LCAO) representation. i, j are the indexes on the atomic position and µ, ν are the indexes on the local atomic basis. The electron-phonon interaction is written in the second quantization form as 16,18

Helph =

X

λ Uiµjν a†iµ ajν

iµjνλ

b†λ

+ bλ

(4)

where λ is an index running on phonon modes and U is the electron-phonon coupling for a given mode. When the electron phonon interaction is not present the current in the system can be calculated by applying the Landauer formula 18,42

5



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2e I= h

Z

Page 6 of 32

+∞

(fL − fR ) T (ε)dε

−∞

(5)

where fL(R) = f (µL(R) , T, ε) is the Fermi-Dirac distribution of the left (right) electrode. The transmission function T (ε) can be calculated with the Caroli formula 43 T (ε) = Tr [ΓL (ε)Gr (ε)ΓR (ε)Ga (ε)] where ΓL(R) (ε) = i

h

ΣrL(R) (ε)

−

ΣaL(R) (ε)

i

(6)

; Σr,a L(R) are the retarded and advanced self energy

of the electrodes and Gr,a are the retarded and advanced equilibrium Green’s function of the device region. The Landauer-Caroli formula only holds for a non-interacting Hamiltonian. When electronphonon interaction is taken into account, we have to resort to the more general Non Equilibrium Green’s Function theory framework. The current flowing through the left (right) electrode is calculated using the Meir-Wingreen formula 44

IL(R)

2e = h

Z

+∞

Tr −∞

h

> Σ< L(R) (ε)G (ε)

−

< Σ> L(R) (ε)G (ε)

i

dε

(7)

Σ are the lesser and greater Self Energy of the lead and are defined from the equilibrium α quantities as follows:

Σ< L(R) (ε) = ifL(R) ΓL(R) (ε) Σ> (ε) = −i 1 − f ΓL(R) (ε) L(R) L(R)

(8)

G represent the lesser and greater Green’s function of the device region and can be calculated via kinetic equation

G

)

=G

r

h

) ΣL

+ 6

) ΣL

+

) Σφ


i

Ga

(9)

Page 7 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Interactions here are represented in terms of additional self energies Σφ which can be calculated by applying many body perturbation theory. Notice however that the above formulation is valid provided the interactions are restrained to the device region only. We will work at the level of Self Consistent Born Approximation (SCBA) 16,17 and retain only the self energies describing the Fock diagram in the electron-phonon interaction: the Hartree diagram does not modify the quasi particle lifetime which is the quantity of interest in the context of dephasing. The SCBA lesser (greater) self energy in matrix representation can be calculated as

Σ(ε)

)

X i Z ) = U λ G) (ε − ε′ )U λ Dλ (ε′ )dε′ 2π λ

(10)

where Dλ is the single phonon mode propagator and U the electron-phonon coupling as previously introduced. The retarded and advanced self energies can be calculated consequently by exploiting the relations between G and Gr,a . 16 This general formulation of the SCBA model has been numerically implemented by other authors in the context of molecular modeling to describe Inelastic Electron Tunneling Spectra (IETS); 33,35 in this case the electron phonon couplings and vibrational modes are commonly calculated by DFT or DFTB by finite difference method. This kind of implementation, even though accurate, is computationally demanding: besides the computational burden of calculating the vibrational spectra, the resulting self energies are dense thus preventing a number of numerical optimization based on the sparsity of the self energies. 20 In fact, in the context of atomistic semiconductor device modeling, a number of ad-hoc approximation are commonly introduced to decrease computational burden and memory consumption. 45 In this work we apply a simple model aimed at describing dephasing effects through a single phenomenological parameter, similar to the work from previous authors 37–39,46 and here generalized beyond the nearest neighbor s-like tight binding model. The model is obtained by introducing some approximations which largely simplify the electron-phonon self energy: we assume that the vibrational modes are represented by Einstein phonons strongly localized on 7



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 32

an atomic site (λ = i) and the electron-phonon coupling is described by a local deformation potential identical on all orbitals. The interaction Hamiltonian reduces to

Helph =

X iµν

Ui a†iν aiµ b†i + bi

(11)

The corresponding retarded (advanced) and lesser (greater) self energies are further simplified if we assume that the phonon energy h ¯ Ωλ is negligible compared to the other energy scales in the system, i.e. room temperature and energy dependent features in the spectral function. Within this approximation the interaction self energy represents elastic scattering events and can be written (see 37 ) as

2 r,a, Σr,a, iµ,iν (ε) = γi Giµ,iν (ε)

(12)

γi is regarded in this work as a phenomenological parameter monotonically increasing with temperature, its physical interpretation and magnitude are discussed at the end of this section. This set of approximations has been referred as electron-phonon dephasing model , 37 or simply dephasing model . 46 The same expression for the self energy can be derived for different local and elastic scattering mechanisms: for example it has been used to describe intermolecular contributions in SAMs. 47 We will hold to the interpretation in terms of electron-phonon scattering as this is likely the process responsible for dephasing in the system under investigation in this paper. It should be however emphasized that in the context of molecular electronics some of the approximations introduced in this work are not valid in general: vibronic states usually have a rich spectrum combining both strongly and weakly localized states and the spectral function exhibits peaks related to localized states weakly bounded with the electrodes or interference features. The purpose of the present work is not to justify the application of this model in substitution of a full electron-phonon treatment, but rather to propose a simple way to introduce dephasing phenomena in a computationally efficient way. The dephasing model shares two key physical assumption with B¨ uttiker probe

8


Page 9 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


based models: they are both local, i.e. the self-energy in eq. 12 is diagonal in a single orbital tight binding basis as in the Green’s function formulation of B¨ uttiker probes 29 and they are both elastic. Note that the choice of local self energies is justified in our case by the assumptions of local oscillators. A similar model but with completely dense self energies has also been proposed and applied to molecular wires. 26,46 The formulation in terms of Keldysh self energy has some peculiarities. It can be integrated flawlessly in NEGF implementations and current conservation is ensured at the level of SCBA approximation also beyond linear response and at finite temperature, while at the level of B¨ uttiker probes a suitable generalization has been formulated recently. 31 The second advantage is that the approximations can be gradually relaxed bridging the phenomenological treatment with a full ab-initio description of electron-phonon interactions, for example to include dissipation. We can use this to our advantage to extend the model, up to now only applied to simple model Hamiltonian or single band effective mass models, to a non orthogonal LCAO basis where the application of the B¨ uttiker model is not straightforward. The self-energy in eq. 12 is diagonal only if we use one basis per atom or alternatively if we consider each mode λ coupled to a single orbital. This assumption is arbitrary: it only allows us to provide a fully diagonal model which is still a valid SCBA self energy fulfilling particle conservation correctly. We will refer to this model as Fully Diagonal (FD). We can alternatively keep the definition in eq. 12 when an arbitrary number of basis orbitals per atom are used: the resulting self energy is non-zero on blocks corresponding to atomic degrees of freedom. We will refer to this model as Block Diagonal (BD). So far we have ignored the overlap between atomic orbitals, i.e. we assumed that the electron-phonon coupling is non zero only on diagonal elements. While this assumption can be justified in an orthogonal basis by considering matrix elements of the local deformation potentials as < i|U |j >= Uii δij , it becomes clearly questionable in a non-orthogonal basis where off-diagonal contribution are expected. In fact at a DFTB level it is easily seen that the electron-phonon coupling should have at least the sparsity pattern of the overlap matrix

9



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 32

S. 35 Even though the electron-phonon coupling is not calculated ab-initio, we can generalize the definition of γi to take into account the overlap in an approximate way. In DFTB the projection of a local potential onto the LCAO basis is based on an expression alike: 20

Viµjν = vi

Siµjν Siµjν + vj 2 2

(13)

where vi is a potential projected on a single atom i. If we consider the operator M diagonal in matrix representation, we can write an expression analogous to the one above as ˜ = M S + SM M 2 2

(14)

The expression above can also be seen as the first order of a series expansion of the inverse Löwdin transformation for the operator M defined in the Löwdin orthogonal basis to the ˜ in the original non-orthogonal basis, i.e. operator M ˜ = S 21 M S 12 = M S + SM + O(S 2 ) M 2 2

(15)

We now consider M i as a matrix representation of the local deformation potential in i an orthogonal basis such that M i is zero everywhere except for Miµiµ = γi and apply the

transformation 14. The total self energy is written as

Σr,a, (ε) =

X

˜ i Gr,a, (ε)M ˜i M

(16)

i

For our particular choice of M i the self energy has the same sparsity pattern as S, thus introducing basis dependent corrections proportional to the overlap. In other words, if we would define a vector of phenomenological couplings γi in a Löwdin orthogonal basis and we then rotate this basis back to the DFTB one (or any non-orthogonal one), the couplings are approximately given by eq. 14. We will refer to this model as Overlap Masked (OM), from the sparsity pattern of the self energy.

10


Page 11 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The three models (FD,BD,OM) are very similar, the only difference being a different arbitrary definition of the relation between the electron-phonon coupling and the coupling parameters γi . However the FD and BD models are computationally much more convenient with respect to the OS model, the reason being that the total self energies in the FD and BD models can be easily calculated at once from the diagonal (or block diagonal) elements of the Green’s functions, while in the OS model an explicit summation over the modes is needed. Furthermore, the FD and BD model do not break the iterative NEGF algorithm implemented in DFTB+ 20 which provides a huge speedup and memory saving in calculations and allows for simulation of systems up to several thousands atoms. All the proposed models are elastic, thus allowing to preserve also the parallelization on energy grid point already implemented for the non-interacting case. The elastic nature of the model also allows us to write the current obtained in a Landauer-like fashion. In fact the energy resolved current is identical (in absolute value) regardless of the lead we choose as a reference, i.e. > < > > < > < IL (ε) = −IR (ε) = Tr [Σ< L (ε)G (ε) − ΣL (ε)G (ε)] = − Tr [ΣR (ε)G (ε) − ΣR (ε)G (ε)] at

any energy. In the non interacting case the Landauer formula includes finite temperature effects via the Fermi distribution of the electrodes. The differential conductance at a given chemical potential µ at any given temperature can be calculated from the zero temperature transmission by applying the convolution 48

G(µ, T ) = G0

Z

+∞ −∞

df (ε, µ, T ) T (ε) − dε

dε

(17)

between the zero temperature transmission and the derivative of the Fermi distribution, or thermal broadening function. When we apply the dephasing model the electrode Fermi distributions are already included in the lesser and greater Green’s function via the kinetic equation. However in the elastic limit it is possible to carry out those contributions (a proof is given in appendix) and to write the energy resolved current in the form 11



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

IL (ε) = (fL − fR )Tef f (ε)

Page 12 of 32

(18)

where Tef f can be thought as an effective transmission and we can apply the thermal broadening as in eq. 17 a posteriori. In this way we can discriminate the broadening in the transmission due to dephasing from the broadening due to finite temperature of the electrodes. This is only possible as the present model does not introduce dissipation in the device region, while it is not in general possible to write the interacting current in terms of an effective transmission. To conclude this section we attempt to define a physical interpretation for the dephasing parameter γ and an indicative value which can be used in practice. We start with a single site model coupled to a single harmonic oscillator with frequency ω. In the high occupation limit, i.e. assuming kB T ≫ h ¯ ω, γ can be approximated by the expression 37 γ 2 ≈ 2U 2 kB T /(¯hω) where U is the electron-phonon coupling matrix element defined in eq. 4. The quantity U 2 /(¯hω) corresponds to the polaron formation energy Vpol , 49 therefore we can write γ ≈ 2kB T Vpol . Alternatively, we recall that the electron-phonon coupling element can be written as U 2 = h ¯ / (2mω) × (∂H/∂Q)2 where m is the reduced oscillator mass, H the hamiltonian and Q the displacement of the ions along the the normal vibrational mode. Therefore we can write the coupling parameter as 1 γ ≈ kB T × 2mω 2

∂H ∂Q

2

= kB T Velph

(19)

Velph represents the intramolecular electron phonon coupling energy, as defined in BCS theory. 50 The same expression was derived by other authors 38 where Velph is substituted by the intramolecular reorganization energy Vreorg . This is not surprising, as all the quantities are closely related to each other. 51 Equation 19 gives an indication of the temperature depen√ dence of the parameter γ, i.e. γ ∝ T and an order of magnitude for the coupling in a single 12


Page 13 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


orbital picture. The electron-phonon coupling energy of aromatic molecules is in the range 0.05 − 0.3eV, 50 corresponding to a parameter γ in a range 0.03 − 0.09eV at room temperature. However the reorganization energy is dependent on the system size and contributed by several modes and the result needs to be adjusted to account for a local basis and for the multi-phonon case. The multi-phonon coupling energy can be written to first approximation as Velph ≈ P

Vi where Vi are contributions given by the coupling of single localized Einstein oscillators

coupled with molecular orbitals. Such a decomposition works generally well in aromatic compounds 50 and it is similar in spirit to the decomposition of the inner sphere reorganization energy used in charge transfer theory. 52 The coupling energy commonly decreases with the size of the system as Velph ≈ 1/N , where N is the atomic extension of the molecular orbital. Therefore to first approximation we have Vi ≈ N Velph . We assume a dephasing parameter proportional to Vi and define γi2 ≈ N kB T Velph . As a consequence the local coupling does not have a strong dependence on the degree of localization of the molecular orbitals and/or on system size, therefore it should be to some degree portable between similar systems. As an example, if we consider aromatic compounds, Velph ≈ 1.8eV/N 50 where N is the number of carbon atoms, leading to γ ≈ 0.2eV at room temperature. In this estimation we ignored many details about the actual calculation of the reorganization energy such as orbital degeneracy or bond order contributions, therefore we can not expect a quantitative fit with more accurate models and we stress that the model should be considered phenomenological in nature. Nevertheless these arguments provide an indication of the order of magnitude for the dephasing coupling and important results about the dependence on temperature and orbital size. In order to demonstrate its validity, we applied the dephasing model on a metabenzene junction, a simple system exhibiting a destructive QI which resembles the behavior of AQ wires. We obtained a perfect agreement in the current behavior between a two sites minimal tight binding model and a H¨ uckel model, provided that we scale properly the dephasing coupling according to the degrees of freedom of the

13



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

system. We also compared our result with available data for the inelastic SCBA model in molecular orbital basis. 25 With an appropriate choice of Velph as fitting parameter, the current-temperature characteristic of the different models get almost indistinguishable at room temperature and generally get closer to each other at increasing temperature above 150K, compatibly with the underlying physical assumptions. The details of the calculations are shown in the supporting information.

Results We apply the model previously described to model the effect of dephasing in a molecular junctions with QI. To this end we consider a wire with a central cross-conjugated anthraquinone (AQ) region and a wire with a linearly conjugated anthracene (AC) region. Molecular wires based on AC and AQ have the same length and similar HOMO-LUMO gap, however their conductance differ by one to almost three order of magnitude at room temperature. 10,12,41 This is attributed to a strong destructive interference effect, the robustness of this effect with respect to temperature has been demonstrated by means of a model Hamiltonian approach 25 and is experimentally confirmed by the fact that the conductance of AQ is systematically some orders of magnitude lower than the reference AC system, at any temperature up to room temperature. In this work we model the molecular wires in a mechanically controlled break junction configuration as shown in in Figure 1. The gold electrodes are oriented in the [111] direction and the extended device region includes three lattice monolayers. The electrodes and the molecular region are described via a DFTB Hamiltonian applying a recently developed parametrization for gold-organic compounds, 53 the dephasing model is only applied on the subspace of the extended molecular region corresponding to the molecular wire. The choice of the electrode configuration is arbitrary: even though the exact value of conductance may depend strongly on the exact morphology and bonding configuration at the

14


Page 14 of 32

Page 15 of 32


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment



Page 16 of 32

Page 17 of 32




1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

several orders of magnitude due to the partial suppression of the QI feature. The effect of the block diagonal model is all over slightly more pronounced. As an additional feature of the transmission spectra, we notice a characteristic transition from a Lorentzian shape of the LUMO resonance to a semi-elliptical one. The transition between a Lorentzian resonance spectrum, typical of weakly coupled MOs and well described in terms of Breit-Wigner theory and semi-elliptical spectrum is characteristic of many body corrections in the SCBA for a strong perturbation. 54 For this effect to appear, as pointed out also by other authors, 38 it is necessary a coupling significantly larger than the coupling between the molecular orbital and the electrode. This may suggest that in near resonance systems in weak coupling regime dephasing may have also a direct impact on charge injection properties. A strong coupling regime is usually not addressed in a full electron phonon inelastic approach: the computational cost usually does not allow for a large number of SCBA iterations and some approximations are usually introduced. 25 The dephasing model on the other hand is computationally very efficient and does not suffer from these problems. We notice that the convergence of the SCBA is extremely sensitive to energy. While close to the antiresonance dip and to the HOMO a few cycles suffice (typically no more than three), close to the LUMO a very large number of iterations can be needed (in the order of twenty cycles), due to the fact that the quasi-particle lifetime is completely determined by the dephasing model. The effective transmission of the AC wire as shown in Figure 4 is for the block diagonal model only. As the transport is simply due to the off resonance transmission in between the HOMO-LUMO gap and the transmission has no strong energy dependent modulation as in the QI case, the dephasing model leaves the transport properties almost unaffected. Some modulation is still present due to broadening of the HOMO level close to the DFTB Fermi energy. However, it should be noticed that the exact position of the Fermi level is hardly captured at a DFTB, DFT level and the experimental value is more likely to be closer to the middle of the HOMO-LUMO gap where almost no modulation is observed. The position

18


Page 18 of 32

Page 19 of 32




1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

In proximity of the dip the thermal broadening can alone explain a conductance increase up to two orders of magnitude at room temperature: yet this effect alone is not large enough to justify differences larger than two orders of magnitude observed experimentally. 12 The dephasing introduces an additional broadening which becomes dominant when γ ≫ kB T . The conductance of the AQ wire is systematically more than two orders of magnitude lower with respect to AC, regardless of thermal broadening and dephasing and consistent with experimental measurements. 41 This confirms numerically that QI effects can survive at room temperature. Our model is obviously missing satellite peaks visible as low temperature features. 12 Those signatures are due to the openings of additional transport channels in energy, the same physical effect exploited in IETS measurements. Being the dephasing model elastic, those signatures can not be captured. However, these effects are completely washed out when temperature is increased and finite lifetime effects begin to become dominant and are experimentally completely indistinguishable already at a temperature around 100K. We can gain better physical insight in the dephasing mechanism by limiting the application of the model to those atomic sites contributing more to the transport properties close to the anti-resonance. In the case of AQ the anti-resonance behavior is supposedly due to destructive interference between the HOMO LUMO and HOMO-1 LUMO paths. The HOMO and HOMO-1 levels are spread all over the molecular wire including the benzene linkers, while the LUMO is more strongly localized in the central AQ molecule. Molecular transport and QI phenomena can be to some extent interpreted in terms of local pathways 55 in a chemical analogous of a mesoscopic interferometer, where the waveguides are substituted by chemical bonds. In this spirit, the dephasing process will not affect the QI features homogeneously all over the system, but we may rather make the hypothesis that decoherence has a more dramatic impact when it occurs on some specific pathways. Therefore, if we apply a finite γ only to some specific atoms, we can identify molecular regions more critical in terms of dephasing process. We have tried to impose a value γ 6= 0 only in some molecular 20


Page 20 of 32

Page 21 of 32




1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

electron-phonon interaction including the Fock diagram in the elastic limit. We proposed three different ways to generalize the model for the application to non orthogonal LCAO basis. The corresponding implementations results in varying degree of localization of the self energy and computational efficiency and deliver close results. We applied the model to an AQ and AC based molecular junctions. By comparing with ab-initio methods from literature we determine that the dephasing model can as well capture room temperature transport characteristics. The conductance exhibit a strong dependence on dephasing which is necessary to justify the large dependence of conductance on temperature observed experimentally. The finite temperature distribution of the injected carriers alone is not sufficient. By applying the model on different regions of the molecular wires we also reveal that the suppression of the quantum interference feature is completely determined by decoherence occurring on the AQ molecular region. The rest of the wires play little or no role in determining transport properties both at zero and finite temperature. This result is probably connected to a strong role of some specific vibrational modes in promoting coherence loss, yet to be verified. It also justifies the fact that the temperature dependent conductance measurements are very well reproduced both in MCBJ and SAM experiments even though the exact morphology of the wires in the two cases is not exactly the same. Both the interference effect and its suppression are dominated by electronic properties of the AQ unit itself and should be therefore robust with respect to different wire lenghts and configurations. To conclude the combination of the dephasing model with DFTB, a method aimed at computational efficiency itself, makes the dephasing model a useful tool for fast modeling of QI molecular junctions at room temperature and can find future applications in demanding simulations, for example for large scale ensemble averaging. The same formulation of the model can also be easily integrated with any LCAO method, including DFT.

22


Page 22 of 32

Page 23 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Ackowledgement Financial support by the Deutsche Forschungsgemeinschaft (FR2833/36-1, FR2833/50-1) is gratefully acknowledged.

Supporting Information Available Numerical comparison between the dephasing model and the inelastic SCBA model for a meta benzene junction. This material is available free of charge via the Internet at http: //pubs.acs.org/.

Appendix A Landauer-like expression for the current In a non interacting system the Landauer-Caroli formula is easily derived from the MeirWingreen. The demonstration is a text-book exercise, we recall here the main steps as we will follow the same path. We can write the MW formula in a more convenient form 42 (we omit energy arguments) 2e IL = i h

Z

+∞ −∞

Tr [ΓL G< − Σ< L A] dε

(20)

where A = i (Gr − Ga ) = Gr (ΓL + ΓR ) Ga is the spectral function. By substituting the < definition of Σ< L from equation 8 and recalling that we can write the non interacting G as

< a G< = Gr (Σ< L + ΣR ) G =

= ifL Gr ΓL Gr + ifR Gr ΓR Gr

23


(21)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 32

and substituting the definitions in equations 8 in the expression above we obtain the identity

r r i Tr [ΓL G< − Σ< L A] = (fL − fR ) Tr [ΓL G ΓR G ]

(22)

. In the interacting case we have an additional self energy which we will refer to as Σr,a, φ The relation above does not hold in general possible because the in general the lesser and greater interaction self energy can not in general be written in the same form. However in the elastic case we can still write an expression in the form 2e I= h

Z

+∞ −∞

(f (ε, µL , T ) − f (ε, µR , T )) Tef f (ε)dε

(23)

as we can write all the contributions to G as sum of products of equilibrium quantities and Fermi distributions at given energy. In the inelastic case this is of course not possible as contributions at different energies are mixed. If we can write the current in such form, then we can apply the thermal broadening a posteriori, via convolution with the Fermi derivative as in the non interacting case. 48 We show in the following an explicit expression resembling the Landauer-Caroli valid with the dephasing model. We introduce a label i on the SCBA cycle such that r,a;i−1 −1 Gr,a;i = εS − H − Σr,a L,R − Σφ a;i ;i−1 G G;i = Gr;i Σ L,R + Σφ

(24) (25)

with ΣL,R = ΣL + ΣR . The self energies Σr,a,;i are calculated from Gr,a,;i according φ to the equation 16. We can now write a Gamma function associated to the dephasing self energy Γφ = i Σrφ − Σaφ at any order in Born Approximation (BA). For example for the

first two orders we have

24


Page 25 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Γ0φ =

X

M λ Gr;0 (ΓL + ΓR ) Ga;0 M λ

(26)

Γ1φ =

X

M λ Gr;1 ΓL + ΓR + Γ0φ Ga;1 M λ

(27)

λ

λ

in the expression for Γnφ we see that the expressions By substituting recursively Γn−1 φ above can be generalized as a function of Gr,a;n , Gr,a;n−1 ...Gr,a;0 , ΓL,R at any BA order as

Γnφ

=

n X

˜ n,n−m ˜ n,n−m + Γ Γ R L

m=0

(28)

˜ n,n−m represents all the terms containing Gr,a;n ...Gr,a;m and is written as where Γ L

˜ n,n−m = Γ L

X λ

M λ Gr;n

"

X λ

"

M λ Gr;n−1 · · ·

"

X λ

#

#

#

M λ Gr;n−m ΓL Ga;n−m M λ · · · Ga;n−1 M λ Ga;n M λ (29)

Similarly, we can apply the definitions in eq. 8 to write the lesser self energy as

Σ

A Self Energy Model of Dephasing in Molecular Junctions - The

Recommend Documents