Is doped BaBiO a conventional superconductor? 3 V. Meregalli and S. Y. Savrasov Max-Planck-Institut fu¨r Festk¨orperforschung, Heisenbergstr. 1, 70569 Stuttgart, Germany. (February 1, 2008) 8 We report density functional calculations based on local density approximation (LDA) of the 9 properties of doped barium bismuthates. Using linear-response approach developed in the frame- 9 work of the linear muffin-tin-orbital method the phonon spectrum of the Ba0.6K0.4BiO3 system is 1 calculated and is compared with the results of the neutron diffraction measurements. The effect of n dopinginthecalculationismodelledbythevirtualcrystalandmassapproximations. Theelectron- a phonon coupling constant λ is then evaluated for a grid of phonon wave-vectors using the change J in the potential due to phonon distortion found self-consistently. A large coupling of the electrons 7 to the bond-stretching oxygen vibrations and especially to the breathing-like vibrations is estab- 2 lished. Also, a strongly anharmonic potential well is found for tilting-like motions of the oxygen octahedra. Thismodeisnotcoupledtotheelectronstolinearorderinthedisplacements,therefore ] n an anharmonic contribution to λ is estimated using frozen–phonon method. Our total (harmonic o plus anharmonic) λ is found to be 0.34. This is too small to explain high-temperature supercon- c ductivity in Ba0.6K0.4BiO3 within the standard mechanism. Finally, based on standard LDA and - LDA+U like calculations, a number of properties of pure BaBiO such as tilting of the octahedra, r 3 p breathing distortion, charge disproportionation and semiconducting energy gap value is evaluated u anddiscussedinconnectionwiththenegativeU extendedHubbardmodelfrequentlyappliedtothis s compound. . t a m - d I. INTRODUCTION. α=0.21±0.03 for BKBO and α=0.22±0.03 for BPBO. n Using their analysis, the authors of Ref. 5 concluded o Since the discovery of superconductivity at T ∼30K c that ”phononic” effects in these materials are only in- c in Ba K BiO (BKBO)1,2, and from earlier studies of [ 1−x x 3 dicative of dressedelectronic excitations. From studying BaPb1−xBixO3 (BPBO)systemwithTc ∼13K,thereisa the imaginary part of optical conductivity in BKBO the 2 fundamental question whether the conventional phonon authorsofRef.6gavethevalueofλ∼0.2. Electronicspe- v mediated pairing mechanism is operative in these high– cific heat measurements7 of Ba K BiO have yielded 1 T superconductors(HTSC).Dopedbariumbismuthates 0.6 0.4 3 5 arce different from the HTSC cuprates3, since no anti- Ns∗(0) ∼0.32 states/[spin×eV×cell] giving a mass en- 2 hancementfactorN∗(0)/Nband(0)∼1.4(Nband(0)∼0.23 ferromagnetic ordering exists for the parent compound s s s 1 states/[spin×eV×cell] for x=0.4). A complicated sit- BaBiO . This seriously doubts that strong electron cor- 0 3 uation exists with the transport measurements. The 8 relationsexistandareresponsibleforthepairing. Simple temperature dependent resistivities for superconducting 9 cubic superconducting phase makes BKBO and BPBO BKBO and BPBO have ranged from metallic to semi- / similar to the isotropic low T superconductors. How- t c conducting and two–channel model of the conductivity a ever, there is a number of features which makes doped m in the bismuthates was discussed8. While good grain– bismuthates similar to the cuprates. Both systems are boundary–free thin films and single crystals of BKBO - perovskite oxide superconductors with surprisingly low d doped well away from the CDW instability seem to ex- density ofstates atthe Fermi level. This canhardlygive n hibitmetallicbehaviour,thevaluesoftheresistivityitself high transition temperatures for the BCS–like supercon- o are (like in the cuprates) unusually high and are of the c ductors. Ashigh–Tccupratesareoriginatedfromantifer- order a few hundred µΩ×cm at room temperature This : romagnetic insulators, the parent BaBiO compound is v 3 couldpointoutthatanadditional(tostandardelectron– a charge–density wave (CDW) insulator in which oxy- i phonon) scattering mechanism is presented. X gen octahedra around the Bi ions exhibit alternating r breathing–inand breathing–outdistortions. The Bi ions The most direct evidence on the importance of a existinthe chargedisproportionatestate whichis chem- electron–phonon interactions in superconductivity of ically interpreted as 2Bi4+ ⇒Bi3++Bi5+. It is therefore the bismuthates has been given by the tunnelling temptingtoconnectthemechanismofsuperconductivity measurements9,10. Although not identical for differ- with the nature of these insulators. ent junctions, the deduced Eliashberg spectral functions Unfortunately,experimentalestimatesoftheelectron– α2F(ω)bearacloseresemblancewiththephonondensity phonon coupling strength do not lead to the firm of states determined by inelastic neutron scattering11. conclusion on the origin of superconductivity in the The estimated values of λ vary from 0.7 to 1.2 which bismuthates. Large isotope effect with α=0.4 has seemtobesufficienttoexplainhighcriticaltemperatures been reported for BKBO4. Other measurements5 give within the standard mechanism. 1 The electron–phonon coupling in doped BaBiO has interest since it may provide an insight on the supercon- 3 been investigated theoretically by several methods. The ductivity mechanism in the bismuthates. Recent calcu- authorsofRef.12studythisproblemusingtight–binding lations using constrained density–functional theory have (TB) fit to the energy bands which are obtained from been carried out to obtain the Coulomb interaction pa- density–functional calculations based on local–density rameters for the Bi 6s orbitals28. No indication for neg- approximation13 (LDA).The computedEliashbergspec- ative U of the electronic origin was reported. tral function α2F(ω) has been found to display promi- In the present work we try to address several prob- nent features in the frequency range corresponding to lems seen from the above introduction by means of the the oxygen stretching modes and the value of λ=1.09 state–of–the-artdensityfunctionalLDAcalculations. As has been reported. Crude calculations based on rigid– a first problem, we study lattice dynamics of the super- muffin–tin approximation (RMTA) also give large λ ∼3 conductingcubicBa K BiO usingrecentlydeveloped 0.6 0.4 3 indicating a strong–coupling regime14. Two estimates of linear-responseapproachimplementedwithintheLMTO λ using total–energyfrozen–phononmethod appearedin method29. The effect of doping is modelled by the vir- the literature15,16. Note that, in contrast to the TB and tual crystal and virtual mass approximations. On the RMTA methods, the frozen–phonon calculations treat basis of this calculation, we estimate electron–phonon screening of the potential due to lattice distortions self– coupling in this compound. The linear–responsemethod consistently. The value of the electron–phonon coupling used by us is advantageous in contrast to the frozen– strength equal to 0.3 for the breathing mode has been phononapproachsinceitallowsthe treatmentofpertur- found15 andtheroughestimateofλ∼0.5wasobtained16 bations with arbitrary wave vectors q. We have demon- using 12 q=0 phonons for ordered cubic Ba K BiO . strateditsaccuracybycalculatinglattice–dynamical,su- 0.5 0.5 3 There was a partial success in predicting structural perconducting and transport properties for a large vari- phase diagram for the parent compound BaBiO within ety of metals30, and we believe that our calculatedvalue 3 the LDA15,17–19. The experimental structure mainly of λ will be a realistic estimate for the electron–phonon consists of combined tilting and breathing distortions of coupling strength in this high–T superconductor. Also, c the oxygen octahedra corresponding to the instable R– a recent publication31 deals with the application of the point phonons of the cubic phase20,21. While rotational linear–responsemethodtostudytheelectron–phononin- instability was found in all calculations, the frozen–in teraction in another high–T superconductor CaCuO . c 2 breathing mode was not described by pseudopotential The second problem, we focus in our work, is study- calculation17 andtwo linear–muffin–tin–orbital(LMTO) ing the effects of anharmonicity in the electron-phonon calculations15,19 give the value for the breathing distor- coupling. It is widely accepted that certain phonon tion about 30% off the experimental one. Less rigor- modes are strongly anharmonic in the high–T mate- c ous potential–induced–breathing (PIB) model obtained rials. Frozen–phonon calculations produce double–well both instabilities22 with similar accuracy. It is not clear potentials for buckling motions of oxygen atoms per- whetherthesediscrepanciesareduetosensitivitytocom- pendicular to the CuO planes in nearly all HTSC31–33, putationaldetails ordue to the localdensity approxima- chain–bucklingdistortionsarefound to be anharmonic33 tion itself. in YBa Cu O , X–point tiltings of the octahedra along 2 3 7 Agreatamountofwork2 hasbeendonetounderstand (110)directions are instable34 in La CuO , andR–point 2 4 the properties of the barium bismuthates on the basis of instabilities corresponding to breathing and tilting exist thenegativeU extendedHubbardmodeloriginallyintro- in the doped barium bismuthates15,19,16. The influence duced in Ref. 23. The valence configuration of semicon- of anharmonicity to high–temperature superconductiv- ductingBaBiO canbeviewedasBa Bi3+Bi5+O which ityhasbeenaddressedinseveralpublications35–37,espe- 3 2 6 representsalatticeofelectronpairscentredateverysec- ciallybecauseofthesmallisotopeeffectfoundforHTSC ond Bi site (Bi3+). The sites occupied with Bi5+ions cuprates. Triple–degenerate pure rotationalmode at the are interpreted as those with no electrons. Rice and co– R point of cubic ordered Ba K BiO was predicted 0.5 0.5 3 workers24 have proposed that such local pairs are sta- to exhibit a double–well potential and some estimates of bilised by polarising the O octahedra and the effective theanharmoniccontributionstoλhavebeengiven16. We on–site U becomes negative due to the large electron– extendthisanalysisbysolvingnumericallySchr¨odinger’s phonon coupling. Recently, Varma25 has pointed out equation for the anharmonic potential well found from thatnegativeU canalsobeoftheelectronicorigindueto frozen–phonon calculations. The anharmonic λ is then theskippingofthevalence”4+”bytheBiion. Thelatter computed along the lines proposed in Refs. 37,38 by es- can provide a possible explanation for a well-separated timating the electron–phonon matrix elements from the optical and transport energy gap in the bismuthates26. energy bands computed for different tilting distortions. The mean–field phase diagram of the negative U Hub- Weconclude,inaccordwiththepreviousfindings16,that bard model exhibits several stable phases involving a this contribution, while not decisive, is not negligible for CDWsemiconductor,andasingletsuperconductor. This the total value of λ. is in qualitative agreement with the experimental phase The third purpose of our work is to study the prop- diagram27. erties of the undoped parent compound BaBiO . We 3 The question on the origin of negative U is of great try to answer the question whether the LDA gives an 2 adequate description of the ground state properties for exchange–correlation formula after Ref. 43 is used. The thischarge–density–waveinsulator. Sincetherewassome valence bands are treated scalar relativistically and the inconsistency reported in previous calculations15–19, we core levels - fully relativistically. A number of k points wanttoruleoutpossiblesensitivityofthefinalresultsto fortheBrillouinzone(BZ)integrationusinganimproved the internal parameters used in our band structure cal- tetrahedronmethod44 istakentobe20per 1 thBZ.The 48 culations with the full–potential LMTO method39. We chargedensity andthe potentialintheinterstitialregion carefullychooseourLMTObasisset,numberofkpoints, are expanded in plane waves with the cutoff correspond- plane waveenergycutoffandotherparametersbyexam- ing to the (28,28,28) fast–Fourier–transform (FFT) grid ining the convergency of the total energy and the cal- in the real space (approximately 10000 plane waves). culated properties with respect to them. Based on the Wefirstsummarisethemainfeaturesofthecalculated well converged data, we come to the conclusion that the electronicstructureindopedBaBiO . Theoccupiedpart 3 breathing distortions are seriously underestimated (ide- ofthebands(seeFig. 1)mainlyconsistsofBi(6s)–O(2p) ally, absent) in the LDA, and, therefore, the insulated hybridisedbandcomplex. Thisisinaccordwiththepre- state is not correctly described. This strongly resembles viouscalculation42. Forthecubicperovskitephase,there the situation with the antiferromagnetic ground state of isonlyonebandcrossingtheFermilevel,whichisanan- thecupratessuperconductorswhichisalsonotdescribed tibonding Bi–O sp(σ) band. Similar situation is found by the LDA40. We perform a number of model calcula- in the cuprate superconductors where Cu–O dp(σ) an- tions in the spirit of the LDA+U method41 in order to tibonding bands dominate at the Fermi energy E . A F clarify this problem. simple tight–binding model involving Bi(6s), O(2p) or- Therestofthepaperisorganisedasfollows: InSec. II, bitals, and two–centre nearest–neighbour sp(σ) interac- our linear–response calculations of the lattice dynamics tion can be used to understand the principal features of and the electron–phonon interactions in doped BaBiO3 these energy bands45. It was early noted45 that for the are described. Sec. III considers anharmonicity correc- case of half–filling (undoped cubic BaBiO ) this model 3 tionstoλforthetiltingmotionsoftheoxygenoctahedra. has a perfectly nested Fermi surface for the wave vector Sec. IV reports our LDA and model calculations of the q corresponding to the R–point. Therefore, it is tempt- ground state properties for pure BaBiO3. In Sec. IV we ing to interpret the appearance of breathing distortions give our conclusions. ascommensuratePeierlsinstability andcubic perovskite II. HARMONIC PHONONS AND λ Ba0.6K0.4BiO3 as a doped Peierls insulator46. To under- stand whether nesting can bring any effect in static sus- This section presents our results on the lattice dy- ceptibility, we have analysed q–dependence of the func- namics and the electron–phonon interaction for the cu- tion bic perovskite superconductor Ba K BiO . We also 0.6 0.4 3 summarise the main features of the calculated electronic structure and discuss our predicted equilibrium lattice δ(Ekj)δ(Ek+qj′) (1) configuration. The band structure calculations are per- kXjj′ formed with the highly precise full–potential LMTO method39. The details of the calculations are the fol- lowing: The effect of potassium doping is taking into account using virtual crystal approximation (VCA) by 1.0 considering a fractional nuclei charge Z=55.6 at the Ba E site. Numerous supercell investigations16,42 of the dop- 0.8 F ing influence on the calculated energy bands justify the applicabilityoftheVCA.Amultiple,three–κLMTOba- Ry0.6 sissetwiththe tailenergiesequalto-0.1,-0.8,and-2Ry. gy, isemployedforrepresentingvalencewavefunctions. The Ener0.4 valencestatesinclude6sand6porbitalsofBi,2porbitals 0.2 of O, and 6s orbitals of Ba. Such semicore states as 5d orbitals of Bi, 2s orbitals of O, and 5p orbitals of Ba are 0.0 treated as bands and are included in the main valence panel using two κ LMTO basis with κ21,2=-0.1,-0.8 Ry. -0.2Γ X M Γ R X The main panel also includes unoccupied 5d orbitals of Ba with the 2κ basis and 4f orbitals of Ba with the 1κ FIG. 1. Calculated LMTO energy bands for cubic basis (κ2=-0.1 Ry). Deeper lying 5s states of Ba are re- Ba0.6K0.4BiO3. The potassium doping is taken into ac- solved in a separated energy panel. All other states are count using virtual crystal approximation. treated as core levels. The muffin–tin sphere radii were takento be: S =3.25a.u., S =2.25a.u., andS =1.80 for the realistic energy bands E (relative E ) us- Ba Bi O kj F a.u. All calculations are performed at the experimen- ing experimental structures. We conclude that the nest- tal lattice constant a=8.10 a.u. The Barth–Hedin–like ing is far from perfect in the case of half–filling and 3 dimerisationof the oxygenoctahedracan hardly be con- phase. According to the neutron diffraction data20, the nected with it. A realistic TB model should also in- average structure is cubic although presence of a weak clude Bi(6p) orbitals and their nearest–neighbour inter- long–range superstructure characterised by the octahe- action with O(2p) states45,28. Upon potassium doping, dra rotations at the angles about 3◦ was also found48,49. the bands hardly change except for a slight lowering of Recent XAFS measurements50 report on the locally dis- E awayfromhalf–filling. ForBa K BiO the Fermi ordered rotations. From their analysis, the authors F 0.6 0.4 3 surfacerepresentsaroundedcube centredatthe Γ point of Ref. 50 conclude that the rotations can either be as shown in Fig. 2. Analysis of the band structure fac- along (1,1,1) or (1,1,0) axe. Previous frozen–phonon tor given by Eq. (1) as a function of q shows featureless calculations16 performed for the ordered Ba K BiO 0.5 0.5 3 behaviour and any effect of the nesting enhancement on investigate (1,0,0) component of the tilting mode which the electron–phonon interaction is not expected for this is found to be unstable with the total energy minimum band dispersion. corresponding to the angle 7◦. Ourowntotal–energycalculationsalsoconfirmtheex- istence of tilting distortions. Fig. 3 shows that the total energy exhibits a double well behaviour as a function of the rotation angle. We choose (1,1,0) axe for the tilting as is the case in the undoped compound. The unit cell in the calculation is doubled according to the R point of the cubic phase. The total energy minimum is found at the angle equal to 5◦. The energy gain compared to the cubicphaseisonly10meV/(1×cell)whichindicatesthat at the temperatures of the order T the rotations can be c dynamic. The double well behaviour at such small en- ergyscaleunambiguouslypointsoutontheimportanceof evaluating anharmonicity contribution in total electron– phonon coupling. This problem will be discussed in the following section. l 30 l e c 1 / 20 FIG. 2. Calculated Fermi surface for cubic Ba0.6K0.4 V BiO3. using the LMTO method. The effect of potassium e dopingis takenintoaccount within thevirtualcrystalap- m proximation. The centre of the cube corresponds to the Γ 10 point of theBrillouin zone. , ε y ω 5 g 4 ε ω We second discuss our results for the calculated equi- r 0 ε4 e ω3 3 librium lattice configuration in Ba K BiO . The theoretical–to–experimental volume0.6rat0i.o4 V/3Vexp is En ω 20 ε~ε 1 1 2 found to be 1.01, and the calculated bulk modulus is -10 equal to 1.25 Mbar. Both neutron diffraction21 and x– ray–absorption–fine–structure47 (XAFS) measurements -10 -5 0 5 10 show that frozen–in breathing distortions are absent in the superconducting phase. We have performed frozen– Tilting angle, degrees phononcalculationsforthedoubledcellcorrespondingto the R–point and for several breathing distortions. The FIG. 3. Frozen–phonon calculation of the total en- total energy minimum shows that the undistorted cubic ergy (meV/1 cell) as a function of the tilting angle in phase is stable in accord with these experiments. The × Ba0.6K0.4BiO3. The levels ǫn are the solutions of the curvature is well fitted with standard parabola, which Schr¨odinger equation for the anharmonic oscillator with showsthatthe breathingmode isharmonicinthe super- thedouble–potentialwellshownonthefigure. Thetransi- conducting phase. tionsωn =ǫn ǫ0 involvingdifferentphononexcitedstates − We further investigate tilting of the octahedra. Ex- are illustrated byarrows. perimentally, for the undoped compound the octahe- dra rotated21 at the angle ∼11.2◦ along (1,1,0) axe. We now report our main results on the calculated lat- Morecomplicatedsituationexistsinthesuperconducting tice dynamical properties of Ba K BiO . The den- 0.6 0.4 3 4 sity functional linear–response approach29 implemented frequency range from 0 to 10 THz. The high–frequency onthebasisofthefull–potentialLMTOmethod39isused modes mainly consist of the oxygen bond–stretching vi- in this calculation. The dynamical matrix is computed brations. The longitudinal branch at the point R cor- at 20 irreducible q points corresponding to the (6,6,6) responds to the breathing mode which in our calcula- reciprocal lattice grid of the cubic BZ. The effect of the tion has a frequency 15.7 THz. From the analysis of potassiumsubstitution onthe phonon spectrumis taken our polarisationvectors,we conclude that oxygenbond– intoaccountbyvirtualmassapproximation. TheLMTO bending vibrations dominate in the frequency interval basissetandothertechniquedetailshavebeendescribed between 6 and 10 THz. The octahedra tilting modes above. One more comment should be said on evaluating are at the low–frequency interval. They exhibit signifi- BZ integrals in the linear–response calculation. Here, cant softening near the q–point R=(1,1,1)π/a. Due to one can essentially improve the accuracy of the integra- symmetry, one can talk about pure tilting at the line tion by using a multigrid technique29. A (6,6,6) grid (20 between q–point M=(1,1,0)π/a and the R point. Ex- irreducible k–points) is used for finding linear–response actlyatthe M pointnearestoctahedratilt in–phaseand functions while the effects of the energy bands and the they tilt out–of–phase at the R–point. A nearly–zero– Fermi surface are taken into account using a (30,30,30) frequencytriple–degeneratemodeexistsatq=(1,1,1)π/a grid (816 irreducible k–points). whichcorrespondsto thepurerotationalT phonon. In 2u The calculated phonon spectrum along major symme- fact, for T=0 this mode should have a slightly imagi- try directions of the cubic BZ is plotted in Fig. 4 . Solid nary frequency for the cubic structure according to our circles denote the calculated points and the lines result frozen–phonon analysis illustrated in Fig. 3 . But, due from interpolation between the circles. Several features tonumericalinaccuraciesthelinear–responsecalculation can be seen from these phonon dispersions. Three high– givesverysmallpositiveω=0.5THz. Nosignificantsoft- frequency optical branches around the ω ∼17 THz are ening of the tilting modes near the point M is predicted well separated from the other modes distributed in the by our calculation. λ (q):0.64 0.41 0.25 0.25 0.25 0.44 0.44 0.17 0.00 0.11 0.46 0.62 0.16 0.21 0.25 λep(q):0.01 0.27 0.13 0.18 0.20 0.51 0.60 0.04 0.00 0.06 0.77 1.02 0.20 0.20 0.13 20 tr 0.06 18 16 0.05 0.170.15 0.14 0.15 0.34 0.24 0.05 0.23 0.12 0.130.13 0.32 z H 14 T , 12 y c n ue 10 0.130.01 0.01 0.01 0.05 0.07 0.01 q 0.03 0.06 0.02 0.06 re 8 0.43 0.04*20.09 F 0.06 0.03 0.09*2 6 0.07 0.03 0.05 4 2 0 Γ (ξ00) X (0.5ξ0) M (ξξ0) Γ (ξξξ) R (0.5ξξ) X DOS2, .5st./THz FIG. 4. Calculated phonon spectrum of Ba0.6K0.4BiO3 using density–functional linear–response method. The potassium dopingis taken into account using virtual crystal and virtualmass approximation. The calculated points areshownbysymbols. Thelinesresultfrominterpolationbetweenthepoints. Horizontalbarsindicatethemeasured49 phononfrequencies. Numbersforeveryphononmodeindicatethecalculatedelectron–phononcouplingconstantsλqν. (Onlythevalueslargerthan0.01areemphasised.) Ontopofthefigureshownare(i)thevaluesofλsummedoverall branchesfor given q,(ii) thevaluesof q–dependenttransport constant λtr. The calculated phonondensityof states F(ω) is shown on the right. 5 The phonondispersioncurvesalongΓX, ΓM,and ΓR the energies ∼70meV (or 17 THz) in undoped BaBiO . 3 symmetry directions for Ba K BiO have been very Thereforeitistemptingtoconnectpossiblesourceforthe 0.6 0.4 3 recently measured by inelastic neutron scattering48,49. discrepancies with our poor treatment of doping. The Horizontal lines in Fig. 4 indicate measured phonon fre- authors of Ref. 11 discuss doping induced appearance of quencies at the symmetry points Γ, X,M, andRe as we localised holes on the oxygen 2p orbitals which screen were able to deduce them from Fig. 1 of Ref. 49. The the charge on the oxygen anions. This charge reduction existence ofsoft rotationalmodes nearR canbe directly willlowertheenergyofthesemodes. Ifthelocalisedhole seenfromthemeasuredphonondispersions. Theauthors picture is correct, the VCA will not capture this since it of Ref. 49 have reported that their samples still have a removes electrons from the conduction band by(uniform weeklong–rangesuperstructurecharacterisedbythetilt- distributingtheholesbetweenO(2p)andBi(6s)orbitals. ing of the octahedra. An extremely low frequency ∼0.9 We now report our results for the calculatedelectron– THz of these modes was measured. This is in agreement phonon interaction. Based on our screened potentials with our calculations. whichareinducedbynucleidisplacementsandarefound Twoothercommentsshouldbesaidonthecomparison self–consistently, we evaluate matrix elements of the between our theory and the experiment. One comment electron–phonon interaction, gqν . The standard k+qj′kj concerns frequency interval from 0 to 10 THz. Here, our expression51 for the electron–phonon matrix elements calculationisseentoreproducemeasuredphonondisper- reads as sions with the accuracy of the order 10%. In particular, the lowest mode in Γ has ωcalc=3.79 THz which can be gqν = k+qj′ QR(qαν) δ+ V kj (2) vcoomlvpesarBeda(Kw)ithanωdexBpi ∼vi3b.r5atTioHnzs.. TThheisnmexotdmeomdaeininlyΓinis- k+qj′kj * (cid:12)(cid:12)XRα 2MRωqν Rα (cid:12)(cid:12) + (cid:12) (cid:12) the oxygen out–of–phase mode. Here, ωcalc=4.88 THz whereQ(qν) aretheorth(cid:12)(cid:12)onorpmalisedpolarisat(cid:12)(cid:12)ionvectors Rα and the measured frequency is less than 6 THz. The associatedwiththemodeqν,M aretheatomicmasses, R so–calledferroelectric mode has a frequency 5.59THz in R runs over basis atoms in the unit cell and α runs over our calculation which is close to ωexp found near 6 THz. directions x,y,z; δ+ V denotes self–consistent change in Rα This mode is bond–bending longitudinal and it has the the potential associated with the q–wave displacements strongestpolar character. Usually it exhibits large split- ofatomsRalongαaxe. Inpracticalcalculationswehave ting from the TO mode at Γ in cubic perovskites. The also added so–called incomplete basis set corrections to presence of free chargecarriersscreensCoulomb interac- the matrix elements (2) according to the method devel- tion at long distances, and therefore, the LO–TO split- oped in Ref. 30 ting is absent in our calculation. The dispersion of the The coupling strength λ for the electrons with the qν ferroelectric mode as a function of q is also seen to be phonon of wave vector q and branch ν is given by the correctly reproduced. following integral The second comment concerns our comparisonfor the hlaitgiho–nfsreaqrueefnocuyndinttoerbvael,lewsshaecrceutrhaeteraesnudlttshoefotvheeracllaldcius-- λqν = Ns2(0)kjj′δ(Ekj)δ(Ek+qj′)|gkq+νqj′kj|2 (3) crepancy consists about 20%. The highest mode at Γ is X Bi–Obond–stretchingmode. Here,wereportthevalueof whereN (0)isthedensityofstatesatE =0percelland s F ω equalto17.91THzandtheω valueonlyslightlylarger per one spin. Indexes j and j′ numerate the bands (not than 15 THz was found experimentally. The authors of spins) and spin degenerate case is assumed throughout Ref. 48 discuss an anomalous dispersion for the longi- thepaper. Thetotalcouplingconstantλresultsbysum- tudinal optical branch of the one–dimensional breathing ming λ over ν and by averaging over BZ. Two delta qν modealongΓX withitspronouncedfrequencyrenormal- functions in (3) impose integration over the space curve isation. Our calculation, on the other hand, gives much resulting from the crossing of two Fermi surfaces sepa- lessdispersiveopticalbranchesalongthisdirectionascan rated by q. For this integral we have used as many as be seen from Fig. 4. It is not clear whether this result 816 k–points per irreducible BZ. is due to inhomogeneity of the potassium distribution or The calculated values of λ at the symmetry direc- qν otherimperfectnessofthesamplesusedintheexperiment tions of the BZ are indicated in Fig. 4 along with the or due to drawbacks in our calculation connected with calculated phonon dispersions. (We emphasise only the the virtual crystal approximation. In fact, it is clearly values larger than 0.01.) On top of the figure shown seen that all our high–frequency branches are overesti- are the values of λ which are summed over allbranches q matedby∼20%incontrasttotheexperimentalones(ex- for given q. It is seen that the electron–phonon cou- cept eventually the breathing vibrations near the point plingislargeforthehigh–frequencybond–stretchinglon- R). This result also follows from the comparison of our gitudinal branch. This result is expectable from band calculated and the measured11 phonon density of state structure calculations45 since bond–stretching and espe- F(ω)[see Fig.5(a)]. It was found experimentally11 that cially breathing vibrations produce modest changes in theoxygenbond–stretchingmodesexhibitsofteningwith the energy bands near E . Coupling is strongly en- F the substitution of Ba by K. These modes are located at hanced near the points M and R, where it reaches the 6 values ∼0.3. Here, the mode corresponds to two– or On the basis of our evaluated phonon dispersion and three–dimensional breathing. The value of λ =0.3 for q–dependent electron–phonon interaction we calculate b the breathing mode at R is in accord with the previous the Eliashberg spectral function α2F(ω). This is plot- frozen–phonon calculation15. The authors of Ref. 16, on ted in Fig. 5(b) by full lines. There, we also show the other hand, give much lower value for λ equal to by symbols two α2F(ω) which were deduced from the b 0.04. For other bond–stretching vibrations, we find λ tunnelling measurements9. Comparing the experiment qν of the order 0.1. with our calculations, it is first seen that the intensities From Fig. 4 we conclude that the electron–phonon of high–energy peaks are approximately the same which coupling is not small for the bond–bending oxygen means that we reproducethe coupling for these phonons modes. It is seen that λ is enhanced for the wave vec- sufficiently accurate. It is also seen both from the the- tors near R and also along ΓX direction. In the latter oryandthe experimentthatwhileinthephonondensity case, the vibrations correspondto the ferroelectricmode of states [Fig. 5(a)] mainly the TO phonons contribute and the value of λ as high as 0.43 is found for the q– to the high–energy structure, for the α2F(ω) these are qν point (1/3,0,0)π/a. Strongly anharmonic tilting modes, the LO phonon modes. Our tendency to overestimate on the other hand, do not exhibit noticeable electron– the phonon frequencies at high energies is again clearly phonon coupling in the linear order with respect to the distinguishable. displacements. ExactlyatthepointR,thesemodeshave The most prominent feature seen from our calculated electron–phonon matrix element equal to zero by sym- α2F(ω) is the absence of any structure for the low– metry, and,therefore,smallvalues ofω donot bringany frequency interval below 40 meV. This strongly con- effect on enhancing the coupling. We refer to the follow- tradicts with the experimental α2F(ω) which exhibits ing section on our evaluated anharmonicity corrections. two intensive peaks in this region centred at 15 and 30 meV. According to our analysis of partial F(ω), the ori- gin of the first peak could be due to low energy Ba(K) 1.2 and Bi vibrations together with the tilting modes, and TO V 1.0 the second peak can result from the bond–bending oxy- e modes m gen modes. It is not clear however why the calculation 0.8 / seriously underestimates the electron–phonon coupling . t s 0.6 for these phonons, while it correctly describes the cou- , ω) 0.4 pling for the bond–stretching modes. From the band F( structure arguments45 one can expect that only bond– 0.2 stretching modes will have large interaction with elec- 0.0 trons. The bond–bending modes cannot produce any 0.6 λ =0.29 LO significant changes in the bands near EF since sp(σ) in- 0.5 ep modes teraction for these kind of distortions is not changed in linearorder. Thesameistrueforthetiltingmodes. Since ω) 0.4 thereisnopartialweightofthe Ba(K)orbitalsatE we ( F F 0.3 also do not expect strong electron–phonon coupling for 2 α 0.2 the low–frequency interval. 0.1 Another possible explanation for the observed peaks 0.0 is due to the contributions connected with anharmonic- ity. For the anharmonic phonons, not only one–phonon 0.6 λ =0.32 virtual transitions must be seen in the α2F(ω), but also tr 0.5 higher order virtual states. The detailed discussion on ) ω 0.4 this subject will be presented in the next section, here ( F 2 0.3 we only give the value 0.04 as our final answer for the αtr anharmonic contribution to λ resulting from the tilting 0.2 motions. Whilenotnegligible,thisvaluealoneagaindoes 0.1 not explain the intensity of the experimental α2F(ω) at 0.0 low–energies. 0 10 20 30 40 50 60 70 80 Our calculatedtotal value of λ resulting fromthe har- Frequency, meV monic phonons and linear electron–phonon coupling is foundtobe 0.29. This is too small to account for the su- FIG. 5. Results for doped BaBiO . (a) Comparison 3 between calculated and experimental11(symbols) phonon perconductivity at 30K in the compound Ba0.6K0.4BiO3. In fact, only the high–energy phonons contribute to our density of states. (b) Calculated Eliashberg spec- tral function α2F(ω) and the results of the tunnelling coupling. Asaresult,ourestimatedvalueofωlog ashigh measurements9 (symbols). (c) Calculated transport spec- as550Kisfound. UsingMcMillan’sTc expression52 with tral function α2trF(ω). µ∗=0 we however find the critical temperature with our 7 set of parameters as low as 4.5K. One can try to esti- response calculations give nearly zero–frequency vibra- mate the error in our λ value due to the overestimation tions forthis mode attheR–pointof cubic BZ.Our own of the frequencies for the bond–stretching modes. Using (andprevious16)frozen–phononcalculationspredicthere the expression52: λ≃NI2/Mω¯2, where NI2 is the elec- a double–potential–wellbehaviour with very shallow en- tronic prefactor and Mω¯2 is an average force constant, ergy minimum at 5◦ as illustrated in Fig. 3. Despite the one sees that λ should increase with lowering ω¯2. Our reportedaveragestructureforthesuperconductingphase 20% error in ω¯ results in 30% error in ω¯2, and this can is cubic21, some experiments48 discuss the existence of lead to the actual λ values which are 30% higher than long–rangesuperstructurecharacterisedbythe rotations we calculate. However, λ ∼0.4 is also not sufficient to of the octahedra. Recent XAFS measurements50 report explain the value of T . on the presence of locally disordered rotations. c Our calculated electron–phonon contribution to the transportpropertiesisafinalsubjectofthissection. The Unfortunately, though formulated54, the problem of quantitiesresponsiblefortheelectronictransportareeas- the influence of anharmonicity to superconductivity is ily deduced from the linear–response calculations30. By not tractable numerically in a full volume. Our sim- inserting the electron velocity factor (vkj −vk+qj′)2 to plified treatment is based on the expression introduced the expression (3), we calculate transport constant λ . by Hui and Allen38 which generalises zero–temperature tr Itsq–dependenceisshownatthetopofFig. 4alongwith electron-phonon coupling to the anharmonic case by in- theq–dependenceoftheelectron–phononλ. Weseethat cluding matrix elements over all phonon excited states. both functions exhibit very similar behaviour in the BZ. The phononstates|niandtheir energiesǫn areobtained The transport spectral function α2 F(ω) is the central by solving the Schr¨odinger equation for an oscillatory tr quantity needed for evaluating temperature–dependent mode characterised by q and ν. (We will not label the electricalandthermalresistivityaslow–ordervariational states|ni,ǫn with(qν) forsimplicity). Forharmonicpo- solutions of the Boltzmann equation53. Our calculated tentialwellsǫn−ǫ0 is just nωqν (in atomicunits), where α2 F(ω) is shown in Fig. 5(c). It is seen that this ωqν isthephononfrequency. Thisleadsonlytothe one– tr function behaves closely to the superconducting α2F(ω) phonon virtual states (n=1) which are involved in the which is usually the case in metals30. The total average matrix elements of the electron–phonon interaction. For λ isfoundtobe0.32. Basedonthesedata,weevaluate anharmonicpotentialwellsthe spectrumisgenerallydif- tr electron–phononlimitedelectricalresistivityρtobe14.3 ferent from the set of equidistant levels. One example is µΩ×cmatT=273K.Thisisatleastoneorderofmagni- a double well of the tilting mode which is shown in Fig. tudelowerthanthevaluesofρreportedintheliterature8. 3 The real spectrum ǫn obtained as the solution of the Itisclearthatthesourceforthisdiscrepancyisthesame Schr¨odinger equation is plotted in Fig. 3 by horizontal asinourdescribingsuperconductingproperties. Itisun- lines. Therefore, our first purpose is to examine what likelythatthestrongelectroniccorrelationsarepresented effect in our λ would bring the proper treatment of all inthesematerialsbecausetheparentcompoundsaredia- virtual phonon states. The second problem is connected magnetic (not antiferromagnetic) insulators. Therefore, withthe modificationofthe electron–phononmatrixele- it is unlikely that other (than electron–phonon) scatter- ments due to higher–orderterms in the expansionof the ing mechanisms take place, as spin fluctuations, for ex- change in the ground state potential with respect to the ample, in HTSC cuprates. Taking into account anhar- atomic displacements. Since we wish to examine these monic phonons, polarons, or bipolarons may be decisive effects only for the tilting mode at the point R of the for describing the superconductivity and transport phe- cubic BZ, we make an essential approximation by as- nomena here. Our basic conclusion is that the conven- suming that the tilting mode is not coupled to the other tional ideas on the electron–phonon mechanism are not modes of either this wave vector or other wave vectors operative in the bismuthates. which is generally not true when anharmonic terms are included into the lattice dynamicalproblem. We assume III. TILTING AND ANHARMONIC λ. thatthepolarisationvectorsforthismodeareknownand are given by exactly the out–of–phase rotations of the Towards further understanding of superconductivity nearest octahedra along the (1,1,0) axe. Any processes in the bismuthates, we try to evaluate anharmonicity of phonon–phonon interactions will be neglected in this corrections in the electron–phonon interaction. As we treatment. We alsoneglectby the finite–temperature ef- have mentioned in the introduction, several experiments fectsusingtheargumentsgiveninRef.37Inthedouble– and numerous frozen–phonon calculations point on pos- well problem one expects that the modifications due to sible importance of these effects in the theory of HTSC frequency and electron–phonon–matrix–elements renor- cuprates33–37. In the bismuthates, the first candidate malisation are quite dramatic35,37 and bring the largest to study anharmonicity is the tilting mode: Our linear– effect in the values of λ. We start from a general zero–temperature expression for the electron–phonon coupling in anharmonic case, which can be written as follows38 8 1 2 λqν = N (0) (fkj −fk′j′)δ(Ekj −Ek′j′ +ωn) G[kn′]j′kj /[ωn]2, (4) s kX′j′kjXn (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) where ω =ǫ −ǫ . For the moment we will not assume as in Ref. 38 that ω are small at the electron energy scale. n n 0 n (The latter reduces the integral with the Fermi step functions fkj −fk′j′ to the integral with the delta functions as givenbyEq. 3.) WeintroduceageneralisedmatrixelementG[n] forthe virtualtransitiontothen’sphononstate: k′j′kj G[kn′]j′kj =hn|∆gk′j′kj|0i. (5) The electron–phonon matrix element ∆gk′j′kj in (5) involves the transitions between the states |kji and |k′j′i near the Fermi surface ∆gk′j′kj =hk′j′|∆V|kji, (6) where∆V isthetotalchangeinthegroundstatepotentialinducedbythelatticedistortionassociatedwiththemode qν. We haveespecially included ∆ in the notation∆gk′j′kj since ∆V is not a derivative ofthe potentialwith respect to the displacements but it is the difference between the self–consistent potential V(r,{t +∆t }) for the distorted R R crystalandthe potentialV(r,{t })forthe undistortedcrystal. The atomic positionsatthe equilibriumaregivenby R the vectors t = t+R, where t denote the translations and R are the basis vectors. The displacements associated R with the mode qν are describedby the vectorfield ∆tR. Therefore ∆V is proportionalto ∆tR and so does ∆gk′j′kj. By introducing complex (infinitesimal) polarisation vectors δQ of the mode55 qν, the displacements in any atomic R cell t can be found using the formula: ∆t =δQ eiqt+c.c. (7) Rα Rα where α runs over directions x,y,z, and c.c. stands for the complex conjugate. (The quantities ∆V, ∆t , and δQ R R should, in principle, be labelled with qν but we omit this for simplicity.) The phonon states |ni are the functions of the displacements ∆t or δQ . In order to compute the matrix Rα Rα element hn|∆gk′j′kj|0i over the phonon states we should expand ∆V with respect to the displacements. Keeping the terms up to second order, this expansion reads as δV ∆V = δQ eiqt + Rα δt Rα Rα t X X +12 δQRαδQR′α′ eiq(t+t′)δt δ(2δ)tV′ + RR′αα′ tt′ Rα R′α′ X X +12 δQRα(δQR′α′)∗ eiq(t−t′)δt δ(2δ)tV′ +c.c. (8) RR′αα′ tt′ Rα R′α′ X X Here, δV/δt is associated with the first–order derivative of the potential when a single nucleus centred at t+R Rα experiences an infinitesimal displacement along α–th direction, and δ(2)V/δt δt′ stands for the second–order Rα R′α′ derivative. [Notationt is shorthandfor (t+R) .] Obviously,both these responsefunctions haveno dependence on Rα α the mode qν. If V(r) has a periodicity of the original lattice, the change δV/δt is a function of general type. One Rα expects that δV/δt is only not zero in the vicinity of the displaced atom and it goes to zero when r departs from Rα the site t+R. However, because of the translational invariance of the original crystal, considering the response at the point r due to the movement of atom in t+R must be equivalent to considering the response at the point r−t due to the movement of atom at R (when t=0). Therefore we can write that δV(r)/δt =δV(r−t)/δR . We now Rα α introduce the lattice sum δV δ+ V = eiqt (9) Rα δt Rα t X which represent a variation of the potential per unit displacement induced by to the movements of atoms R along α–th axeby infinitesimalamountδt proportionalto exp(iqt)in everycellt. Itis easyto provethatthe expression Rα (9) translates like a Bloch wave with wave vector q in the original lattice, i.e. δ+V(r+R) = eiqRδ+V(r). (We will sometimes omit indexes Rα) Notation δ+V refers to the travelling wave of vector +q while complex conjugated quantity δ−V would refer to the travelling wave of vector −q. 9 Onecananalogouslydefine latticesums associatedwiththe second–orderchangesofthe potential. Theseenterthe secondandthirdcontributionsin(9). Consider,forexample,thelatticesumassociatedwiththesecondcontribution: δ+ δ+ V = eiq(t+t′) δ(2)V . (10) Rα R′α′ δt δt′ tt′ Rα R′α′ X This expression translates like a Bloch wave of vector 2q because eiq(t+t′)δ(2)V(r+t′′) = δt δt′ tt′ Rα R′α′ X eiq(t+t′) δ(2)V(r) = tt′ δ(t−t′′)Rαδ(t′−t′′)R′α′ X e2iqt′′ eiq(t+t′) δ(2)V(r) . (11) δt δt′ tt′ Rα R′α′ X Analogously, the lattice sum associated with the third contribution in (9) can be denoted as δ+δ−V. It represents a travelling wave of wave vector 0, i.e. it is periodical at original lattice. Using the notations (9), and (10), the change in the potential ∆V given by the formula (8) now has the form ∆V = δQ ×δ+ V + Rα Rα Rα X 1 +2 δQRαδQR′α′ ×δR+αδR+′α′V + RR′αα′ X 1 +2 δQRα(δQR′α′)∗×δR+αδR−′α′V +c.c. (12) RR′αα′ X It is clear that when this expansion is used in the matrix element (6), certain selection rule will occur for the wave vectors k′ and k. Namely, the matrix element hk′j′|δ+V|kji is equal to zero unless k′ = k+q, hk′j′|δ+δ+V|kji=0 unless k′ =k+2q, and hk′j′|δ+δ−V|kji=0 unless k′ =k. Let us now introduce the electron–phonon matrix element associated with the first–order change in the potential G[kn′]j{′1k}j =δk′k+q k+qj′ hn|δQRα|0i×δR+αV kj . (13) * (cid:12) (cid:12) + (cid:12)XRα (cid:12) (cid:12) (cid:12) The electron–phonon matrix elements associated with(cid:12) the second–order cha(cid:12)nges in the potential have two forms (cid:12) (cid:12) according to the second and third contributions in (12) 1 G[kn′]j{′2k}j = 2δk′k+2q*k+2qj′(cid:12)(cid:12)RRX′αα′hn|δQRαδQR′α′|0i×δR+αδR+′α′V(cid:12)(cid:12)kj+, (14) (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) Gk[n′]j{′2k′j} = 12δk′k*kj′(cid:12)(cid:12)RRX′αα′hn|δQRα(δQR′α′)∗|0i×δR+αδR−′α′V(cid:12)(cid:12)kj+. (15) (cid:12) (cid:12) Then, the expression (4) for λ splits into(cid:12)three contributions (cid:12) qν (cid:12) (cid:12) λ =λ{1}+λ{2}+λ{2′} (16) qν qν qν qν associatedwithoneelectron–phononmatrixelementfromthe firstorder,Eq. (13),andtwomatrixelements fromthe second–order,Eqs. (14), and (15), i.e. 2 2 λ{q1ν} = N (0) (fkj −fk+qj′)δ(Ekj −Ek+qj′ +ωn) G[kn+]{q1j}′kj /[ωn]2, (17) s kXjj′Xn (cid:12) (cid:12) (cid:12) (cid:12) λ{q2ν} = N2(0) (fkj −fk+2qj′)δ(Ekj −Ek+2qj′ +ωn(cid:12)) G[kn+]{22q}j(cid:12)′kj 2/[ωn]2, (18) s kXjj′Xn (cid:12) (cid:12) (cid:12) (cid:12) λ{q2ν′} = N2(0) (fkj −fkj′)δ(Ekj −Ekj′ +ωn) Gk[nj]′{k2j′(cid:12)} 2/[ωn]2.(cid:12) (19) s kXjj′Xn (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) (cid:12) 10