Table Of Content

Dielectricbreakdown inspinpolarizedMottinsulator Zala Lenarcˇicˇ1 and Peter Prelovsˇek1,2 1J. Stefan Institute, SI-1000 Ljubljana, Slovenia and 2Facultyof Mathematics and Physics, University of Ljubljana, SI-1000 Ljubljana, Slovenia NonlinearresponseofaMottinsulatortoexternalelectricfield,correspondingtodielectricbreakdownphe- nomenon,isstudiedwithinofaone-dimensionalhalf-filledHubbardmodel.Itisshownthatinthelimitofnearly spinpolarizedinsulatorthedecayrateofthegroundstateintoexcitedholon-doublonpairscanbeevaluatednu- mericallyaswelltohighaccuracyanalytically.Resultsshowthatthethresholdfielddependsonthechargegap asF ∝∆3/2. Numericalresultsonsmallsystemsindicateonthepersistenceofasimilarmechanismforthe th 2 breakdownfordecreasingmagnetizationdowntounpolarisedsystem. 1 0 PACSnumbers:71.27.+a,71.30.+h,77.22.Jp 2 n Thenonlinearresponsetoexternalfieldsandmoregeneral Inthefollowingwestudytheprototype1DHubbardmodel, a J nonequilibriumpropertiesofstronglycorrelatedelectronsand 8 Mottinsulatorsinparticular[1]aregettingmoreattentionin H = t (eiφc† c +H.c)+U n n , (1) 1 recent years, also in connectionwith powerfulnovelexperi- − i+1,σ iσ i↑ i↓ iσ i mentaltechniques,e.g. thepump-probeexperimentsonMott X X ] l insulators [2], as well as novel systems, the prominent ex- withperiodicboundaryconditions(p.b.c.) wherec† ,c are e iσ iσ ample being the drivenultracoldatoms within the insulating - creation(annihilation)operatorsforelectronsatsiteiandspin r phase [3]. In this connection, one of the basic phenomena t σ = , . TheactionofanexternalelectricfieldF isinduced s to be understood is the dielectric breakdown in Mott insu- ↑ ↓ viathe Peierlsphaseφ (vectorpotential)andits timedepen- . at lators, studied experimentallyin effectivelyone-dimensional dence, i.e. φ˙(τ) = e0F(τ)a0/~. Furtheron we use units m (1D) systems more than a decade ago [4]. The concept of ~ = e0 = a0 = 1, as well as we put t = 1 defining the Landau-Zener(LZ)single-electrontunneling[5,6]asastan- unitofenergy. Insuchamodelweinvestigatefinitesystems - d dardapproachtodielectricbreakdownofbandinsulators[7] of length L and at half-filling N +N = L but in general u d n isnotstraightforwardtogeneralizetocorrelatedelectrons[8– at finite total spin, Sz = (N N )/2 and magnetization o 10]. Theoreticaleffortshavebeensofarrestrictedtothepro- m=Sz/L. u − d c totype Hubbard model at half-filling. In 1D numerical ap- [ Letusfirstconsidertheproblemofasingleoverturnedspin, proaches have given some support to analytical approxima- i.e. ∆Sz =L/2 Sz =1. Here,basiswavefunctions ϕ 1 tionsforthemostinterestingquantitybeingthethresholdfield − | jmi v correspondtoanemptysite(holon)atsitej andadoublyoc- F anditsdependenceonthechargegap∆[9],typicallyre- 7 th cupiedsite(doublon)atsitem. Takingintoaccountthetrans- vealing a LZ type dependence F ∆2. Differentdepen- 1 th ∝ lationalsymmetryofthemodel(1)withp.b.c.(evenwithtime 8 denceisfoundnumericallywithinthedynamical-mean-field- dependent φ(τ)) at given (total) momentum q = 2πm /L 3 theoryapproach[11]asrelevantforhighdimensionsD 1. q 1. ≫ the relevant basis is |Ψlqi = (1/√L) jeiqj|ϕj,j+li,l ∈ [0,L 1]. At fixed φ adiabatic eigenfunctions can be then 0 search−edintheform ψ = d Ψj Pleadingtothe eigen- 2 InthisLetterweapproachtheproblemofadielectricbreak- | i j j| qi 1 valueequation, downfrom a partially spin polarizedMott insulator. We use P : v thefactthatthegroundstate(g.s.) ofthe1DHubbardmodel 1 1 1 Xi isinsulatingatanyspinpolarizationwiththechargegapmod- − U = L E U +2(cos(q′ φ)+cos(q′ φ q)). estlydependentonthemagnetizationm. Inparticular,asin- q′ − − − − r X (2) a gle spin excitation in fully polarized system m 1/2, i.e. ∆S = 1 state, canbe studiedexactlynumerically∼aswellas In the limit L the g.s. energy E0 representing the → ∞ holon-doublon(HD) boundstate can be expressed explicitly tohighaccuracyanalytically.Therelevantmechanismforthe decay of the g.s. under constant external field F is the cre- as E0 = U (U2 + 16cos2(q/2))1/2. We note that (in − spite of the q-dependence) g.s. states for all q are noncon- ationof holon-doublon(HD) pairs. We show thatdueto the ducting since from Eq. (2 ) it follows that the charge stiff- dispersion-lessg.s. thesimilaritytotheLZtunnelingisonly partialandleadstoadifferentscalingF ∆3/2. Furtheron ness 0 ∂2E0/∂φ2 0 for L . On the other th D ∝ → → ∞ ∝ hand, excited states form a continuum with lower edge at we study numerically on small systems also the model with ∆S > 1, m < 1/2 in a finite field F. Results indicate that E1 =U 4cos(q/2). − thedecaymechanismremainsqualitativelyandevenquantita- Since φ(τ) conserves total q we furtheron consider only tivelysimilaratpolarizationsm<1/2,inparticularforlarger solutions within the q = 0 subspace representing the ab- ∆wherebythemostinterestingcaseisclearlytheunpolarized solute g.s. wavefunction |0i with d0j = Ae−κ|j|eiφj and m=0system. A = √tanhκ. Here, the charge gap ∆ = E1 E0 and − 2 therelatedg.s. localizationparameterκaregivenby Here,wecanalreadyrealizesomeessentialdifferencestothe usual concept of of LZ tunneling, i.e., Φ and Eq. (5) do ∆= 4+ U2+16=4(coshκ 1). (3) k − − not favor transitions to lowest lying excited state but rather When we consider theptime-dependentφ(τ) we have to deal to k κ/√3, hencethe reductionto a two-levelproblemis ∼ at finite L with adiabatic states E (φ) as, e.g., shown in notappropriate. n Fig. 1 for finite L. At finite L 1/κ E0 is essentially φ- Therateofak(τ)followingfromEqs.(5),(7)isnotsteady. ≫ independent but the same holds as well for lowest excited SinceweareinterestedinlowF weaverageitovertheBloch states E ,n & 1 which makes an usual application of two- periodτ =2π/F togeta¯=a (τ )whichisapproximately n B k B level LZ approach not straightforward to apply. If we sup- thesameformajorityofk(fixingherek =π), pose that due to field F > 0 the transition probability be- tweenneighbouringstatesishigh(neglectingfinitesizegaps AU π 1 ′ i ξ a¯= dξ exp dξ′ω (ξ′) (8) π between them) excited states are well represented by ’free’ −√L Z−π (cid:18)ωπ(ξ)(cid:19) F Z−π ! HDpairstateswith iAU π i ξ 1 2π dξexp dξ′ω (ξ′) , (9) dkj = √Leikj, ǫk =U−4cos(φ−k), k = L mk. (4) ∼ F√LZ−π F Z−π π ! As shown further relevant transitions due to time-dependent afterper partesintegrationofEq. (8) andneglectingthe first φ(τ) happento effectivestates k with m 1 since the fastoscillatingpart,smalleralsoduetoanadditionalprefactor k | i | | ≫ g.s. 0 iswelllocalized. F. FinalsimplificationforsmallFcanbemadebyreplacing | i coshκ cosξ ξ2/2+ ∆/4 and consequently extending 8 − ∼ integrationsinEq.(9)toξ = . Thisleadstoananalytical expressionforthedecayrateΓ±,∞definedby a0 2 exp( Γτ) 6 whereΓ=La¯2/τ , | | ∼ − B | | 4 ∆3/2B(∆) √2∆3/2 B(∆) (2∆)3/2 En Γ= K21 exp 2 3πF 3 3F !∼ √8 (cid:18)− 3F (cid:19) (10) 0 whereK1/3(x)isthemodifiedBesselfunctionandB(∆) = ∆(∆+8)3/2/(∆+4), and the last exponentialapproxima- -2 tion is valid for small enough Γ. The main conclusion of 0.0 0.2 0.4 0.6 0.8 1.0 the analysis is that Γ in Eq. (10) depends on ∆3/2/F un- Φ(cid:144)2Π like usual LZ theory applications [8, 9] yielding ∆2/F. As thethresholdfieldisusuallydefinedwiththeexpressionΓ Figure 1. Energy levels En (in units of t) vs. phase φ for holon- exp( πF /F),Eq.(10)directlyleadsF =(2∆)3/2/(3π∝). doublon pair states in the system with L = 21 sites and U = 4. − th th It is straightforward to verify the validity of approxima- Thicklinerepresentstheg.s. andtheeffectiveHDpairstatedisper- sion. tionsforNd = 1viaadirectnumericalsolutionofthetime- dependentSchro¨dingerequation(TDSE)withφ=Fτ within Letusnowconsiderthedecayoftheg.s. 0 afterswitch- thefullbasisatq = 0andfinitebutlargeL >100. Timede- | i ing constant field F(τ > 0) = F,φ = Fτ. We present an pendenceoftheg.s. weight a0(τ)2 ispresentedinFig.2for | | analysisfortheinitialdecaywheremostweightisstillwithin typicalcaseU =4anddifferentfieldsF =0.2 0.5.Results − the g.s., i.e. a0(τ) an6=0(τ). In such case the excited for the case of an instantaneous switching F(τ > 0) = F | | ≫ | | stateamplitudetime-dependencean(τ)isgivenby (shown for F = 0.5) reveal some oscillations (with the fre- τ τ′ quencyproportionalto the gap ∆) butotherwise clear expo- a (τ)= F dτ′Φ (τ′)exp(i ω (τ′′)dτ′′), (5) nential decay with well defined Γ. In order to minimize the n n n − Z0 Z0 fast-switching effect we use in Fig. 2 and furtheron mostly whereΦn = n∂/∂φ0 andωn(τ)=En(φ) E0. smooth transient [11], i.e., field increases as F(τ < 0) = h | | i − Analytically progress can be made by using effective HD F exp(3τ/τ )toitsfinalvalueF(τ >0)=F. B states k as approximate excited states with ω (ξ) = ǫ InFig.3wecompareresultsforΓasobtainedviathreedif- k k | i − E0 =4(coshκ cosξ),ξ =Fτ k. Byusingtherelation ferentmethods:a)directnumericalsolutionofTDSE,b)ana- − − lyticalapproximationwithanaveragedecayrateintofreeHD hk|0iωk =hk|H0+U −H|0i=Uhk|n0↓|0i=UA/√L, states,numericallyintegratingEq.(8),andc)theexplicitex- (6) pression(10)whereadditionalsimplificationoftheparabolic whereH0denotesonlykineticterminEq.(1),onecanexpress dispersion of excited states is used. The agreementbetween Φ inEq.(5)as k differentmethodsissatisfactoryessentiallywithinthewhole ∂ ∂ UA∂ω−1 regimeofsmallΓ anddeviationsbetweenanalyticalandnu- Φ = k 0 = k 0 = k . (7) k h |∂φ| i ∂φh | i √L ∂φ mericalresultsbecomevisibleonlyforlargeΓ 0.1. More- ∼ 3 0.4 F=0.5 ofmagnetization1/2> m≥0(relevant1≤Nd ≤L/2)for two casesofU = 4,10,respectively,andthespanofappro- priatefieldsF. Examplesarechosensuchtorepresentcharge 0.3 2 gap (for a single HD pair) being small ∆ 1.3 < W and 0 ∼ a large∆ 6.5> W,respectively,relativetothenoninteract- F=0.4 n 0.2 ingband∼widthW =4. l - 0.1 F=0.3 m»1(cid:144)2 F=0.8 0.6 m=3(cid:144)8 F=0.2 0.0 Nd m=1(cid:144)4 -0.2 0.0 0.2 0.4 0.6 (cid:144) m=0 2 Τ(cid:144)ΤB a0 0.4 F=0.6 Figure 2. (Color online) a) Ground state weight ln|a0|2 vs. time -ln 0.2 τ/τB forU = 4anddifferentfieldsF = 0.2−0.5. ForF = 0.5 F=0.3 thecomparisonofresultsforsmoothlyandinstantaneouslyswitched F(τ)ispresentedwhileforF <0.5onlysmoothswitchingisused. 0. 0. 0.2 0.4 Τ(cid:144)Τ B over, results confirm the expected variation lnΓ 1/F es- ∝ sentiallyinthewholeinvestigatedrangeofF. Figure4. (Coloronline)Normalizedg.s. weight(1/Nd)ln|a0|2vs. timeτ/τ forU = 4andfieldsF = 0.3,0.6,0.8,forvariousspin B states1≤N ≤L/2. d 0.1 0.04 m»1(cid:144)2 F=2.6 G m=3(cid:144)8 0.01 m=1(cid:144)4 d TDSE N m=0 (cid:144) HDnumerical 2 HDanalytical a0 0.02 F=2.2 0.001 n 2. 3. 4. -l F=1.6 1(cid:144)F 0 0. 0.2 0.4 Figure 3. (Color online) Ground state decay rate Γ (log scale) vs. Τ(cid:144)Τ 1/F forU = 4asevaluatedbydirectnumericalsolutionofTDSE B (fullline),decayintofreeHDstates,numericallyintegratingEq.(8) (dottedline),andanalyticalexpression,Eq.(10)(dashedline). Figure 5. (Color online) The same as in Fig. 4 for U = 10 and F =1.6,2.2,2.6. Onecanassumethatasimilarmechanismofthedielectric breakdownvia the decayinto freeHD pairsremainsvalidat The main conclusion following from Figs. 4,5 is that the finite deviationsNd > 1 andm < 1/2. Inorderto test this g.s.weight a0 2indeeddecaysproportionaltoNdconfirming | | scenariowe performthe numericalsolutionofTDSE forthe the basic mechanismof the field-inducedcreationof (nearly model,Eq.(1), withthefinitefieldF(τ). Calculationforall independent)HDpairs. ThedecayrateΓ definedas a0 2 | | ∝ Sz sectorscoveringthewholeregime0 m < 1/2areper- exp( ΓN τ) is only moderately dependent on N and m. d d ≤ − formedonfiniteHubbardchainswithuptoL=16sitesusing ResultsconfirmthatΓisessentiallyindependentofN inwell d theLanczosprocedurebothforthedeterminationoftheinitial polarized systems with m 1/4, which is compatible with ≥ g.s. wavefunction 0 aswellasforthetimeintegrationofthe independent decay into low concentration of HD pairs. For | i TDSE [14] within the full basis for given quantum numbers larger U = 10 in Fig. 5 the invariance of Γ extendseven to N ,N ,q =0reachinguptoN 107basisstates. Weuse unpolarizedsituationm = 0(N /L = 1/2)forintermediate d u st d ∼ everywheresmoothtransientforthefieldF(τ). Sincethede- fieldsF 2.2. cayrateoftheg.s. weight a0 2 isexpectedtoscalewiththe There≥are some visible deviations at m 1/4 for weak- | | ≤ numberofoverturnedspinsN therelevantquantitytofollow est fields both in Fig. 5 for F = 1.6 and even more for d andcompareis(1/Nd)ln a0 2(τ). smaller U = 4 and F = 0.3 in Fig. 4, indicating on larger | | In Figs. 4,5 we present numerical results for time depen- Γ and correspondinglyfaster decay of unpolarizedg.s. with denceofnormalizedg.s.weightln a0 2/Ndasobtainedviaa m = 0 relative to nearly saturated m 1/2. Part of | | ∼ directsolutionoftheTDSEforL = 16withthewholerange this enhancement of Γ can be attributed to the dependence 4 inparticularofanunpolarizedg.s. [10,11]. 8 HD LZ ìæ Thecase of a nearlypolarizedstate N = 1describesthe d m»1(cid:144)2 ìæ à mechanism of the field-induced decay of the g.s. into sin- 6 m=1(cid:144)4 ìæ à gle HD pair. Here one can follow differences to usual LZ- Fth4 m=0 ìæ à type approaches: a) the g.s. is localized and dispersionless æ ìæ à withintheinsulator,b)the transitionisnotbetweentwo iso- ì à 02 ææææàæìææàìæ àìæ àìæ à lftairootenmdstlEoeqvlose.wls(e5bs)tu,(et8xr)caittthheeadrtstmtaotaeatsr,icxco)neiltneinsmtueeuanmdts,ofdmeooxranecoottveefxarcviiottersfotsrltalaontwessis-, 0 2 4 6 8 onecanwelluseeffectivefreeHDstates,d)dispersionofef- D fectiveHDstatesisunlikeinLZapplicationsnothyperbolic, e.g., ω (k2 +κ2)1/2 butratherparabolicω = k2 +κ2 k k ∝ Figure6.(Coloronline)ThresholdfieldF vs.chargegap∆fordif- which is presumablythe main origin for qualitatively differ- th ferentmagnetizationsm∼1/2(givenbyNd =1)andm=1/4,0 entbehaviorofthethresholdfieldFth ∆3/2whichisafinal asobtainednumericallyforL = 16. Fullcurve(HD)representthe manifestationofthedistinctiontousu∝alLZapplications. On analyticalapproximation,Eq.(10),whilethedashedcurveistheLZ the other hand there are some similarities. In particular the approachresultfromRef.[9]. analytical expression for the average transition rate, Eq. (9), wherematrixelementisintegratedout,appearsanalogousto two-levelproblemandreadyforphase-integraltransformation of the charge gap on the magnetization ∆(m). The ther- intoimaginaryplaneasusedoriginallybyLandau[5]andthen modynamic (L ) value ∆0 = ∆(m = 0) is known generalized[17, 18] and appliedaswell to breakdownprob- → ∞ via the Bethe Ansatz solution given by the equation ∆0 = lem [10, 13]. Still it is straightforwardto verify that for the (16/U) ∞dx√x2 1/sinh(2πx/U) [15, 16]. Values for 1 − levelsunderconsiderationωkdonotsatisfycriteriaforitsap- ∆(m 1/2) as given by Eq. (3) are somewhat larger than plication,buttheanalogyratheremergesthroughtheapplica- ∼R ∆0 with the relative difference becoming more pronounced tionofthesteepestdescentapproximationtoEqs.(9). for U < 4. Still taking into accountactual ∆(m) some en- The picture of the decay of the driven Mott insulator into hancement seems to remain at m 0 at least for weaker ∼ HD pairs remains attractive for magnetization approaching fields F and smaller U. This could indicate that the decay the unpolarized g.s. There seem to be two characteristic into HD pairsare notindependentprocessesbutcorrelations lengthscalescontrollingthemechanism,theHDpairlocaliza- duetofiniteconcentrationofN /Lenhancedecay. d tionlengthζ = 1/κandtheStark (Bloch)localizationscale FinallyletusconsiderthethresholdfieldforthedecayF th L =8/F.Ourresultsindicatethatforlarger∆(smallζ)and asdefinedagainbyΓ exp( πF /F). Wepresentresults S ∝ − th well localized HD pairs the mechanism of decay into nearly in Fig. 6 for F as function of the gap ∆. To extract F th th independentHDpairsremainsatleastqualitativelyvalid. On vs. ∆weusenumericaldataforΓ(F)obtainedfromnumer- the other hand, we find indications that for smaller ∆ and cichaalrg|ae0g|2a(pτ∆) a(sm,)e.wg.e,suhsoewfonrimnFigs1./22,4a,n5d.mFo=rth1e/4reEfeqr.e(n3c)e, weaker F (larger LS), the decay is enhanced, i.e., pointing ∼ intothedirectionofmorecollectivedrivenexcitationsfavored while form = 0we use exact∆0. Somedeviationbetween also in the interpretationof experiments[4]. It should be as m 1/2 and m = 1/4 results can be still attributed to ac- ∼ wellpointedoutthatthephenomenonofHDpairgeneration tually slightly smaller gap for the latter magnetization. For is not particularly specific to 1D systems discussed here but comparisonwe plot also the analyticalresult emergingfrom cangeneralisedtohigherdimensionalMottinsulatorsaswell. Eq. (10), F ∆3/2, as well as the dependence following th from the LZ a∝pproach [9] with F = ∆2/8. From Fig. 6 This work has been supported by the Program P1-0044 th and the project J1-4244 of the Slovenian Research Agency we concludethatthegeneraltrendF (∆) isquite wellrep- th (ARRS). resented by the single HD pair result which deviates signif- icantly from the LZ dependence at least for larger ∆ > 6. At the same time, we should note that our numerical results in the range1 < ∆ < 2.1agree also well with data analyz- ingnumericallytheg.s. decayusingthet-DMRGmethod(at [1] M. Imada, A. Fujimori, and Y. Tokura, Rev. Mod. Phys. 70, m=0)forthesamemodelbutbiggerL 50[9]. 1039(1998). ∼ In conclusion, we have presented an analysis of the di- [2] H.Okamoto,H.Matsuzaki,T.Wakabayashi,Y.Takahashi,and electric breakdown within the Mott-Hubbard insulator start- T.Hasegawa,Phys.Rev.Lett.98,037401(2007). [3] N.Strohmaieretal.,Phys.Rev.Lett.104,080401(2010). ingfromaspinpolarizedgroundstate. Suchanapproachhas [4] Y. Taguchi, T. Matsumoto, and Y. Tokura, Phys. Rev. B 62, clearlyanadvantagesincetheproblemcanbesolveduptode- 7015(2000). sired accuracynumericallybut as well captured analytically. [5] L.Landau,Sov.Phys.1,89(1932). Assuchthesituationcanserveatleastaswellcontrolledtest [6] C.Zener,Proc.R.Soc.A137,696(1932). formoredemandingsituationsofanarbitrarymagnetization, [7] C.Zener,Proc.R.Soc.A143,523(1934). 5 [8] T. Oka, R. Arita, and H. Aoki, Phys. Rev. Lett. 91, 066406 [12] Q.NiuandM.G.Raizen,Phys.Rev.Lett.80,3491(1998). (2003). [13] T.Oka,arXiv11056.3143. [9] T.OkaandH.Aoki,Phys.Rev.Lett.95,137601(2005). [14] forareviewseee.g.P.PrelovsˇekandJ.Boncˇa,arXiv1111.5931. [10] T.OkaandH.Aoki,Phys.Rev.B81,033103(2010). [15] E.H.LiebandF.Y.Wu,Phys.Rev.Lett.21,192(1968). [11] M. Eckstein, T. Oka, and P. Werner, Phys. Rev. Lett. 105, [16] A.A.Ovchinnikov,Zh.Eksp.Teor.Fiz57,2137(1969). 146404(2010). [17] A.M.Dykhne,Sov.Phys.JETP14,941(1962). [18] J.P.DavisandP.Pechukas,J.Chem.Phys.64,3129(1976).