4 0 0 2 Quantum theory of shuttling instability in a n movable quantum dot array a J 0 Andrea Donarini†, Tom´aˇs Novotny´†,‡ and Antti–Pekka 2 Jauho† ] †MikroelektronikCentret, Technical UniversityofDenmark,ØrstedsPlads, ll Bldg.345east,DK-2800Kgs.Lyngby,Denmark a ‡DepartmentofElectronicStructures,FacultyofMathematics andPhysics, h CharlesUniversity,KeKarlovu5,12116Prague,CzechRepublic - s E-mail: [email protected],[email protected] e m Abstract. We study the shuttling instability in an array of three quantum t. dots the central one of which is movable. We extend the results by Armour a and MacKinnon on this problem to a broader parameter regime. The results m obtained by an efficient numerical method are interpreted directly using the Wigner distributions. We emphasize that the instability should be viewed as - d acrossoverphenomenonratherthanaclear-cuttransition. n o c [ 1. Introduction 1 v Since the first proposal by Gorelik et al. [1] of the shuttling instability in a 7 genericnanoelectromechanicalsystem(NEMS)consistinginamovablesingleelectron 5 transistor the phenomenon of the shuttling transport attracted much attention. 3 However, until recently the fully quantum theory of the phenomenon was not 1 developed. Thefirstquantummechanicalstudyonamodifiedsetupwasaccomplished 0 4 by Armour and MacKinnon [2], closely followed by the work of the present authors 0 onthe originalGorelik’ssetup[3]. Inthis paperweextendthe resultsbyArmourand / MacKinnon. The phase space analysis in terms of the Wigner functions introduced t a in [3] reveals directly the nature of the transport (incoherent tunnelling versus m shuttling) in different regions in contrast to the indirect evidence used by Armour - and MacKinnon. Also the new numerical scheme that we use enables us to access a d n wider range of parameters fast and reliably. o c 2. Model and method of solution : v i Let us consider a simple NEMS consisting of an array of three quantum dots X (device) connected to two leads [2]. The central dot is assumed to be movable in r a a parabolic potential. The Hamiltonian of the device consists of the mechanical and the electronic parts H = H +H where H = h¯ωaˆ†aˆ is the Hamiltonian of mech el mech the harmonic oscillator and H = |αiǫ (xˆ)hβ|. We assume strong el Pα,β∈{0,L,C,R} αβ Coulomb blockade regime with no double occupancies so that the vectors |αi with 2 α = 0,L,C,R span the entire electronic Hilbert space of the device. Each matrix element ǫ (xˆ) is still a full matrix in the oscillator space. ǫ(xˆ) reads explicitly αβ 0 0 0 0  0 ∆V t (xˆ) 0  ǫ(xˆ)= 2 L (1)  0 tL(xˆ) −∆2xV0xˆ tR(xˆ)   0 0 t (xˆ) −∆V  R 2 wherexˆ= h¯ (aˆ+aˆ†)isthepositionoperator,∆V,thedevicebias,isthedifference q2mω betweentheenergyoftheleftandtherightdot. x ishalfthedistancebetweenthetwo 0 outer dots and represents the maximum amplitude of the inner dot oscillation. The three dots areelectrically connectedonly via a tunnelling mechanism. The tunnelling length is given by 1/α and the tunnelling strengths depend on the position operator xˆ of the inner grain as t (xˆ)=V e−α(x0+xˆ), t (xˆ)=V e−α(x0−xˆ). L 0 R 0 The dynamics of the device is described by the generalized master equation (GME) formalism [2]. The density matrix of the device evolves according to ρ˙ = −i[H,ρ]+Ξρ+ρ˙ . The first term represents the coherent evolution of the isolated d device. The coupling to the leads responsible for the electronic transfer from/to the device is introduced in the wide band approximation following Gurvitz et al. [4] and yields the second term in the equation (in the block notation of [2]) ρ −ρ 0 0 0 RR 00   0 ρ 0 −ρ /2 Ξρ=Γ 00 LR (2)  0 0 0 −ρCR/2     0 −ρRL/2 −ρRC/2 −ρRR  withΓ being the injectionrateto the leads. The thirdtermdescribesthe effectofthe environmentontheoscillator,consistinginmechanicaldampingandrandomquantum and thermal excitation (Langevin force). It reads [2] γ γ ρ˙ =− n¯(aa†ρ−2a†ρa+ρaa†)− (n¯+1)(a†aρ−2aρa†+ρa†a) (3) d 2 2 where γ is the damping rate and n¯ is the mean occupation number of the oscillator at temperature T. The stationary version of the above GME was solved numerically after the oscillator Hilbert space was truncated at sufficiently large N so that all dynamically excited states were contained within the basis. Utilizing the decoupling properties of the GME in the block notation resulted in the problem of finding a unique null space ofa(super)matrix withthe lineardimension10N2 (with uptoN =40). The Arnoldi iteration [5] with preconditioning was found to be superior to direct methods, such as singular value decomposition or inverse iteration. The preconditioning involving the inversionof the Sylvester part of the problem (a fast procedure) is crucial for the convergence of the method which is very fast and low memory consuming (the whole (super)matrix does not have to be stored in the memory). 3. Results and discussion Theresultingstationarydensitymatrixwasusedtoevaluatethemeanvalueofcurrent I = ΓTr ρ [2] and the phase space distribution of the charged central dot using osc RR the transformation into Wigner coordinates [3]. These quantities as a function of the device bias ∆V and for three different values of the injection rate Γ are plotted 3 0.05 Γ = 0.05 0.045 Γ = 0.1 Γ = 0.2 0.04 0.035 0.03 0.025 I 0.02 0.015 0.01 0.005 0 0 0.5 1 1.5 2 2.5 3 3.5 ∆ V Figure 1. Plotofthecurrentthroughthetripledotsystemasafunctionofthe device bias for different injection rates Γ. The other parameters are V0 = 0.5, α=0.2,x0=5,γ=0.0125,ω=1. in the figures 1 and 2, respectively. It was found in [2] that the triple dot system exhibits differentregimesoftransportatdifferentdevice biases. The currentpeaksat ∆V ≈nω(seefigure1)wereidentifiedaseffectsofelectromechanicalresonanceswithin the device. Yet, the different peaks may correspond to different physical mechanisms — while the peak around ∆V ≈ω is mainly due to the incoherent oscillator-assisted tunnelling the peak at ∆V ≈2ω reveals a clear shuttling component. This finding by ArmourandMacKinnonbasedonindirectevidenceofparametricdependenciesofthe currentcurves (e.g.the dependence of the currentcurve onthe tunneling length 1/α) isconfirmedbythedirectinspectionofthephasespacedistributions(seethefirstrow offigure 2). The half-moon-likeshape characteristicforshuttling transport[3]is only present around ∆V ≈ 2ω while all other plots show just the fuzzy spot indicative of incoherenttunnelling. However,ourdirectcriterionfordetectingtheshuttling regime reveals a close similarity between the resonances. For increasing injection rate we can see that the shuttling regime gradually sets in also in the vicinity of the first resonance peak. This reveals the crossover character of the onset of the shuttling instability found also previously [3]. The sharp transition into the shuttling regime reportedinsemiclassicalstudiesissmearedintothecrossoverduetothe noisepresent in the system and properly accounted for in our approach. 4 Figure2. Wignerdistributionsforthecentraldotinthechargedstate. Different devicebiasesarethepointsofminimumormaximumcurrent(seefigure1). The authors thank T. Eirola for indispensable advice concerning the numerical methods. Support of the grant 202/01/D099of the Czech grant agency for one of us (T.N.) is also gratefully acknowledged. [1] L.Y.Gorelik,A.Isacsson,M.V.Voinova,B.Kasemo,R.I.Shekhter,andM.Jonson 1998Phys. Rev. Lett.804526 [2] A.D.ArmourandA.MacKinnon. 2002 Phys. Rev. B66035333 [3] T.Novotny´, A.Donarini,andA.-P.Jauho. 2003 Phys. Rev. Lett.90256801 [4] S.A.GurvitzandYa.S.Prager. 1996Phys. Rev. B5315932 [5] G. H. Golub and C. F. Van Loan 1996 Matrix Computations (The John Hopkins University Press)