Charge-memory effect in a polaron model: equation-of-motion method for Green functions

We analyze a single-level quantum system placed between metallic leads and strongly coupled to a localized vibrational mode, which models a singlemolecule junction or an STM setup. We consider a polaron model describing the interaction between electronic and vibronic degrees of freedom and develop and examine different truncation schemes in the equation-of-motion method within the framework of non-equilibrium Green functions. We show that upon applying gate or bias voltage, it is possible to observe charge-bistability and hysteretic behavior which can be the basis of a charge-memory element. We further perform a systematic analysis of the bistability behaviour of the system for different internal parameters such as the electron-vibron and the lead-molecule coupling strength.


Introduction
Within the field of single-molecule electronics [1,2,3,4], beside experimental progress with regard to vibrational properties and their signatures in transport [5,6,7,8,9], related phenomena such as switching, memory effects and hysteretic behavior in molecular junctions have gained increasing importance and attention.Random and controlled switching of single molecules [10,11,12], as well as conformational memory effects [13,14,15] have been recently explored.Related to these effects, there is the so-called charge-memory effect, that is basically a hysteretic behaviour of the chargevoltage, respectively, current-voltage characteristics arising from the interplay between the polaron shift and Franck-Condon blockade [16] in the presence of electron-vibron interaction.Several works have recently addressed this interesting feature of molecular junctions, both experimentally [17,18] and theoretically [19,20,21,22,23].
The charge-memory effect can be explained in the framework of a simple singlelevel polaron model [20,23,24,25], where the electronic state is coupled to a vibronic mode with frequency ω 0 , see the sketch in Fig. 1.If the energy of the unoccupied electron level without electron-vibron interaction is ǫ 0 , the occupied (charged) state of the interacting system will have the energy ǫ 1 = ǫ 0 − ǫ p , where ǫ p is so-called polaron shift (called also recombination energy).Neutral and charged (polaron) states correspond to local minima of the potential energy surface and get metastable, if the electron-vibron interaction is strong enough.Upon applying an external voltage, one can change the state of this bistable system, an effect that is accompanied by hysteretic charge-voltage and current-voltage curves.In this model it is not necessary to include Coulomb interaction explicitly, though one can straightforwardly incorporate charging effects.
It was alternatively suggested in Refs.[21,22] that quantum switching between bistable states rather results in telegraph noise at finite voltage than in a memory effect.In a recent paper [23], considering the problem in the weak molecule-to-lead coupling limit, i.e. for the level width Γ ≪ ω 0 , ǫ p , we considered the crossover between these two pictures, if one takes into account the time scale of the switching process.Indeed, the switching time τ between the two states of interest should be compared with the characteristic time of the external voltage sweeping, τ s ∼ V (t)/(dV (t)/dt).For τ ≫ τ s , quantum switching can be neglected, and hysteresis can be observed, while in the opposite limit, τ ≪ τ s , the averaging removes the hysteresis.In Ref. [23] we found that at large enough electron-vibron interaction strenght λ/ω 0 (where λ is introduced below in Eq. ( 1)), the switching time τ is exponentially large compared to the inverse bare tunneling rate 1/Γ.
In this paper we consider the polaron problem in another limiting case, ω 0 ≪ Γ < ǫ p , corresponding to the regime where the Born-Oppenheimer approximation holds true.In this situation of intermediate molecule-to-lead coupling we address the case of time-independent applied voltages, considering then the stationary problem and focusing on the properties of the two states of interest.In this paper we consider the parameter ranges for which the fluctuations between the two charge states of interest are negligible.It should be noted that the mean-field solutions we consider are metastable, but they are physical in the case of very long switching time τ .Of course, the exact equilibrium solution is a superposition of these two states.We assume that all possible relaxation processes, with the exception of the switching between metastable states, are much faster than τ .Within the framework of nonequilibrium Green functions [26,27,28,29,30] we use an equation-of-motion (EOM) approach (Color) Schematic representation of a single-level model system interacting with a vibronic mode and coupled to left and right leads.
which allows for studying the appearance of the charge-memory effect addressed at different levels of approximation, starting with the self-consistent Hartree level.The work is partly built on and further develops ideas introduced in Refs.[24,25,20] and complements the results based on a master-equation approach in Ref. [23].
The paper is organized as follows: In Sec. 2 we introduce the model and the relevant quantities necessary to describe transport through the molecular junction.In Sec. 3 we outline the equation-of-motion method, used for the calculation of nonequilibrium Green functions, and discuss different truncation schemes.In the final Sec. 4 we present and discuss the results of our approach.

The single-level electron-vibron Hamiltonian
We consider a three-terminal device with the central part given by a single level interacting with a vibronic mode (see Fig. 1).Possible experimental implementations include metal-molecule-metal junctions or STM spectroscopy of a single molecule on a conducting substrate.To describe such a system we use the standard electron-vibron Here the terms in the first line describe the central system including one electronic state with energy ǫ 0 , one vibronic state with frequency ω 0 and their mutual interaction with a coupling strength λ.The second line in Eq. (1) contains the Hamiltonian of the two leads with independent-particle states and the tunneling between the leads and the central region via the couplings V i,k .The index i denotes the left and right leads, while k labels the electronic states of electrons in the leads.Below it proves furthermore convenient to employ the vibronic "position" and "momentum" operators The energy level ǫ 0 in Eq. ( 1) can be shifted through the gate voltage V G .We choose as reference energy ǫ 0 = 0 for V G = 0 and assume a linear capacitive coupling, ǫ 0 = αeV G putting α = 1.Note that the presence of a bias voltage V B = ϕ L − ϕ R can also change the energy of the electronic level via the parameter ϕ 0 = ϕ R + ηV B , where 0 < η < 1 describes the symmetry of the voltage drop across the junction: η = 0 corresponds to the completely asymmetric case, while η = 0.5 stands for the symmetric case.Thus in all approximations the bias and gate voltages are taken into account twofold, through the potentials of the leads, which can be chosen e.g. as ϕ L = V B /2, ϕ R = −V B /2, and the effective energy of the level, which is correspondingly ǫ0 = ǫ 0 + eϕ 0 = eV G + e(η − 0.5)V B .From this expression for ǫ0 it follows that in the case of asymmetric bias-voltage drop across the junction, η = 0, the energy of the unoccupied electron level will be centered around the electrochemical potential of the right lead and moved away from this value through the gate voltage.This ingredient will be crucial for the effect addressed in this paper because the additional presence of the polaron shift will then fix the energy of the occupied (charged) state below the electrochemical potential of the right lead resulting in a blocked charged state under appropriate parameter conditions.In the case of symmetric voltage drop, η = 0.5, the energy of the unoccupied electron level will be centered around the zero of the energy resulting in a different scenario for what concerns the memory effect.In this paper we consider the case of an asymmetric junction for which we show that the memory effect occurs for small bias voltage in a wide range of the parameters entering the model Hamiltonian (1).The symmetric situation can also give rise to a hysteretic behaviour but only at finite bias voltages, making this case less interesting for the memory effect addressed here.

Spectral function, average charge and current
To obtain physical information from the Hamiltonian (1) we use the Green function method within the equation of motion (EOM) formalism.This method is an alternative to the Green function techniques earlier applied to the considered problem in Refs.[20,23,24,25].We start with the retarded Green function for two generic operators a and b which is defined as where , denotes a thermal average and {a, b} = ab + ba is the anti-commutator.For the stationary case, the general expression (3) reduces to an object with only one time argument, For a single-level system the retarded Green function in Eq. ( 4) is obtained by replacing the generic operators a and b by the electronic operators.For the Hamiltonian (1), this reads From Eq. ( 5) one obtains the spectral function A(ǫ) of the system through the expression where G r (ǫ) is the Fourier transform of G r (t).The spectral function is the basic ingredient for obtaining the transport properties of the system such as average current and charge on the molecule.The expression for the current through the molecule is given by where f 0 i is the equilibrium Fermi function in the i-th lead.The tunneling couplings to the right (Γ R ) and left (Γ L ) leads are where the matrix elements V i,k are assumed to be energy-independent (wide-band limit).The full level broadening is given by the sum Γ = Γ L + Γ R .Below Γ R and Γ L are assumed to be the same.
The average charge (number of electrons), n = d † d , is given by where f (ǫ) is the distribution function of electrons inside the molecule.For the approximations used in this paper, we employ the same distribution function as in the non-interacting case, because we are focusing on the case of intermediate molecule-lead coupling.Fast tunneling into and out of to the molecule makes plausible the assumption that the electrons are in strong a non-equilibrium situation and can then be decribed via Eq.( 10) that is obtained assuming only elastic processes.Moreover, within the first approximation that we discuss in the paper, it is possible to obtain the distribution function (10) analytically by calculating the lesser Green function of the problem or applying the Hartree approximation directly to the Hamiltonian as in Ref. [20].
For a more complete explanation of the basic formulas introduced here the reader is referred to Refs.[28,29,30].

General formalism
The equations of motion for non-equilibrium Green functions are obtained from the Heisenberg equation for a Heisenberg operator a(t), Here and below all Hamiltonians are assumed to be time-independent because we consider the stationary problem, the applied voltages enter the problem only as timeindependent parameters, changing the position of the molecular energy level and the electrochemical potentials in the left and the right lead.With Eq. ( 11) the time derivative of the Green function, Eq. (4), reads After performing a Fourier transform of Eq. ( 12) we obtain where , indicates the Fourier transform of a given Green function.Equation ( 13) is the starting point for the EOM method.By applying successively a time derivative there are new high-order Green functions appearing.The idea of the EOM method is thereby to truncate this iterative process at some point making the mean-field like approximation of the highest-order Green function through lower-order functions, in order to obtain a closed set of equations.In our case we start from the Hamiltonian (1) and consider the first-oder equation for the function d, d † and second-order equations for the functions xd, d † and pd, d † .

EOM method for the single-level electron-vibron Hamiltonian
The EOM method for Hamiltonian (1) generates the expressions The equation for c i,k , d † is closed (including only the function d, d † ).By substituting Eq. ( 15) into Eq.( 14) and introducing the self-energy Σ of the leads through we obtain eventually The last term, describing the interaction between electron and vibron, has to be truncated at this level or found from higher-order equations and then truncated at a higher level of approximation.The lead self-energy will be used below in the wide-band approximation Σ(ǫ) = −iΓ.

Self-consistent Hartree approximation
The simplest way to close Eq.( 17) is to perform the truncation by approximating Then we obtain immediately for the Green function Here the quantity x remains to be calculated.To this end we compute, respectively, the time derivatives of the x-operator, and the p-operator, Upon combining Eqs. ( 20) and ( 21) we get In the stationary case addressed, Eq. ( 22) yields a direct connection between the "position" of the vibron and the particle number in the dot: In view of Eq. ( 19), we finally obtain for the spectral function ( 6) the following selfconsistent expression , This result is equivalent to the one obtained earlier in Refs.[24,25,20] using alternative approaches.The same spectral function can be found if one takes the self-energy in Hartree approximation.

Second approximation
In the first approximation above, the self-consistent Hartree treatment, fluctuations of the particle number n and the vibron coordinate x are completely neglected.In order to go one step further and estimate possible corrections, we start from the generated equations for the second-order Green functions, xd, d † and pd, d † , The second approximation that we consider here is based on the factorization The corresponding set of equations reads from which we obtain the second approximation for the Green function, where we introduced ∆ = ǫ + iη − ǫ0 , Ω = and G r H is given by the Green function obtained in the Hartree approximation, Eq.( 19).After inserting the expression (23) for the level population x , Eq. ( 29) reads We then calculate the spectral function (6) and the average number of electrons, (9).The self-consistent calculation is performed following the chain G r (ǫ) → A(ǫ) → n → G r (ǫ).Before entering into the discussion of the calculated quantities, we consider the structure of the Green functions obtained in the two different approximations: • In the limit where Ω → 0 the second approximation reduces to the first one.
Although the EOM method is not a systematic expansion, it tells us that the second approximation consistently extends the first one and reproduces it in a limiting case.• The second term on the right site of Eq. ( 30) represents an additional shift with respect to Eq. ( 19).In the case of very small frequencies ω 0 , Eq. (30 , involving the first term of a series expansion in ω 0 .

Results, discussion and conclusions
Starting from the expressions derived above within the first (self-consistent Hartree) and second approximation we have performed numerical simulations for the average population and current in the molecular junction as a function of gate and bias voltage for the different parameters λ, Γ and T entering the model Hamiltonian.Furthermore we compare the two approximations.The underlying nonlinear equations give rise to bistability in the level population thereby enabling memory effects and affecting the current.Below we analyze in detail the parameter ranges and conditions for memory (Color) Level population n in the self-consistent Hartree approximation (left) and the second approximation (right) as a function of normalized bias voltage for different electron-vibron interaction strength λ = 3ω 0 (black), λ = 4ω 0 (red), λ = 5ω 0 (green), the other parameters are Γ = 5ω 0 , effects to occur.We focus here on the case of an asymmetric junction, η = 0, while keeping the coupling to the leads symmetric, Γ L = Γ R = Γ, to reduce the parameter space.
First we investigate the gate-voltage dependence of the level population.In Fig. 2 it is clearly seen that bistability takes place only at larger values of the electron-vibron coupling λ.The critical values at which bistability occurs and disappears depend on the coupling Γ to the leads and on temperature T ; we discuss this parameter dependence below.Note that at large values of λ the level population of the stable states is close to 0 and 1; thus these two memory states are well distinguishable in charge.
A comparison between the left and right panel of Fig. 2, representing the two different levels of approximations in the EOM method, shows that the self-consistent Hartree treatment underestimates the parameter range where bistability occurs: the critical value of λ for the occurence of bistable behavior is close to 3ω 0 (see left panel).At this value a bistable regime has already developed for the second approximation.This can be partially understood taking into account the additional term appearing in the Green function for the second approximation, Eq. ( 30).This term increases the polaron shift thereby enhancing the bistable behaviour in the second approximation.Fig. 3 shows the bias-voltage dependence of the level population.In this case a qualitative difference arises between electrostatically symmetric (η = 0.5) and asymmetric (η = 0) junctions.For asymmetric junctions both states are stable at zero voltage, and both charge states are easily accessible.As we showed recently [23], asymmetric junctions are thus favourable, since they exhibit memory effects and hysteretic behavior at zero bias, enabling controlled switching upon ramping the bias voltage.For symmetric junctions hysteresis is expected only at finite bias voltage (nonequilibrium bistability [20]), and hence only a single stable state exists at zero bias.Furthermore, at finite voltage the level is only partially occupied by tunneling electrons.These two features render the symmetric system less suited for a memory setup.for different electron-vibron interaction strength λ = 3ω 0 (black), λ = 4ω 0 (red), λ = 5ω 0 (green), for the same parameters as in Fig. 3.
We further display in Fig. 4 the current-voltage characteristics, which reflects the switching behavior of the system.The characteristic feature is a current jump at the bias voltage value where recharging sets in.This behavior can be used to test the state of the system and as readout.
Finally, we depict in Fig. 5 "phase diagrams" showing the boundaries between the parameter regions where bistable memory states (below the boundaries in the two panels of Fig. 5) and single-valued states exist.The left panel of Fig. 5 shows the Γ − T parameter plane.The curve separating single-valued and bistable states can be roughly approximated by the condition Γ + T = c(λ), where c(λ) is extracted to be c(λ) ≈ λ 1.7 .This means that either thermal or quantum tunneling broadening suppresses the hysteresis.The disappearance of the memory effect at high temperature is due to enhanced electron tunneling into the higher-energy state.We note that in the analysis of the temperature dependence the effect of vibronic coordinate fluctuations is neglected that can be relevant close to the threshold between single-valued and bistable regime.The right panel in Fig. 5, displaying the λ 2 − Γ parameter plane, shows on the one hand that bistable behaviour requires, for growing Γ, increasing electron-vibron coupling λ.The boundary between the two regimes is approximately a straight line in the λ 2 − Γ plane.Hence the condition for finding the memory effect is given by Γ ≤ 0.63λ 2 /ω 0 .This clearly shows that for the appearence of the memory effect at low temperature, the two energy scales to be compared are the level broadening Γ and the polaron shift λ 2 /ω 0 .A further important conclusion is that the memory effect is suppressed for large coupling to the leads.Since, on the other hand, larger coupling favours fast information writing and reading and also can additionally suppress effects from quantum tunneling between states, the problem arises to find optimal parameters for utilizing the memory effect.This will be the subject of future work.
To conclude, we considered a charge-memory and switching phenomena within a single-level polaron model of a molecular junction in the framework of the equation-of-motion approach to the nonequilibrium Green function theory at different levels of approximation.Electrostatically symmetric and asymmetric junctions show qualitatively different bistability behavior.In the latter case, controlled switching of the molecule is achieved by applying finite voltage pulses.We showed that bistability takes place for sufficiently large elelctron-vibron coupling for a wide range of further parameters such as molecule-lead coupling and temperature and performed a systematic analysis of this parameter dependence by computing phase diagrams for the memory effect.We focussed on the investigation of bistability for stationary states, leaving the problem of time-dependent fluctuations and fluctuation-induced switching as an open question to be investigated in future research.