Positron surface state as a spectroscopic probe for characterizing surfaces of topological insulator materials
Vincent Callewaert, K. Shastry, Rolando Saniz, Ilja Makkonen, Bernardo Barbiellini, Badih A. Assaf, Donald Heiman, Jagadeesh S. Moodera, Bart Partoens, Arun Bansil, A. H. Weiss
PPositron surface state as a new spectroscopic probe for characterizing surfaces oftopological insulator materials
Vincent Callewaert, ∗ K. Shastry, Rolando Saniz, Ilja Makkonen, Bernardo Barbiellini, BadihA. Assaf,
4, 5
Donald Heiman, Jagadeesh S. Moodera,
6, 7
Bart Partoens, Arun Bansil, and A. H. Weiss Department of Physics, Universiteit Antwerpen, Antwerpen 2020, Belgium Department of Physics, University of Texas at Arlington, Arlington, Texas 76019, USA Department of Applied Physics, Aalto University School of Science,P.O. Box 15100, FI-00076 Aalto, Espoo, Finland Department of Physics, Northeastern University, Boston, Massachusetts 02115, USA D´epartement de physique, Ecole Normale Sup´erieure, CNRS,PSL Research University, 24 rue Lhomond, 75005 Paris, France Department of Physics, MIT, Cambridge, Massachusetts 02139, USA Francis Bitter Magnet Laboratory, MIT, Cambridge, Massachusetts, 02139, USA (Dated: October 16, 2018)Topological insulators are attracting considerable interest due to their potential for technologicalapplications and as platforms for exploring wide-ranging fundamental science questions. In order toexploit, fine-tune, control and manipulate the topological surface states, spectroscopic tools whichcan effectively probe their properties are of key importance. Here, we demonstrate that positronsprovide a sensitive probe for topological states, and that the associated annihilation spectrum pro-vides a new technique for characterizing these states. Firm experimental evidence for the existenceof a positron surface state near Bi Te Se with a binding energy of E b = 2 . ± . I. INTRODUCTION
Quickly after their initial discovery, topological insu-lators (TIs) were recognized to hold significant potentialfor new technological applications and as playground forfundamental physics . An intrinsic challenge with TIs,which arises due to the fact that their interesting prop-erties originate from Dirac states located in a nanoscopiclayer near the surface, remains to separate the fingerprintof the topological surface states from the bulk behaviourof the sample. Highly surface sensitive techniques suchas angle resolved photoemission spectroscopy and scan-ning tunnelling microscopy have thus proven to be anindispensable tool to establish the existence of the gap-less states in several systems and to confirm various ofthe predicted quasi-particle properties .In this article, we demonstrate that positrons provide ahighly surface sensitive probe for the topological Diracstates. Since positron annihilation spectroscopy (PAS)techniques, with measurements of the 2D angular correla-tion of the annihilation radiation (2D-ACAR) in particu-lar, are well suited to measure both the low and high mo-mentum components of the annihilating electronic stateswithout complication of matrix element effects, they canprovide useful information on the Dirac state orbitals.Our calculations show that spin-polarized positron beamscan additionally resolve the spin-textures associated withthe topological states, owing to the predominant annihi-lation between particles with opposite spins .In section II, we present the experimental evidence forthe existence of a bound positron state at the surface ofthe TI Bi Te Se and the measured binding energy . Sec- tion III contains a discussion of the theory and compu-tational details used in our first principles investigation.In section IV, we show that the theory confirms the ex-perimental interpretation and predicts a significant over-lap between the positron and the topological states. Wealso demonstrate that spin-polarized positron measure-ments can reveal the spin-structure at the surface. Insection V we summarize the results and discuss possibleapplications and advantages of PAS over other spectro-scopic techniques. II. EXPERIMENTAL RESULTS
Our Bi Te Se films are grown by molecular beam epi-taxy on Si (111). The substrates are etched in hydroflu-oric acid prior to loading in vacuum. A stoichiometric2:2:1 Bi:Te:Se flux ratio is used. The substrate temper-ature is fixed at 200 ◦ C during the growth. The filmsused in this study are typically 40 nm thick. A 100 nmSe cap is then deposited, in-situ, on the sample surfaceafter cooling down the substrate to room temperature.The capping layer protects the film surface from oxida-tion and atmospheric contaminants.X-ray diffraction (XRD) is systematically used to char-acterize the samples, as briefly discussed in ref. 5. Thec-axis lattice constant for the film used in this work isfound to be equal to 30 . ± .
03 ˚A. Energy dispersiveX-ray spectroscopy confirmed stoichiometry within a 5%error on samples resulting from an identical growth.The samples are then transferred to the experimentalpositron chamber. In order to decap the samples, the a r X i v : . [ c ond - m a t . m t r l - s c i ] S e p protective Se layer is evaporated under UHV conditions,prior to the positron annihilation experiment. A heaterbutton is placed behind the sample in a holder and asuitable current was passed to heat the sample for 20minutes at 200 ◦ C. This procedure is similar to the de-capping sequence used in ref. 6. The technical detailsconcerning the setup of the positron experiments can befound in ref. 4 and references therein.Positrons annihilate predominantly with the valence elec-trons but the small fraction that annihilates with coreelectrons produces highly unstable core holes, which arefilled by the Auger process. Therefore, if positrons anni-hilate in a surface state (SS), positron-induced Auger-electron spectroscopy (PAES) provides a particularlyclean method to determine the composition of the sur-face, free from a secondary electron background . Aschematic picture of the process is drawn in figure 1(a).Results of PAES experiments from the TI Bi Te Se sur-face are shown in figure 1(b), where signals from Bi, Te,Se, C and O can be identified; the latter two are caused bythe presence of a small concentration of contaminants ad-sorbed on the surface . These results reveal the presenceof a bound positron SS. Were this not the case, positronswould either get trapped between the blocks of quintuplelayers (QL) of the material or would be re-emitted be-fore they annihilate. Since the extent of one QL block isabout 10 ˚A, which corresponds roughly to the mean freepath of a 60 eV electron, any Auger signal coming frombelow the first QL is too weak to be detected. Thus, thefact that the annihilation induced Auger peak intensitiesare observable is clear evidence that the positron is in astate localized at the surface at the time it annihilates.Auger Mediated Positron Sticking (AMPS) experimentsprovide an independent proof for the existence of the SSand allow us to determine its binding energy . In theAMPS mechanism, the excess energy from a positrondropping into the image potential well is transferred toa valence electron. This can result in the emission ofan Auger electron if the energy difference between thepositron SS and the initial state, determined by the in-cident positron’s kinetic energy, is larger than the elec-tron workfunction . The maximum kinetic energy of theAuger electrons is then given by E max = E p + E b − φ − ,where E p is the energy of the incident positron, E b isthe binding energy of the positron surface state, and φ − is the electron workfunction. Figure 2(a) illustrates theAMPS mechanism schematically. The observed increasein amplitude of the Auger signal at low energies as theenergy of the incident positrons is increased, is shown infigure 2(b), and it confirms the presence of the SS. Know-ing the electron workfunction, the binding energy of theSS can be determined from the positron energy thresholdvalue for Auger electron emission: E T h = E p for which E max = 0. The linear fit shown in figure 2(c), yields E T h = 1 . E a = 0 . , one can eliminate the electron work-function using the expression , E a = E b + φ − − .
80 eV, which gives a binding energy of E b = 2 . ± . III. THEORY AND COMPUTATIONALDETAILS
Our first-principles calculations are carried outin the zero-positron-density limit of the two-component electron-positron density functional theory(2CDFT) . In this limit, which is exact in the caseof a delocalized positron in a perfect crystal or at asurface, the electron density remains unperturbed bythe presence of the positron. The computations thusconsists of an electronic and positronic groundstatecalculation which are performed subsequently.
A. Electronic structure
The electronic ground state is obtained using theprojector-augmented wave (PAW) method as imple-mented in the VASP software package . Electronexchange-correlation effects are treated using thePerdew-Burke-Ernzerhof (PBE) functional , and spin-orbit coupling is included in the computations. Thekinetic energy cutoff for the plane-wave expansion of thewavefunctions is set to 275 eV. For the bulk calculations,we use the rhombohedral unit cell with a Γ-centered11 × × k -grid in combination with a Gaussiansmearing of width 0 . k -gridwith 11 × × . . B. Positron state
The effective potential for the positron in the zero-density limit of the 2CDFT is determined by theCoulomb interaction with the nuclei, the Hartree interac-tion with the electron density and the electron-positroncorrelation potential. The latter is usually described withlocal density (LDA) or generalized gradient approxima-tions (GGA), which give reliable results for bulk systems.A fundamental limitation of these semi-local approxima-tions is that they always describe the formation of Ps − in the limit of a dilute electron gas. In the case of asurface, however, the correct limit is given by the imagepotential − / z − z ), where z denotes the distanceto the surface and z represents the image potential ref-erence plane. We impose this limit in the vacuum regionby considering the corrugated mirror model , in whichthe image potential is constructed to follow the same iso-surfaces as the electron density. In the vacuum region z > z , we take the least negative of the LDA potential and the image potential. The strength of the image po-tential is given by : V im ( r ) = − z eff ( n − ( r )) − z ) , (1)where n − ( r ) is the electron density and the effective dis-tance to the surface is determined by: z eff ( n − ( r )) = ∞ (cid:90) z dz (cid:48) z (cid:48) δ ( n − ( r ) − (cid:104) n − (cid:105) ( z (cid:48) )) . (2)Here, (cid:104) n − (cid:105) is the electron density averaged over theplanes parallel to the surface and δ denotes the Diracdelta function. We approximate the image potential ref-erence plane z by the background edge position, whichis determined by the position outside the surface wherethe electron density starts decaying exponentially.We used the MIKA/doppler package to obtain thepositron ground state. These calculations are performedin an all-electron way in the sense that a superposition offree atomic core quantities, e.g. density and Hartree po-tential, are added to the self-consistent valence electronproperties. The Kohn-Sham equations for the positronare solved on a real space grid using a Rayleigh multigridimplementation . C. Electron-positron momentum density
The goal of the present paper is to investigate whetherPAS can be used to measure the properties of the TI’sDirac states. We thus need to calculate the electron-positron momentum density, which contains informationabout a sample’s electronic structure, and determine if itcontains a clear fingerprint of the topological states.To the best of our knowledge, electron-positron momen-tum density calculations in which the electronic wave-functions are not collinear, have not been studied in lit-erature. Hence, we present in some detail a generalizationof the theory to the non-collinear case.Spin-polarized positron annihilation measurements ex-ploit the fact that the two gamma annihilation only oc-curs for electron-positron pairs in a singlet state. If onespecifies the initial spin of the positron, this translatesto saying that the positron will only annihilate with elec-trons of the opposite spin. The magnetization of theelectron-positron momentum density along a specifiedaxis can thus be obtained by taking the difference be-tween spectra obtained by aligning the positrons paralleland anti-parallel to that axis. As long as the electronand positron spins are good quantum numbers, i.e. theyare position independent, the effect of the spin is easilytaken into account by realizing that the positron will bein a singlet state with exactly half of the electron states with the opposite spin. In systems where the spin can-not be considered a good quantum number, however, amore careful examination is required. In general, we canwrite the momentum density of the annihilating electron-positron pairs as : ρ ( p ) = 4 πr e c (cid:88) j g j (cid:88) s e ,s p (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (cid:90) d r e − i p · r ˆ S s α j ( r , s e ; r , s p ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) (3)where | α j (cid:105) are the natural geminals which diagonalizethe reduced two-body density matrix, sometimes also re-ferred to as electron-positron pairing wavefunctions, andthe g j are their occupation numbers. The spin of theelectron and positron in the geminal are denoted by s e and s p , respectively, and j represents a set of quantumnumbers (excluding the spin of the particles). The fac-tor 4 πr e c , with r e the classical electron radius and c thespeed of light, is the annihilation rate constant . Theoperator ˆ S s = 1 − ˆ S , where ˆ S is the total spin opera-tor for the electron-positron pair, projects on the singletstate. For the purpose of notation as well as practicalcalculations, it is convenient to define: A j,s e ,s p ( p ) = (cid:90) d r e − i p · r α j ( r , s e ; r , s p ) (4)as well as the matrix:Γ j ( p ) = (cid:18) | A j, ↑↓ ( p ) | A j, ↑↓ ( p ) A ∗ j, ↓↑ ( p ) A j, ↓↑ ( p ) A ∗ j, ↑↓ ( p ) | A j, ↓↑ ( p ) | (cid:19) . (5)In measurements with unpolarized positron beams, thepositron has statistically a 50% chance to be either in thespin-up or spin-down state. In this case, upon evaluationof eq. (3), the off-diagonal terms of Γ j ( p ) drop since thegeminals with opposite spin orientations, e.g. α j ( r , ↑ ; r , ↓ )and α j ( r , ↓ ; r , ↑ ), are not simultaneously occupied. Theresult for the momentum density then becomes: ρ ( p ) = πr e c (cid:88) j g j Tr[Γ j ( p )] , (6)where Tr[ . . . ] denotes taking the trace. In case thepositron beam is perfectly polarized parallel or anti-parallel to the z -axis, we obtain: ρ ↑ z ( p ) = 2 πr e c (cid:88) j g j | A j, ↑↓ ( p ) | , (7)and: ρ ↓ z ( p ) = 2 πr e c (cid:88) j g j | A j, ↓↑ ( p ) | , (8)respectively. The magnetization along the z -axis is ob-tained by taking the difference between these two spectra,and can conveniently be written as: ρ z ( p ) = 2 πr e c (cid:88) j g j Tr[ σ z Γ j ( p )] , (9)where σ z denotes the Pauli matrix. Analogous obser-vations can be made for a positron polarized along thedifferent axes, thus we can write in general: ρ i ( p ) = 2 πr e c (cid:88) j g j Tr[ σ i Γ j ( p )] , (10)where i = { x, y, z } and the σ i are the Pauli matrices. Adetailed derivation of the above formulas can be found inthe appendix.In electron-positron momentum density calculationsbased on the 2CDFT, one assumes that the natural gem-inals can be written in terms of a product of the electronand positron single particle Kohn-Sham orbitals ψ − j,s e and ψ + s p , where the positron is assumed to reside in its ground-state, and the occupation numbers of the electronic or-bitals replace those of the natural geminals g j . Electron-positron correlation effects are included by introducing amultiplicative term γ , i.e. the enhancement factor, whichcan be state and/or space dependent. We thus have: α j ( r , s e ; r , s p ) = (cid:113) γ j,s e ,s p ( r ) ψ − j,s e ( r ) ψ + s p ( r ) . (11)Note that, in general, it is justified to consider thepositron wavefunction to be collinear even though theelectronic states are not. Indeed, electron-positron spin-spin interactions are small and generally neglected in PASstudies and positrons stay too far away from the nucleito experience any significant spin-orbit interaction. Wethus assume that the orbital part of the positron wave-function is independent of the chosen spin-polarization: ψ + s p ( r ) = ψ + ( r ) χ s p , where χ s p denotes a two-componentspinor for the positron. Note that for the calculation ofthe momentum density from eqs. (6) and (10), we haveto set ψ + ↑ ( r ) = ψ + ↓ ( r ) instead of explicitly setting a po-larization.In our calculations, we consider the state-dependent en-hancement factors : γ j,s e ,s p = λ LDAj,s e ,s p /λ IP Mj,s e ,s p . The λ ’sdenote the partial annihilation rates in the LDA and in-dependent particle model (IPM), respectively, and theformer is calculated as: λ LDAj,s e ,s p = πr e c (cid:90) d r | ψ − j,s e ( r ) | | ψ + s p ( r ) | γ ( n − ( r )) , (12)with γ ( n − ( r )) the LDA enhancement factorparametrized by Drummond . The IPM annihila-tion rates are obtained by setting γ ( n − ( r )) = 1.The high-momentum components of the wavefunctionsare important to accurately calculate the electron-positron momentum density. It is thus necessary to usethe all-electron wavefunctions in the above formulae in-stead of the soft pseudo wavefunctions, i.e. we explicitlyperform the PAW transformation : | ψ − (cid:105) = | ˜ ψ − (cid:105) + (cid:88) i (cid:16) | φ − i (cid:105) − | ˜ φ − i (cid:105) (cid:17) (cid:104) ˜ p i | ˜ ψ − (cid:105) . (13)Here, | ˜ ψ − j (cid:105) are the soft pseudowavefunctions, (cid:104) ˜ p i | arethe projectors and | φ − i (cid:105) and | ˜ φ − i (cid:105) are the localized all-electron and soft pseudo partial waves of the ionsrespectively. The details on how we performed thistransformation can be found in refs. 21 and 24. D. Positronium model
We can theoretically determine the activation energyfor Ps desorption from a Bi Te Se, of which the experi-mental results are described in ref. 4, by calculating theparticle’s binding energy to the surface. In order to modelthe Ps state, we consider the Schr¨odinger equation for aneutral particle in an effective potential well . Here, theeffective potential outside the surface is determined byan attractive and a repulsive contribution. The repulsivecontribution, due to the overlap of the electron of the Pswith electrons of the material, is given by V R ( z ) = | φ P s | e − ( z − z ) /λ , (14)where φ P s is the Ps workfunction, z the backgroundedge position and λ the characteristic length of the elec-tron density decay outside the surface. The Ps workfunc-tion can be calculated by taking the sum of the work-functions of the constituent particles minus their bindingenergy: φ P s = φ + + φ − − .
25 Ha. The attractive part ofthe interaction is given by the Van der Waals interactionand can be written as V vdW ( z ) = − C ( z − z (cid:48) ) F (( z − z (cid:48) ) /λ ) , (15)where the strength of the interaction is given by the ex-pression : C = ¯ h π ∞ (cid:90) dω α ( iω ) (cid:18) (cid:15) ( iω ) − (cid:15) ( iω ) + 1 (cid:19) . (16)The bulk dielectric function (cid:15) at imaginary frequenciescan be obtained by first evaluating the dielectric functionat real frequencies, which is readily calculated from first-principles in the RPA approximation, and then applyinganalytic continuation. The Ps polarizability α can beobtained from the analytic expression for H-like atoms,given in ref. 31, by rescaling. Indeed, the Ps problemcan be solved by going to the center of mass coordinates,which then yield the same equations as for the H atom.The only differences are that the Bohr radius is twice aslarge and the ionization energy is half the value of that ofH. The analytic damping function F , for which we takeexpression (17) of ref. 32, describes the saturation of theVan der Waals interaction as the particle draws closer tothe surface and regularizes the divergence at the refer-ence plane position z = z (cid:48) . The reference plane positioncan in principle take another value than the backgroundedge position but since they are both, in the case of anelementary metal with lattice parameter a , located closeto a/
2, we make the approximation z (cid:48) = z . For z < z ,we extend the repulsive interaction, and add V vdW ( z ) toensure the continuity of the potential, with a cutoff setby the Ps workfunction: V ( z ) = min { φ P s , V R ( z ) + V vdW ( z ) } Θ( z < z )+ { V R ( z ) + V vdW ( z ) } Θ( z ≥ z ) . (17)The different contributions to the potential are show infigure 3. The Ps state and its energy are obtained bysolving the resulting Schr¨odinger equation − ψ (cid:48)(cid:48) V ( z ) ψ = Eψ. (18)
IV. COMPUTATIONAL RESULTS
We start our discussion of the computations by show-ing that the measured Ps activation energy E a = 0 . is consistent with the theoretical predictions. We takethe activation energy to be equal to the groundstate en-ergy predicted by the Ps model discussed in the previ-ous section. For the parameters in the model, we findthat the Van der Waals interaction strength evaluates to C = 2 .
306 eV · ˚A and from the electronic and positronicworkfunctions φ − = 4 .
904 eV and φ + = 2 .
392 eV, we ob-tain φ P s = 0 .
493 eV. The values for the background edgeposition and the characteristic length of the electron den-sity decay in the vacuum region are given by z = 1 .
250 ˚Aand λ = 0 .
365 ˚A. Using these values, the model predictsthat the Ps forms a delocalized state in the bulk of thematerial. We note, though, that the experimental valuefor the electronic workfunction φ − = 4 . φ − ∈ [4 . , .
9] eV. Over therange φ − = 4 .
90 eV to φ − = 4 .
72 eV, we find that thegroundstate gradually shifts from the bulk to the sur-face. To determine when we have a surface state, we setthe criterion that the Ps density should decay below 1%of its maximum value beyond the first QL block insidethe material. In the range φ − ∈ [4 . , .
72] eV, the Psmodel predicts a surface state with a binding energy of E P s = 0 . ± .
05 eV, in good agreement with the ex-perimental results.Next, we investigate the predictions of the 2CDFT calcu-lations to determine whether they support the proposedinterpretation of the PAES and AMPS experiments. Ourfirst observation is that the positron in its groundstate in-deed resides in the surface’s image potential well ratherthan the gaps in between the QLs, which also act asstrong positron traps. We obtain the binding energy ofthe positron by taking the difference between the vacuumlevel and the positron’s chemical potential. The vacuumlevel is determined in the usual way by the taking thevalue of the Hartree potential in middle of the vacuumregion. We find that the positron SS has a binding en-ergy of E b = 2 .
69 eV, in excellent agreement with themeasured value. We find that the lifetime evaluates to τ = 309 .
25 ps. This value seems reasonable compared with the lifetime of 340 −
380 ps measured for positronstrapped at the surface of colloidal PbSe quantum dots .On the other hand, a lifetime of 580 ps has been deter-mined for positrons trapped at an Al surface , whichcan not be reproduced within the LDA approximation .One workaround suggested in literature is to set the en-hancement factor to zero for z > z , i.e. assume that thepositron will not annihilate in the vacuum region . Wefind, though, that this operation makes the result for thelifetime depend sensitively on the value for the image po-tential reference plane z . For this reason, as well as thescarcity of experimental data that show this operation isjustified, the rest of our calculations have been carriedout without modifying the LDA enhancement factor.Now that the calculations confirmed the existence of thebound positron SS, we turn to the important question ofthe extent to which this SS overlaps with the Dirac coneelectrons. This overlap is of central importance because itdetermines the annihilation rate of the positron with theelectrons occupying the topological states and thus thesensitivity with which positron annihilation spectroscopycan probe the Dirac states. This can be seen from eq.(12), where the partial annihilation rate is determinedby the sum over λ j where j denotes the states on thecone.The computed densities of the positron SS, ρ + , and thetopological Dirac states, ρ − Dirac , are shown in figure 4.The density of the topological states is obtained by sum-ming the one-particle densities for all states on the conebetween the Dirac point and a specific value for the elec-tron chemical potential µ − . Although the positron is seento probe only the topmost atomic layers of the material,it still penetrates the material sufficiently to have a sig-nificant overlap with the Dirac states. Moreover, the leftpanel of figure 4 shows that the overlap with the Diracstates changes sensitively depending on the population ofthe Dirac states near the Fermi-level. Our calculations ofthe momentum density, discussed below, further demon-strate that this underlying overlap translates into a clearsignal coming from the annihilation of the positron withthe Dirac fermions.A partially filled energy band when it crosses the Fermienergy gives rise to a break in the electron momentumdensity, which is the basis of the measurement of Fermisurfaces in materials via 2D-ACAR experiments. A stan-dard procedure for enhancing the Fermi surface signal inthe spectrum is the Lock-Crisp-West (LCW) map ob-tained by folding all the higher momentum (Umklapp)contributions into the first Brillouin zone . Figure 5shows the calculated LCW map together with a cut alongΓ − M over a range of values of the electron chemicalpotential, which simulates different doping levels of theDirac cone. The evolution of the plateau around theΓ-point clearly indicates the sensitivity of the positronto the Dirac cone states. The relative drop in intensitybetween 5% −
7% at the Fermi momentum compares fa-vorably with, for example, the 1% drop found for theNd − x Ce x CuO − δ high T c superconductor in which 2D-ACAR experiments have been shown previously to beviable in detecting Fermi surface sheets due to Cu-Oplanes .A topic which has drawn considerable interest in the caseof topological insulators, is the spin-momentum lock-ing of the topological states. Measurements using spin-polarised positron beams exploit the fact that a two pho-ton decay is only possible between electrons and positronswith opposite spins . In recent work, spin-effects in theelectronic structure of simple ferromagnets were observedusing differences between the doppler broadening of theannihilation radiation (DBAR) measured with positronaligned parallel and anti-parallel to a polarizing mag-netic field. . In a similar ACAR experiment, Weber etal. successfully resolved the spin-dependent Fermi sur-face of the ferromagnetic Heusler compound Cu MnAl.This motivates us to investigate whether spin-polarisedpositrons can be used to detect the spin-structure of thetopological states at the surface. The signal from theFermi-surface can be extracted from the LCW map bytaking the difference between the signal obtained at dif-ferent doping levels. In figure 6, we show the results ob-tained by taking the difference between the LCW mapsobtained with µ − = E F + 0 . µ − = E F in thevicinity of the Γ-point. As expected, we see the plateaudue to the extra occupation of the cone in the total am-plitude. Our results for the magnetization along the x -and y -directions, agree well with the results obtained inseveral studies of various tetradymite TIs , which allpredict a clockwise rotation of the spin. We see thatthe z -component of the magnetization increases gradu-ally away from the Γ-point. This out of plane compo-nent develops due to the hexagonal warping of the Diraccone, as pointed out by Fu . We note that the differencein amplitude for the magnetic components is quite pro-nounced w.r.t. to the Fermi-surface signal. Indeed, wefind that the signal from the magnetization about halfthat of the Fermi-surface signal obtainable with an unpo-larized beam. This means that the magnetization signalstill constitutes a promising 2% −
4% of the total sig-nal. We note, though, that in real experiments, positronbeams are not perfectly polarized, as we have assumed inour calculations. Thus, in experiment, a proper weight-ing has to be performed which will lead to a smaller sig-nal.
V. CONCLUSION AND OUTLOOK
Our study establishes the existence of a positronsurface state near the topological insulator Bi Te Se.The results of our calculations show that this surfacestate can be exploited as a spectroscopic characterizationtool for probing surfaces of topological materials. Since asignificant fraction of positrons annihilate with electronsoccupying Dirac cone states, 2D-ACAR experimentsshould be able to measure their momentum distributionwith high precision , and thus obtain information concerning the nature of the Dirac states which iscomplementary to that accessed through angle-resolvedphotoemission, scanning tunnelling and other surface-sensitive spectroscopies without complications of relatedmatrix element effects . PAES and Doppler broadeningof the annihilation radiation measurements can, inturn, be used to characterize the chemical compositionof surfaces. In combination with 2D-ACAR experi-ments, these positron spectroscopies could be exploitedto determine effects of various surface impurities onthe topological states, in addition to the role of bulkdefects . Now our study identified a positron surfacestate, positron spectroscopies can prove valuable for thecharacterization of nano-structured topological insula-tors. Indeed, positrons have shown to act as effectiveself-seeking probes for nano-crystal surfaces withoutrequiring the preparation of single crystal specimens ,whereas the applicability of conventional spectroscopictechniques is limited. Finally, our calculations show thatthe spin-textures of the Dirac states should be accessiblethrough 2D-ACAR measurements using a spin-polarizedpositron beam since positrons predominantly annihilatewith electrons of the opposite spin . VI. ACKNOWLEDGEMENTS
I. M. acknowledges discussions with M. Ervasti andA. Harju.V. C. and R. S. were supported by the FWO-Vlaanderenthrough Project No. G. 0224.14N. The computationalresources and services used in this work were in partprovided by the VSC (Flemish Supercomputer Center)and the HPC infrastructure of the University of Antwerp(CalcUA), both funded by the Hercules Foundation andthe Flemish Government (EWI Department). I. M. ac-knowledges financial support from the Academy of Fin-land (projects 285809 and 293932). The work at North-eastern University was supported by the US Depart-ment of Energy (DOE), Office of Science, Basic EnergySciences grant number DE-FG02-07ER46352, and bene-fited from Northeastern University’s Advanced ScientificComputation Center (ASCC) and the NERSC supercom-puting center through DOE grant number DE-AC02-05CH11231. K. S. and A. W. acknowledge financialsupport from the National Science Foundation throughgrants DMR-MRI-1338130 and DMR-1508719. D. H. re-ceived financial support of the National Science Founda-tion (grant ECCS-1402738). J. S. M. was supported bythe STC Center for Integrated Quantum Materials underNSF grant DMR-1231319, NSF DMR-1207469 and ONRN00014-13-1-0301. B. A. A. also acknowledges supportfrom the LabEx ENS-ICFP: ANR-10-LABX-0010/ANR-10-IDEX-0001-02 PSL.
20 30 50 100 200 500
Outgoing electron energy (eV) . . . . . . . . C o un t s p e r s ec o nd
200 150 125 100 90 80
Time of flight (ns) vac ∆E∆E EE vac Annihilation Auger emission
Te Se Bi Bi/C O (a) (b)
FIG. 1. (a) Schematic representation of the PAES mecha-nism. In the first step, a positron (blue) annihilates with anelectron (red) occupying a core level and creates a highly un-stable hole. In the second step, an electron from a higher levelfills this hole and transfers the energy difference between thetwo levels to a second electron. If the energy difference is suf-ficiently large, and the second electron is close enough to thesurface, it can traverse the surface dipole and escape from thesample. The measured outgoing electron energy correspondswith the transferred energy in the Auger process minus theenergy difference between the second electron’s state and thevacuum level. (b) Results of the PAES measurements on theBi Te Se sample in which Auger signals from the differentelements are indicated.
VII. REFERENCES ∗ [email protected] M. Z. Hasan and C. L. Kane, Rev. Mod. Phys. , 3045(2010); A. Bansil, H. Lin, and T. Das, ibid . , 021004(2016); J. E. Moore, Nature , 194 (2010); A. R.Mellnik, J. S. Lee, A. Richardella, J. L. Grab, P. J.Mintun, M. H. Fischer, A. Vaezi, A. Manchon, E.-A. Kim,N. Samarth, and D. C. Ralph, ibid . , 449 (2014);V. Mourik, K. Zuo, S. M. Frolov, S. R. Plissard, E. P.A. M. Bakkers, and L. P. Kouwenhoven, Science ,1003 (2012). D. Hsieh, D. Qian, L. Wray, Y. Xia, Y. S. Hor, R. J.Cava, and M. Z. Hasan, Nature , 970 (2008); D. Hsieh,Y. Xia, D. Qian, L. Wray, J. H. Dil, F. Meier, J. Oster-walder, L. Patthey, J. G. Checkelsky, N. P. Ong, A. V. Fe-dorov, H. Lin, A. Bansil, D. Grauer, Y. S. Hor, R. J. Cava,and M. Z. Hasan, ibid . , 1101 (2009); P. Roushan,J. Seo, C. V. Parker, Y. S. Hor, D. Hsieh, D. Qian,A. Richardella, M. Z. Hasan, R. J. Cava, and A. Yaz-dani, ibid . , 1106 (2009); Y. Xia, D. Qian, D. Hsieh,L. Wray, A. Pal, H. Lin, A. Bansil, D. Grauer, Y. S. Hor,R. J. Cava, and M. Z. Hasan, Nat. Phys. , 398 (2009). S. Berko and J. Zuckerman, Phys. Rev. Lett. , 339(1964). K. Shastry, A. H. Weiss, B. Barbiellini, B. A. Assaf, Z. H.Lim, P. V. Joglekar, and D. Heiman, J. Phys.: Conf. Ser. , 012006 (2015). B. A. Assaf, T. Cardinal, P. Wei, F. Katmis, J. S. Moodera,and D. Heiman, Appl. Phys. Lett. , 012102 (2013). D. Zhang, A. Richardella, D. W. Rench, S.-Y. Xu, A. Kan-dala, T. C. Flanagan, H. Beidenkopf, A. L. Yeats, B. B. Buckley, P. V. Klimov, D. D. Awschalom, A. Yazdani,P. Schiffer, M. Z. Hasan, and N. Samarth, Phys. Rev.B , 205127 (2012). A. Weiss, R. Mayer, M. Jibaly, C. Lei, D. Mehl, and K. G.Lynn, Phys. Rev. Lett. , 2245 (1988). S. Mukherjee, M. P. Nadesalingam, P. Guagliardo, A. D.Sergeant, B. Barbiellini, J. F. Williams, N. G. Fazleev, andA. H. Weiss, Phys. Rev. Lett. , 247403 (2010). S. Chu, A. P. Mills, and C. A. Murray, Phys. Rev. B ,2060 (1981). B. Chakraborty and R. W. Siegel, Phys. Rev. B , 4535(1983). E. Boronski and R. M. Nieminen, Phys. Rev. B , 3820(1986). P. E. Bl¨ochl, Phys. Rev. B , 17953 (1994). G. Kresse and J. Furthm¨uller, Comput. Mater. Sci. , 15(1996). G. Kresse and J. Furthm¨uller, Phys. Rev. B , 11169(1996). G. Kresse and D. Joubert, Phys. Rev. B , 1758 (1999). J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev.Lett , 3865 (1996). The distance between the QL blocks is severely overesti-mated when using the PBE functional. As positrons arestrongly repelled by the ions, the separation between theQL strongly influences the value of the positron workfunc-tion and in order to obtain reliable results, we deem itappropriate to work with the experimental lattice param-eters instead. The lattice parameters only slightly affectthe electronic structure as the results of our bandstructure . . . . . . . . . Incident Positron Energy (eV) . . . . . I n t e g r a t e d I n t e n s i t y S E p e a k FitExperiment0 2 4 6 8 10
Outgoing electron energy (eV) . . . . . . . . . N o r m a li ze d i n t e n s i t y ( a r b . un i t s ) E p = 2 . E p = 3 . E p = 3 . E p = 4 . E p = 5 . E p = 6 . EE vac ф - -E b E p ΔESample VacuumPositronpotential (a)(b) (c)
FIG. 2. (a) Schematic representation of the AMPS mecha-nism. The left part of the diagram shows the incident positron(blue) that drops in the image potential well. In this process,the positron transfers an energy ∆ E , determined by the inci-dent kinetic energy E p and the binding energy of the SS E b ,to an electron of the system through a virtual photon, as indi-cated in the right part of the figure. If the energy difference islarger than the electronic workfunction φ − , the electron canescape to the vacuum. (b) The measured low-energy Augersignals for the Bi Te Se sample. The outgoing electron en-ergy is determined by the transferred energy ∆ E minus therequired energy to escape from the sample. The different linesshow the result for varying energies of the incident positrons.(c) The integrated peak amplitudes of the low-energy Augersignal associated with the AMPS mechanism as a function ofthe incident positron energy.calculations agree very well with the previously reportedfirst-principles results and those of ARPES mea-surements . N. D. Lang and W. Kohn, Phys. Rev. B , 3541 (1973). R. M. Nieminen and M. J. Puska, Phys. Rev. Lett. , 281(1983). We are updating the standard corrugated mirror modelfor the potential at the surface where GGAcorrections are traditionally not included. I. Makkonen, M. Hakala, and M. J. Puska, Phys. Rev. B , 035103 (2006). T. Torsti, M. Heiskanen, M. J. Puska, and R. M. Niemi-nen, Int. J. Quantum Chem. , 171 (2003). T. Torsti, T. Eirola, J. Enkovaara, T. Hakala, P. Havu,V. Havu, T. H¨oyn¨al¨anmaa, J. Ignatius, M. Lyly, I. Makko-nen, T. T. Rantala, J. Ruokolainen, K. Ruotsalainen,E. R¨as¨anen, H. Saarikoski, and M. J. Puska, Phys. StatusSolidi B , 1016 (2006). I. Makkonen, M. Hakala, and M. J. Puska, J. Phys. Chem.Solids , 1128 (2005). A. Zubiaga, M. M. Ervasti, I. Makkonen, A. Harju,F. Tuomisto, and M. J. Puska, J. Phys. B. , 064005(2016). R. A. Ferrell, Rev. Mod. Phys. , 308 (1956). Distance to topmost atomic layer ( ˚ A)-1 0 1 2 3 4 5 E n e r g y ( e V ) -1.5-1-0.50 φ Ps V vdW ( z ) V rep ( z ) V rep ( z ) + V vdW ( z ) V rep ( z ) + V vdW ( z ) FIG. 3. Potential obtained for the Ps model with the valuesfor C , z , λ and φ + mentioned in the discussion. For theelectronic workfunction, we took φ − = 4 .
612 eV, which givesa Ps workfunction of φ Ps = 0 . D i s t a n ce t o t o p m o s t a t o m i c l a y e r ( ˚ A ) -8-6-4-20246 Density (arb. units) ρ + ρ − Dirac : µ − = E F ρ − Dirac : µ − = E F + 0 . eV (a) (b) FIG. 4. Overlap of the positron SS with the Diracstates. (a) Planar average of the positron (blue) and electron(red/yellow) density associated with the Dirac states belowthe Fermi-energy for two different values of the chemical po-tential µ − . (b) Density of the topological surface state andthe positron in the same spatial region as panel (a). The pro-gressively lighter blue isosurfaces show the positron densityat 80%, 20% and 2% of the maximum value, respectively, andthe red isosurfaces show the electronic charge density asso-ciated with the electron states on the Dirac cone below theFermi-level at 10% of the maximum value. The Bi, Te and Seatoms are shown in purple, brown and green colors, respec-tively. µ − = E F µ − = E F + 0 . eV . . . . . . . Γ M . . . . . . . R e l a t i v e I n t e n s i t y µ − = E F µ − = E F + 50 meV µ − = E F + 100 meV µ − = E F + 150 meV µ − = E F + 200 meV − . . . . . E n e r g y ( e V ) (a) (b) (c) KMΓ
FIG. 5. Theoretical momentum densities. (a) LCW map with the chemical potential located at the Fermi-level. (b) LCWmap with the chemical potential raised by 0 . M direction for different values of thechemical potential. The inset shows the band structure near the Fermi-level ( E F = 0 . − . − . . . . k y ( ˚ A − ) ρ tot ρ x − . − .
05 0 .
00 0 .
05 0 . k x ( ˚A − ) − . − . . . . k y ( ˚ A − ) ρ y − . − .
05 0 .
00 0 .
05 0 . k x ( ˚A − ) ρ z . . . . . × − − . − . . . . × − − − × − − − − × − FIG. 6. Difference between the LCW maps obtained withdifferent doping levels of the Dirac cone: µ − = E F + 0 . µ − = E F . The top left pane of the figure show the to-tal amplitude of the LCW map. The top right, bottom leftand bottom right figures show the magnetization componentsalong the x -, y - and z -axes respectively. We only show the re-sult zoomed in around the Γ-point as the difference betweenthe LCW maps is exactly zero in the rest of the Brillouinzone. The inner- and othermost edges of the non-zero partin the plots correspond with the dashed lines shown in fig-ures 5(a) and 5(b), respectively. The length of the recipro-cal axes is | b | = 1 .
688 ˚A − and the amplitudes are given inps − ˚A . (It is readidly seen that the units of the LCW mapare in ps − ˚A by realizing that the integral over the LCWmap yields the positron’s annihilation rate, or, in the case ofthe magnetic LCW maps, the difference in annihilation ratebetween two measurements with opposite spin polarizationsfor the positron.) M. Alatalo, B. Barbiellini, M. Hakala, H. Kauppinen,T. Korhonen, M. J. Puska, K. Saarinen, P. Hautoj¨arvi,and R. M. Nieminen, Phys. Rev. B , 2397 (1996);B. Barbiellini, M. J. Puska, M. Alatalo, M. Hakala,A. Harju, T. Korhonen, S. Siljam¨aki, T. Torsti, and R. M.Nieminen, Appl. Surf. Sci. , 283 (1997). N. D. Drummond, P. L´opez R´ıos, R. J. Needs, and C. J.Pickard, Phys. Rev. Lett. , 207402 (2011). R. Saniz, B. Barbiellini, P. M. Platzman, and A. J. Free-man, Phys. Rev. Lett. , 096101 (2007). E. Zaremba and W. Kohn, Phys. Rev. B , 2270 (1976). R. Szmytkowski, Phys. Rev. A , 012503 (2001). S. H. Patil, K. T. Tang, and J. P. Toennies, J. Chem.Phys. , 8118 (2002). L. Chai, W. Al-Sawai, Y. Gao, A. J. Houtepen, P. E. Mi-jnarends, B. Barbiellini, H. Schut, L. C. van Schaarenburg,M. A. van Huis, L. Ravelli, W. Egger, S. Kaprzyk, A. Ban-sil, and S. W. H. Eijt, APL Materials , 022111 (2013). K. G. Lynn, W. E. Frieze, and P. J. Schultz, Phys. Rev.Lett. , 1137 (1984). R. M. Nieminen, M. J. Puska, and M. Manninen, Phys.Rev. Lett. , 1298 (1984). D. G. Lock, V. H. C. Crisp, and R. N. West, J. Phys. F:Metal Phys. , 561 (1973). B. Barbiellini, M. J. Puska, A. Harju, and R. M. Niem-inen, J. Phys. Chem. Solids , 1693 (1995); A. Shukla,B. Barbiellini, L. Hoffmann, A. A. Manuel, W. Sadowski,E. Walker, and M. Peter, Phys. Rev. B , 3613 (1996). A. Kawasuso, M. Maekawa, Y. Fukaya, A. Yabuuchi, andI. Mochizuki, Phys. Rev. B , 024417 (2012). J. A. Weber, A. Bauer, P. B¨oni, H. Ceeh, S. B. Dug-dale, D. Ernsting, W. Kreuzpaintner, M. Leitner, C. Pflei-derer, and C. Hugenschmidt, Phys. Rev. Lett. , 206404(2015). J. Henk, A. Ernst, S. V. Eremeev, E. V. Chulkov, I. V.Maznichenko, and I. Mertig, Phys. Rev. Lett. , 206801(2012). S. Basak, H. Lin, L. A. Wray, S.-Y. Xu, L. Fu, M. Z. Hasan,and A. Bansil, Phys. Rev. B , 121401 (2011). H. Lin, T. Das, L. A. Wray, S.-Y. Xu, M. Z. Hasan, andA. Bansil, New J. Phys. , 095005 (2011). L.-L. Wang and D. D. Johnson, Phys. Rev. B , 241309(2011). L. Fu, Phys. Rev. Lett. , 266801 (2009). S. B. Dugdale, J. Laverock, C. Utfeld, M. A. Alam, T. D.Haynes, D. Billington, and D. Ernsting, J. Phys.: Conf.Ser. , 012083 (2013). A. Bansil and M. Lindroos, Phys. Rev. Lett. , 5154(1999); J. C. Campuzano, L. C. Smedskjaer, R. Benedek,G. Jennings, and A. Bansil, Phys. Rev. B , 2788 (1991);J. Nieminen, H. Lin, R. S. Markiewicz, and A. Bansil,Phys. Rev. Lett. , 037001 (2009). F. Tuomisto and I. Makkonen, Rev. Mod. Phys. , 1583(2013). T. R. Devidas, E. P. Amaladass, S. Sharma, R. Rajaraman,D. Sornadurai, N. Subramanian, A. Mani, C. S. Sundar,and A. Bharathi, Europhys. Lett. , 67008 (2014). S. W. H. Eijt, A. Van Veen, H. Schut, P. E. Mijnarends,A. B. Denison, B. Barbiellini, and A. Bansil, Nat. Mat. ,23 (2006). X.-Q. Dai, B. Zhao, J.-H. Zhao, Y.-H. Li, Y.-N. Tang, andN. Li, J. Phys.: Condens. Matter , 035502 (2012). J. Chang, L. F. Register, S. K. Banerjee, and B. Sahu,Phys. Rev. B , 235108 (2011). M. Neupane, S.-Y. Xu, L. A. Wray, A. Petersen,R. Shankar, N. Alidoust, C. Liu, A. Fedorov, H. Ji, J. M.Allred, Y. S. Hor, T.-R. Chang, H.-T. Jeng, H. Lin, A. Ban-sil, R. J. Cava, and M. Z. Hasan, Phys. Rev. B , 235406(2012). T. Arakane, T. Sato, S. Souma, K. Kosaka, K. Nakayama,M. Komatsu, T. Takahashi, Z. Ren, K. Segawa, andY. Ando, Nat. Commun. , 636 (2012). N. G. Fazleev, J. L. Fry, K. H. Kuttler, A. R. Koymen,and A. H. Weiss, Phys. Rev. B , 5351 (1995). N. G. Fazleev, J. L. Fry, and A. H. Weiss, Phys. Rev. B , 165309 (2004). B. Barbiellini and J. Kuriplach, Phys. Rev. Lett. ,147401 (2015). G. G. Ryzhikh and J. Mitroy, J. Phys. B. , 4051 (1999). Appendix A: Detailed derivation spin-resolved momentum density
In this appendix, we give a detailed derivation of the momentum density formulas given in the main manuscript.We start from the general expression given in ref. 57, which defines the rate to start with the ground-state with N electron and a single positron Ψ and end up in the final state Φ, with N − p : Λ ν ( p ) = 4 πr e c N (cid:88) i =1 (cid:88) s i ,s p (cid:12)(cid:12)(cid:12)(cid:12) (cid:90) dτ (cid:90) d r i (cid:90) d r p Φ ∗ ν ( x , . . . , x i − , x i +1 , . . . , x N ) e − i p · r p × ˆ S si δ ( r i − r p ) Ψ( x , . . . , x N ; x p ) (cid:12)(cid:12)(cid:12)(cid:12) . (A1)To keep notation in check, we define x i = { r i , s i } to denote the particle position and spin, and dτ = (cid:81) Nj =1 ,j (cid:54) = i d x j to represent integration over all the non-annihilating electron coordinates, thus including a sum over the possiblespin directions of the particle. The contraint that only electron-positron pairs in a singlet state contribute to the 2 γ annihilation is taken into account through the operator ˆ S si = 1 − ˆ S /
2, with ˆ S i the total spin operator of electron i and the positron. As we will show further on, this operator projects the electron-positron pair on the singlet statesof the respective pairs.Thanks to the anti-symmetry of the wavefunction, we can swap the electron indices around such that the annihilatingelectron always has the label N . We make use of the delta function to get perform the integration over the positroncoordinate to get:Λ ν ( p ) = 4 πr e c N (cid:88) s N ,s p (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) dτ (cid:90) d r N Φ ∗ ν ( x , . . . , x N − ) e − i p · r N ˆ S sN Ψ( r , s , . . . , r N , s N ; r N , s p ) (cid:12)(cid:12)(cid:12)(cid:12) . (A2)If we are not concerned with the precise final state Φ ∗ ν we end up with, we can define the 2 γ transition rate ρ ( p ) = (cid:80) ν Λ ν ( p ) from Ψ to some final state with a pair of photons of total momentum p . The completeness relation: (cid:88) ν Φ ∗ ν ( x , . . . , x N − )Φ ν ( x (cid:48) , . . . , x (cid:48) N − ) = δ ( x − x (cid:48) ) . . . δ ( x N − − x (cid:48) N − ) , (A3)then allows us to write: ρ ( p ) = 4 πr e cN (cid:88) s e ,s p (cid:90) d x . . . (cid:90) d x N − (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) d r N e i p · x N ˆ S sN Ψ( x , . . . , r N , s N ; r N , s p ) (cid:12)(cid:12)(cid:12)(cid:12) . (A4)We now recognize the electron-positron two-body reduced density matrix, defined as :Γ ep ( r e , s e , r p , s p ; r (cid:48) e , s (cid:48) e , r (cid:48) p , s (cid:48) p ) = N (cid:90) d x . . . (cid:90) d x N − Ψ ∗ ( x , . . . , x N − , r e , s e ; r p , s p ) × Ψ( x , . . . , x N − , r (cid:48) e , s (cid:48) e ; r (cid:48) p , s (cid:48) p ) . (A5)It is convenient to introduce the natural geminals, also called pairing-wavefunctions, which diagonalize the abovedensity matrix: Γ ep ( r e , s e , r p , s p ; r (cid:48) e , s (cid:48) e , r (cid:48) p , s (cid:48) p ) = (cid:88) j g j α ∗ j ( r e , s e ; r p , s p ) α j ( r (cid:48) e , s (cid:48) e ; r (cid:48) p , s (cid:48) p ) . (A6)We thus arrive at the general expression for the 2 γ momentum density: ρ ( p ) = 4 πr e c (cid:88) j g j (cid:88) s e ,s p (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) d r e − i p · r ˆ S s α j ( r , s e ; r , s p ) (cid:12)(cid:12)(cid:12)(cid:12) , (A7)where we dropped the now unnecessary label on the singlet operator.Let us now examine the effect of the singlet projection operator. We have:ˆ S s = 1 − ˆ S
2= 1 − (cid:16) ˆ S s e + ˆ S s p + 2 ˆ S s e · ˆ S s p (cid:17) = 1 − (cid:16) ˆ S s e + ˆ S s p + 2 (cid:104) ˆ S s e ,x ˆ S s p ,x + ˆ S s e ,y ˆ S s p ,y + ˆ S s e ,z ˆ S s p ,z (cid:105)(cid:17) (A8)2If we make use of: ˆ S s e | s e , s p (cid:105) = ˆ S s p | s e , s p (cid:105) = 34 | s e , s p (cid:105) ˆ S s e ,x ˆ S s p ,x | s e , s p (cid:105) = | s e s p ||− s e , − s p (cid:105) = 14 |− s e , − s p (cid:105) ˆ S s e ,y ˆ S s p ,y | s e , s p (cid:105) = − s e s p |− s e , − s p (cid:105) ˆ S s e ,z ˆ S s p ,z | s e , s p (cid:105) = s e s p | s e , s p (cid:105) , (A9)then we find: ˆ S s | s e , s p (cid:105) = (cid:18) − s e s p (cid:19) ( | s e , s p (cid:105) − |− s e , − s p (cid:105) ) . (A10)It is obvious that this gives zero if the electron and positron spins have the same value and the prefactor becomes 1 / ± ( |↑↓(cid:105) − |↓↑(cid:105) ).After performing the sum over all electron and positron spins, we find for the momentum density expression: ρ ( p ) = 4 πr e c (cid:88) j g j (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) d r e − i p · r √ α j ( r , ↑ ; r , ↓ ) − α j ( r , ↓ ; r , ↑ )) (cid:12)(cid:12)(cid:12)(cid:12) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) (cid:104) α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) − α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) − α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) (cid:105) . (A11)In the rest of the derivation, it is more convenient to work with a spinor representation for the geminals, which wedefine as: α j ( r e ; r p ) = α j ( r e , ↑ ; r p , ↑ ) α j ( r e , ↑ ; r p , ↓ ) α j ( r e , ↓ ; r p , ↑ ) α j ( r e , ↓ ; r p , ↓ ) . (A12)Equation (A7) then becomes: ρ ( p ) = 4 πr e c (cid:88) j g j (cid:12)(cid:12)(cid:12)(cid:12) (cid:90) d r e − i p · r ˆ S s α j ( r ; r ) (cid:12)(cid:12)(cid:12)(cid:12) = 4 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) (cid:104) ˆ S s α † j ( r (cid:48) ; r (cid:48) ) (cid:105) (cid:104) ˆ S s α j ( r (cid:48) ; r (cid:48) ) (cid:105) . (A13)Indeed, making use of: ˆ S s α j ( r e ; r p ) = 12 α j ( r e , ↑ ; r p , ↓ ) − α j ( r e , ↓ ; r p , ↑ ) α j ( r e , ↓ ; r p , ↑ ) − α j ( r e , ↑ ; r p , ↓ )0 , (A14)it is straightforward to check that one obtains the same result asin equation (A11).Next, if we assume that the geminals are collinear in the positron spin, we can write: α j ( r e ; r p ) = (cid:18) α j ( r e , ↑ ; r p ) α j ( r e , ↓ ; r p ) (cid:19) ⊗ χ p , (A15)where ⊗ denotes a direct product, the remaining arrow indicates the electron spin, and χ p is the (position-independent)spinor for the positron. For a positron fully polarized along the positive and negative z -axis, respectively, we havethe single particle spinors: χ z + p = (cid:18) (cid:19) and χ z − p = (cid:18) (cid:19) . (A16)3For the geminal spinor this gives: α z + j ( r e ; r p ) = α j ( r e , ↑ ; r p , ↑ )0 α j ( r e , ↓ ; r p , ↑ )0 , and α z − j ( r e ; r p ) = α j ( r e , ↑ ; r p , ↓ )0 α j ( r e , ↓ ; r p , ↓ ) . (A17)and after applying the singlet operator to them:ˆ S s α z + j ( r e ; r p ) = 12 α j ( r e , ↓ ; r p , ↑ )0 , and ˆ S s α z − j ( r e ; r p ) = 12 α j ( r e , ↑ ; r p , ↓ )00 . (A18)From (A13) we obtain: ρ z + ( p ) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) ρ z − ( p ) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) . (A19)The magnetization is obtained as the difference between these two spectra and gives: ρ z ( p ) = ρ z − ( p ) − ρ z + ( p ) = 2 πr e c (cid:88) j g j (cid:34)(cid:12)(cid:12)(cid:12)(cid:12) (cid:90) d r e − i p · r α j ( r , ↑ ; r , ↓ ) (cid:12)(cid:12)(cid:12)(cid:12) − (cid:12)(cid:12)(cid:12)(cid:12) (cid:90) d r e − i p · r α j ( r , ↓ ; r , ↑ ) (cid:12)(cid:12)(cid:12)(cid:12) (cid:35) (A20)A positron polarized along the x -axis is represented by the single particle spinors: χ x + p = 1 √ (cid:18) (cid:19) and χ x − p = 1 √ (cid:18) − (cid:19) . (A21)Thus we get: α x + j ( r ; r ) = 1 √ α j ( r , ↑ ; r , ↑ ) α j ( r , ↑ ; r , ↓ ) α j ( r , ↓ ; r , ↑ ) α j ( r , ↓ ; r , ↓ ) , α x − j ( r ; r ) = 1 √ α j ( r , ↑ ; r , ↑ ) − α j ( r , ↑ ; r , ↓ ) α j ( r , ↓ ; r , ↑ ) − α j ( r , ↓ ; r , ↓ ) , (A22)and:ˆ S s α x + j ( r ; r ) = 12 √ α j ( r , ↑ ; r , ↓ ) − α j ( r , ↓ ; r , ↑ ) α j ( r , ↓ ; r , ↑ ) − α j ( r , ↑ ; r , ↓ )0 , ˆ S s α x − j ( r ; r ) = 12 √ − α j ( r , ↑ ; r , ↓ ) − α j ( r , ↓ ; r , ↑ ) α j ( r , ↓ ; r , ↑ ) + α j ( r , ↑ ; r , ↓ )0 , (A23)which result in the momentum densities: ρ x + ( p ) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) (cid:104) α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) − α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) − α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) (cid:105) .ρ x − ( p ) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) (cid:104) α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ )+ α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) (cid:105) . (A24)So the magnetization in this case is: ρ x ( p ) = ρ x − ( p ) − ρ x + ( p ) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) × (cid:104) α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) (cid:105) (A25)4Finally, a positron polarized along the y -axis is represented by the single particle spinors: χ y + p = 1 √ (cid:18) i (cid:19) and χ y − p = 1 √ (cid:18) − i (cid:19) . (A26)The geminal spinors become: α y + j ( r ; r ) = 1 √ α j ( r , ↑ ; r , ↑ ) iα j ( r , ↑ ; r , ↓ ) α j ( r , ↓ ; r , ↑ ) iα j ( r , ↓ ; r , ↓ ) , α y − j ( r ; r ) = 1 √ α j ( r , ↑ ; r , ↑ ) − iα j ( r , ↑ ; r , ↓ ) α j ( r , ↓ ; r , ↑ ) − iα j ( r , ↓ ; r , ↓ ) , (A27)ˆ S s α y + j ( r ; r ) = 12 √ iα j ( r , ↑ ; r , ↓ ) − α j ( r , ↓ ; r , ↑ ) α j ( r , ↓ ; r , ↑ ) − iα j ( r , ↑ ; r , ↓ )0 , ˆ S s α y − j ( r ; r ) = 12 √ − iα j ( r , ↑ ; r , ↓ ) − α j ( r , ↓ ; r , ↑ ) α j ( r , ↓ ; r , ↑ ) + iα j ( r , ↓ ; r , ↑ )0 , (A28)from which we find: ρ y + ( p ) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) (cid:104) α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ )+ iα j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) − iα j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) (cid:105) .ρ y − ( p ) = 2 πr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) (cid:104) α j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + α j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) − iα j ( r , ↓ ; r , ↑ ) α ∗ j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) + iα j ( r , ↑ ; r , ↓ ) α ∗ j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) (cid:105) . (A29)This gives the final component of the magnetization: ρ y ( p ) = ρ y − ( p ) − ρ y + ( p ) = 2 iπr e c (cid:88) j g j (cid:90) d r (cid:90) d r (cid:48) e − i p · ( r − r (cid:48) ) × (cid:104) α j ( r (cid:48) , ↑ ; r (cid:48) , ↓ ) α ∗ j ( r , ↓ ; r , ↑ ) − α j ( r (cid:48) , ↓ ; r (cid:48) , ↑ ) α ∗ j ( r , ↑ ; r , ↓ ) (cid:105) . (A30)If we introduce the notation: A j,s e ,s p ( p ) = (cid:90) d r e − i p · r α j ( r , s e ; r , s p ) (A31)and the matrix: Γ j ( p ) = (cid:18) | A j, ↑↓ ( p ) | A j, ↑↓ ( p ) A ∗ j, ↓↑ ( p ) A j, ↓↑ ( p ) A ∗ j, ↑↓ ( p ) | A j, ↓↑ ( p ) | (cid:19) , (A32)then the above results can be written as: ρ i ( p ) = 2 πr e c (cid:88) j g j Tr [ σ i Γ j ( p )] , (A33)where i = { x, y, z } and σ i are the Pauli matrices.Let us now derive the momentum density as measured in experiments with unpolarized positron beams. In this case,there is a statistical uncertainty on the direction of the positron spin. We can assume that the 50% of the positronsare in the up state and 50% in the down state, w.r.t. whatever direction of the quantization axis. This means wemeasure: ρ tot ( p ) = 12 ( ρ z + ( p ) + ρ z − ( p )) = πr e c (cid:88) j g j Tr [Γ j ( p )] ,,