text:
# Energy dependent counting statistics in diffusive superconducting tunnel junctions.

## I Model and theory

We consider a diffusive normal-insulating-normal-insulating-superconducting ($`NI_1N^{}I_2S`$) junction (see Fig. 1), where both the resistances of the insulating tunnel barriers ($`I_1`$ and $`I_2`$) are much larger than the resistance of the normal diffusive region ($`N^{}`$) between the barriers. The short diffusive region $`N^{}`$, which does not contribute to the junction resistance, leads to complete isotropization of the system parameters between the barriers. Nazcirc  This makes the junction effectively zero dimensional dotcom  and the exact shape and size of the diffusive region is thus not important.

The conductances $`g_1`$ and $`g_2`$ of the two tunnel barriers can be written as $`g_{1(2)}=(2e^2/h)_{n=1}^{N_{1(2)}}\mathrm{\Gamma }_{1(2),n}`$, where $`N_{1(2)}`$ is the number of transport modes at the barrier $`1(2)`$ and $`\mathrm{\Gamma }_{1(2),n}`$ is the transparency of mode $`n`$. We consider the limit of large number of modes, $`N_{1(2)}1`$ and small transparency of each mode, $`\mathrm{\Gamma }_{1(2),n}1`$. However, the conductances $`g_1`$ and $`g_2`$ are much larger than the conductance quanta $`2e^2/h`$, i.e. $`_{n=1}^{N_{1(2)}}\mathrm{\Gamma }_{1(2),n}1`$. Under these conditions we can neglect Coulomb blockade effects as well as weak localization like corrections to the cumulants of the current.

It is further assumed that the dwell time of the particles in the the normal region $`N^{}`$ is much shorter than the inelastic scattering- and phase breaking times. We are also only interested in zero frequency properties.spin  With these assumptions, we can analyze the full counting statistics of the junction within the framework of the circuit theory, recently developed by NazarovNazcirc ; Nazgenfcn  and Belzig and Nazarov.Belzig01a ; Belzig01b  Very recently, a detailed discussion of the circuit theory was given. Belzig02  We note that the zero energy counting statistics was studied for the same junction in Ref. \[Boerlin02, \], here we study the full energy dependence.

In the circuit theory, the mesoscopic conductor is described in terms of nodes connected via various circuit elements, just as in an ordinary electric circuit theory. Due to the effective zero-dimensionality of the junction, the circuit theory representation (see Fig. 1) only contains three nodes, the two reservoirs ($`N`$ and $`S`$) and the isotropic normal region ($`N^{}`$) between the tunnel barriers, and two circuit elements, the tunnel barriers ($`I_1`$ and $`I_2`$).

However, instead of characterizing each node with an electric potential, as is done in an ordinary electrical circuit theory, they are characterized with a $`4\times 4`$ matrix Green’s function (in Keldysh-Nambu space). Here we have the three node Green’s functions $`\stackrel{ˇ}{G}_N,\stackrel{ˇ}{G}_S`$ and $`\stackrel{ˇ}{G}`$, as shown in Fig. 1. They obey the normalization condition $`\stackrel{ˇ}{G}_N^2=\stackrel{ˇ}{G}_S^2=\stackrel{ˇ}{G}^2=\stackrel{ˇ}{1}`$. Throughout the paper we use the convention with “check”, $`\stackrel{ˇ}{.}`$ denoting $`4\times 4`$ matrices and “hat” $`\widehat{.}`$ denoting $`2\times 2`$ matrices.

The matrices of the reservoir nodes are known. In the normal reservoir we have

$`\stackrel{ˇ}{G}_N`$ $`=`$ $`e^{i\chi _N\stackrel{ˇ}{\tau }_K/2}\stackrel{ˇ}{G}_N^0e^{i\chi _N\stackrel{ˇ}{\tau }_K/2},\stackrel{ˇ}{\tau }_K=\left(\begin{array}{cc}0& \widehat{\sigma }_z\\ \widehat{\sigma }_z& 0\end{array}\right)`$ (3)
$`\stackrel{ˇ}{G}_N^0`$ $`=`$ $`\left(\begin{array}{cc}\widehat{\sigma }_z& \widehat{K}_N\\ 0& \widehat{\sigma }_z\end{array}\right),\widehat{K}_N=2\left(\begin{array}{cc}12f_+& 0\\ 0& 2f_{}1\end{array}\right),`$ (8)
where $`f_\pm (E)=f(E\pm eV)`$, with $`f(E)=(\text{exp}[E/kT]+1)^1`$ the Fermi distribution function. Here $`\chi _N`$ is a counting field, “counting” the electrons passing in and out of the normal reservoir and $`\widehat{\sigma }_z`$ is a Pauli matrix. In the superconducting reservoir we have

$`\stackrel{ˇ}{G}_S`$ $`=`$ $`e^{i\chi _S\stackrel{ˇ}{\tau }_K/2}\stackrel{ˇ}{G}_S^0e^{i\chi _S\stackrel{ˇ}{\tau }_K/2}`$
$`\stackrel{ˇ}{G}_S^0`$ $`=`$ $`\left(\begin{array}{cc}\widehat{R}_S& \widehat{K}_S\\ 0& \widehat{A}_S\end{array}\right),\widehat{R}(\widehat{A})_S=\left(\begin{array}{cc}g_{R(A)}& f_{R(A)}\\ f_{R(A)}& g_{R(A)}\end{array}\right),`$ (13)
$`\widehat{K}_S`$ $`=`$ $`(\widehat{R}_S\widehat{A}_S)[f(E)f(E)],`$ (14)
where $`\mathrm{\Delta }`$ is the superconducting gap and $`f_{R(A)}=i\mathrm{\Delta }[(E\pm iϵ)^2\mathrm{\Delta }^2]^{1/2}`$ and $`g_{R(A)}=[1f_{R(A)}^2]^{1/2}`$ are Green’s functions characterizing the superconducting reservoir ($`ϵ`$ is an infinitesimal positive number). The counting field $`\chi _S`$ of the superconducting reservoir is introduced for notational clarity, due to current conservation it is in principle not needed.

The matrix $`\stackrel{ˇ}{G}`$ of the node between the tunnel barriers is to be determined. It is found from a conservation law for matrix currents (a matrix Kirchoff’s rule), as

$`\stackrel{ˇ}{I}_1+\stackrel{ˇ}{I}_2+\stackrel{ˇ}{I}_E=0`$ (15)
where the matrix currents are given by

$`\stackrel{ˇ}{I}_1`$ $`=`$ $`{\displaystyle \frac{g_1}{2}}[\stackrel{ˇ}{G}_N,\stackrel{ˇ}{G}],\stackrel{ˇ}{I}_2={\displaystyle \frac{g_2}{2}}[\stackrel{ˇ}{G}_S,\stackrel{ˇ}{G}]`$
$`\stackrel{ˇ}{I}_E`$ $`=`$ $`{\displaystyle \frac{2iE}{\delta }}[\stackrel{ˇ}{G}_E,\stackrel{ˇ}{G}],\stackrel{ˇ}{G}_E=\left(\begin{array}{cc}\widehat{\sigma }_z& 0\\ 0& \widehat{\sigma }_z\end{array}\right)`$ (18)
where $`[A,B]=ABBA`$ and $`\delta `$ is the mean level spacing in the $`N^{}`$-region divided by the conductance quanta $`2e^2/h`$. The Thouless energy is given by $`E_{Th}=(g_1+g_2)\delta /4`$, i.e. it depends only on the escape time through the tunnel barriers, not on the diffusion time through the $`N^{}`$-region. The term $`\stackrel{ˇ}{I}_E`$ leads to a modification of the proximity effect at finite energies, and eventually a suppression for $`EE_{Th}`$.

From Eq. (15) and the condition $`\stackrel{ˇ}{G}^2=1`$, we can determine the missing node matrix $`\stackrel{ˇ}{G}(\chi _N,\chi _S)`$ as a function of the counting fields. Following the lines of Ref. \[Boerlin02, \], noting that Eq. (15) can be written in the form $`[\stackrel{ˇ}{G},\stackrel{ˇ}{G}_x]=0`$, we find the physically relevant solutionrelevant

$`\stackrel{ˇ}{G}`$ $`=`$ $`\stackrel{ˇ}{G}_x(\stackrel{ˇ}{G}_x^2)^{1/2},`$
$`\stackrel{ˇ}{G}_x`$ $`=`$ $`{\displaystyle \frac{g_1}{2}}\stackrel{ˇ}{G}_N+{\displaystyle \frac{g_2}{2}}\stackrel{ˇ}{G}_Si{\displaystyle \frac{E}{2E_{Th}}}(g_1+g_2)\stackrel{ˇ}{G}_E.`$ (19)
Knowing $`\stackrel{ˇ}{G}(\chi _N,\chi _S)`$, we can find the cumulant generating function $`F(\chi _N,\chi _S)`$. It is given from the general relation between the matrix current and the generating function,Nazgenfcn

$`{\displaystyle \frac{F(\chi _N,\chi _S)}{\chi _N}}`$ $`=`$ $`{\displaystyle \frac{it_0}{8e^2}}{\displaystyle 𝑑E\text{tr}\left[\stackrel{ˇ}{\tau }_K\stackrel{ˇ}{I}_1(\chi _N,\chi _S)\right]}`$
$`{\displaystyle \frac{F(\chi _N,\chi _S)}{\chi _S}}`$ $`=`$ $`{\displaystyle \frac{it_0}{8e^2}}{\displaystyle 𝑑E\text{tr}\left[\stackrel{ˇ}{\tau }_K\stackrel{ˇ}{I}_2(\chi _N,\chi _S)\right]}`$ (20)
where $`t_0`$ is the measurement time. Using that under the trace,

$`{\displaystyle \frac{(\stackrel{ˇ}{G}_x^2)^{1/2}}{\chi _{N(S)}}}`$ $`=`$ $`{\displaystyle \frac{\stackrel{ˇ}{G}_x}{\chi _{N(S)}}}\stackrel{ˇ}{G}_x(\stackrel{ˇ}{G}_x^2)^{1/2}={\displaystyle \frac{i}{4}}g_{1(2)}[\stackrel{ˇ}{\tau }_K,\stackrel{ˇ}{G}_{N(S)}]\stackrel{ˇ}{G}`$ (21)
$`=`$ $`{\displaystyle \frac{i}{4}}g_{1(2)}\stackrel{ˇ}{\tau }_K[\stackrel{ˇ}{G}_{N(S)},\stackrel{ˇ}{G}]={\displaystyle \frac{i}{2}}\stackrel{ˇ}{\tau }_KI_{1(2)},`$
we get the counting field dependent part of the cumulant generating function

$`F(\chi _N,\chi _S)`$ $`=`$ $`{\displaystyle \frac{t_0}{4e^2}}{\displaystyle 𝑑E\text{tr}\left[(\stackrel{ˇ}{G}_x^2)^{1/2}\right]}.`$ (22)
We note that the standard definition of a generating function demands $`F(\chi _N=0,\chi _S=0)=0`$. This implies that a counting field independent constant $`\text{tr}[G_x^2(\chi _N=0,\chi _S=0)]^{1/2}`$ should in principle be added to the integrand in Eq. (22). However, such a constant does not alter the cumulants, and below we will for simplicity omit it.

The probability $`P(Q)`$ that $`Q`$ electron charges have been transported through the junction in time $`t_0`$ is given by

$`P(Q)`$ $`=`$ $`{\displaystyle _\pi ^\pi }{\displaystyle \frac{d\chi _N}{2\pi }}e^{i\chi _NQF(\chi _N,\chi _S=0)}`$ (23)
Eqs. (8) - (23) completely determine the full charge counting statistics for arbitrary temperatures and voltages. For our purposes, it is however more convenient to focus on the properties of the cumulant generating function in Eq. (22). From the cumulant generating function we can derive all cumulants of the current by repeatedly taking the derivative with respect to the counting field $`\chi _N`$ (or $`\chi _S`$). For the two first cumulants, current and noise, we have (evaluated at the $`N`$-contact),

$`I`$ $`=`$ $`i{\displaystyle \frac{e}{t_0}}{\displaystyle \frac{F(\chi _N,\chi _S)}{\chi _N}}|_{\chi _N=\chi _S=0},`$ (24)
$`S`$ $`=`$ $`{\displaystyle \frac{2e^2}{t_0}}{\displaystyle \frac{^2F(\chi _N,\chi _S)}{\chi _N^2}}|_{\chi _N=\chi _S=0}.`$ (25)
and similarly for higher order cumulants.

## II Full counting statistics.

In what follows we focus on the limit $`kT,eV,E_{Th}\mathrm{\Delta }`$. In this limit, the physics is dominated by the induced proximity effect in the diffusive region $`N^{}`$ between the two barriers. proximity  The proximity effect affects the transport properties on the energy scale $`E_{Th}`$. It can on a microscopic level be described as resulting from coherence between the two quasiparticles in the Andreev reflection, e.g. an incident electron and a retroreflected hole. An induced proximity effect leads in general to large modifications of the transport properties, i.e. modifications scaling with the system size. We note that this is in strong contrast to the situation in normal systems, where phase coherence only gives rise to small quantum corrections. Buttikerrew

However, being a coherent phenomenon, the proximity effect is also sensitive to modifications of the electron and hole phases. At finite energies, the electron and hole pick up different phases and at $`EE_{Th}`$, the electron-hole coherence is destroyed and consequently, the proximity effect is suppressed. Moreover, a magnetic flux in the $`N^{}`$-region large enough to break time reversal symmetry, also suppresses the proximity effect, since the electron and the hole pick up different phases in a finite magnetic field. Here we however only consider the case with negligibly small magnetic fields.

Formally, at energies well below the superconducting gap $`\mathrm{\Delta }`$, the matrix Greens function $`\stackrel{ˇ}{G}_S^0`$ in Eq. (14) simplifies considerably, since $`\widehat{R}_S=\widehat{A}_S=\widehat{\sigma }_x`$ and consequently $`\widehat{K}_S=0`$. We then find $`F(\chi _N,\chi _S)`$ in Eq. (22) by diagonalizing $`(\stackrel{ˇ}{G}_x^2)^{1/2}`$, giving

$`F(\chi _N,\chi _S)={\displaystyle \frac{t_0}{e^2}}{\displaystyle \frac{dE}{2\sqrt{2}}\sqrt{\alpha +\sqrt{\alpha ^2+\beta ^24g_1^2g_2^2\mathrm{\Lambda }}}},`$
(26)
where we have introduced

$`\alpha `$ $`=`$ $`g_1^2+g_2^2{\displaystyle \frac{E^2}{E_{Th}^2}}(g_1+g_2)^2,`$
$`\beta `$ $`=`$ $`2g_1(g_1+g_2){\displaystyle \frac{E}{E_{Th}}},`$
$`\mathrm{\Lambda }`$ $`=`$ $`\left(1e^{2i(\chi _N\chi _S)}\right)f_+(1f_{})`$ (27)
$`+`$ $`\left(1e^{2i(\chi _N\chi _S)}\right)f_{}(1f_+).`$
We note that at $`E=0`$, we recover the generating function of Ref. \[Boerlin02, \]. The terms $`\text{exp}(\pm 2i[\chi _N\chi _S])`$ show that electrons are transported in pairs out of ($`+`$) and into ($``$) the superconductor and into ($`+`$) and out of ($``$) the normal reservoir. This pair transport follows from the fact that the charge transfer mechanism across the $`NS`$-interface is Andreev reflection. The form of the cumulant generating function, with the double square roots, tells that the pair transfer is highly correlated. There are however several limits where the generating function is simply proportional to $`\mathrm{\Lambda }`$, describing an uncorrelated, (generalized) Poissonian, transfer of pairs of electrons (as follows from Eq. 23).

In the limit of high energy, $`EE_{Th}`$, for arbitrary $`g_1`$ and $`g_2`$, the integrand of the generating function in Eq. (26) takes on the form

$`{\displaystyle \frac{F(\chi _N,\chi _S)}{E}}={\displaystyle \frac{t_0}{4e^2}}{\displaystyle \frac{g_1g_2^2}{(g_1+g_2)^2}}\left({\displaystyle \frac{E_{Th}}{E}}\right)^2\mathrm{\Lambda }.`$ (28)
For these energies, the probability of pair tunneling is strongly suppressed and the pairs are thus emitted in a Poisson process. This is a consequence of the different phases picked up by the electrons and holes, $`\pm E/E_{Th}`$, at finite energies. At $`EE_{Th}`$ the electrons and holes loose their coherence and the pair tunneling probability is suppressed.

In the limit of a dominating coupling of the diffusive region $`N^{}`$ to the normal reservoir, $`g_1g_2`$, the integrand of the generating function becomes

$`{\displaystyle \frac{F(\chi _N,\chi _S)}{E}}={\displaystyle \frac{t_0}{4e^2}}{\displaystyle \frac{g_2^2}{g_1}}{\displaystyle \frac{1}{1+(E/E_{Th})^2}}\mathrm{\Lambda }.`$ (29)
Thus, in this limit, the pair transfer is also uncorrelated. In this case, the barrier $`I_1`$ has a negligible resistance compared to $`I_2`$ and the junction is effectively a $`NI_2S`$-junction. As discussed above, in this case the pair transfer is uncorrelated. We note that for high energies $`EE_{Th}`$, it coincides with Eq. (28), (for $`g_1g_2`$).

We also noteBoerlin02  that for the limit $`g_1g_2`$, strong coupling of $`N^{}`$ to the superconductor and for the additional limit of low energies, $`EE_{Th}`$, the pair transfer is uncorrelated. In this case, the integrand of the generating function is simply

$`{\displaystyle \frac{F(\chi _N,\chi _S)}{E}}={\displaystyle \frac{t_0}{4e^2}}{\displaystyle \frac{g_1^2}{g_2}}\mathrm{\Lambda }`$ (30)
The reason for the uncorrelated pair transfer in this case is the following: In the limit $`g_1/g_20`$, a proximity gap opens up in $`N^{}`$ and the junction becomes an effective $`NI_1S^{}`$ junction ($`S^{}`$ denoting the gapped normal region between the barriers). For energies below the induced gap, the pair transfer is thus Poissonian, as in a standard $`NIS`$-junctions. Khlus87 ; deJong94 ; Muz94

Making a connection to the proximity effect, it can be noted that the in two cases, Eq. (28) and Eq. (29), the strongly suppressed proximity effect in $`N^{}`$ is responsible for the uncorrelated pair transport. However, in the last case, Eq. (30), the Poissonian pair transfer is due to the induced proximity gap in $`N^{}`$.

It is instructive to also study qualitatively the full counting statistics, starting from the scattering point of view in Ref. \[Muz94, \]. The cumulant generating function, generalized to many modes, is

$`F(\chi _N,\chi _S)={\displaystyle \frac{t_0}{h}}{\displaystyle 𝑑E\text{ln}(1+R_{eh}(E)\mathrm{\Lambda })}`$ (31)
where the $`R_{eh}(E)`$’s are the Andreev reflection eigenvalues, i.e. the real and energy dependent eigenvalues of the hermitian electron-hole scattering matrix product $`S_{eh}^{}S_{eh}`$ (see e.g. Ref. \[deJong97, \]). Here $`..`$ denotes ensemble average over impurity configurations and $`\mathrm{\Lambda }`$ is the same as in Eq. (27). For a general function $`h(R_{eh})`$, the ensemble average can be written as an integration $`𝑑R_{eh}\rho (R_{eh})h(R_{eh})`$, where $`\rho (R_{eh})`$ is the distribution of Andreev reflection eigenvalues (the corresponding quantity in NS-systems to the transmission eigenvalue distribution in normal systems).

Expanding the logarithm in Eq. (31), makes it is clear that the limits where the charge transfer is uncorrelated, Eqs. (28) to (30), corresponds to the limit where $`R_{eh}R_{eh}^2`$ (and consequently $`R_{eh}R_{eh}^n`$ for all higher moments $`n>2`$ since $`0R_{eh}1`$).

From this we can draw the following conclusions. The limits of uncorrelated pair transfer corresponds to the cases where there are no or very few “open” Andreev channels, i.e with Andreev reflection probability close to unity (leading to $`R_{eh}R_{eh}^2`$). In the case where there are more open channels, the various moments $`R_{eh}^n`$ are of the same order of magnitude, resulting e.g. in a noise below the Poissonian noise of uncorrelated pair transfer across the NS-interface, further discussed below. As is also further discussed below, this Andreev channel picture provides a simple explanation for the voltage and temperature dependence of the noise as well as the current.

## III Noise.

The rest of the paper is devoted to a detailed analysis of the properties of the the second cumulant, the noise, which is within reach of existing experimental techniques. We note that the cumulant generating function $`F(\chi _N,\chi _S)`$ in Eq. (26) can be written as $`g_1`$ (or $`g_2`$) times a function which only contains the ratio $`g_1/g_2`$, i.e. all transport properties studied below, (normalized) depend only on the ratio of the conductances.

For the benefit of the following discussion, we first study the current, given from Eq. (24) and (26),

$`I={\displaystyle \frac{1}{e}}{\displaystyle 𝑑E\frac{f_{}f_+}{\sqrt{2}}\frac{g_1^2g_2^2}{\sqrt{\alpha ^2+\beta ^2}\sqrt{\alpha +\sqrt{\alpha ^2+\beta ^2}}}}`$ (32)
The current has been studied both theoretically Volkov ; Nazcirc  and experimentally, Lefloch  here we briefly summarize the findings. The differential conductance $`G=dI/dV`$ at zero temperature is plotted as a function of $`eV/E_{Th}`$ in Fig. 2, for different ratios $`g_1/g_2`$.

The conductance for $`eVE_{Th}`$ is given by

$`G(eVE_{Th})G_{NS}={\displaystyle \frac{g_1^2g_2^2}{(g_1^2+g_2^2)^{3/2}}}.`$ (33)
In the opposite limit, $`eVE_{Th}`$, the conductance decays as a power law with applied voltage, as $`(E_{th}/eV)^2g_2^2g_1/(g_1+g_2)^2`$.

The behavior of the conductance for voltages of the order of $`E_{Th}`$ depends strongly on the ratio $`g_1/g_2`$. For $`g_1g_2`$, i.e. for a strong coupling of the normal region $`N^{}`$ to the normal reservoir, the conductance decreases monotonically with voltage as $`(g_2^2/g_1)/[1+(eV/E_{Th})^2]`$. In the opposite regime, $`g_2g_1`$, the normal region is strongly coupled to the superconducting reservoir. The conductance has a strong peak at voltages $`eVE_{Th}`$, which can be seen as an ensemble averaged electron-hole resonance. In other words, for energies $`EE_{Th}`$ the distribution of Andreev reflection eigenvalues, $`\rho (R_{eh})`$, is shifted towards more “open” channels.

As pointed out above, in the limit of $`g_1/g_20`$, i.e. when decoupling the normal reservoir, a proximity induced gap $`E_{Th}`$ opens up in the spectrum in the normal region $`N^{}`$. In this limit, the conductance shows a singularity at $`eV=E_{Th}`$, similar to the standard NIS-tunneling conductance, but with the induced gap $`E_{Th}`$ instead of the superconducting gap $`\mathrm{\Delta }`$.

The temperature dependence of the zero voltage conductance, $`G(T)`$, is shown in Fig. 3. The temperature dependence is qualitatively similar to the voltage dependence, shown in Fig. 2, due to the form of the spectral current density in Eq. (32). We note that for $`kTE_{Th}`$ we have $`G=G_{NS}`$ from Eq. (33). In the opposite limit, $`kTE_{Th}`$ the conductance decays as $`1/kT`$. In the intermediate regime, the behavior of the conductance depends strongly on the ratio $`g_1/g_2`$. For $`g_1g_2`$, it decays monotonically with temperature, while in the opposite limit, $`g_2g_1`$, it has a peak at temperatures $`kT`$ of order $`E_{Th}`$.

We then turn to the noise, found from Eq. (25) and (26) to be

$`S`$ $`=`$ $`{\displaystyle 𝑑E\frac{2\sqrt{2}g_1^2g_2^2}{\sqrt{\alpha ^2+\beta ^2}\left(\alpha +\sqrt{\alpha ^2+\beta ^2}\right)^{3/2}}}`$ (34)
$`\times `$ $`[[f_+(1f_{})+f_{}(1f_+)](\alpha +\sqrt{\alpha ^2+\beta ^2})`$
$`(f_{}f_+)^2{\displaystyle \frac{g_1^2g_2^2}{\alpha ^2+\beta ^2}}(2\alpha +3\sqrt{\alpha ^2+\beta ^2})].`$
First we consider the equilibrium noise, for $`eVE_{Th},kT`$. In this limit we have $`f_{}f_+=0`$ and the noise is

$`S_{eq}`$ $`=`$ $`{\displaystyle 𝑑E\frac{4\sqrt{2}f(1f)g_1^2g_2^2}{\sqrt{\alpha ^2+\beta ^2}\sqrt{\alpha +\sqrt{\alpha ^2+\beta ^2}}}}.`$ (35)
Noting that we have $`2f(1f)=kT[d(f_{}f_+)/dV]_{eV=0}=2kT(df/dE)`$ we see from Eq. (32) and (35) that

$`S_{eq}(T)`$ $`=`$ $`4kTG(eV=0,T),`$ (36)
for arbitrary temperature $`T`$. This is just what we expect from the fluctuation-dissipation theorem. The temperature dependence of the noise, normalized with $`4kT`$, is thus just the temperature dependence of the zero voltage conductance, plotted in Fig. (3).

At $`kTE_{Th}`$ we find just $`S_{eq}(kTE_{Th})=4kTG_{NS}`$. In the opposite limit, $`kTE_{Th}`$ the equilibrium noise tends towards a constant, temperature independent value

$`S_{eq}(kTE_{Th})`$ $`=`$ $`E_{Th}{\displaystyle \frac{\pi g_1g_2^2}{(g_1+g_2)\sqrt{g_1^2+g_2^2}}}.`$ (37)
In the intermediate temperature regime, the noise, as the conductance in Fig. 3, depends strongly on the ratio $`g_1/g_2`$.

The voltage dependence of the noise is qualitatively different from the temperature dependence. In the limit $`kTeV,E_{Th}`$, the shot noise regime, we have

$`S`$ $`=`$ $`{\displaystyle _{eV}^{eV}}𝑑E{\displaystyle \frac{\sqrt{2}g_1^2g_2^2}{\left(\alpha +\sqrt{\alpha ^2+\beta ^2}\right)^{3/2}}}`$ (38)
$`\times `$ $`\left(1+{\displaystyle \frac{3g_1^2g_2^2}{\alpha ^2+\beta ^2}}{\displaystyle \frac{\alpha (\alpha ^2+\beta ^22g_1^2g_2^2)}{(\alpha ^2+\beta ^2)^{3/2}}}\right).`$
Here we focus on the differential noise, $`dS/dV`$, which is plotted as a function of voltage for different ratios $`g_1/g_2`$ in Fig. 4.

In the low voltage limit, $`eVE_{Th}`$, the differential noise is given by deJong97 ; Boerlin02

$`{\displaystyle \frac{dS}{dV}}|_{eVE_{Th}}P_{NS}=2eG_{NS}\left(2{\displaystyle \frac{5g_1^2g_2^2}{(g_1^2+g_2^2)^2}}\right).`$ (39)
The Fano factor, $`P_{NS}/2eG_{NS}`$, in this limit, is thus $`(25g_1^2g_2^2/[g_1^2+g_2^2]^2)`$. It varies from $`2`$ in the limits $`g_1g_2`$ and $`g_2g_1`$, to $`3/4`$ forFauchere98  $`g_1=g_2`$. A Fano factor $`2`$ is an indication of uncorrelated emission of pairs of electrons into or out of the superconductor. As pointed out in connection to Eq. (29) and (30), in these limits, for $`eVE_{Th}`$ and $`g_1g_2`$ or $`g_2g_1`$, the generating function $`F(\chi _N,\chi _S)`$ in Eq. (26) is proportional to Boerlin02  $`\text{exp}(2i[\chi _N\chi _S])`$, showing that the charge transfer is indeed an uncorrelated emission (for negative voltages) of pairs of electrons from the superconductor.

In the high voltage limit, $`eVE_{Th}`$, the differential noise decreases with voltage as

$`{\displaystyle \frac{1}{e}}{\displaystyle \frac{dS}{dV}}=2\left({\displaystyle \frac{E_{th}}{eV}}\right)^2{\displaystyle \frac{g_2^2g_1}{(g_1+g_2)^2}}`$ (40)
i.e. the same power law behavior as the conductance, which is clear directly from the generating function in Eq. (28).

At intermediate voltages, the noise depends strongly on the ratio $`g_1/g_2`$. For a dominating coupling of the normal region to the normal reservoir, $`g_1g_2`$, the differential noise decreases monotonically with applied voltage, as

$`{\displaystyle \frac{1}{e}}{\displaystyle \frac{dS}{dV}}=2{\displaystyle \frac{g_2^2}{g_1}}{\displaystyle \frac{1}{1+(eV/E_{Th})^2}}.`$ (41)
This can also be obtained directly from the generating function in Eq. (29). We note that this behavior is very similar to what was found for a diffusive $`NIS`$-junction,Hekkila02  a diffusive normal region of non-negligible resistance, connected to a superconducting reservoir via a tunnel barrier.

In the opposite limit, $`g_1g_2`$, the noise shows a double peak structure around $`eVE_{Th}`$, the voltage where the conductance shows a single peak (see Fig. 2). Such a double peak behavior has been discussed before in single mode junctions with sharp Andreev resonances. Fauchere98  In single mode junctions, the double peak behavior follows directly from the fact that if the Andreev reflection probability $`R_{eh}`$ has a resonance at some energy $`E_0`$, the noise, $`R_{eh}(1R_{eh})`$, shows a double peak around $`E_0`$. Interestingly, as is clear from Fig. 4, in the $`NI_1N^{}I_2S`$, this double-peak behavior survives the ensemble average. This indicates that the whole distribution $`\rho (R_{eh})`$ of Andreev reflection eigenvalues is shifted from being dominated by “closed” towards being dominated by “open” channels at resonance $`EE_{Th}`$.

We note that such a double peak behavior can also be found in normal single-mode double-barrier junctions, it is not an effect of the superconductivity (see Ref. \[Buttikerrew, \]). However, in normal, many-mode systems, such a behaviour is, to leading order in number of modes, washed out when performing an ensemble average.

Additional insight can be obtained by studying the differential Fano factor, $`[dS(V)/dI(V)]/2e`$. In the differential Fano factor, the voltage dependence of the conductance is canceled out from the noise, and we obtain the “bare” voltage-dependence of the noise. We have

$`{\displaystyle \frac{1}{2e}}{\displaystyle \frac{dS(V)}{dI(V)}}`$ $`=`$ $`2{\displaystyle \frac{2g_1^2g_2^2\left(2\alpha +3\sqrt{\alpha ^2+\beta ^2}\right)}{(\alpha ^2+\beta ^2)\left(\alpha +\sqrt{\alpha ^2+\beta ^2}\right)}}.`$ (42)
This quantity is plotted in Fig. 5. We note that the differential Fano factor shows a much stronger dependence of the applied voltage, compared to earlier studied Belzig01a ; Reulet02 ; Hekkila02 ; Belzig02b ; Reulet02b  diffusive $`NS`$ and $`NIS`$ junctions. The resonant double-peak in the noise, for $`g_2g_1`$ is not visible in the differential Fano factor, instead it shows a dip at $`eVE_{Th}`$. This is in accordance with our explanation of the double peak-behavior of the noise above, since, again making the comparison to a single mode junction, the differential Fano factor is the differential shot noise divided by the differential conductance, $`[R_{eh}(1R_{eh})]/R_{eh}=1R_{eh}`$, which thus has a dip at the resonant energy.

Moreover, in the high voltage limit, $`eVE_{Th}`$, the differential Fano factor saturates at a constant value $`2`$, i.e. independent on the relation between the conductances $`g_1`$ and $`g_2`$. This is a manifestation of the uncorrelated charge transfer at these energies, as discussed in connection to Eq. (28). Also in the limit $`g_2g_1`$, the Fano factor approaches $`2`$, independent on voltage. This is again a manifestation of the uncorrelated pair transfer, and can be seen directly from the generating functions in Eq. (29). From this we conclude that the differential Fano factor is thus the relevant quantity for studying the energy dependence of the charge transfer mechanism.

## IV Conclusions

In conclusion, we have investigated the energy dependent full counting statistics of a diffusive normal-insulating-normal-insulating-superconducting ($`NI_1N^{}I_2S`$) junction. We have used the recently developed circuit theory of full counting statistics, allowing us to access the full temperature and voltage dependence of the statistics. In general, the charge is transported via correlated transfer of pairs of electrons. However, in the case of strongly asymmetric contacts or energies much larger than the Thouless energy, the pair transfer is uncorrelated. The second cumulant, the noise, was studied in detail. It was found to depend strongly on voltage and temperature. For low temperatures and a junction resistance dominated by the tunnel barrier to the normal reservoir ($`I_1`$), the noise shows a double-peak behavior at voltages $`eVE_{Th}`$, a signature of an ensemble averaged electron-hole resonance.

We acknowledge discussions with F. Lefloch, W. Belzig and M. Büttiker. This work was supported by the Swiss National Science Foundation and the project for Materials with Novel Electronic Properties.


format: text

text:
# Intrinsic accuracy in 3-dimensional photoemission band mapping

## I Introduction: Why intrinsic accuracy?

The quasiparticle bandstructure $`E(𝐤)`$, which reflects the peaks of the spectral function $`A(E,𝐤)`$ depending on energy and wavevector $`𝐤`$, is the key property of the crystalline solids. Angle-resolved photoemission (PE) spectroscopy (for reviews see, e.g., ) is the main experimental method to map $`E(𝐤)`$ with resolution in energy $`E`$ and, in principle, in the 3-dimensional $`𝐤`$ (3D band mapping). This is based on conservation of $`𝐤`$ in the process of photoexcitation from the occupied initial state to unoccupied final state in the bulk of the crystal. However, the PE experiment is inherently performed on the crystal terminated by the surface. By virtue of the remaining surface-parallel periodicity, the surface-parallel wavevector $`𝐤_{}`$ is conserved in the whole system and can be determined from the surface-parallel vacuum wavevector $`𝐊_{}`$ with the accuracy limited only by instrumental resolution. As the surface-perpendicular periodicity is however broken, the surface-perpendicular wavevector $`k_{}`$ is not conserved in the whole system and can not be directly determined. Control over $`k_{}`$ is therefore a fundamental problem of 3D band mapping.

The control over $`k_{}`$ has two aspects. First, $`k_{}`$ is distorted when the photoelectrons exit from the bulk of the solid into vacuum. As $`k_{}`$ is conserved in the very photoexcitation process in the bulk, it can nevertheless be recovered if the $`k_{}`$-dispersion of the final state $`E^f(k_{})`$ or that of the initial state $`E^i(k_{})`$ back in the bulk is known. Commonly an empirical free-electron-like approximation is used here for the final states, though the accuracy of this can often be insufficient. Nevertheless, the true quasiparticle $`E^f(k_{})`$ including the non-free-electron and self-energy effects in the final state can be determined using very-low-energy electron diffraction (VLEED) spectroscopy (see  and references therein) based on that the PE final states, as proved by the one-step theory of PE, are the time-reversed LEED states . It is a common point of view that knowledge of the final-state $`E(𝐤)`$ exhaust the whole problem of 3D band mapping. Indeed, variation of the final-state energy $`E^f`$ in the PE experiment (commonly achieved by variation of the photon energy $`h\nu `$) would through known $`E^f(k_{})`$ be translated into variations of $`k_{}`$, and accurate 3D mapping of the occupied $`E(𝐤)`$ would be performed.

There is another aspect however. Due to inelastic absorption and elastic reflection from the crystal potential the PE final state is intrinsically damped towards interior of the solid . Such a confinement in the surface-perpendicular coordinate is equivalent, by the uncertainty principle, to certain intrinsic broadening in $`k_{}`$. The PE signal via this final state is then formed by all initial states within the $`k_{}`$-broadening interval, and the PE peak reflects an average of the quasiparticle $`E^i(k_{})`$ in the valence band. A difference of the measured averaged $`k_{}`$-dispersion from the true $`E^i(k_{})`$ limits intrinsic accuracy of 3D band mapping. It is borne by fundamental physics of the PE process, and can not be improved instrumentally.

The intrinsic accuracy appears thus as an extrinsic factor to the excited-state self-energy effects as deviations of the quasiparticle excited-state $`E(𝐤)`$, measured in the ideal PE experiment with undamped final states, from the ground-state bandstructure.

The intrinsic accuracy problem has in pieces been tackled in quite a few papers. Here, my aim is to delineate an entire concept, including essential physics and mechanisms of the intrinsic accuracy, and unveil its non-trivial effects encountered in the PE experiment. The concept is built upon a simple picture involving interplay of the final- and initial-state spectral functions, which gives a clearer physical insight compared to the explicit one-step PE theory with its heavy computational machinery. Finally, I discuss optimization of the PE experiment with respect to the intrinsic accuracy of 3D band mapping.

## II Basic physics of the PE process

In the following analysis of the PE process it is implied that $`E(𝐤)`$ is the quasiparticle bandstructure, already incorporating the excited-state self-energy effects. We will neglect the surface effects and concentrate on the bulk derived states with significant dispersion in $`k_{}`$ (so-called 3D states). The ideal $`E`$\- and $`𝐤_{}`$-resolution of the PE experiment is assumed.

### A Initial and final states of PE

First we will recall physics of the PE initial and final states based on the ideas from . For simplicity we will assume that their wavefunctions have only one Bloch wave constituent.

The final-state wavefunction is the time-reversed LEED wavefunction. As sketched in Fig.1 (the upper panels represent the final state), the wavefunction is described by a Bloch wave which is damped in the surface-perpendicular direction $`r_{}`$ towards the crystal interior . There are two sources of such damping: 1) inelastic absorption, expressed by absorption potential $`V_i`$ connected with the electron lifetime as $`2V_i=`$ $`\mathrm{}/\tau _e`$, and 2) elastic reflection from the crystal potential in the final-state band gaps. The damping is described by complex $`k_{}^f`$ component of the final-state wavevector $`𝐤^f`$, with $`Imk_{}^f`$ connected to the photoelectron escape length $`\lambda `$ as $`2Imk_{}^f=\lambda ^1`$ (the factor of 2 comes from squaring the wavefunction amplitude for intensity).

The corresponding final-state dispersion $`E^f(k_{})`$, shown in Fig.1 by black lines, is radically different from this one if the Bloch waves were undamped (propagating), as shown in gray: the dispersions as a function of $`Rek_{}^f`$ are smooth and gapless, passing continuously through the band gaps of the undamped $`E(𝐤)`$ . Such regions will in the following be somewhat loosely referred to as final-state band gaps. Note that due to surface-parallel invariance of the PE process the Bloch waves are undamped in the surface-parallel direction, retaining real $`𝐤_{}`$ and gapped $`𝐤_{}`$-dispersions .

The final-state spectral function $`A^f(E,𝐤)`$, due to the finite photoelectron lifetime, is in principle characterized by certain distribution in $`E`$. However, in the PE experiment this is reduced to the $`\delta `$-function at the energy given by the PE analyser. Due to damping of the final-state wavefunction $`A^f(E,𝐤)`$ is characterized, on the contrary, by Lorentzian distribution in real $`k_{}`$, as shown in Fig.1, centered at $`k_{}^0=Rek_{}^f`$ of the final-state Bloch wave and having the fullwidth $`\delta k_{}=2Imk_{}^f`$ (in the angle-resolved PE experiment $`𝐤_{}`$ is fixed and can be omitted). Therefore, $`A^f(E,𝐤)`$ is given by

$$A^f(k_{})\frac{\delta k_{}}{\left(k_{}k_{}^0\right)^2+(\delta k_{}/2)^2}$$
(1)
where the nominator stands for the normalization $`\underset{\mathrm{}}{\overset{+\mathrm{}}{}}A^f(k_{})𝑑k_{}=1`$. The final state of PE is thus characterized by fixed $`E`$ and broadening in $`k_{}`$ whose width is determined by $`Imk_{}^f`$ of the damped final-state Bloch wave .

Nature of the PE initial state is different and complementary. The initial-state wavefunction is basically undamped: Removal of holes by the finite hole lifetime is compensated by their simultaneous generation by the electromagnetic field, which is almost homogeneous, extending into the crystal interior over a light absorption length of $``$100 Å, very large compared to the final-state damping length. As sketched in Fig.1 (the lower panels represent the initial state), initial-state wavefunction is therefore described by an almost undamped Bloch wave, having well-defined real $`k_{}`$.

The corresponding initial-state dispersion $`E^i(k_{})`$ shows up sharp and gapped $`k_{}`$-dispersions typical of the undamped states, radically different from the final-state.

The initial-state spectral function $`A^i(E,𝐤)`$ is characterized, owing to well-defined $`k_{}`$, by distribution in real $`k_{}`$ reduced to the $`\delta `$-function. Due to finite hole lifetime $`\tau _h`$ it is characterized, on the contrary, by Lorentzian distribution in energy, as shown in Fig.1, centered at the band energy $`E^i(k_{})`$ and having the fullwidth $`\delta E=\mathrm{}/\tau _h`$. Therefore, $`A^i(E,𝐤)`$ is given by

$$A^i(E)\frac{\delta E}{\left(EE^i(k_{})\right)^2+\left(\delta E/2\right)^2}$$
(2)
where the nominator ensures $`\underset{\mathrm{}}{\overset{+\mathrm{}}{}}A^i(E)𝑑E=1`$. The initial state of PE is thus characterized, complementarily to the final state, by fixed $`k_{}`$ and broadening in $`E`$ whose width is determined by the hole lifetime .

Note that some implementations of the one-step PE theory  describe the hole lifetime effects by damping of the initial states. However, such an approach is in clear contradiction with the PE experiment, which unambiguously finds in the valence band the gapped $`E^i(k_{})`$ dispersions typical of undamped states. The same fundamental flaw is implied by an PE lineshape analysis (see, e.g., ) where the $`E`$-broadening of the initial state is replaced by an equivalent $`k_{}`$-broadening $`\delta E\left(\frac{E^i}{k_{}}\right)^1`$. Although such a replacement does not alter the lineshape while the $`k_{}`$-dispersions remain linear, it leaves out some of the intrinsic effects discussed here.

Inverse photoemission (IPE) is commonly treated as the time-reversed PE process  in which the IPE final state in the conduction band is equivalent to the PE initial state in the valence band. However, the above arguments suggest a fundamental difference between the two: as the IPE final state is excited by the damped LEED wavefunction, it is, in contrast to the PE initial state, also damped. Correspondingly, the conduction band $`E(𝐤)`$ measured in the IPE experiment should be characterized by gapless $`k_{}`$-dispersions similar to the PE final-state dispersions in Fig.1. Surprisingly, such a fundamental difference between the IPE and PE processes, tracing back to fundamental principles of quantum mechanics, has so far been addressed nor theoretically neither experimentally.

### B Development of the PE spectrum: Intrinsic $`k_{}`$-resolution

We will now apply the above principles to analyze the PE response of 3D states. We will use a simple formalism based on interplay of the spectral functions , which is extremely simple and gives a clear physical insight. Its accuracy, despite principal limitations compared to the explicit one-step PE theory (see III E), is nevertheless sufficient for evaluation of the PE peak profiles .

As illustrated in Fig.2, the photoemission current $`I(E^f,E^i)`$ at the final-state energy $`E^f`$ and initial-state energy $`E^i=E^fh\nu `$ is formed by adding up the elementary photocurrents $`dI(E^f,E^i)`$ from all $`dk_{}`$ intervals, in which direct transitions take place, with $`k_{}`$ spanning the whole Brillouin zone (BZ). Each $`dI(E^f,E^i)`$ is weighted by the above spectral functions, $`A^f(k_{})`$ to express closeness of $`k_{}`$ to $`k_{}^0`$ of the final-state Bloch wave, and $`A^i(E^i)`$ to express closeness of $`E^i`$ to the initial-state band energy $`E^i(k_{})`$ at this $`k_{}`$ with $`\delta E`$ varying along $`E^i(k_{})`$. Moreover, each $`dI(E^f,E^i)`$ is multiplied by amplitude factors of the final-state surface transmission $`|T^f|^2`$ and photoexcitation matrix element $`\left|M_{fi}\right|^2`$. This yields $`dI(E^f,E^i)=dk_{}\left|T^f\right|^2\left|M_{fi}(k_{})\right|^2A^f\left(k_{}\right)A^i\left(E^i\right)`$. To obtain the whole photocurrent, this expression should be integrated in $`k_{}`$:

$`I(E^f,E^i)`$ $``$ $`{\displaystyle \underset{\mathrm{}}{\overset{+\mathrm{}}{}}}dk_{}\left|T^f|^2\right|M_{fi}(k_{})|^2`$ (5)
$`{\displaystyle \frac{\delta k_{}}{\left(k_{}k_{}^0\right)^2+(\delta k_{}/2)^2}}`$
$`{\displaystyle \frac{\delta E}{\left(E^iE^i(k_{})\right)^2+(\delta E/2)^2}}`$
This integral can be evaluated at negligible computational cost (it can even be performed analytically using Taylor expansion of $`E^i(k_{})`$ ).

Cuts of $`I(E^f,E^i)`$ along the $`E^f=const`$ line correspond to the constant-final-state (CFS) measurement mode, along the $`h\nu =E^fE^i=const`$ line to the energy-distribution-curve (EDC) mode, and along $`E^i=const`$ to the constant-initial-state (CIS) mode. Fig.2 illustrates, for simplicity, the PE spectrum taken using the CFS mode in which $`k_{}^0`$ is constant throughout the spectrum.

The PE spectrum shows up, whichever the measurement mode, a peak centered roughly at the energy $`E^i(k_{}^0)`$ dictated by the direct transition at $`k_{}^0`$ between the final and initial bands. Broadening of the peak results from the final-state $`k_{}`$-broadening combined with the initial-state $`E`$-broadening. In this context the $`k_{}`$-broadening appears as an intrinsic $`k_{}`$-resolution of the PE experiment (in the following we will use these terms interchangeably depending on the context). Note that this resolution, in contrast to the instrumentally limited resolution in $`𝐤_{}`$, is limited by intrinsic mechanisms involved in the PE process and can not be improved instrumentally.

### C Regimes of the PE experiment

If the $`k_{}`$-resolution is at its best extreme $`\delta k_{}0`$, then $`A^f\left(k_{}\right)`$ becomes $`\delta `$-function $`\delta (k_{}k_{}^f)`$ and the integral (5) for $`I(E^f,E^i)`$ is reduced to $`A^i\left(E\right)`$ whose maximum points to the initial band energy $`E^i(k_{}^0)`$ in the $`𝐤`$-point dictated by the direct transition. The PE peaks will then exactly follow the 3D quasiparticle dispersions in the valence band. In practice, such a situation is realized under the condition

$$\delta k_{}k_{}^{BZ}$$
(6)
where $`k_{}^{BZ}`$ is the surface-perpendicular dimension of the BZ. This condition identifies the so-called bandstructure regime (BS-regime) of the PE experiment which ensures accurate 3D band mapping . In the real space, it reads as $`\lambda c_{}`$, where $`c_{}`$ is the surface-perpendicular dimension of the unit cell.

If the PE final state comprise a few Bloch waves, the condition (6) should be generalized as $`\delta k_{}\frac{k_{}^{BZ}}{N}`$, where $`N`$ is the number of the corresponding final bands. Often in the literature (see, e.g., ) this number is confused with the number of all unoccupied bands available for given $`E^f`$ and $`𝐤_{}`$; as the latter becomes with energy progressively immense, an erroneous conjecture is made that at high energies any 3D band mapping is impossible. However, in the multitude of the available unoccupied bands only those are effective as the final bands, whose Bloch waves provide effective transport of the photoelectrons out of the solid by significant coupling to the outgoing plane wave . Typically there are only a few ($`N=12`$) final bands in this sense (see examples in III A).

If the $`k_{}`$-resolution is at the opposite extreme $`\delta k_{}\mathrm{}`$, any resolution in $`k_{}`$ is lost, $`A^f\left(k_{}\right)`$ becomes constant and the integral (5) yields, assuming that $`|T^f|^2`$ and $`\left|M_{fi}\right|^2`$ are independent of $`k_{}`$ and $`E^i(k_{})`$ is linear, the 1-dimensional DOS (1DOS) $`\frac{k_{}}{E^i}`$. The PE peaks will then manifest the 1DOS maxima independently of $`k_{}^0`$. In practice, such a situation is achieved under the condition $`\delta k_{}k_{}^{BZ}`$, which identifies the so-called 1DOS-regime .

Interpretation of the PE data between the bandstructure and 1DOS regimes is doubtful : the PE peaks retain some $`h\nu `$-dispersion, although damped by strong averaging in $`k_{}`$, but they can be far from the true quasiparticle $`E^i(k_{})`$.

Evolution of the regimes of the PE experiment with the $`k_{}`$-resolution will now be illustrated by a PE simulation. The model employs the initial state whose parameters are chosen to resemble roughly the sp-band of Cu : It is described by the free-electron dispersion $`E^i(k_{})=\frac{\mathrm{}^2}{2m}k_{}^2+V_{000}`$ and the energy broadening $`\delta E`$ varying linearly as a function of energy to match zero at the Fermi level $`E_F`$ and an experimental value of 1.7 eV at the bottom of the band. The Fermi-Dirac cut is introduced as the ideal step function. The photocurrent $`I(E^f,E^i)`$ is evaluated by integration of the spectral functions in $`k_{}`$ (5), assuming constant amplitude factors $`|T^f|^2`$ and $`\left|M_{fi}\right|^2`$. The integration also extends above the Fermi vector $`k_{}^F`$ to include certain contribution from the tails of $`A^i(E)`$ in the unoccupied region of the valence band (see III B). For clarity, the PE spectra were evaluated for the CFS mode, in which $`k_{}^0`$ is constant (deviations from the EDC mode are insignificant, because normally the final-state $`\frac{E^f}{k_{}}`$ is large enough to ensure small variations of $`k_{}^0`$ through the peak profile).

The simulation was performed for four $`\delta k_{}`$ values: zero for the ideal $`k_{}`$-resolution, $`0.2`$$`|`$$`\mathrm{\Gamma }X`$$`|`$ as a representative for the BS-regime in a real PE experiment, $`0.6`$$`|`$$`\mathrm{\Gamma }X`$$`|`$ for degraded $`k_{}`$-resolution, and $`2`$$`|`$$`\mathrm{\Gamma }X`$$`|`$ for the 1DOS-regime. The results are shown in Fig.3 as series of the PE spectra with $`k_{}^0`$ scanning along $`\mathrm{\Gamma }X`$ (upper panels, each series in different intensity scale) and positions of the spectral peaks mapped on top of the true valence band (lower panels).

In the ideal $`k_{}`$-resolution limit the PE peak reflects the initial-state $`A^i\left(E\right)`$ and ideally follows the true $`E^i(k_{}^0)`$. Interestingly, upon passing $`k_{}^F`$ there appears an extremely small peak just below $`E_F`$; it replicates the low-energy tail of the spectral function $`A^i(E)`$ whose maximum is already above $`E_F`$ (see III B).

In the BS-regime ($`\delta k_{}=0.2`$$`|`$$`\mathrm{\Gamma }X`$$`|`$ in our simulation) the $`k_{}`$-resolution remains sufficient for the PE peaks to closely follow $`E^i(k_{})`$. This regime can be used for accurate 3D band mapping. Two peculiarities should be noted here: 1) Despite $`A^i(E)`$ becomes the $`\delta `$-function upon approaching $`E_F`$, the PE peak does not sharpen up to singularity; 2) When $`k_{}^0`$ passes $`k_{}^F`$ and enters into the unoccupied region, the peaks do not disappear but become ’ghost’ peaks which are highly asymmetric and have their maximum stationary in energy just below $`E_F`$. Both peculiarities are appearances of the $`k_{}`$-broadening and will be analyzed in III B and III C.

As the $`k_{}`$-resolution degrades ($`\delta k_{}=0.6`$$`|`$$`\mathrm{\Gamma }X`$$`|`$), the concomitant averaging of $`E^i(k_{})`$ causes asymmetry of the peaks and strong intrinsic shifts of their maxima from the true $`E^i(k_{})`$. Although the PE peaks deceptively keep dispersion on $`k_{}^0`$, their use for 3D band mapping will return a distorted valence band dispersion. Note appearance of an additional shoulder-like structure, which originates from the large 1DOS in the bottom of the band and becomes notable as the main peak shifts upwards.

In the 1DOS-regime ($`\delta k_{}=2`$$`|`$$`\mathrm{\Gamma }X`$$`|`$) the $`k_{}`$-resolution is completely lost. The maxima of the PE peaks become stationary in energy reflecting the maximum of the valence band 1DOS (and a tiny replica of the $`A^i(E)`$ singularity at $`E_F`$, omitted from the plot). The $`k_{}^0`$ changes show up only in the intensity modulation. Note however that the PE peaks are notably shifted from the true 1DOS maximum.

The regimes of the PE experiment in general follow the well-known ”universal curve” $`\lambda (E)`$  which gives the $`k_{}`$-resolution as $`\delta k_{}=\lambda ^1`$. The BS-regime holds while $`\lambda `$ remains large, which takes place at low $`E^f`$ below or not far above a plasmon excitation energy of 20-30 eV (see examples in III A). Upon further increase of energy $`\lambda `$ reaches its minimum around 50-100 eV, and the PE experiment enters into the 1DOS-regime, although in practice the $`k_{}`$-dispersion signatures are never completely suppressed. The BS-regime recovers in the high-energy region above $`300`$ eV where $`\lambda `$ rises again (despite at these energies the unoccupied bands become immense in number, only one of them remains effective in the PE final state in the sense of effective coupling to the outgoing plane wave ).

The high-energy region is characterized by low PE intensity resulting from small valence band crossection and from the Debye-Waller factor, in which the increase of $`h\nu `$ reduces the PE intensity on equal footing with temperature. However, the advent of new high-brilliance synchrotron radiation sources and multidetection PE analysers has allowed for PE experiments in this region with reasonable count rates and resolution in $`E`$ and $`𝐤_{}`$ . $`𝐤`$-resolved experiments were reported, for example, in  where the $`sp`$-band of Al was found to be mirrored in the $`h\nu `$-dispersion of the PE signal up to 750 eV. Interestingly, this study has found an increase of the surface to bulk signal intensity ratio with $`h\nu `$, opposite to the tendency expected from the ”universal curve”, and attributed this surprising phenomenon to intrinsic phonon excitation effects involved in the PE process. These results call for further studies on the role of the phonon excitation in high-energy PE spectroscopy.

## III Mechanisms limiting the intrinsic accuracy

Degradation of the $`k_{}`$-resolution, as we observed above, results in intrinsic shifts of PE peaks from the true valence band $`k_{}`$-dispersions, limiting the intrinsic accuracy of 3D band mapping. These shifts can originate through various mechanisms, which will now be analyzed.

### A Non-linearity mechanism

Degradation of the $`k_{}`$-resolution can crawl into the position of the PE peaks, firstly, by non-linearity of the valence band $`k_{}`$-dispersion. As evident from Fig.2, in this case the integral number of the $`dk_{}`$-intervals within the $`k_{}`$-broadening will be different for the initial states whose $`E^i`$ is above the direct-transition energy $`E^i(k_{}^0)`$ and for the states whose $`E^i`$ is below $`E^i(k_{}^0)`$, being larger where the 1DOS is larger. Correspondingly, the PE peak becomes asymmetric. Its maximum somewhat shifts towards larger integral number of $`dk_{}`$, deviating from the true $`E^i(k_{}^0)`$. In the band interior the PE peaks shift towards the band edges, where the 1DOS is larger. Near the band edges the peaks, on the opposite, shift towards the band interior (so-called in-band shifting ) because beyond the band edge there are no states and the 1DOS is zero. We will refer to this mechanism to limit the intrinsic accuracy as the non-linearity mechanism. It was behind the intrinsic shifts observed in the PE simulation in Fig.3, where the PE peaks clearly followed the characteristic pattern of shifting towards larger 1DOS. Note that the intrinsic shifts due to the non-linearity mechanism depend, apart from $`\delta k_{}`$ and non-linearity of $`E^i(k_{})`$, on the initial-state $`\delta E`$.

The non-linearity mechanism can be illustrated by experimental data for VSe<sub>2</sub> (for a detailed discussion see ). Information about the final-state dispersions and lifetimes, required for determination of $`\delta k_{}`$, was achieved in this case by VLEED (see IV). The experimental energy dependence of $`V_i`$ is shown in Fig. 4 (a). It is characterized by a sharp increases of $`V_i`$ at the plasmon excitation energy $`\mathrm{}\omega _p`$. The experimental $`E^f(k_{})`$ is shown in Fig. 4 (b). Note that these final band can hardly be described by free-electron-like dispersions: there are two bands, and each can be free-electron fitted only locally with the inner potential strongly depending on $`E^f`$ and $`𝐤_{}`$ (such a complicated structure of the final states is in fact typical of many materials with large unit cell, containing more than one atom, due to backfolding of many free-electron bands into the first BZ and their mutual hybridization). An estimate of $`\delta k_{}`$ was obtained, neglecting some increase in the final-state band gaps, as $`2V_i\left(\frac{E^f}{k_{}}\right)^1`$ (see IV). The obtained $`\delta k_{}`$ values are shown in Fig. 4 (b) superimposed on the final bands. Following the increase of $`V_i`$, the $`k_{}`$-resolution degrades above $`\mathrm{}\omega _p`$. The results of PE band mapping, the binding energy of the PE peaks as a function of $`k_{}^0`$, are shown in Fig. 4 (c) with distinction for the final states below and above $`\mathrm{}\omega _p`$. Near the band edges of the 3D bands the experimental points above $`\mathrm{}\omega _p`$ are systematically shifted towards the band interior compared to those below $`\mathrm{}\omega _p`$. It should be stressed that here we have an explicit example how the same initial states yield different energies of the PE peaks depending on the $`k_{}`$-resolution. The observed shifts demonstrate the intrinsic effect of in-band shifting, occurring mainly by the non-linearity mechanism. Interestingly, increase of $`E^f`$ towards 50 eV results in further increase of $`\delta k_{}`$ and suppression of any dispersion of the PE peaks .

Our analysis demonstrates that for VSe<sub>2</sub> the energy region of accurate 3D band mapping is exhausted by low $`E^f`$ below or not far above $`\mathrm{}\omega _p`$. Such a situation is in fact fairly general, because a pronounced plasmon threshold is typical of many materials with some exceptions such as Cu and Ni. The intrinsic accuracy effects similar to VSe<sub>2</sub> were also observed for TiS<sub>2</sub> , WSe<sub>2</sub> , Ge , etc.

In the final-state band gaps the Bloch waves experience additional damping due to elastic scattering from the crystal potential . Although this is normally much weaker than the damping due to $`V_i`$, as in the above example, for materials with exceptionally large final-state band gaps the additional damping can cause notable degradation of the $`k_{}`$-resolution. Graphite is a typical example of this (for a detailed discussion see ). Its final bands were simulated using the empirical pseudopotential method based on the VLEED data for their energies in the $`\mathrm{\Gamma }`$-point (the simulation is therefore less accurate in the BZ interior) and $`V_i`$. The resulting $`E^f(k_{})`$ is shown in Fig.5 (a) and (b) as a function of $`Imk_{}^f`$ and $`Rek_{}^f`$, respectively, using the double-BZ representation to reflect the dipole selection rules in graphite. The final states are characterized by a strong increase of $`\delta k_{}=2Imk_{}^f`$ in the two band gaps in the $`\mathrm{\Gamma }`$-point. Based on these data, we simulated the PE response of the valence $`\pi `$-band, having the 3D character, using the same formalism (5). The results are shown in Fig.5 (b). When the final-state $`k_{}^0`$ moves towards the $`\mathrm{\Gamma }`$-point and enters the band gaps, progressive degradation of the $`k_{}`$-resolution results, by the non-linearity mechanism, in increasing of the in-band shifts. In the lower band gap, which has the maximal width, the increase of $`\delta k_{}`$ is particularly strong; the concomitant in-band shifting is so dramatic that it overcomes the trend dictated by the band dispersions themselves and reverses the dispersion of the PE peaks compared to the true $`\pi `$-band. Such an unusual reversed dispersion has indeed been found in the PE experiment .

Typically, intrinsic shifts near the band extrema due to the non-linearity mechanism are of the order of a few tenths of eV. The experiment often suggests however larger shifts . As degradation of the $`k_{}`$-resolution goes along with reduction of the escape length $`\lambda `$, this reveals certain contribution of surface effects such as suppression of the 1DOS in the band extrema (see III D). It should be noted that the condition (6) identifying the BS-regime is particularly restrictive for layered materials owing to their small $`k_{}^{BZ}`$ ( see  for a specialized survey).

### B ’Ghost’ peaks at $`E_F`$

An interesting phenomenon occurs when a 3D band crosses the Fermi surface. This situation is described by the PE simulation in Fig.3: When the direct-transition $`k_{}^0`$ moves outside the Fermi surface and starts sampling the unoccupied part of the valence band, the PE peaks, mysteriously, do not disappear. Instead, they become asymmetric and stationary in energy just below the Fermi level $`E_F`$, and only gradually reduce in amplitude. Experimentally, such a phenomenon has been observed, for example, in a famous PE study on Na, where the $`sp`$-band crosses $`E_F`$ . We will refer to this phenomenon as ’ghost’ peaks. It has in fact two contributions due to different mechanisms.

The first contribution is due to a general mechanism taking place for 3D as well as 2D bands. It comes from the unoccupied states next to $`E_F`$: The spectral function $`A^i(E)`$, centered at the band energy above $`E_F`$, protrudes its low-energy tail, however small the remaining amplitude, below $`E_F`$ and gives therefore some PE intensity (such an effect for 2D bands is discussed, for example, in ). Such a contribution can be isolated in our PE simulation by restricting the integration (5) to $`k_{}>k_{}^F`$. It survives even in the limit of the ideal $`k_{}`$-resolution, see Fig.3. In any case, the ’ghost’ peak intensity due this mechanism is very small.

The second contribution, much predominating for the 3D bands, comes from the occupied states: As illustrated in Fig.6, even if $`k_{}^0`$ moves outside the Fermi surface, the $`dk_{}`$-states at the initial-state $`E^i`$ below $`E_F`$ give certain ’ghost’ intensity by virtue of being accessible via the $`k_{}`$-broadening. Upon moving $`k_{}^0`$ further from $`k_{}^F`$ the $`A^f\left(k_{}\right)`$ profile shifts, reducing its amplitude at the initial band and, concomitanly, the ’ghost’ peak amplitude. Note that this mechanism gives a significant ’ghost’ peak intensity even if the initial-state $`E`$-broadening is vanishing which takes place upon approaching $`E_F`$.

The ’ghost’ peak phenomenon complicates identification of the $`E_F`$ crossings in mapping of the Fermi surface formed by 3D bands. A well-known example is NbSe<sub>2</sub>, whose 4p<sub>z</sub> band forms a 3D pocket of the Fermi surface with extension in $`k_{}`$ small compared to the $`k_{}`$-broadening; direct mapping of this pocket is precluded by non-dispersive ’ghost’ intensity . Even if the mapping is performed as a function of $`𝐤_{}`$, the ’ghost’ peaks stay at $`E_F`$ through the entire interval of $`𝐤_{}`$ where at least some part of the 3D band remains occupied. It should be noted that often the PE data, especially for quasi-2D materials, is interpreted neglecting the 3D effects; in this case such ’ghost’ peaks in the $`𝐤_{}`$-dispersion are misinterpreted as a sign of unusual self-energy effects. A way to circumvent the ’ghost’ peak problem in the Fermi surface mapping was suggested in : one varies $`𝐤_{}`$, simultaneously changing $`E^f`$ to remain in the $`E^i\left(k_{}\right)`$ minimum.

### C Linewidths near $`E_F`$

Similar in origin to the ’ghost’ peaks is an interesting PE linewidth behavior upon approaching $`E_F`$: As seen in our simulations in Fig.3 for the BS-regime, the PE peak from the 3D band, mysteriously, retains a finite linewidth up to $`E_F`$ without sharpening up as might be expected from the behavior of the initial-state $`A^i(E)`$, which becomes singular at $`E_F`$ because of the vanishing $`E`$-broadening. The origin of such a counter-intuitive linewidth behavior is again due to the 3D character of the band and the $`k_{}`$-broadening : As one can conjecture from Fig.6, even if the $`E`$-broadening is zero, the initial states around the direct-transition energy $`E^i(k_{}^0)`$ are accessible via the $`k_{}`$-broadening and contribute to the PE peak width. In the ideal $`k_{}`$-resolution limit, as demonstrates our simulation in Fig.3, the PE peak well becomes singular towards $`E_F`$.

An experimental example of such an effect of the $`k_{}`$-broadening on the linewidths can be seen, for example, in the PE peak from the $`sp`$-band of Cu whose linewidth increases upon approaching $`E_F`$ . In  it was shown that the same effect can explain an unusual linewidth energy dependence near $`E_F`$ observed in Bi<sub>2</sub>Sr<sub>2</sub>CaCu<sub>2</sub>O<sub>8</sub>, which otherwise would indicate deviations from the Fermi liquid picture.

Significant PE intensity next to $`E_F`$, seen in all regimes of the simulation in Fig.3, also originates from existence of a 3D band crossing $`E_F`$ and the $`k_{}`$-broadening. Such a peculiarity is typical of the PE experiments on 3D materials such as Cu.

### D Other mechanisms

Intrinsic shifts in the 3D bands can also appear due to a matrix element mechanism: If $`M_{fi}`$ as a function of energy and $`k_{}`$ undergoes sharp variations within the $`k_{}`$-broadening, a contribution of the states above and below $`E^i(k_{}^0)`$ is different (see Fig.2). The peak becomes asymmetric and it maximum experiences an intrinsic shift from the true band energy. Such a mechanism has been identified, for example, in TiS<sub>2</sub> . Similar effects due to $`T^f`$ should be less pronounced, because the energy variations of $`T^f`$ are normally smooth and become sharper only at extremely low $`E^f`$, as can be seen in VLEED data (see, e.g., )

The above mechanisms were all based on the bulk bandstructure picture, in which the bandstructure near the surface is considered identical to the bulk one. This implies that the crystal potential, and thus the wavefunctions, ideally repeat the bulk ones up to the surface. In fact, within a few atomic layers of the PE escape depth they can experience notable modification. Such surface effects also give rise to intrinsic shifts. For example, the 1DOS singularities in the band extrema, characteristic of the bulk $`E(𝐤)`$, are smeared in the local DOS near the surface (see, e.g., ), which results in an increase of the in-band shifting. The surface effects are harder to control compared to the $`k_{}`$-resolution effects. Qualitatively, their contribution reduces with increase of $`\lambda `$ simultaneously with reduction of the intrinsic shifts of the bulk origin due to the concomitant improvement of the $`k_{}`$-resolution.

Surface photoelectric effect, generated by abrupt change in the dielectric response at the surface barrier and described by the $`𝐀`$ part of $`M_{fi}`$, can influence the PE signal from the 3D bands through interference of the surface and bulk PE components (see, e.g., ). Normally this effect only causes some asymmetry of the PE peaks. However, near the surface and bulk plasmon excitation energies, especially for free-electron metals such as Na, the lineshape changes can be drastic and return large intrinsic shifts . In  it was demonstrated that the surface and bulk components can be experimentally separated by changing the light incidence angle.

### E Intrinsic accuracy within the one-step PE theory

The intrinsic effects are naturally embodied in the one-step PE theory (see, e.g., ) which finds the PE intensity as a matrix element $`\left|\mathrm{\Phi }^f\left|\widehat{𝐀}\widehat{𝐩}\right|\mathrm{\Phi }^i\right|`$ between the final- and initial-state wavefunctions $`\mathrm{\Phi }^f`$ and $`\mathrm{\Phi }^i`$. Interference between different Bloch wave constituents in the final and initial states, matrix element effects and surface effects, all ignored in the simplified formalism (5) used above, are naturally included here. However, a heavy computational machinery often obscures the physical mechanisms.

An example of the intrinsic accuracy analysis using the one-step PE theory can be found in  in application to TiSe<sub>2</sub>. Similarly to the above results for VSe<sub>2</sub>, for the final states above $`h\omega _p`$ the $`k_{}`$-resolution was found to degrade, causing strong broadening and intrinsic shifts of the PE peaks. Surface effects due to smearing of the bulk 1DOS singularities have also been identified. Another example is a calculation on Na , which described the ’ghost’ peaks at $`E_F`$ and in-band shifting in the bottom of the valence band.

## IV Optimization of the PE experiment using VLEED

3D band mapping remains relevant if the intrinsic shifts are small. Of the factors affecting these shifts, the non-linearity of the $`k_{}`$-dispersion, initial-state $`E`$-broadening and, partly, variations of $`M_{fi}`$ are all due to inherent properties of the valence bands under study and remain beyond control by the experimentator. The $`k_{}`$-resolution, on the contrary, is an exclusively final-state property, and by tuning $`E^f`$ can taken to values which ensure negligible intrinsic shifts. This is achieved in the BS-regime $`\delta k_{}k_{}^{BZ}`$ (see II C). The exact energy ranges of this regime are however much material dependent. Knowledge of the $`k_{}`$-resolution energy dependence is therefore required.

VLEED, covering the energy range below $``$40 eV, has recently been established as the experimental method giving the most direct access to the PE final states with resolution in the 3D wavevector (see  and references therein). In the VLEED spectra of elastic electron reflection, the energies of the spectral structures give the critical points in the final bands. $`E^f(k_{})`$, obtained from this experimental data, can further be used to control $`k_{}`$ in the PE experiment. Broadening of the VLEED structures, moreover, gives an estimate of the corresponding $`V_i`$ values . By performing calculations of the complex bandstructure with the experimental $`V_i`$ , one obtains the $`k_{}`$-resolution as $`\delta k_{}=2Imk_{}^f`$ of the generated Bloch waves (and simultaneously the energy dependence of $`\lambda `$, as it was recently demonstrated for graphite in comparison with determination of $`\lambda `$ by the conventional overlayer method ). As such calculations can however be time-consuming, $`\delta k_{}`$ can be estimated in a simplified manner directly from the $`V_i`$ values and the final-state group velocity $`\frac{E^f}{k_{}}`$ as

$$\delta k_{}2V_i\left(\frac{E^f}{k_{}}\right)^1$$
(7)
This estimate is accurate however only outside the final-state band gaps. Nevertheless, increase of$`\delta k_{}`$ associated with additional damping in the gaps (see Fig.5) is normally less significant, except for a few materials with exceptionally wide band gaps such as graphite.

Use of the VLEED predictions on the $`k_{}`$-resolution allows the experimentator to optimize the PE experiment on the intrinsic accuracy by selecting the $`E^f`$ intervals which ensure the BS-regime (see the examples in III A).

As a concluding remark, it should be noted that one can significantly improve accuracy of 3D band mapping by computational modelling of the intrinsic accuracy effects, as described above, and correcting the experimental data correspondingly. This is often critical for proper evaluation of the true excited-state self-energy effects, for example in graphite  and Na .

## V Conclusion

Fundamental principles of PE have been analysed, focussing on emission from 3D states in the valence band. Damping of the final state due to inelastic and elastic scattering in the crystal results in equivalent broadening in $`k_{}`$, which acts as the intrinsic $`k_{}`$-resolution of the PE experiment. The PE peak is then formed as a matrix-element weighted average of the quasiparticle $`k_{}`$-dispersion in valence band, which can be shifted from the dispersion itself. Such shifts intrinsically limit the accuracy of 3D band mapping, which remains accurate only in the bandstructure regime, characterized by the intrinsic $`k_{}`$-broadening small compared to the BZ extension in $`k_{}`$. Upon degradation of the $`k_{}`$-resolution the PE experiment yields increasing intrinsic shifts and, finally, enters into the 1DOS regime. The intrinsic $`k_{}`$-resolution can be controlled using VLEED, which gives the final-state dispersions and lifetimes. The mechanisms forming the intrinsic shifts have been exposed, including non-linearity of the $`k_{}`$-dispersion, variations of the matrix element and surface effects. Appearance of the intrinsic shifts and transformation of the regimes of the PE experiment has been analysed using a physically transparent picture of interplay of the final- and initial-state spectral functions, and illustrated by PE simulations and experimental examples. The $`k_{}`$-broadening gives rise to other surprising effects typical of 3D valence states such as ’ghost’ peaks whose direct-transition $`k_{}`$ falls outside the Fermi surface, and finite width of the PE peaks at the Fermi level.

I thank R. Claessen, P.O. Nilsson and R. Feder for promoting discussions. Financial support by the Deutsche Forschungsgemeinschaft (CL124/5-1) is gratefully acknowledged.


text:
# The high-pressure 𝜶/𝜷 phase transition in lead sulphide (PbS)

## 1 Introduction

Among the IV-VI narrow gap semiconductors lead sulphide, PbS, has been a subject of interest during the last three decades, owing to its potential technological relevance, e.g. in long-wavelength imaging, as diode lasers Philipps2000  or thermovoltaic energy converters. Many of the experimental and theoretical investigations have been focused on the electronic structure Santoni92 ; Wei97 . Currently, PbS attracts new interest because of its importance for the understanding of the complex misfit layer compounds (MX)1+x(TX2)m (M=Sn, Pb, Sb, Bi, rare earth; T=Ti, V, Cr, Nb, Ta; X=S, Se; 0.08$`<`$x$`<`$0.28; m=1-3) Rouxel95 , of which PbS can be regarded as a component structure. These composite structures are single-crystal multi-layer systems, exhibiting strong anisotropic material properties.

At ambient conditions $`\alpha `$-PbS crystallises in the rock salt structure type (B1) Wasserstein51  with the lattice parameter $`a`$=5.9240(4) Å and space group $`Fm\overline{3}m`$. The elastic constants of $`\alpha `$-PbS were determined by ultrasonic velocity measurements Padaki81 . Resistivity measurements under high pressure on several lead chalcogenides revealed the existence of a phase transition Samara62  which occurs in PbS at 2.2 GPa. First structural studies under high pressure were performed using X-ray film methods Wakabayashi68 . Later, the compression behaviour of PbS was studied by means of high-pressure energy dispersive X-ray powder diffraction (HP-EDX) up to 35 GPa Peresada76 ; Chattopadhyay84 ; Chattopadhyay86 ; Quadri96 ; Jiang00 .

The phase transition at 2.2 GPa Samara62  to a $`\beta `$-phase having orthorhombic symmetry was characterised as a first order transition Wakabayashi68 ; Chattopadhyay84 ; Chattopadhyay86 . However, the structure of this phase is still a matter of debate. Two different structure types were proposed, viz. the GeS structure type (B16) with the space group Pbnm Samara62 ; Wakabayashi68 ; Quadri96 ; Jiang00  and the closely related, but higher symmetric CrB structure type (B33) with the symmetry Cmcm Chattopadhyay84 ; Chattopadhyay86 . A definitive decision between these two possibilities has not been possible until now on the basis of the experimental data available.

A second phase transition at about 25 GPa was detected by HP-EDX and assigned to the transformation into cubic $`\gamma `$-PbS with the CsCl (B2) type structure Chattopadhyay84 ; Chattopadhyay86 .

Since the equation-of-state parameters for $`\alpha `$-PbS, the symmetry of the $`\beta `$ phase and its compression behaviour are hitherto unknown it was found worth reexamining the PbS system under pressure. To elucidate the structure of the intermediate $`\beta `$ phase, angle dispersive X-ray powder diffraction and density functional theory (DFT) calculations were performed. Here the results of the study of the crystal structure of PbS under high pressure are presented. The mechanism of the phase transition is discussed using the symmetry relations between the space groups of the two phases.

## 2 Experimental methods

Crystals were grown by chemical vapour transport (CVT) in evacuated quartz tubes in a gradient furnace (1053 to 923 K) with stoichiometrically mixed elements as starting material and iodine as transport agent. In the low-temperature zone cube shaped crystals were obtained after several days with edge lengths up to 3 mm. The chemical composition of the crystals was checked by electron micro-probe analysis and the phase purity was studied by conventional X-ray powder diffraction.

The high-pressure experiments were performed using diamond anvil cells Merrill74 . NaCl was used as internal standard for the pressure determination by the Decker equation-of-state Decker71 . Several runs were performed using a 4:1 methanol/ethanol mixture for maintaining hydrostatic conditions. The highest pressure obtained was 6.8 GPa. Angular dispersive diffraction data were collected with Mo-$`K\alpha `$ radiation on a Mar2000 image plate diffractometer. The geometry correction for the radial integration of the two-dimensional data and the transformation into standard one-dimensional powder patterns were performed using Fit2d Hammersley96 . Whole powder pattern fitting employing Fullprof Rodriguez93  was performed for the precise determination of the unit cell parameters.

## 3 Computational aspects

Geometry optimisations were performed for $`\alpha `$-PbS and for both proposed structure types of $`\beta `$-PbS. The ab initio plane-wave CASTEP code Payne92 ; Milman2000  was used which is based on the density functional theory (DFT). The exchange correlation interaction was accounted for by using the Perdew, Burke and Ernzerhof version of the generalised gradient approximation (PBE-GGA) of the exchange-correlation functional Perdew92 ; White94 . Ultrasoft pseudopotentials Vanderbilt90 ; Kresse94  with a maximum cutoff energy of the plane waves of 290 eV were used to describe the electron-ion interaction. The integration over the Brillouin zone was performed employing the Monkhorst-Pack Monkhorst76  sampling scheme of reciprocal space, with a distance of 0.04 Å-1 between the sampling points. All structural parameters not constrained by the space group symmetry were relaxed.

The calculation of the elastic constants is based on applying small strains to the already relaxed structure, followed by relaxation of the atomic positions. In order to avoid anharmonic effects the maximum applied distortion was about 0.2 %. Several amplitudes for each strain component were used. Elastic constants were determined from a linear fit of stress to the applied strain.

## 4 Results and discussion

Some integrated powder patterns are presented in figure 1. The onset of a phase transition was detected at 2.5 GPa by the appearance of additional reflections. The new phase could be indexed with an orthorhombic unit cell and cell parameters $`a`$=3.98(2) Å, $`b`$=11.11(7) Å and $`c`$=4.16(4) Å. Coexistence of the two phases was observed up to a pressure of 6 GPa.

The powder patterns obtained above 2.5 GPa could be indexed with the same unit cell for both, B16 and B33, structure types. A structure refinement of the $`\beta `$-phase in either case turned out to be impossible because of the strong peak overlap with reflections stemming from the gasket material (Inconel steel). Thus, a discrimination between both possible structure types was not possible by the diffraction method.

### 4.1 Lattice parameters and elastic properties of $`𝜶`$\- and $`𝜷`$-PbS

The cell parameters of the orthorhombic $`\beta `$-phase are related to those of the cubic $`\alpha `$-phase by $`a_o=\sqrt{2}a_c/2,b_o=2a_c,c_o=\sqrt{2}a_c/2,`$ with subscripts $`o`$ and $`c`$ referring to the orthorhombic and the conventional cubic lattices, respectively. Figure 2a shows the dependence of (pseudo-) cubic lattice parameters on pressure for both phases.

The normalised volume changes are given in figure 2b. Because of the narrow pressure range observed in the experiment the volume at zero pressure V0 of $`\beta `$-PbS could not be obtained unambiguously from the equation-of-state (eos) fit. Instead, the volume at zero pressure from the ab initio calculations was used to normalise the experimental data of the $`\beta `$-phase.

The non-zero components of the calculated elastic stiffness tensor $`(c_{ij})`$ (in Voigt’s notation, units GPa) are $`c_{11}=141(3)`$, $`c_{12}=34(2)`$ and $`c_{44}=19.97(5)`$ for $`\alpha `$-PbS at ambient pressure. They agree well with the available experimental data ($`c_{11}=149`$, $`c_{12}=35`$ and $`c_{44}=29`$ GPa) determined at 5 K Padaki81 . For $`\beta `$-PbS no experimental data are available, the computed values at 4 GPa are

$$\left(\begin{array}{cccccc}\hfill 105(5)& \hfill 11(6)& \hfill 67(5)& & & \\ & \hfill 89(8)& \hfill 30(3)& & & \\ & & \hfill 103(7)& & & \\ & & & \hfill 21(1)& & \\ & & & & \hfill 71.02(6)& \\ & & & & & \hfill 12(2)\end{array}\right).$$

Linear bulk moduli determined from elastic constants are 70(2) GPa for $`\alpha `$-PbS and 53(3) GPa for $`\beta `$-PbS. The value of the linear volume bulk modulus of the $`\alpha `$-phase is higher than that of the high-pressure $`\beta `$-phase. This reflects the different compression behaviour in the two phases (section 4.3). The bulk moduli $`b_0`$ and their pressure dependencies $`b^{}`$ were determined for $`\alpha `$\- and $`\beta `$-PbS by fitting Birch Murnaghan equations-of-state Birch78  to the normalised volume data. The fit results are given in table 1 together with the linear compressibility data $`k`$.

The bulk modulus $`b_0`$, as well as the pressure dependence of the bulk modulus $`b^{}`$, for $`\alpha `$-PbS show excellent agreement between experimental and calculated values. Also, the linear compressibilities derived from a linear fit to experimental data and from calculated compliances agree excellently. The limited pressure range of the experimental data prevented us from determining $`b^{}`$ for $`\beta `$-PbS. Hence, the bulk modulus was determined with a second order Birch Murnaghan equation-of-state ($`b^{}:=4`$). It is well known that $`b_0`$ and $`b^{}`$ are highly correlated variables in eos-fitting Angel2000 . Consequently, $`b_0`$ is higher in the experiment, since $`b^{}`$ is smaller comparing the numbers between experiment and theory. The resulting bulk modulus is smaller than data reported recently Chattopadhyay86 . However, direct comparison is difficult since the authors of Chattopadhyay86  used a different equation of state.

### 4.2 Stability of the $`𝜷`$ phase

As stated earlier, the quality of our diffraction data did not allow the unambiguous assignment of the symmetry and, thus, discrimination between both possible structure types B16 and B33 for $`\beta `$-PbS. By contrast, the quantum mechanical calculations made the distinction possible, because only for the CrB-type structure (B33) convergence was achieved. The enthalpy for this phase with respect to the $`\alpha `$-phase is shown in figure 3 as a function of pressure. A transition pressure of 9.6(1) GPa for the $`\alpha /\beta `$-transformation follows from the crossover of the fitted curves. This transition pressure is higher than the experimentally determined values (onset of the transition at 2.6 GPa, this study; 2.2 GPa Samara62 ). However, DFT-GGA calculations such as those presented here suffer from systematic ’underbinding’ and the neglect of temperature effects. Furthermore, they refer to perfect ideal crystals. Hence, a discrepancy of 7 GPa between observed and calculated transition pressures can be tolerated. Also, the main aim of such calculations is not to reproduce numbers, but instead provide insight into the transition mechanism, and for this purpose such calculations are well suited.

### 4.3 Compression mechanism

In $`\alpha `$-PbS the lead atoms are sixfold coordinated by sulphur atoms in a regular octahedral arrangement. The only structural parameter in the B1 structure type not constrained by the space group symmetry is the magnitude of the cell parameter. The Pb–S bond length is $`a_c/2`$ and, of course, the compression mechanism results in shortening of the bond length.

No major change in the linear compressibility between cubic and orthorhombic PbS follows from the experimentally determined cell parameters (Fig. 2b). However, Pb–S bond lengths derived from the ab initio calculations (Fig. 4a) clearly show that the compression behaviour is anisotropic. Layered chalcogenides show a strongly anisotropic compression with the strongest compression perpendicular to the layers (e.g. TiS2 Allan1998 , SnS2 Knorr2001-1 , NbS2 Ehm2002-2 ). The main compression mechanism in these compounds is narrowing of the van-der-Waals gaps, whereas the layers themselves are bonded predominately ionic. The structure of $`\beta `$-PbS can be thought of as a layer structure as well, consisting of rock salt like bilayers, stacked along the $`b`$-axis of the orthorhombic unit cell. This is shown in figure 4b where two simple coordination polyhedra can be distinguished. Within the bilayers Pb is found in a fivefold coordination for which a square pyramid is an idealised geometry. The basal pinacoid face of the square pyramid forms together with two sulphur atoms from the next bilayer an almost trigonal prism. Figure 4a shows the strongest compression for this particular coordination polyhedron. The Pb–S(6,7) bonds connect the layers. The two different bonds within the bilayers are less compressible. While the pyramids deform almost isotropically, the main compression mechanism is a distortion of the trigonal prism, due to a change of the Pb–S(6,7) bond lengths. This mechanism effectively influences the length of the $`b`$-axis. Changes in the $`a`$ and $`c`$ axes are related to lateral distortions within the bilayers, i.e. changes of the Pb–S(2-5) bond lengths, and of the angles between them.

Admittedly, a description of the coordination in $`\beta `$-PbS by two different polyhedra is not fully appropriate since this approach doesn’t account for the increase of the coordination number with pressure. Alternatively, the coordination of the lead atoms can be described as sevenfold with the geometry of a single-capped trigonal prism Kepert1982 . The capped trigonal prism in $`\beta `$-PbS is shown in Figure 4c. The ideal geometry of a single-capped trigonal prism requires that all apices are located on the surface of a sphere and the edges $`d_i`$ fulfil the conditions given in table 2 Kepert1982 . The deviation of this coordination polyhedron from ideal can be quantified by applying e.g. the Pinsky and Avnir continuous-symmetry-measure (CSM) Pinsky1998 . Here a simplified CSM formulation was used.

From the five different edges $`d_i`$ radii $`r_i`$ were determined according to table 2 and from the weighted average radius

$$\overline{r}_w=\underset{i=1}{\overset{5}{}}n_ir_i/\underset{i=1}{\overset{5}{}}n_i,$$
(1)
the edge lengths of the corresponding ideal polyhedron $`d_i^{}`$ could be calculated. The average distortion measure ADM is then given by

$$\text{ADM}=\frac{_{i=1}^5n_i\left(d_id_i^{}\right)^2}{_{i=1}^5\left(n_id_i^{}\right)^2}\times 100$$
(2)
and the mean edge length distortion is

$$\overline{\mathrm{\Delta }}_d=\sqrt{\underset{i=1}{\overset{5}{}}n_i\left(d_id_i^{}\right)^2}/\underset{i=1}{\overset{5}{}}n_i.$$
(3)
The ADM is normalised with respect to the edge lengths of the ideal polyhedron in order to compensate for the effect of the compression. Hence, it is a dimensionless number. The compression of the coordination polyhedron is represented by the change of of the circumsphere’s radius $`\overline{r}_w`$ with pressure. The mean edge length distortion $`\overline{\mathrm{\Delta }}_d`$ is a measure of the variance of the edges. The resulting values for selected pressures are given in Table 3. At low pressures the conformation of the coordination polyhedron deviates most strongly from that of an ideal single-capped trigonal prism Kepert1982  (ADM=0). With increasing pressure ADM and $`\overline{\mathrm{\Delta }}_d`$ indicate increasing regularity of the coordination polyhedron. Above approximately 10 GPa some residual deformation remains virtually constant. The position of the central lead atom does not coincide with the conjoint centroid of the seven sulphur ligands (the difference $`\delta _y`$ is along the orthorhombic $`b`$ direction and given in table 3). With increasing pressure the position of Pb tends towards that of the sulphur’s centroid. A similar behaviour has recently been reported for $`\alpha `$-PbF2 Ehm2002-3 , where Pb is elevenfold coordinated by F in a tri-capped trigonal prism.

### 4.4 Characterisation of the $`𝜶/𝜷`$ phase transition

Coexistence of the two phases between 2.5 and 6 GPa was observed. Since the peak profile remained sharp throughout the whole pressure range covered by the experiment, pressure inhomogeneities can be excluded as the origin of this phase coexistence. Instead, it reflects the first order character of the transition (e.g. Chattopadhyay84 ; Chattopadhyay86 ). This is further supported by the observed discontinuity in the molar volume of 4.7 % at the phase transition.

The transition is reconstructive. This is justified by i) the change of the coordination number from six to seven, and ii.) because there is no direct group/subgroup relationship between the two space groups $`Fm\overline{3}m`$ and $`Cmcm`$. Figure 5 shows the symmetry relations in terms of klassengleiche (index $`k`$) and translationengleiche (index $`t`$) subgroups together with the changes of the unit cell dimensions. The crucial step here is the transition to space group $`Pbnm`$ realised by a displacement $`0<\mathrm{\Delta }a<1/2`$ of every second (001)-bilayer along the orthorhombic $`a`$-, i.e. cubic $`1\overline{1}0`$, direction. The space group $`Cmcm`$ is obtained as a limit state, when $`\mathrm{\Delta }a`$ takes the value of $`1/2`$. It is interesting to note that the displacement is the same as for the preferential glide system $`1\overline{1}0\{001\}`$ of NaCl. In the $`1\overline{1}0`$ directions easy displacements of $`\{001\}`$ layers are possible because movement of the ions occurs without having contact of equally charged ions.

The $`\alpha /\beta `$ transition in PbS is accompanied by a quadrupling of the unit cell. The wave vector of the transition is $`(0,0,\pi /a)`$ in the midpoint of the $`\mathrm{\Delta }`$-line connecting the $`\mathrm{\Gamma }`$-point at the zone centre with the $`X`$-point on the boundary of the Brillouin zone of the cubic face-centred lattice.

## 5 Summary and conclusion

The elastic behaviour of $`\alpha `$\- and $`\beta `$-PbS was measured up to 6.8 GPa. Bulk moduli, elastic constants and linear compressibilities were determined by DFT calculations and compared to experiment. The symmetry of $`\beta `$-PbS was found to be $`Cmcm`$. The differences in the compression behaviour between $`\alpha `$ and $`\beta `$-PbS can be understood on the basis of the different crystal structures and consequently, different compression mechanisms. PbS bilayers as found in the $`\beta `$-phase occur as MX subsystem of the chalcogenide misfit layer compounds Rouxel95 . Structural investigations under high pressure are not known for the misfit layer compounds. Therefore, our findings for the compression mechanism in $`\beta `$-PbS may serve at least as model for the MX subsystem in misfit layer compounds. Of course, for a better understanding of the elastic properties of the MX-TX2 multi-layer systems further investigations are required.

## Acknowledgement

This research was performed in the framework of the Kieler Forschergruppe ”Wachstum und Grenzflächeneigenschaften von Sulfid- und Selenid-Schichtstrukturen” funded by the German Science Foundation DFG (De412/21-1).


text:
# TWO PHOTON WIDTH OF HEAVY PSEUDOSCALAR MESONS

## 1 Introduction

We discuss the partial width of the pseudoscalar charmonium state $`\eta _c`$ and bottomonium state $`\eta _b`$ into two photons. Predictions from potential models are examined and compared with experimental values for the $`\eta _c`$ case. Through the NRQCD factorisation procedure results for $`\eta _c`$ are also compared with those from $`J/\psi `$ data, and results for $`\eta _b`$ to the $`\mathrm{{\rm Y}}`$ decay data.

## 2 Experimental values and relation to vector electromagnetic width

In this work we revisit the calculation of the two photon width of $`\eta _c`$, highlighting newest experimental results and updating the potential model calculation. Unlike the $`\eta _c`$, the $`\eta _b`$ state hasn’t been observed yet. We will examine various theoretical predictions for the electromagnetic decay of the lightest $`b\overline{b}`$ bound state. We start with the two photon decay width of a pseudoscalar quark–antiquark bound state with first order QCD corrections. The ratio of the pseudoscalar decay rate to the vector one is given by

$$\frac{\mathrm{\Gamma }(P\gamma \gamma )}{\mathrm{\Gamma }(Ve^+e^{})}3Q^2\frac{(13.38\alpha _s/\pi )}{(15.34\alpha _s/\pi )}=3Q^2\left[1+1.96\frac{\alpha _s}{\pi }+𝒪(\alpha _s^2)\right].$$
(1)
This expression can be used to estimate the radiative width of pseudoscalar state from the measured values of the leptonic decay width of the vector state.

## 3 Potential models predictions for $`\eta _c`$ and $`\eta _b`$ $`\gamma \gamma `$ decay width

We present the results one can obtain for the absolute width, through the extraction of the wave function at the origin from potential models. The “prototype” potential is given by the Cornell potential model

$$V(r)=\frac{k}{r}+\frac{r}{a^2}$$
(2)
which allows us to compute the Born decay width. The full expression of the pseudoscalar decay width is given by:

$$\mathrm{\Gamma }(P\gamma \gamma )=\mathrm{\Gamma }_{Born}^P\left[1+\frac{\alpha _s}{\pi }\left(\frac{\pi ^220}{3}\right)\right]$$
(3)
$`\mathrm{\Gamma }_{Born}^P`$ depends on $`\psi (0)`$ of the particular potential model taken into account.

## 4 Octet component procedure

We present after another procedure which admits other components to the meson decay beyond the one from the colour singlet picture (Bodwin, Braaten and Lepage). The decay width expression is given by means of NRQCD from the expansion

$$\mathrm{\Gamma }=\underset{n=1}{\overset{2}{}}\frac{2\mathrm{}f_n(\alpha _s)}{M^{d_n}}𝒪_n$$
(4)
which is a sum of terms, each of which factors into a short–distance perturbative coefficient $`\mathrm{}f_n`$ and a long–distance nonperturbative matrix element $`𝒪_n`$. From the experimental values of electromagnetic decay and light hadron of the vector state we obtain the two photon decay width of the pseudoscalar state. We compute also the decay width from the lattice calculation of the nonperturbative long–distance terms.

## 5 Conclusions

The $`\mathrm{\Gamma }(\eta _c\gamma \gamma )`$ prediction gives the value $`7.5\pm 1.6`$ keV. The $`\mathrm{\Gamma }(\eta _b\gamma \gamma )`$ prediction gives the value $`466\pm 101`$ eV. Prediction of the BBL procedure and other theoretical results are in good agreement with each other.