The Spectral Coupled Cluster Method with $k$ dependency for strongly correlated lattices

0.92 MB·

by  Alessandro Mirone

TheSpectralCoupledClusterMethodwithk dependencyforstronglycorrelatedlattices.
Alessandro Mirone
European Synchrotron Radiation Facility, BP 220, F-38043 Grenoble Cedex, France
(Dated:February15,2017)
WeadapttheCoupledClusterMethodtosolidstatestronglycorrelatedlatticeHamiltoniansextendingthe
CoupledClusterlinearresponsemethodtothecalculationofelectronicspectraandobtainingthespace-time
Fourier transforms of generic Green’s functions. We apply our method to the MnO2 plane with orbital and
magneticordering,tointerpretelectronenergylossexperimentaldata,andtotheHubbardmodel,whereweget
insightintothepairingmechanism.
7
1
0 1. INTRODUCTION whereSgivestheidealexactsolutionandSN isasum, trun-
2 catedtoNterms,ofproductsofelectron-holepairexcitations:

b The recent developments of the Linear Response Coupled
e Cluster1havestretchedtheapplicabilityoftheCoupledClus- S =∑N t Symm ∏ni cˆ† cˆ† (2)
F ter Method2 (CCM), originally formulated as a ground state N i=1 i (k=1 αi,k ai,k)

5 approximation, to excited states. In particular Crawford and
1 Ruudhavecalculatedvibrationaleigenstatescontributionsto Eachterminthesumistheproductofasetofelectron(hole)-
Raman optical activity3 while Govind et al. have calculated creation operators cˆ†, and is determined by a choice of in-

] excitonic states in potassium bromide4. More recently the dexesαi,k(ai,k),withthegreek(latin)letterα(a)rangingover
l
e Spectral Coupled Cluster Method (SCCM) has been derived empty(occupied)orbitalsofthereferencestateΦ0 . Onecan
| i
- fromthelinearresponsetocalculateelectronicspectra. This considerthereferencestateasthevacuumstateandthateach
r
t isdonebyiterativelyrefiningaresolventoperatorwhichde- term i in the sum SN creates, from the vacuum, an excited
s
. scribestheresponseoftheCoupledClustersolutiontoprobe state which is populated by ni holes and ni electrons. The
at operatorsandthenbyobtainingthespectralfunctionsbythe SymmoperatormakestheansatzsymmetricfortheHamilto-
m diagrammaticalexpansionoftheirexpectationvalue5. niansymmetrysubgroupwhichtransforms,uptoafactor,the
reference state Φ into itself. The t’s are free coefficients
- InthisworkwefurtherexpandSCCMforsolidstatemodel | 0i i
d thatareobtainedfromtheSchrodinger’seigen-equation:
Hamiltonians, byadaptingittothecalculationofspace-time
n
Fourier transforms of generic Green’s functions. In the next
o H Ψ E Ψ =0 (3)
section(2)werecallfirsttheCCformalismandthenintroduce | i− | i
c
[ thenewformulationofSCCMwithkdependency. by setting the eigen-equation residue to zero in the space of

3 With the goal of simulating the energy loss spectra of the excitedstateswhichentertheSN sum:
orbitally ordered MnO plane in La Sr MnO we detail,
42v ilnatesetchteioMn n(3O)2, tphleanme,oda2enldHaalsmoilttwonoia0sn.i5mtph1lai.t5fiewdefuo4sremtsoosfimthuis- 0=hΦ0| k∏=ni1cˆai,kcˆαi,k!e−SNHeSN|Φ0i ∀i∈[1,N] (4)
5 Hamiltonian which describe the Hubbard model and a sim-
8 ple linear chain that we can solve exactly for a preliminary This equation enforces the orthogonality between
0 validationofourSCCMimplementation. Wediscussinsec- e−SNHeSN Φ0 and the basis of excited states. The
| i
. tion (4) the results for the simple linear chain, for the Hub- overlap with the reference state is non-zero and gives
1
0 bardmodeland,finally,fortheMnO2 plane. FortheHubbard anapproximationoftheenergy:
17 mthoed4el w4eHcuobmbpaardresyosutremca,lcgueltatitniognintosigthheteinxtaoctthreesroolluetioofnthoef E’EN =hΦ0|e−SNHeSN|Φ0i (5)
: collec×tivespinonicexcitationsinthegroundstateandapossi- (6)
v
i blepairingmechanism.FortheorbitallyorderedMnO2plane The Coupled Cluster Method expands these equations by
X wesimulate, bycalculatingthedynamicstructurefactor, the means of the Baker-Campell-Hausdorff expansion formula6
r lowenergyresonancesobservedinanelectronenergylossex-
whichfortwoarbitraryoperatorsAandBstatesthat:
a
periment.
e−ABeA=B+[B,A]+1/2[[B,A],A]+..1/n![[...[B,A]...],A]+...
(7)
The numerical applicability of CCM relies on the fact that,
2. METHOD
when A is substituted with S which contains creation oper-
N
ators only, and B is replaced by the Hamiltonian H with the
The Coupled Cluster method2 uses an exponential ansatz propertiesdescribedbelow,thenonlythefirstfivetermsinthe
topreserveextensiveproperties.Givenareferencestate Φ0 , seriescanbenonzero. Thiscanbedemonstratedconsidering
| i
thesolution Ψ isrepresented,inthisansatz,by: that each term of H is formed by a limited number of anni-
| i
hilationoperators. TheCoulombinteraction,infact,contains
Ψ =eS Φ eSN Φ (1) productsoffourone-particleoperators.
0 0
| i | i’ | i
2
TheaccuracyoftheCCMsolutionincreaseswithN.Toin- where c† denotes the operator which creates an electron on
i
creaseN,ateachiteration,theexcitationnotyetcontributing orbital i, not to be confused with cˆ† which is defined with
to SN and giving the largest scalar product in equation (4) is respect to the reference state taken as vacuum. In this work
added,forthefollowingstep,tothelistofexcitationsconcur- we will consider also, numerically, the partial contributions
ringtoSN+1. OncewehaveobtainedtheCCMgroundstate S↑↑,S↑↓,S↓↑,S↓↓toS(q,ω)obtainedbysplittingNq=Nq↑+Nq↓.
anditsgroundenergyE,weareinterestedintheGreen’sfunc- Another example is the spectral function A(q,ω), which is
tionsthatareobtainedapplyingallsortsofprobeoperatorson theimaginarypartoftheone-particleretardedGreen’sfunc-
thegroundstate. tion, related to angle resolved photo-emission spectroscopy
We are interested here in solid states model Hamiltonians (ARPES)andtoinversephoto-emissionspectroscopy(IPES):
appliedonsystemshavingperiodicboundaryconditions,and
being composed of a Bravais lattice of identical unit cells. 1
A(q,ω)= (G(C ,C ,ω,q,γ) G(C ,C ,ω,q, γ)) (13)
Eachunitcelliscomposedofoneormoreatoms. Thebasis 2iω q q − q q −
statesoftheFockspaceareformedasantisymmetrisedtensor
whereC is given, for ω>µ (with µ the chemical potential,
products of atomic localised orbitals. Considering Hamilto- q
wherewetakethegroundenergyasthezeroenergy),by
nians that are invariant under translation by a Bravais lattice
vector,andwhicharealsotimeindependent,itisconvenient,
Carpes= ∑exp(iqr)c† (14)
forcalculationpurposes,tostudytheGreen’sfunctionsspace- q i i
i uc
time Fourier transform. All kinds of Green’s functions (ad- ∈
vanced, retarded, time-ordered, correlationfunctions)canbe and,forω<µ,by
obtained from the following function (takingh¯=1 for nota-
tionalsimplification): Cipes= ∑exp( iqr)c (15)
q − i i
i uc
∈
Φ eS†(P (B))†(H ω iγ) 1P (A)eS Φ
G(B,A,ω,q,γ)=h 0| q − − − q | 0i In the general case, to calculate equation (8), we have to
Φ0 eS†eS Φ0 solvetwoproblems. Thefirstoneistofindanapproximated
h | | i
(8) solutionRfortheresolventequation:
In this expression A and B represent, each, a product or a
sum of products of cˆ† and cˆ operators, and Pq is an opera- (H ω iγ) R =Pq(A)eS Φ0 (16)
− − | i | i
torwhichsymmetrizesitsargumentundertheBravaislattice
vectortranslationsoperators,herenotedasT : andthesecondoneistocalculatethescalarproduct. Theop-
R
erator P (A) can contain, for a generic A, both annihilation
q
P(D)=∑exp(ikR)T (D). (9) andcreationoperators. Inordertobetteraccommodateequa-
k R
R tion (16) into the framework of the coupled cluster method,
whichdealswithexcitationswhicharewrittenasaproductof
As an example of the different choices of A and B, we can
creation operators only, we rewrite A in a commutated form
considertheadaptationofequation8tothedynamicstructure
A whereannihilationoperatorsdisappear:
S
factor calculation. The dynamic structure factor is observed
experimentally in different kind of spectroscopies like elec- (H ω iγ) R =eSP (A ) Φ (17)
q S 0
tronenergyloss(EELS)andInelasticX-RayScatteringSpec- − − | i | i
troscopy(IXSS).Thedynamicstructurefactorisgivenbythe whereA isdefinedastheHausdorffexpansion(7)ofe SAeS
S −
space and time Fourier transform of the density-density cor- discardedofallthosetermshaving,withrespecttotherefer-
relation function7,8 which is defined, writing n(q,t) for the encevacuum Φ ,non-contractedannihilationoperators:
0
spatialFouriercomponentqoftheelectrondensityattimet, | i
as: eSPq(AS) Φ0 =eSPq(e−SAeS) Φ0 . (18)
| i | i
1 ∞
S(q,ω)= dtexp(iωt) n†(q,t)n(q,0) (10) The same transformation can be applied to the B operators
2ωZ−∞ D E whichappearintheGreen’sfunctionequation(8). Notethat
intheaboveexpressionwehavemadeuseofthetranslational
andcorresponds, inequation8, tothespectraldensityofthe
invarianceoftheSoperator.
Green’sfunction:
We represent an approximated solution for the resolvent,
1 introducingtheapproximatingoperatorRω,γandthefollowing
S( q,ω)= (G(Nq,Nq,ω,q,γ) G(Nq,Nq,ω,q, γ)) ansatz, which is refined at each increment of the number Nr
− 2iω − −
(11) byaddinganewexcitationoftheform∏n cˆ† :
k=1 jk
where A and B have been replaced by N , written summing
q
overalltheorbitalsibelongingtoanunitcelluc(atvariance R = Rω,γeS Φ0 (19)
| i | i
watipthosniqtiownhricihinissidsuemucm:edoverthewholelattice),andlocated Rω,γ=rω0,γPq(AS)+∑Ni=r1rωi,γ(cid:16)Pq(cid:16)∏nk=ri 1cˆ†ji,k(cid:17)(cid:17)⊥i (20)
Nq= ∑exp(iqri)c†ici (12) Isntatnhdissfeoxrprtehsesioornthroigaorneaflriesaetipoanraomfeittesrasragnudmtehnet(w.)i⊥thi sryemspbeoctl
i uc
∈
3
to Pq(AS) and to all the operators Pq ∏nk=rl 1cˆ†jl,k with l <i. skiemeppilnegMthneOscahmaeini.mpWleempeenrtfaotrimon,thtehseescarmosesmvaoldidela,tiaonnds bbyy
We build our spectral CCMequations(cid:16)(SCCM e(cid:17)quations) by
multiplyingequation(17)attheleftwithe S,andbysetting simplyadaptingthesystemparameterstothetestcases.
−
theresiduetozero:
nr † 3.1. ModelandparametersfortheMnO2plane
i
0∀i∈=[1,Nr]hΦ0|Pq k∏=1cˆ†ji,k!⊥i e−S(H−ω−iγ)Rω,γeS|Φ0i theWperekveioepusinpatpheisr9waoprakrtthfoersiatms seizmeoadnedl fHoarmailmtoondiaifincaatsioinn
  (21)
tothespindependenttermattheMnsiteswhichisaimedto
Theequation(21)isexpandedbytheHausdorff(7)expansion implement the spin ordering. For sake of clarity we report
formulasubstitutingAwithSN andBwith(H ω iγ)Rω,γ. hereallthedetailsofthemodelHamiltonian.
The Hausdorff expansion contains a finite nu−mbe−r of non- TheMnsitesareplacedatintegercoordinates(2i,2j)with
zerotermsbecauseRω,γcontainsonlycreationoperatorswhile j∈[0,Nx−1]andi∈[0,Ny−1].Theoxygenatomsareatpo-
eachtermofH containsamaximumoffourannihilationop- sitions(2i+1,2j)and(2i,2j+1). Periodicboundarycondi-
erators. tionsareapplied,withperiods2Ny,2Nx inthetwodirections.
Theexpansiongivesasetoflinearequationsforther pa- Inordertolimitthecomputationalcost,werestrictthedegrees
rameters. Theaboveequation(21)istheequivalentofequa- offreedomtothoseorbitalswhicharethemostimportantfor
tion15ofCrawfordandRuud3transposedtothecaseofape- the electronic properties of the MnO2 plane in the mangan-
riodiccrystalandforexcitationswhichdiagonalisethetrans- ites. These are the eg 3d orbitals of Mn, namely the x2 y2
−
lationalsymmetryoperators. and 3z2 r2 orbitals, and the p oxygen orbitals which point
−
Theresolventequationaccuracyisimprovedbysystemati- towardMnsites. Thex,yaxisareintheMnO2 planeandzis
callyincreasingNr,selecting,ateachiteration,thesetofnum- outofplane.TheMnt2gorbitals(xy,xz,yz)ofagivenspindi-
bers rection(upordowndependingontheatomicsite)have,each,
occupancyfixedto1. Theothert ofoppositespindirection
2g
nr ,j ,....,j (22) havetheiroccupancyfreezedto0.
Nr+1 Nr+1,0 Nr+1,nrNr+1 Wefreezethet degreesoffreedomwiththeaimofreduc-
n(cid:16) (cid:17)o 2g
notyetcontributingtotheresolventandcorrespondingtothe ingthecomputationalcostofthemodel. Thejustificationfor
largestscalarproductintheSCCMequation(21). Whenwe thisfreezingisthat,neglectingtransitionto(much)higheren-
calculatetheresiduewefixω=ω atthecenterofthespectral ergyorbitals,theelectron-electroninteractionactingbetween
r
regionofinterest. Theparametersraregiven,overthespec- two t2g electrons remains within the t2g sector, the one be-
tral region of interest, by a linear algebra operations of the tweentwoeg’s: withintheeg one. Thisistrueforsymmetry
kindr=(M1)−1(M2+ωM3)wheretheM’sarematricesob- reasons. For the same reasons the interaction between a t2g
tainedfromSCCMexpansion. OnceweknowtheRoperator electron and a eg one does not change the occupation num-
wecancalculatethespectrawiththefollowingequation: beroft2g ore2g orbitals. Moreover,inthecaseofaparallel-
spint -e pair,theCoulombinteractiondoesnotchangethe
2g g
Φ0 eS†(Pq(BS))†Rω,γeS Φ0 spinvalueoftheorbitals. Theexchangeinteractionbetween
G(A,B,q,ω,γ)= h | | i (23) at electronandae oneofoppositespinistheonlyinterac-
Φ eS†eS Φ 2g g
h 0| | 0i tionwhichcanflipt2g spins. Ifweneglectthisinteractionwe
This expression is calculated following the procedure al- canfreezethet2g electronvariablesandsubstitutethemwith
readydescribedinSCCM5,usingtheWick’stheoremandthe an effective mean-field exchange interaction. This approxi-
linked-clustertheoremtoobtainasumofproductsofGreen’s mation should not alter too much the physics of the MnO2
function of different orders, and a hierarchical set of Dyson plane model which is dominated by eg electrons. The full
equations for the Greens functions, that we truncate at the many-bodytreatmentcanbeextendedtot2g byrestoringthe
Hartree-Fockapproximation. off-diagonal part of the t2g e2g exchange interaction at the
−
expense of a higher computational cost, but this will not be
doneinthispaper.
3. MODEL For oxygens we restrict to py for the (2i+1,2j) sites and
to p at the (2i,2j+1) sites. These are the oxygen orbitals
x
Inthispaperweextendouroriginalapplication5ofSCCM, whichplaythemajorroleinbridgingtheMnsitesalongthe
xandydirections. Forsymmetryreasonstheyarecoupledto
which was tested on a small 2x2 MnO lattice with periodic
2
thee sectoronly.
boundaryconditions,toalargerlattice. Consideringalarger g
The Hamiltonian acting on the Hilbert space spanned by
latticeallowsustofinelyresolvethekdependencyoftheex-
theseorbitalsiscomposedofthefollowingparts
citations and to accommodate in the model different choices
of charge, spin, and orbital ordering. We detail below the
model Hamiltonian used for the MnO2 plane. We also val- H=H +H +HMn+HMn+HO (24)
bare hop U J U
idate our implementation, deriving from it simpler systems
that can be solved exactly: the Hubbard 4 4 model and a namelyH whichcontainstheone-particleenergiesofthe
bare
×
4
orbitals,thehoppingHamiltonianH whichmoveselectron FinallytheoxygenHubbardtermis
hop
between neighboring sites, the Hubbard correlations HMn,
U
HMn and HO for Manganese and Oxygen. The bare Hamil-
J U HO=∑U n n +
tonian,usingletters panddforthesecondquantisationoper- U i,j p (cid:18) px,s=+12,2i+1,2j px,s=−12,2i+1,2j
atorsonoxygenandmanganeseorbitalsrespectively,is:
n n (30)
Hbare=i,j∑,gd,s(εd+(1/2−smi,j)h)dg†d,s,2i,2jdgd,s,2i,2j+ py,s=+12,2i,2j+1 py,s=−12,2i,2j+1(cid:19)
i∑,j,sCapid3†z2−r2,s,2i,2jd3z2−r2,s,2i,2j+ weTuosfiexktnhoewfrleedegpearfarommeteorusropfrtehveiomuosdwelofrokrothnemMannOga2npitleasn9e.
(ε U )∑ p† p +p† p (25) Parameters are given in eV units. The effective Slater inte-
p− p x,s,2i+1,2j x,s,2i+1,2j y,s,2i,2j+1 y,s,2i,2j+1 grals used in that work correspond, in the present model, to
i,j,s
(cid:16) (cid:17) U =6.88,U =5.049,J =+0.917. Theexchangesplitting
d d0 d
where the mi,j is a function of the site position and takes h induced by occupied polarisedt2g orbitals is h’2eV. We
the values 1 depending on the Mn t magnetisation. The use a hopping t =1.8 taken from our previous work9. For
2g
±
parameter h is the exchange splitting induced on the eg or- oxygen we choseUp =5eV. The parameter εp controls the
bitals by the t electrons. The index g takes the val- amountofchargeback-donationfromoxygentomanganese.
2g d
ues g = x2 y2,3z2 r2. The Mn one-particle energies The predominant O 2p character of doped holes found in
d
(εd+(1/2 −s mi,j)h)−are spin (s= 1/2) and position de- manganites10 correspondstoavalueεp whichraisesthebare
− ±
pendent and take into account the mean-field exchange with oxygenorbitalsenergiesabovethebareMnones.Thevalueof
theMnt2goccupiedorbitals(xy,xz,yz)(whosedegreesoffree- εpinfluencestheaverageoccupationoftheegorbitals. These
domarediscardedfromthemodel).Theenergyofthe3z2 r2 occupancies match the ones found in the previous works for
orbitalisraisedbyavalueCapiabovethebareenergytosi−mu- avalueεp’2. WehavefixedtheCapi effectiveapicalligand
latetheligand-fieldsplittinginducedbythehoppingtoapical field splitting in the following way. We have considered an
oxygens which are not otherwise treated in the model. The isolatedplaquetteformedbyaMnsiteandthefourneighbour-
oxygen orbitals one-particle energy term takes into account ingoxygenatoms. Usingthepresentmodelisation, withone
dthoeuHbluebabnadrdsicnogeleffiocciecnutpa−tiUonpstooncoomxypgeennsas.teHUO andfavoring eegteresl.ecFtroornl,otwwovaeluleecstroofnCsappeirthoexyoguet-no,fa-pnldantheeegabiosvoeccpuapraiemd-.
Thehoppingtermis We have increased Capi till the in-plane orbital is occupied.
Thisoccursat1.5eV.Thenwearbitrarilyscaleitdownto1eV
in order to realise a system which, without the in-plane ki-
H =t ∑ ∑ l f p† d +
hop gd,x x,s,2i l,2j gd,s,2i,2j netic effect, which are enforced when the complete plane is
i,j,gd,sl=±1 (cid:16) − considered,wouldhaveout-of-planeoccupiede orbital.
g
f p† d +c.c. (26)
gd,y y,s,2i,2j l gd,s,2i,2j
−
(cid:17)
where
3.2. ParametersfortheHubbardmodel
f = f =1/2
3z2 r2,x 3z2 r2,y
− −
f = f =√3/2 (27) In order to cross-check our implementation on the
− x2−y2,y x2−y2,x Hubbard-Modelwecompareittotheexactnumericalcalcula-
The Coulombintra-site repulsiveinteraction forMn is given tionofa4 4systemathalf-fillingwithU/t=8(parameters
×
byUd foranelectronpaironthesameorbital,andbyUd0 for of the Hubbard model). To obtain the Hubbard model from
twodifferentorbitals: our MnO plane model, we set to zero all the model param-
2
etersexceptHMn. Thenweaddanhoppingtermactingonly
U
HUMn=i,∑j,gdUd ngd,s=+21,2i,2jngd,s=−21,2i,2j+ betweenx2−y2orbitalsofneighbouringsites:
i,j∑,s1,s2Ud0 n3z2−r2,s1,2i,2jnx2−y2,s2,2i,2j (28) Hhop=−t ∑ dx†2 y2,s,2i+2,2jdx2 y2,s,2i,2j+
i,j,s − −
(cid:16)
whereletterndenotesthenumberofelectronsonagivenor- d† d +c.c. (31)
bital. TheCoulombexchangefore orbitalsis x2 y2,s,2i,2j+2 z2,s,2i,2j
g −
(cid:17)
TospanthesameHilbertastheHubbardmodelweinitialise
HMn=J ∑ d† d† d d the CC method with a ground state where the only occupied
J d i,j,s1,s2 3z2−r2,s2,2i,2j x2−y2,s1,2i,2j 3z2−r2,s1,2i,2j x2−y2,s2,2i,2j orbitalsarethedx2 y2 ones. TheHamiltonian,withthesetof
(29) parametersherede−scribed,doesnotconnecttootherorbitals
whilethee -t exchangeisincludedasamean-fieldtermin- than the d ’s and thus the effective degrees of freedoms
g 2g x2 y2
sideH . coincidewit−hthoseoftheHubbardmodel.
bare
5
3.3. Parametersforthesimplechain Necc 2nd quantisationform Asciiart
Φ0
| i | i
ssNypysiWn=te-emd1oc,wowamnnhpdiceahlureescwetoretuohrnoesmbftoeaoltinlhnoloywid,nintatoghoeepnxafeaorcaldtlmoismweotlieeunnrtgssio:iwontnsa=ayfl.o1rcW.h8aa,esiCnicmaosppnielsti=tifiidnee1gdr, 123 d3†z2−dpr3††x2z11201−ddxxr2220−−pyyx22801|||ΦΦΦ000iii ||| iii
εfepre=nt2vaanludeεsdo,fhU,Up.,TUhde,Unujm=b0e.rsCoaflcsuitleastioinntahreexdodnireecattiodnifs-, 4 p†x7d3†z2−r26px5dx2−y26|Φ0i | i
d0
Nx is kept small enough so that an exact numerical solution, ··· d† d† p† d† ··· ···
forthegroundstateandforthespectra,canbefound. 300 3z2−r26 3z2−r24 x3 3z2−r22
dx2−y22dx2−y24dx2−y26dx2−y210|Φ0i | i
4. DISCUSSION TABLEII.SecondquantisationandAscii-artrepresentationofsome
selected excitations, appearing in the S operator which gives the
4.1. Thesimplechain groundstateofthesimplechainmodel,forNx=6andUd0 =3. The
chainispopulatedwithspin-downelectronsThepositionswheresec-
ondquantisationoperatorsactonthereferencestatearemarkedwith
Weshowintable(II)thereferencestate Φ0 andsomese- a black background. The Mn sites are indexed with even indexes
| i
lectedexcitationsforthesimplestofourthreestudiedmodels: from0to10,theOsiteswithoddindexesfrom 1to11,with 1
− −
thesimplechain,forU =3. Inthefirstcolumnwemarkthe beingequivalentto11giventheperiodicboundaryconditions.
d0
iterationnumberatwhichtheexcitationappearsintheexpo-
nentS (seeequation(1))whenwecalculatethegroundstate.
N
ThethirdcolumnshowsanAscii-artviewofthechainwhich
canbereadfollowingthelegendoftable(I). Thechaincon-
sistsofsixMnOunitsandispopulatedwithninespin-down
a) b) c)
d) e) f)
TABLEI.GuidetoAsciiartrepresentationofelectronicstates. a)a
spin-downelectronoccupyingax2 y2orbital.b)aspin-upelectron
inz2 orbital. c)OnthesameMns−ite: aspin-downinx2 y2 anda
spin-x2 y2upinz2. d)thexsymbolisusedfordoubleo−ccupancy:
here x2−y2 is doubly occupied. e) on the right of the Mn site an
−
oxygenorbitalisdoublyoccupied. f)ablackbackgroundhighlights
the orbital where an excitation is occurring. Here, with respect to
areferencestategivenbye),aspin-downelectronhasbeenmoved
fromtheoxygenorbitaltotheMnz2orbital. Emptyorbitalsarenot
drawn.
electrons. The occupied Mn orbitals are represented by
x2 y2
acrosscomposedofminussig−nstosignifyaspin-downelec- FIG.1.Simplechaingroundstateenergyasafunctionofthenumber
tron(seelegendintable(I)). Inthereferencestate|Φ0i,they ofexcitationinSoperatorforUd0 =1.
are occupied at all six sites. The remaining three electrons
populate one every two oxygen orbitals. These orbitals are
represented,whentheyareoccupiedbyaspin-downelectron, particle caseU =0 (not shown) the exact ground energy is
d0
byaminussign. Thespin-upsectorisdiscardedinthissim- recoveredwith13one-particlesymmetrisedexcitations. We
ple model. The positions where excitations are created are showinfigure(3)thespectralfunctionA(q,ω)(seeequation
highlightedbyablackbackground. Theground-stateenergy, (13))foranARPESoperatoractingonMn3z2 r2orbitals
−
versusthenumberofexcitations,isshowninfigures(1,2),for
twochoicesofthecorrelationenergy:Ud0 =1andUd0 =3.We Carpes= ∑ exp i2πkl d† (32)
have let the algorithm to iterate long enough so that the per k l<Nx (cid:18) Nx (cid:19) 3z2−r2,s=−21,2l
Mn-site energy converges to within few meV from the exact
energy.Infigures(1,2)weplotthegroundstateenergyversus The SCCM spectral response , calculated with 1.5 104
×
thenumberofexcitations,forU =1andU =3. Wereach spectralexcitationoperators,isshownforU =1,upperrow,
d0 d0 d0
a convergence error of the order of a hundredth of eV after andU =3,lowerrow,andiscomparedtotheexactcalcula-
d0
250iterationsforU =3whileforU =1theconvergenceis tionandtotheHartree-Fockmeanfieldmethod. ForU =1,
d0 d0 d0
fasterasexpected,duetotheweakercorrelation. Intheone- at k=0 the spectra consists in a unique peak at ω=1.75.
6
Atk=1thispeakisfoundatlowerenergyandnewspectral
featuresappeararoundω=3.35. Athigherk onlytheselat-
terspectralfeaturessurviveandshifttohigherenergy. Given
thesizeofthechain,andtimereversalsymmetry,thespectra
fork=4(5), notplotted, isidenticaltotheonefork=2(1).
ForU =3weshowthemainspectralfeatureswhichaccount
d0
for more than 90% of the spectral density. The main spec-
tral feature is found around ω=2.9 for k=0 and moves to
ω=3.5 for k=1. Still for k=1 we note the emergence of
new spectral features at the slightly higher energy ω=3.9.
At variance with the U =1 case these features keep mov-
d0
ing at higher energy but, for k =2 and k =3, they remain
marginally small and are not plotted. For these latter values
ofkwenoteinsteadthepersistenceofthelowestenergypeak
whichmovestolowerenergies.Wecanseefromthesenumer-
ical experiments that the SCCM reproduces well the shapes,
the positions, and the strength of all the spectral features for
alltheconsideredvaluesofU , whilethediscrepancyofthe
d0
Hartre-Fockmethodincreases,asexpected,increasingU .
d0
FIG.2.Simplechaingroundstateenergyasafunctionofthenumber 4.2. TheHubbardmodelathalffillingandU/t=8
ofexcitationinSoperatorforUd0 =3(right).
.
4.2.1. GroundstatewiththeCCmethodandbyexactnumerical
solution.
#8 #7
#4
#2
#5
#6
#1 #3
FIG.4. Hubbard8 8modelathalffillingandU/t=8. Ascii-art
×
representationofthefirsteightexcitationsappearingintheSoperator
forthegroundstate.
FIG. 3. The spectral function for the simple chain atUd0 =1 and Weshowinfigures(4)and(5),usingasciiart,aselectionof
fifoerldCskaorlpuestio=n;∑dlo<tNtexde:xtph(cid:16)eie2Nxπxakclt(cid:17)ndu3†zm2−er2r,iσc=a−l12s,2oll.utSioonli;ddalisnhee:d(trheedmcoelaonr ethxeciHtautibobnasrdap8pea8rinmgoidnetlhaetShaelxfpfiolnlienngt,fUor/tth=eg8roaunnddpsetraitoedoicf
online):theSCCMsolution. boundary cond×itions. At each iteration step we introduce a
new excitation and add it to the operator S together with all
7
#50 #75
#80
#70
#40
FIG.5.Hubbard8 8modelathalffillingandU/t=8.Aselection
×
ofexcitationsappearingatfurtherstagesoftheconvergenceprocess.
FIG.7.Convergenceforasmaller4 4latticeandcomparisonwith
×
exactdiagonalisation.
theexcitationsthatcanbeobtainedbyapplyingallthesystem
translation,rotationandspinreversalsymmetries.
We show in figure (6) the ground state energy versus the thefirstvectoroftheKrylovspace,wetakethesamereference
number of excitations for the 8 8 lattice and in figure (7) statethatwehaveusedintheCCmethod.Ifwerestraintothe
×
the comparison with exact diagonalisation for a 4 4 lat- sector of even spin-reversal symmetry we still get the same
×
tice, whose smaller size allows for an exact numerical solu- groundstate.TherestrictionisrealisedbystartingtheKrylov
tion. The dimensionality of the Hubbard 4 4 at half filling spacewithavectorwhichisthesumoftheCCreferencestate
is 1.7 108 whichisstillaffordablewith×thecurrentcapa- plus the state obtained by swapping the spin-plus and spin-
∼ ×
bilityofatypicalworkstation. Thenumericalsolutionofthe minus checker-boards. The Hamiltonian commutes with the
spin-reversaloperatorandgeneratesthespannedspacewhich
isthereforeaneigen-spaceofthespin-reversaloperator. An-
othergroundstateisfound,instead,ifwerestraintotheodd
spin-reversal symmetry sector, by starting the Krylov space
with the difference between the two checker-board-swapped
states. Thenonnegligibleenergydifferencebetweenthetwo
ground states indicates that collective spinonic fluctuations
playanimportantroleintheHubbardmodel.Theimportance
ofthesefluctuationsisalsovisibleintheasymptoticbehaviour
oftheCCgroundstatesenergy. Afterreachingtheenergyre-
gionwhichcorrespondstothecollectivespinonicexcitationof
thegroundstate, anyfurtherimprovementoftheCCground
energy is very slow. This is due to the fact that the collec-
tivespinonicfluctuations,appearingatallthepossiblewave-
lengths,requirealargenumberofexcitationoperators,forall
thesub-partsofthelatticeconcernedbyaspin-reversalfluc-
tuation, and each excitation concerns all the electrons inside
these regions. We see instead, in figure (5), that the excita-
tions appearing in S are localised to small sub-parts of the
lattice and thus cannot represent collective spinonic fluctua-
tions. Thishappensdespitetheelevatednumberofiterations
that we have considered and, in particular, in the region be-
FIG.6. groundstateenergyversusthenumberofexcitationsforthe tween the 40th and 80th term where the CC energy seems to
8 8lattice. bestabilised.
×
Weshowinfigure(8)thecomparisonbetweenthespectra
groundstateisobtainedbyLanczosdiagonalisationwhere,for calculated in the framework of the exact numerical solution
8
andthespectracalculatedbySCCM,fora4 4lattice. The
×
exact numerical spectra was calculated applying the ARPES
operatorCkarpes=∑ix,iyexp(i(Qxix+Qyiy))c2ix,2iy,s=−12 onthe
ground state obtained by Lanczos diagonalisation, and then
obtainingtheChebychevpolynomialcoefficientsofthespec-
tra by iteratively applying the Hamiltonian on the resulting
state11. IntheHubbardmodelthereisjustonetypeofactive
orbital,thatinourimplementationisd ,andwerepresent
x2 y2
it,inthiscontext,withthesymbolcfor−notationalsimplifica-
tion.TheSCCMspectrahasbeenobtainedwith16excitations
( N =16 in formula (2)) and 3 104 terms in the resolvent
×
( formula (20)). We show in the upper row the spectra for
Q=(π,π). Ontherightgraph,upperrow,theSCCMspectra
is compared to the exact one. We can see that apart for the
mainlineposition,manyfeaturesintheSCCMspectraarein
disagreementwiththeexactcalculation. Ontheleftside,the
sameSCCMspectraiscomparedtotheexactcalculationfor
amodifiedHubbardmodelwhereanextratermH hasbeen
pin
addedtotheHamiltonian, inordertopincollectivespinonic
fluctuation.Thistermconsistsinaweakmagneticfieldacting
onhalfofthesites: thechecker-boardblacksitesdefinedby
(i+j) mod2=0,tofavorspin-upoccupation:
H = ∑ hc† c (33)
pin − 2i,2j,1/2 2i,2j,1/2
(i+j) mod2=0
FIG.8. Thespectralfunctionforthe4 4HubbardmodelatU =
8,t=1andforCkarpes=∑ix,iyexp i(Qxix×+Qyiy) c2ix,2iy,σ=−12.,cal- samndalwlceohmavpeargeidvetontthheishtoeprmpinthgeasntdrecnogrtrheloaftihon=te0r.m1swbhuitchcains
culatedwithN=16inSN and3×(cid:0)104termsinth(cid:1)eresolventRω,γ. becomenonnegligibleforlong-rangefluctuationswhichswap
largenumberofspin-upandspin-downsites.Theeffectofin-
troducing this term is that, as expected, the modified system
can be represented with improved fidelity by SCCM which
starts, in our implementation, from a defined reference state
with no fluctuations. The SCCM calculation at ARPES mo-
mentumQ=(0,π)andQ=(0,0), shownonthelowerline,
reproduce the pinned Hubbard model with a similar level of
fidelity. AlargernumberN ofexcitationsintheCCexponent
S couldimprovethefidelityoftheSCCMcalculation.How-
N
ever, despite the extensive properties of CC method which
render it capable of treating larger system than the exact di-
agonalisation,forthisspecificexamplewhichisstilltreatable
byexactdiagonalisation,theCCrequiresimportantcomputa-
tionalresources. Eachexcitation intheS exponentis sym-
N
metrised by the system symmetry group, which for the half-
filled4 4modelconsistof128elements,sotheoperatorS
N
×
is the sum of thousands of terms already with N as small as
16. TheSCCMproceedsbykeepingtrackoftheresiduesin
equation(21)whicharegeneratedbye−S(H ω iγ)Rω,γeS
− −
wheretheresolventR consistsinourcaseof3 104terms.
ω,γ
×
All this, traduced into our implementation, almost saturates
the 250Gb memory of one node of the ESRF cluster. The
situation can be alleviated by distributing the memory over
FIG. 9. The spectral function for Carpes = severalnodes. Byanappropriatedistributionandmessaging
k
∑miox,diyeelxwpitih(QUxi=x+8Q,tyi=y)1.c2iTx,2hiey,σS=C−C12M, fcoarlcutlhaetion8×wit8h NHu=bb8arids swcihthem4enothdiess,pwrohfiitcshaalmsoouthnetsetxoe1cu1t2iocnput-imcoer.es,Htohweeevxeercuetvioenn
(cid:0) (cid:1)
comparedtotheSCCMcalculationforthe4 4modelwithN=16. timeremainoftheorderof24hours.Forcomparisontheexact
Bothspectrahavebeencalculatedwith3 10×4termsintheresolvent numericalcomputationofthespectrawithLanczosdiagonal-
×
Rω,γ isation and spectra decomposition in Chebychev polynomia
completeson28coreswithinawalltimeof2hours.
9
Considering the application to other systems, all this is zosdiagonalisationcancopewiththehalf-fillingcase,itisa-
nonetheless encouraging for the following reasons. First of fortiorisuitedtosolvetheproblemwhereoneormoreholes,
allthesizeofthesystemwehaveconsideredisalreadyvery or electrons are introduced in the system. It is therefore in-
close to the limit of what can be done with the exact diag- teresting, as a side-study of our work on the Hubbard sys-
onalisation, but not of what we can do with the SCCM. As tem,tocalculatethegroundstateenergyofthe4 4Hubbard
×
an example a 6 6 Hubbard system has a dimensionality of modelbyadding, tothehalf-filledsystem, oneortwoholes,
8.2 1019which×exceedsbyfarthememorycapabilitiesofthe with a determined total kinetic moment. As in the case of
×
mostpowerfulsuper-computernowadays,butwecanshowin spin-reversal parity treated in the preceding sub-section, the
figure(9)theSCCMspectraforthe8 8Hubbardmodelcal- totalmomentumofthesystemisdeterminedbythefirstvec-
culatedwithN=8inS and3 104×termsintheresolvent. toroftheKrylovspace. Thegroundenergyofthehalffilled
N
×
ThenumberofexcitationshasbeenlimitedtoN=8,forthe ground state is E = 8.468875. Adding one hole we get
0
−
ground state, to limit the memory usage, because the string the ground energy E = E +∆ with ∆ = 1.678365,
1h 0 1h 1h
−
which represent the residues in second-quantisation are now while by adding two hole to the half-filled system we get
four times longer. Despite the poor level of approximation E =E +∆ with∆ = 3.399965. Thegroundstatefor
2h 0 2h 2h
−
ofthegroundstate,whichcouldbeimprovedusingmorere- the two-holes doped system is realised introducing a couple
sourcesonasuper-computer,thespectrastillreproduceswith formed by a spin-up and a spin-down hole with (π,π) to-
goodagreementthetwomainprincipalfeatureofthespectra. tal momentum. These energies satisfy the strict inequality
Moreover, the Hubbard model maybe considered as an ex- ∆ <2∆ . We can now perform a conceptual experiment
2h 1h
byconsideringanensembleof4 4Hubbardsystemsathalf-
×
filling. By adding two holes, each in a different Hubbard
system, we increase the energy by 2∆ . The resulting total
1h
energy is therefore higher than what we obtain placing both
#3144
holesonthesameHubbardsystem.
Tocharacterisethefluctuations,giventhecalculatedground
stateofthe4 4system,wecalculatetheexpectationvalues
×
ofthefollowingoperator:
S = ∑ S† S
AB A,α B,α
α=x,y,z
with
SA,α= ∑ ∑σαs1s2c†2ix,2iy,s1c2ix,2iy,s2
(ix+iy) mod2=0 s1s2
SB,α= ∑ ∑σαs1s2c†2ix,2iy,s1c2ix,2iy,s2 (34)
(ix+iy) mod2=1 s1s2
where σ are the Pauli matrices, and we are taking the spin-
FIG. 10. An excitation occurring in the ARPES resolvent for the spin scalar product between the black and white sites of the
Hubbardsystemhighlightingtheroleofcollectivespinonicfluctua- checkerboard. For the nominal reference state, with spin-up
tionindressingchargecarriers. blacksitesandspin-downwhitessites, theexpectationvalue
is < S >= 64. This is also the maximal value of the
AB
−
tremecasefortheapplicationoftheCCmethod,becauselong operator. The ground state of the half-filled system, calcu-
rangefluctuationsbringthesystemveryfarfromthereference lated by exact diagonalisation is <SAB >= 60.4. For the
−
state.Inthenextsectionwearegoingtostudytheenergy-loss exact ground state of the one and two-holes systems we get
spectroscopyoforbital-orderinginmanganitesforwhichwe <SAB>1h= 42.5and<SAB>2h= 30.8respectively. By
− −
expectthattheselongrangefluctuationbelesssevere.Inthese namingrtheratiobetweentheSAB expectationvaluesonthe
systems, infact, thephysicsisdominatedbythemobilityof groundstatesandthelargestmagnitudevalues,whichare64,
thee chargeinthemeanexchangefieldofthet electrons. 56 and 49 for zero, one and two holes respectively, we get
2g 2g
As the long range magnetic ordering is well defined in this r=0.944 , r1h =0.76 and r2h =0.63 for zero, one and two
case,webuildtheCCmodelingstartingfromareferencestate holes,whileatthesametimethetotalspinofthewhole4 4
×
whichreproducessuchordering. latticeisfoundtobeasinglet,adoublet,andagainasinglet.
These numerical evidences highlight the following sce-
nario:
4.2.2. Observationofapairingmechanism The collective spinonic fluctuations occurring in the
•
ground state have a relatively long range so that, for
The dimensionality of the Hubbard system is maximal at the considered case, short range antiferromagnetic or-
half-filling.Byintroducingoneormoreholesorelectronsthe dering,onthescaleofthe4 4lattice,ispreserved(r=
×
dimensionality decreases. Our implementation of the Lanc- 0.944).
10
#22 #10 #21
#12 #11
A
#8
#7
#17 #20
B #6 #13 #23
#9
#4 #18 #19
FIG.11.ReferencestatefortheMnO2planeathalf-dopingwithor- #3 #14
bitalordering. Alongaferromagneticchaindouble-occupancyoxy-
gen sites (called A sites) alternate with single occupancy sites (B #2 #5 #15 #16
sites).TwocircleshighlightoneAsiteandaBsite.
#1
Adding one hole to the half-filled system brings down
•
r to r =0.76. This is caused by the hole mobility
1h
whichdisturbstheantiferromagneticordering.Thehole FIG. 12. First 23 excitations concurring to the MnO2 plane
groundstate
is dressed by a collective spinonic halo. Figure (10)
shows one of the excitation concurring to the ARPES
#133
resolvent,fortheHubbardmodel,whichhighlightsthe
#86
roleofcollectivespinonicfluctuationindressingcharge
carriers.
#128
Addingonemoreholebringsdownrtor =0.63,but
• 2h #137 #60 #198
thetwoholesmoveinacorrelatedwayformingasin-
#119
glet.
4.3. ModelingthehalfdopedMnO2planeofLa0.5Sr1.5MnO4 #77 #35
#122
Figure(11)shows,inasciiart(seelegendintable(I)),the #134
unit cell of the reference CC state used for the orbitally or-
deredgroundstateofaMnO plane. TheorbitalsatMnsites #116
2
are represented by symbols arranged in a cross pattern, for
x2 y2 orbitals, or a vertical line pattern for 3z2 r2 ones
− −
(notpresentinthereferencestate).Theorbitalsattheoxygen
FIG.13.Aselectionofexcitationsforfurtherrefinementsofthe
sitesarerepresentedbyasinglesymbolfortheorbitalwhich
MnO2ground-state.
bridges between two Mn sites. The symbols can be a minus
sign,plussignsortheXsymbol,forspindown,spinup,and
doubleoccupancyrespectively(seelegendintable(I)). This
referencestatecorrespondstotheorbitalandmagneticorder- reversalaccompaniedbya(1, 1)translation, theswapofx
ingobservedintheMnO planeinLa Sr MnO 9.Thehalf andyaxisaccompaniedbya(1−,1)translation.Theenergy,as
2 0.5 1.5 4
doping is realised, in this reference state, by removing one afunctionofthenumberofresolventterms,isshowninfigure
electron, every two oxygen sites, from the otherwise doubly (14). Wecanidentifythreedifferentregions. Thefistone,for
occupied oxygen orbitals. In particular, we remove, in each N upto(about)25, isarapiddescentoftheenergyandcor-
chain, only electrons with majority spin. This is done in or- respondsmainlytoshortrangeone-electronexcitations. The
dertostabilisetheferromagneticalignmentalongthezig-zag secondone,forN between25and90,correspondsmainlyto
chains. shortrangetwo-particlesandthree-particlescorrelationterms,
We show in figures (12,13), using ascii art, a selection of ascanbecheckedinfigure(13).Thelastpartistheonewhere
the excitations appearing in the S exponent for the ground the energy has the lesser variations, and is increasingly con-
state of the 8 8 MnO model at half doping and periodic tributed by many-body correlations of longer range, like ex-
2
×
boundary conditions. At each iteration step we introduce a citations134and198(figure(13))whichcorrelatesthreeMn
new excitation and add it to the operator S together with all and one O atoms. The further we push the refinement pro-
the excitations that can be obtained by applying all transla- cess, the closer to the true ground state the CC solution is
tion, rotation and spin reversal symmetries of the reference supposed to be. In particular the CC solution is supposed to
state. These symmetry operations generators are: the trans- restorethesymmetryofthegroundstate,whichishigherthan
lations generated by the (2,2) and (2, 2) vectors, the spin- for the reference state. A good indicator of this is the occu-
−