text:
# Slow temperature equilibration behind the shock front of SN 1006

## 1 Introduction

It has been suspected for many years that the high Mach number, collisionless shocks of young supernova remnants (SNRs) do not produce electron-ion temperature equilibration. The Rankine-Hugoniot relations require that for very high Mach numbers the temperature for each particle species $`i`$:

$$kT_i=\frac{2(\gamma 1)}{(\gamma +1)^2}m_iv_s^2=\frac{3}{16}m_iv_s^2,$$
(1)
with $`m_i`$ the particle mass, $`v_s`$ the shock velocity, and $`\gamma `$ the ratio of specific heats, usually taken to be $`5/3`$ (e.g. McKee & Hollenbach, 1980). The heating process in collisionless shocks is not well understood, but the Coulomb collision times are too long to provide the required heating. So other, collective, processes should be responsible for the heating. This raises the question whether the heating process leads to temperature equilibration or not, i.e. is the electron temperature very low compared to the proton temperature, which, according to (1), should be lower than the oxygen or iron temperature? If temperatures are not equilibrated at the shock front, and subsequent equilibration proceeds through Coulomb interactions, full equilibration takes $`10^{12}/n_\mathrm{e}`$ s (see Itoh, 1984).

A clear hint for non-equilibration is the low electron temperature in young SNRs, which in no object seem to exceed 5 keV, whereas a typical shock velocity of $`4000`$ km/s should give rise to a mean plasma temperature of 19 keV (e.g. McKee & Hollenbach, 1980; Hughes et al., 2000). X-ray observations usually allow only the electron temperature to be determined from the continuum shape and line ratios. The ion temperature is difficult to measure, as it does not alter the continuum shape, and hardly influences the ionization and excitation balance.

Some measurements of temperature non-equilibration based on optical and UV spectroscopy have been reported. The shock velocity and amount of electron-proton equilibration can be determined from the line widths and ratio of H$`\alpha `$ and H$`\beta `$ emission from non-radiative shocks. For the northwestern shock front in SN 1006, which is also the subject of this letter, Ghavamian et al. (2002) measured a shock velocity from the H$`\alpha `$ width of $`v_s=2890`$ km/s and inferred an electron to proton temperature ratio $`T_\mathrm{e}/T_\mathrm{p}<0.07`$. Together with the measured proper motion Winkler et al. (2003) were able to determine a distance to SN 1006 of $`2.17\pm 0.08`$ kpc. UV spectra obtained with the Hopkins Ultraviolet Telescope showed broad C IV, N V and O VI lines, indicating temperature non-equilibration of these elements (Raymond et al., 1995), and a comparison of these line intensities with He II allowed a measurement of $`T_\mathrm{e}/T_\mathrm{p}<0.05`$ (Laming et al., 1996). X-ray observations of SN 1006 indicate a low value of typically $`\mathrm{log}(n_\mathrm{e}t)9.5`$ (Vink et al., 2000; Dyer et al., 2001), implying that if the equilibration is governed by Coulomb interactions, the plasma did not have sufficient time to equilibrate.

Here we report on a direct X-ray measurement of the O VII temperature in the northwest of SN 1006. It confirms with a high statistical confidence the slow equilibration of electron and ion temperature, but for a position further downstream from the shock than for the optical and UV measurements.

## 2 Observations and methods

We observed SN 1006 (G327.6+14.6) with the XMM-Newton X-ray observatory (Jansen et al., 2001) with the goal of measuring the ionization and equilibration in SN 1006. The two observations are part of the guest observer program and were made on August 10 and 11, 2001. One observation pointed on the bright knot in the northwest, the main topic of this letter, the other on the narrow filaments in the east of SN 1006. The exposure times varied per instrument, but were close to 60 ks for both observations. In order to produce Fig. 1 we also used additional archival XMM-Newton and Chandra data.

Our main result is obtained with the Reflective Grating Spectrometer (RGS, den Herder et al., 2001). The two RGS instruments, RGS1 and RGS2 cover the wavelength range of $`640`$ Å, with a first order resolution of $`\lambda /\mathrm{\Delta }\lambda 300`$ at 20 Å. Emissions from different orders are separated using the intrinsic energy resolution of the CCDs. Unfortunately one CCD chip on the RGS1 and one on the RGS2 no longer function. As a result the RGS1 spectrum does not cover the wavelength range encompassing Ne line emission, and the RGS2 does not cover the wavelength range around 22 Å, which contains the O VII He$`\alpha `$ line emission. We limited the extraction region to an 1′ wide strip across the bright northwestern knot.

The extended emission from SN 1006 means that the standard response matrices had to be convolved with the spatial emission profile displayed in Fig. 1, with an attenuation due to off-axis telescope vignetting. The bright knot has a spatial width of 0.4′ (FWHM), which gives an apparent spectral broadening of $`\mathrm{\Delta }\lambda =0.05\AA `$ (FWHM), or $`\lambda /\mathrm{\Delta }\lambda =430`$ at 21.6 Å. The intrinsic resolution of the RGS is $`\mathrm{\Delta }\lambda =0.06`$ Å(den Herder et al., 2001). As the knot is compact and close to the edge of SN 1006 it is unlikely that any measurable broadening is due to bulk motions along the line of sight.

XMM-Newton also has CCD cameras behind each of its three mirrors, called the European Photon Imaging Camera (EPIC). We use the spectra extracted from the EPIC data to measure the electron temperature and abundances. All standard data reductions and response matrix calculations were done with SAS 5.3.3. Background spectra for the EPIC data were taken from a region outside SN 1006, but close to the center of the field. As RGS background spectra we used archival data of targets with no apparent line emission, such as gamma ray bursts.

## 3 The electron and ion temperatures in the northwest of SN 1006

The measurement of the non-equilibration of the electron and ion temperatures requires the measurement of both the ion temperature, here oxygen temperature, and the electron temperature.

The electron temperature is in this case most accurately measured from the spectral continuum shaped observed with the EPIC CCD spectra, as the low ionization time scales in SN 1006 make the available line ratios only weakly temperature dependent. The EPIC spectra from the northwestern knot were fitted with the SPEX non-equilibration ionization (NEI) code (Kaastra et al., 1996). As the knot is relatively compact, temperature and ionization gradients are probably of minor importance. The plasma parameters obtained by fitting the spectra indicate $`\mathrm{log}(n_\mathrm{e}t)=9.4`$, and $`kT_\mathrm{e}`$$`1.31.7`$ keV, higher than the $`kT_\mathrm{e}`$$``$0.7 keV reported by Long et al. (2003) based on Chandra data. Note, however, that the Chandra spectra have a lower spectral resolution, and have currently more calibration problems. Moreover, as there are two different kind of EPIC instruments, we were able to verify the consistency of the results (Fig. 2). Nevertheless, the measured $`kT_\mathrm{e}`$ is substantially lower than the $`kT_\mathrm{e}20`$ keV expected for a fully equilibrated plasma, but still higher than that expected from Coulomb equilibration alone behind a 3000 km/s shock. Fig. 3 shows the predicted $`kT_\mathrm{e}`$ against $`n_\mathrm{e}t`$ for current observations of the SN 1006 knot, as well as the various ion temperatures. Our fitted $`n_\mathrm{e}t`$ corresponds to plasma shocked 200-300 years ago, which is predicted to have $`kT_\mathrm{e}=600`$ eV. Hence we infer a small degree of collisionless electron heating (around 5% of the shock energy, i.e. $`T_\mathrm{e}/T_\mathrm{p}0.1`$, see Ghavamian et al. (2001)) consistent with optical and UV observations (Ghavamian et al., 2002; Laming et al., 1996) of SN 1006 and with observations of high Mach number shocks in other supernova remnants such as Tycho (Ghavamian et al., 2001) and SN 1987A (Michael et al., 2002).

The spectra are dominated by line emission from O VII, the other line complexes are, however, not from helium-like stages of Ne, Mg and Si, as identified by Long et al. (2003), but from lower ionization stages. For instance the EPIC spectra show consistently that the Mg line centroid is $`1.33\pm 0.01`$ keV, whereas the Mg XI line centroid is $`1.35`$ keV. The Si line centroid is $`1.80\pm 0.01`$ keV, which differs significantly from the Si XIII centroid of 1.85 keV. Instead these centroids indicate ionization stages around Mg IX and Si IX. The Ne centroid, as determined from the RGS2 spectrum, indicates $`0.900\pm 0.004`$ keV, consistent with Ne VII. These centroids provide clear evidence for extreme NEI conditions for the northwestern knot and corroborate the measured $`n_\mathrm{e}t`$ value.

In agreement with Long et al. (2003) we find that the knot emission indicates an overabundance of Si, but contrary to the Chandra data we do find O, Ne and Mg to be somewhat underabundant (Table 1). This may indicate that instead of shocked ISM the knot is actually ejecta or ISM mixed with ejecta. This does not diminish the evidence for temperature non-equilibration, but it makes it harder to combine the results reported here with the H$`\alpha `$ shock velocity measurements.

The RGS1 spectrum (Fig. 2) shows that the emission around 0.66 keV is dominated by O VII He$`\beta `$ emission (O VIII Ly$`\alpha :`$ O VII He$`\beta 1:1.6`$). Ne is detected in the RGS2 spectrum, but Mg and Si are too weak for the RGS instruments. No evidence for Fe XVII line emission is seen at 15.01 Å, 16.78 Å, 17.0 Å, or 17.10 Å, presumably because Fe has not yet reached the Fe XVII charge state. The O VII emission is well described by $`n=i1`$, line ratios that scale with the ratios of the oscillator strengths, as is expected for $`kT_\mathrm{e}`$ $`1`$ keV.

In order to measure the ion temperature through the thermal Doppler broadening we fitted the RGS1 spectrum in the range from $`21.022.3`$ Å, dominated by oxygen line emission, with six absorbed gaussians and a bremsstrahlung continuum with $`kT_\mathrm{e}`$$`=1.5`$ keV fixed to the continuum outside the fitted range ($`N_\mathrm{H}=6.8\times 10^{20}`$ cm<sup>-2</sup>, Dubner et al., 2002). The six gaussian components had centroids fixed at the energies of the bright O V, O VI, and O VII lines. The ratios of those lines were fixed according to calculations with the FAC atomic code (Gu, 2002), for $`kT_\mathrm{e}`$$`=1.5`$ keV and a grid of $`n_\mathrm{e}t`$ values between $`\mathrm{log}(n_\mathrm{e}t)=9.09.4`$. The line broadening was taken to be proportional to the line energies. Spectral fitting was done with spectral fitting package xspec (Arnaud, 1996). This allowed us to use the C-statistic, which is the maximum likelihood statistic appropriate for Poisson noise (Cash, 1979).

The best fits corresponded to $`\mathrm{log}(n_\mathrm{e}t)=9.18`$ with a maximum likelihood statistic of $`C=93.6`$ for 104 data bins, with a possible range of $`\mathrm{log}(n_\mathrm{e}t)=9.049.30`$. This is somewhat lower, but arguably more accurate, than the value derived from the CCD spectra. It is also consistent with the measured oxygen Ly$`\alpha `$ to He$`\beta `$ ratio with the RGS. For higher value of $`\mathrm{log}(n_\mathrm{e}t)`$ the emission would be much more dominated by the O VII resonance line at 21.6 Å, whereas the RGS1 spectrum indicates a substantial contribution of forbidden line emission at 22.1 Å due to O VI inner shell ionization and excitation.

The line broadening needed to obtain an acceptable fit is $`\sigma _E=3.4\pm 0.5`$ eV (68% confidence), if fitted within the wavelength range 21.4-22.6 Å. As Fig. 2 shows the broadening underpredicts the wing of the resonance line between 21.0-21.4 Å, which is a possible indication that the oxygen ions may not be completely thermalized (i.e. are not fully described by a Maxwellian distribution). The statistical confidence of the detection of line broadening is at the $`9.8\sigma `$ level ($`\mathrm{\Delta }C=96.1`$). The 99.73% (3$`\sigma `$) confidence range is $`2.26.5`$ eV. Even if we do not fix the line ratios and optimize each individual line, we obtain $`\sigma _E=2.5\pm 0.4`$ eV with a statistical confidence still at the $`6.5\sigma `$ level. During our analysis we have compared different methods of analyzing the data, with various wavelength fitting ranges. This has given us some feeling of the systematic uncertainties involved. We estimate that the systematic error on $`\sigma _E`$ is about 0.7 eV.

## 4 Discussion

Our analysis of CCD and grating data of a bright knot in SN 1006 obtained by XMM-Newton shows that $`kT_\mathrm{e}`$$`1.5`$ keV and $`\mathrm{log}(n_\mathrm{e}t)=9.049.30`$ The Doppler broadening at 574 eV is $`\sigma _E=3.4\pm 0.5`$ eV, which corresponds to an oxygen temperature of $`kT_{oxygen}=530\pm 150`$ keV (e.g. eq. 10.69 Rybicki & Lightman, 1979). This indicates a shock velocity of $`>4000`$ km/s, if no significant temperature equilibration has taken place, but allowing for the possibility of some adiabatic cooling (eq. 1). The implied shock velocity is higher than the $`3000`$ km/s indicated by recent H$`\alpha `$ measurements (Ghavamian et al., 2002), but given the systematic errors both measurements agree.

The excess emission in the wings of the O VII line emission, although possibly due to calibration uncertainties concerning the wings of the intrinsic RGS line profiles, may also be real and caused by a lack of complete oxygen thermalization. This is not unexpected as the oxygen self equilibration time at high temperatures is comparable to the oxygen-proton equilibration time.

To conclude, the measured high oxygen ion temperature is a clear indication that the shock heating processes resulted in only a small degree (5%) of electron-ion equilibration at the shock front, and that the subsequent equilibration process is slow. This also has some bearing on the acceleration of cosmic rays in high Mach number shocks, which are thought to be injected into the shock acceleration process from the pool of the thermal gas behind the shock front. A low electron-ion equilibration will make it relatively more difficult for electrons than for ions to be accelerated. However, the southwestern and northeastern limbs of SN 1006 are a prime examples of efficient electron acceleration as the X-ray emission from these parts is dominated by synchrotron emission (Koyama et al., 1995). Although we could not directly measure the electron temperature at those sites, the electron temperature close to the northwestern limb seems to be similar to that of northeastern knot. However, the notion of a thermal pool from which particles are accelerated is likely to be too simple, as recent simulations seem to indicate (Schmitz et al., 2002).

We thank John Peterson for making to us available his Monte Carlo Code which helped us verify the validity of our results. JV and MFG were supported for this work by the NASA through Chandra Postdoctoral Fellowship Award Number PF0-10011 & PF0-10014 issued by the Chandra X-ray Observatory Center, which is operated by the Smithsonian Astrophysical Observatory for NASA under contract NAS8-39073. JML was supported by the NASA Contract S 92540F from the XMM GI Program and by basic research funds of the Office of Naval Research.


format: text

text:
# Transport equations for a two-dimensional electron gas with spin-orbit interaction

## Abstract

The transport equations for a two-dimensional electron gas with spin-orbit interaction are presented. The distribution function is a $`2\times 2`$-matrix in the spin space. Particle and energy conservation laws determine the expressions for the electric current and the energy flow. The derived transport equations are applied to the spin-splitting of a wave packet and to the calculation of the structure factor and the dynamic conductivity.

Spin injection and coherent control of spins in various nanostructures represent two principal challenges for the field of spintronics. Recently, the amount of spintronics research has grown up extensively with the ultimate goal of applications to the quantum computing and information processing ALS . A number of spin-based devices have been designed and studied DD ; WAB ; KNT ; GBZ . Spin manipulation in such devices can be achieved by optical KA  or electric FKR ; OYB ; ZRK ; RPS  methods or by ferromagnetic gating fer . A controlled coupling between spin and orbital degrees of freedom is considered to be a particularly promising tool of efficient spin manipulation dating back to the seminal proposal by Datta and Das DD .

Spin-orbit interaction in two-dimensional electron gas (2DEG) confined at GaAs/AlGaAs, GaN/AlGaN or similar heterojunctions arises because of the quantum well asymmetry in the perpendicular \[$`z`$\] direction. The resulting perpendicular electric field leads to the coupling of spin to the electron momentum BR . The strength of this coupling can be experimentally tuned by a gate voltage NAT ; Koga .

Experimental advances in spin manipulation present a certain challenge to develop a proper theoretical description for various phenomena related to the spin-orbit interaction. In particular, modification of universal conductance fluctuations and weak localization has been studied in quantum dots H ; AF ; M . The phenomenon of weak localization has been considered in 2DEG as well Koga ; MZM . The Friedel oscillations in the presence of spin-orbit interactionCR  and the ac conductivity and the plasmon attenuationMCE  are calculated.

The principal goal of the present paper is to derive general transport equations for the spin-dependent distribution function of 2DEG including the effects of spin-orbital coupling. We assume that the spin-orbit interaction in a two-dimensional electron gas has the form,

$$H_{so}=\alpha (\widehat{\sigma }_xp_y\widehat{\sigma }_yp_x)+\beta (\widehat{\sigma }_xp_x\widehat{\sigma }_yp_y),$$
(1)
where the first term is the Bychkov-Rashba term BR  and the second term is the linear Dresselhaus (or anisotropy) term present in semiconductors with no bulk inversion symmetry D . The expression of Eq. (1) corresponds to the confinement along the (001) growth direction DK . Hereinafter we neglect cubic Dresselhaus terms. The free particle Hamiltonian can therefore be written in compact notations as,

$$H=[\frac{𝐩^2}{2m}\mu ]\widehat{\sigma }_0+\alpha _{ik}\widehat{\sigma }_ip_k,\alpha _{ik}=\left(\begin{array}{cc}\beta & \alpha \\ \alpha & \beta \end{array}\right).$$
(2)
Here, as usual, $`\widehat{\sigma }_0,\widehat{\sigma }_x,\widehat{\sigma }_y,\widehat{\sigma }_z`$ constitute the set of Pauli matrices. We use Latin subscripts for spatial coordinates and reserve Greek subscripts for the spin indexes.

Spectral properties. Before deriving the transport equation we describe briefly the spectral properties of the Hamiltonian (1). Its diagonalization is straightforward and reveals the existence of two spin-split subbands in the electron spectrum,

$$ϵ_{1𝐩}=\xi _p+\mathrm{\Delta }_𝐩,ϵ_{2𝐩}=\xi _p\mathrm{\Delta }_𝐩,$$
(3)
where

$$\xi _p=\frac{p^2}{2m}\mu ,\mathrm{\Delta }_𝐩=\sqrt{p^2(\alpha ^2+\beta ^2)+4\alpha \beta p_xp_y},$$
(4)
and $`p^2=p_x^2+p_y^2`$ is the total electron momentum. The eigenfunctions corresponding to the eigenstates (3) are,

$$\psi _{1,2}(𝐱)=\frac{1}{\sqrt{2}}\left(\begin{array}{c}e^{i\chi _𝐩/2}\\ \pm e^{i\chi _𝐩/2}\end{array}\right)e^{i\mathrm{𝐩𝐱}},$$
(5)
with the upper (lower) sign corresponding to the $`\psi _1`$ ($`\psi _2`$) state. The phase factor $`\chi _𝐩`$ depends on the direction of the electron momentum,

$$\mathrm{tan}\chi _𝐩=\frac{\alpha p_x+\beta p_y}{\alpha p_y+\beta p_x}.$$
(6)
For the isotropic 2DEG ($`\beta =0`$) the phase $`\chi _𝐩`$ coincides with the angle between the electron momentum and the $`y`$-axis.

Further and more convenient description of the spectral properties can be obtained by considering the spin-dependent retarded Green function, defined as usual by,

$$iG_{\alpha \beta }^R(x,x^{})=\theta (tt^{})\psi _\alpha (x)\psi _\beta ^{}(x^{})+\psi _\beta ^{}(x^{})\psi _\alpha (x).$$
(7)
Here we have used the shorthand notation for the space and time variables, $`x=(𝐱,t)`$. In a homogeneous system the correlation functions depend on the relative coordinates $`xx^{}`$ only. Using the above expressions (3-5) we can write the expression for the retarded Green function of free electrons in the momentum representation, $`\widehat{G}^R(ϵ,𝐩)=𝑑t𝑑𝐱\widehat{G}^R(x)e^{iϵti\mathrm{𝐩𝐱}}`$, which after simple transformations takes the form,

$`\widehat{G}^R(ϵ,𝐩)`$ $`=`$ $`{\displaystyle \underset{\mu =1,2}{}}\widehat{G}_\mu ^R(ϵ,𝐩)={\displaystyle \frac{1}{2}}{\displaystyle \frac{\widehat{\sigma }_0+\mathrm{cos}\chi _𝐩\widehat{\sigma }_x\mathrm{sin}\chi _𝐩\widehat{\sigma }_y}{ϵϵ_{1𝐩}+i\eta }}`$ (8)
$`+{\displaystyle \frac{1}{2}}{\displaystyle \frac{\widehat{\sigma }_0\mathrm{cos}\chi _𝐩\widehat{\sigma }_x+\mathrm{sin}\chi _𝐩\widehat{\sigma }_y}{ϵϵ_{2𝐩}+i\eta }},`$
with the indexes $`\mu =1,2`$ corresponding to the first and second terms, respectively, in the last expression of Eq. (8).

The central quantity in the transport theory is the density matrix,

$$f_{\alpha \beta }(x,x^{})=\psi _\beta ^{}(x^{})\psi _\alpha (x).$$
(9)
Its value in the thermal equilibrium is related to the imaginary part of the retarded Green function via the fluctuation-dissipation theorem,

$$\widehat{f}(ϵ,𝐩)=\underset{\mu =1,2}{}n_{\mu 𝐩}[\widehat{G}_\mu ^R(ϵ,𝐩)G_\mu ^R(ϵ,𝐩)],$$
(10)
here $`n_{\mu 𝐩}`$ is the Fermi-Dirac distribution function for the $`\mu `$-th state (3). It is convenient to expand the density matrix over the complete set of the Pauli matrices. In particular, for the equal-time density matrix,

$$\widehat{f}_𝐩=\frac{dϵ}{2\pi }\widehat{f}(ϵ,𝐩)=\frac{1}{2}f_𝐩\widehat{\sigma }_0+\frac{1}{2}𝐠_𝐩\widehat{𝝈},$$
(11)
we observe according to Eq. (10) that in the thermal equilibrium,

$$\begin{array}{cc}f_𝐩=(n_{1𝐩}+n_{2𝐩}),\hfill & g_{𝐩x}=\mathrm{cos}\chi _𝐩(n_{1𝐩}n_{2𝐩}),\hfill \\ g_{𝐩z}=0,\hfill & g_{𝐩y}=\mathrm{sin}\chi _𝐩(n_{1𝐩}n_{2𝐩}),\hfill \end{array}$$
(12)
Transport equations. In a generic nonequilibrium state, the density matrix (9) obeys a set of conjugated equations that can be obtained from the equations of motion for the electron operators $`\psi (x)`$ and $`\psi ^{}(x)`$ determined by the Hamiltonian (2),

$`\left[i_t+{\displaystyle \frac{^2}{2m}}+\mu e\varphi _x\right]\widehat{f}(x,x^{})+i\alpha _{ik}\widehat{\sigma }_i_k\widehat{f}(x,x^{})=0,`$
$`\left[i_t^{}{\displaystyle \frac{^2}{2m}}\mu +e\varphi _x^{}\right]\widehat{f}(x,x^{})+i\alpha _{ik}_k^{}\widehat{f}(x,x^{})\widehat{\sigma }_i=0.`$
The equations (Transport equations for a two-dimensional electron gas with spin-orbit interaction) neglect impurity scattering. This is justified for ballistic systems when the mean free path exceeds the characteristic system size, e.g. in high-mobility 2DEG in semiconductor heterostructures (we discuss the impurity scattering at the end of the paper). In Eq. (Transport equations for a two-dimensional electron gas with spin-orbit interaction) we allowed for the scalar external field $`\varphi _x=\varphi (𝐱,t)`$. In the absence of electron-impurity or electron-electron collisions no self-energy terms appear in the equations (Transport equations for a two-dimensional electron gas with spin-orbit interaction) which makes it sufficient to consider the equal-time ($`t=t^{}`$) functions only.

Following the known route of deriving kinetic equations LP  we utilize the Wigner transformation for the density matrix,

$$\widehat{f}_𝐩(𝐱,t)=𝑑𝐫e^{i\mathrm{𝐩𝐫}}\widehat{f}(𝐱+\frac{𝐫}{2},𝐱\frac{𝐫}{2},t).$$
(14)
By taking the sum of equations (Transport equations for a two-dimensional electron gas with spin-orbit interaction) in the Wigner representation we obtain,

$`[_t+𝐯]\widehat{f}_𝐩+ie{\displaystyle 𝑑𝐪\varphi _𝐪(\widehat{f}_{𝐩\frac{𝐪}{2}}\widehat{f}_{𝐩+\frac{𝐪}{2}})e^{i\mathrm{𝐪𝐱}}}`$
$`+i\alpha _{ik}p_k[\widehat{\sigma }_i,\widehat{f}_𝐩]+{\displaystyle \frac{1}{2}}\alpha _{ik}_k\{\widehat{\sigma }_i,\widehat{f}_𝐩\}=0,`$ (15)
where $`𝐯=𝐩/m`$. Here we introduced the spatial Fourier transform for the scalar potential $`\varphi _𝐪=𝑑𝐱\varphi (𝐱,t)e^{i\mathrm{𝐪𝐱}}`$, with the shorthand notation for the momentum integration $`d𝐪=d^2q/(2\pi )^2`$.

Finally, to present Eq. (Transport equations for a two-dimensional electron gas with spin-orbit interaction) in a more transparent way we turn to the Pauli matrix representation (11) to write,

$`[_t+𝐯]f_𝐩`$ $`+`$ $`ie{\displaystyle 𝑑𝐪\varphi _𝐪(f_{𝐩\frac{𝐪}{2}}f_{𝐩+\frac{𝐪}{2}})e^{i\mathrm{𝐪𝐱}}}`$ (16)
$`+`$ $`\alpha _{ik}_kg_{𝐩i}=0,`$
$`[_t+𝐯]g_{𝐩i}`$ $`+`$ $`ie{\displaystyle 𝑑𝐪\varphi _𝐪(g_{𝐩\frac{𝐪}{2}i}g_{𝐩+\frac{𝐪}{2}i})e^{i\mathrm{𝐪𝐱}}}`$ (17)
$``$ $`[𝐛_𝐩\times 𝐠_𝐩]_i+\alpha _{ik}_kf_𝐩=0,`$
with the following notation for the precession frequency, $`b_{𝐩i}=2\alpha _{ik}p_k=2\mathrm{\Delta }_𝐩(\mathrm{cos}\chi _𝐩,\mathrm{sin}\chi _𝐩,0)`$.

Conservation laws. The transport equations (16,17) are of the Boltzmann type and therefore fulfill certain particle and energy conservation conditions which will now be obtained. By integrating Eq. (16) with respect to the momentum we find the continuity equation for the particle flow,

$$_t\rho +\frac{1}{e}𝐣=0,$$
(18)
where the electron density and the electric current are given respectively by,

$`\rho ={\displaystyle 𝑑𝐩f_𝐩},j_k=e{\displaystyle 𝑑𝐩[v_kf_𝐩+\alpha _{ik}g_{𝐩i}]}.`$ (19)
The terms containing the external potential $`\varphi _𝐪`$ cancel as is readily seen by the change of integration variables. To obtain the energy continuity condition we multiply Eq. (16) by $`\xi _p`$ and Eq. (17) by $`𝐛_𝐩`$ and add them together. After simple transformations the conservation of energy can be written in the conventional form,

$$_t\rho ^ϵ+𝐣^ϵ=𝐣𝐄,$$
(20)
where the energy density and energy current are,

$`\rho ^ϵ`$ $`=`$ $`{\displaystyle 𝑑𝐩[\xi _pf_𝐩+b_{𝐩i}g_{𝐩i}]},`$
$`j_k^ϵ`$ $`=`$ $`{\displaystyle 𝑑𝐩v_k[\xi _pf_𝐩+b_{𝐩i}g_{𝐩i}]}+\alpha _{ik}{\displaystyle 𝑑𝐩[\xi _pg_{𝐩i}+b_{𝐩i}f_𝐩]}.`$
The equation (20) means that the local energy change is due to the energy flow to the neighboring points in space as well as a result of the local Joule heating (right-hand side).

Wave packet splitting. To give a specific application of the derived equations let us now use them to describe the propagation of a wave packet in 2DEG with a spin-orbit coupling. We neglect a spin-orbit anisotropy $`\beta =0`$ for simplicity. The wave packet propagates along the $`y`$-direction and is uniform along the $`x`$-axis. The transport equations (16,17) are then one-dimensional and (with no external field applied) yield,

$`\left[_t+v_y\right]f=\alpha _yg_x,`$
$`\left[_t+v_y\right]g_x=\alpha _yf,`$
$`\left[_t+v_y\right]g_y=2\mathrm{\Delta }_pg_z`$
$`\left[_t+v_y\right]g_z=2\mathrm{\Delta }_pg_y.`$ (21)
First, we consider a spin-unpolarized Gaussian wave packet injected at the point $`y=0`$ at the time $`t=0`$ and moving with the average momentum $`\overline{p}`$,

$$\widehat{f}_𝐩(𝐱,t=0)=\widehat{\sigma }_0F(y),F(y)=e^{y^2\delta p^2\frac{(p_y\overline{p})^2}{\delta p^2}}.$$
(22)
In this geometry the phase factor $`\chi _𝐩=0`$, which means that the precession vector $`𝐛`$ is directed along the $`x`$-axis. We also observe that $`g_y=g_z=0`$. The remaining two of the equations (Transport equations for a two-dimensional electron gas with spin-orbit interaction) are easily solved by Fourier transforming them into a set of linear algebraic equations. The general solution of Eqs. (Transport equations for a two-dimensional electron gas with spin-orbit interaction) takes the form,

$`f(y,t)=A(yv_+t)+B(yv_{}t),`$
$`g_x(y,t)=A(yv_+t)B(yv_{}t),`$
where we have introduced subband velocities $`v_\pm =v\pm \alpha `$. So far, $`A(x)`$ and $`B(x)`$ are two arbitrary functions which have to be determined from the initial condition (Transport equations for a two-dimensional electron gas with spin-orbit interaction) yielding, $`A(y)=B(y)=F(y)`$. We find that the incident wave packet (22) is decomposed into two independent constituents oppositely polarized along $`x`$-direction and moving with different velocities. The spatial distribution of the electron density is given by the integral over all momenta, i.e.

$`\rho _\pm (y,t)={\displaystyle \frac{1}{2\pi ^{1/2}\delta x(t)}}\mathrm{exp}[{\displaystyle \frac{(x\overline{v}_\pm t)^2}{\delta x^2(t)}}],`$ (23)
with the average velocities $`\overline{v}_\pm =\overline{p}/m\pm \alpha `$, and the Gaussian width at finite times, $`\delta x^2(t)=\delta p^2+\frac{t^2\delta p^2}{m^2}`$. To observe the spin-orbit induced splitting of a wave packet the following conditions should be satisfied,

$$\frac{\delta p}{m}\alpha \frac{t}{\delta p}.$$

The first of the two conditions ensures that the splitting dominates over the wave packet broadening, while the second condition means that enough time has to elapse before the splitting becomes larger than the intrinsic packet width.

Now let us consider an injection of a packet initially polarized along the $`y`$-direction.

$$\widehat{f}_𝐩(𝐱)=(\widehat{\sigma }_0+\widehat{\sigma }_y)F(y),$$
(24)
The equations for the $`f`$ and $`g_x`$ components of the density matrix remain unchanged with the above analysis still valid. The second pair of Eqs. (Transport equations for a two-dimensional electron gas with spin-orbit interaction) is independent of the first pair and have a solution,

$`g_y(y,t)=F(yvt)\mathrm{cos}(2\mathrm{\Delta }_pt),`$
$`g_z(y,t)=F(yvt)\mathrm{sin}(2\mathrm{\Delta }_pt).`$ (25)
According to the expressions (Transport equations for a two-dimensional electron gas with spin-orbit interaction) the initial spin polarization precesses with a frequency $`2\mathrm{\Delta }_p`$ around the axis perpendicular to the propagation direction. Note that the precessing spin propagates with the center-of-mass velocity $`\overline{v}`$ rather than with the subband velocities $`\overline{v}_\pm `$.

The above analysis assumes that a wavepacket is injected with a given momentum $`\overline{p}`$. Such an injection into 2DEG with a spin-orbit coupling is not easy to achieve. For example, injection through an interface with a ’normal’ (with no spin-orbit interaction) 2DEGHM ; MH ; LLF  would not result in a spatial splitting of a wave packet. This is due to the fact that the injection happens with a conservation of energy rather than momentum. As seen from Eqs. (3-4) the two states with the same energy propagate with the same velocityMSB ; LLF  within the approximations of this paper. However, if we take into account the cubic Dresselhaus terms, which have been omitted in our discussion, there can be a splitting of velocities at the same energy. In order to achieve splitting without the cubic terms we need to consider a more complicated setup. As a demonstration of principle, we consider the following example. Let us inject a wave packet propagating along the $`y`$-direction with the spin polarized along the interface ($`x`$-axis), e.g. by injection from a ferromagnetic contact. The states forming the wavepacket belong to the subband 1, with a spin polarization $`\frac{1}{\sqrt{2}}(1,1)`$. Let us now switch on ac magnetic field along the $`y`$-axis rotating the spin direction until it is aligned with the $`z`$-axis, $`(1,0)`$, and then switch the magnetic field off. The resulting state will be an equal mixture of both eigenstates $`\frac{1}{\sqrt{2}}(1,1)`$ and $`\frac{1}{\sqrt{2}}(1,1)`$ without any change of momentum (the energy is no longer conserved). The velocities of these states are different and the packet will split.

The above picture holds not only for the injection of initially polarized packet. If the incident packet is unpolarized and has a given energy, upon entering the interface it will become a mixture of two states: $`\frac{1}{\sqrt{2}}(1,1)`$ with the momentum $`p_0m\alpha `$, and $`\frac{1}{\sqrt{2}}(1,1)`$ with the momentum $`p_0+m\alpha `$. Both velocities remain equal to $`v_0=p_0/m`$. After switching on the ac magnetic field with the frequency $`\omega 2\alpha p_0`$ (which is a resonant frequency for the transition between the two subbands), the first state will evolve into the mixture of the states: $`\frac{1}{\sqrt{2}}(1,1)`$ and $`\frac{1}{\sqrt{2}}(1,1)`$, both with the momentum $`p_0+m\alpha `$, meaning two different velocities $`v_0`$ and $`v_02\alpha `$. The same reasoning shows that the other initial state will develop two velocities $`v_0`$ and $`v_0+2\alpha `$. Therefore, the initially unpolarized packet will split into three parts.

Ballistic spin injection. We envisage a spin injection from ferromagnetic contacts into ballistic 2DEG among the applications for the equations derived above. In this case the injection occurs with conservation of energy, and can be described by the time-independent solution of the equations (16-17) with the appropriate boundary conditions, which require a conservation of the normal components of the electric current, Eq. (19), at the interfaces. A corresponding theory would generalize the existing approach for the ballistic spin-injection based on the ordinary Boltzmann equationRashba . In the latter case the Boltzmann equation method is more convenient for the calculation of spin polarization of current and magnetoresistance than the direct solution of the single-particle Schrödinger equation.

Structure factor. The electron density fluctuations are described by the structure factor PN  defined as the retarded correlation function,

$`\chi (x,x^{})=i\theta (tt^{})\rho (x)\rho (x^{})\rho (x^{})\rho (x),`$ (26)
of the electron density operators $`\rho (x)=\psi _\alpha ^{}(x)\psi _\alpha (x)`$. At equilibrium the structure factor (26) depends on the relative coordinates $`xx^{}`$ only. The imaginary part of the Fourier transform $`\chi (\omega ,𝐪)`$ measures the energy dissipation of the external field at a given frequency $`\omega `$ and a wavevector $`𝐪`$. In the isotropic system the structure factor is related to the ac conductivity by the relation,

$`\sigma (\omega )=\underset{q0}{lim}{\displaystyle \frac{ie^2\omega }{q^2}}\chi (\omega ,𝐪).`$ (27)
The formula (27) is readily checked using the Kubo formula for the conductivity and the continuity equation (18).

According to the fluctuation-dissipation theorem, the structure factor can be determined from calculations of the linear response to an external scalar field. The field-induced modulation of electron density is related to the magnitude of the external perturbation through the structure factor according to PN ,

$`\delta \rho (\omega ,𝐪)=e\chi (\omega ,𝐪)\varphi (\omega ,𝐪).`$ (28)
The electron density modulation is given by the deviation of the function $`f_𝐩(t,𝐱)`$ from its equilibrium value, $`\delta \rho (\omega ,𝐪)=𝑑𝐩\delta f_𝐩(\omega ,𝐪),`$ and can be found from the linearized equations (16,17). In the linear approximation by the external field $`\varphi (\omega ,𝐪)`$, the distribution function is a small deviation

$`f_𝐩=f_𝐩^0+\delta f_𝐩,𝐠_𝐩=𝐠_𝐩^0+\delta 𝐠_𝐩,`$ (29)
from its equilibrium value (12). The linearized transport equations (16,17) take the form,

$`(\omega \mathrm{𝐪𝐯})\delta f_𝐩`$ $``$ $`\alpha _{ik}q_i\delta g_{𝐩k}=e\varphi (\omega ,𝐪)(f_{𝐩+\frac{𝐪}{2}}^0+f_{𝐩\frac{𝐪}{2}}^0),`$
$`(\omega \mathrm{𝐪𝐯})\delta g_{𝐩i}`$ $``$ $`i[𝐛_𝐩\times \delta 𝐠_𝐩]_i\alpha _{ik}q_k\delta f_𝐩=`$ (30)
$``$ $`e\varphi (\omega ,𝐪)(g_{𝐩+\frac{𝐪}{2}i}^0+𝐠_{𝐩\frac{𝐪}{2}i}^0).`$
Solving these equations for the variation of the electron density (19) we obtain the structure factor with the help of the relation (28),

$`\chi (\omega ,𝐪)`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \underset{\mu \mu ^{}}{}}{\displaystyle 𝑑𝐩[1+(1)^{\mu \mu ^{}}\mathrm{cos}(\chi _𝐩\chi _𝐩^{})]}`$ (31)
$`\times {\displaystyle \frac{n_{\mu 𝐩_{}}n_{\mu ^{}𝐩_+}}{\omega ϵ_{\mu ^{}𝐩_+}+ϵ_{\mu 𝐩_{}}}},`$
where $`𝐩_\pm =𝐩\pm 𝐪/2`$. The expression (31) with $`\omega =0`$ corresponds to the previously derived result for the static dielectric functionCR . To simplify further the subsequent discussion we will disregard the anisotropy, $`\beta =0`$, and consider the zero-temperature limit $`T=0`$. The two spin-orbit subbands are axially symmetric, shown on Fig. 1. The subbands are filled up to the same Fermi energy level $`ϵ_F`$ but have two different Fermi momenta, $`p_1`$ and $`p_2`$, determined from the equations $`ϵ_i(p_i)=ϵ_F`$, where $`ϵ_i(p)`$ are given by Eqs. (3,4) with $`\beta =0`$. This leads to the values,

$`p_1=p_0m\alpha +O(m^2\alpha ^2/p_0^2),`$
$`p_2=p_0+m\alpha +O(m^2\alpha ^2/p_0^2),`$ (32)
where $`p_0`$ is determined by $`ϵ_F=p_0^2/2m`$, namely $`p_0`$ is the Fermi-momentum in the absence of spin-orbit interaction. Note that the Fermi velocities for the two subbands,

$$v_i=\frac{ϵ_i(p)}{p}|_{p=p_i}=\frac{p_0}{m}+O(m^2\alpha ^2/p_0^2),$$
(33)
are the same and (up to higher order terms) equal to the Fermi velocity in 2DEG with no spin-orbit coupling $`\alpha =0`$.

The imaginary part of the structure factor $`\chi (\omega ,𝐪)`$ determines the absorption, or Landau damping, of the external field at given frequency and wavevector. The points in the electron momentum space that contribute to the Landau damping correspond to the zeros of the denominators. There are total four determined by the equations,

$$\omega =\mathrm{𝐪𝐯}\pm \alpha p_+\pm \alpha p_{},$$
(34)
with the opposite signs of the last two terms corresponding to the (gapless) transitions within the same subbands \[Eq. (3)\] and equal signs describing the transitions between different subbands.

The terms with $`\mu =\mu ^{}`$ in Eq. (31) represent the effect of intrasubband transitions. Only indirect ($`q0`$) transitions contribute to the imaginary part of the structure factor. For small transferred momenta $`qp_0`$ one can disregard the deviation of the cosine factor from unity and also approximate $`n_{i𝐩_{}}n_{i𝐩_+}𝐪n_i/𝐩`$. Taking the momentum integral we obtain for the contribution of the $`i`$-th subband

$`\mathrm{}\chi _i(\omega ,q)`$ $`=`$ $`\nu _i{\displaystyle \frac{\omega }{\sqrt{q^2v_0^2\omega ^2}}}\theta (q^2v_0^2\omega ^2),`$ (35)
where $`v_0=p_0/m`$ and $`\nu _i`$ stands for the density of states of the $`i`$-th subband at its Fermi surface $`p=p_i`$: $`\nu _1=\frac{m}{2\pi }(1\frac{m\alpha }{p_0})`$, $`\nu _2=\frac{m}{2\pi }(1+\frac{m\alpha }{p_0})`$. Note that the sum of the two contributions (35) is independent of the spin-orbit interaction (up to higher-order terms), a consequence of the fact that the two subbands have the same value of the Fermi-velocity. Spin-orbit interaction results only in a redistribution of the spectral weight between the subbands controlled by the changes in the densities of states.

The terms with $`\mu \mu ^{}`$ in Eq. (31) correspond to the intersubband transitions. Their contribution to the structure factor for $`m\alpha qp_0`$ is negligible compared to the above considered intrasubband transitions by the factor $`q^2/p_0^2`$ (due to the small $`\mathrm{sin}^2`$ prefactor). However, the presence of the two subbands is important as it makes the direct, $`q=0`$, transitions possible. The factor $`n_{1p}n_{2p}`$ then defines the momentum space available for the direct transitions, $`p_1<p<p_2`$ (see Fig. 1), which corresponds to the frequency domain $`2\mathrm{\Delta }_02m\alpha ^2<\omega <2\mathrm{\Delta }_0+2m\alpha ^2`$, where $`\mathrm{\Delta }_0=\mathrm{\Delta }_{p_0}`$,

$$\chi (\omega ,q0)=\frac{\alpha q^2}{4\pi }\underset{p_1}{\overset{p_2}{}}\frac{dp}{(\omega +i0)^24\mathrm{\Delta }_p^2}.$$
(36)
The imaginary part of this expression is

$$\mathrm{}\chi (\omega ,q0)=\frac{q^2\text{sgn}\omega }{32\mathrm{\Delta }_0}\theta [4m^2\alpha ^4(\omega 2\mathrm{\Delta }_0)^2].$$
(37)
The equation (36) corresponds to the previously obtained result MCE  for the optical conductivity $`\sigma (\omega )`$. The expression (36) goes to zero with the wavevector, which is easily understood by noting that the matrix elements for the transitions between $`\psi _1(𝐩_{})`$ and $`\psi _2(𝐩_+)`$ states are suppressed at small transferred momenta since they are orthogonal at $`q=0`$. However, their contribution to the conductivity \[according to Eq. (27)\] remains finite, which is clear since the operator of electron velocity has nonzero matrix elements for the intersubband transitions even at $`q=0`$.

The experimental observation of the direct transitions (37) is feasible in the measurements of the resonant microwave absorption in high-mobility semiconductor heterostructures.

Screened electron-electron interaction and plasmon excitations. So far our analysis has been restricted to the noninteracting electron gas. To incorporate the effects of the electron-electron interaction in the random phase approximation one has to account for the self-consistent electric field induced by the variations of the electron density. The potential for this field $`\varphi _{sc}`$ obeys the Poisson equation. In two dimensions the Fourier transform of the Poisson equation has the form,

$$e\varphi _{sc}(\omega ,𝐪)=V_q\rho (\omega ,𝐪),$$
(38)
where $`V_q=2\pi e^2/q`$ is the bare Coulomb propagator. The random phase approximation (RPA) is then equivalent to the substitution $`\varphi (\omega ,𝐪)\varphi _{sc}(\omega ,𝐪)+\varphi (\omega ,𝐪)`$ in the right-hand side of Eq. (Transport equations for a two-dimensional electron gas with spin-orbit interaction). It is straightforward to see that the structure factor takes the familiar RPA form,

$$\chi _{RPA}(\omega ,𝐪)=\frac{\chi (\omega ,𝐪)}{1V_q\chi (\omega ,𝐪)}.$$
(39)
The pole of this expression determines the plasmon spectrum $`\omega =\omega _q+i\gamma _q`$, where $`\omega _q^2=v^2\kappa q/2`$, with $`\kappa =2\pi e^2\nu `$ standing for the static screening radius. The plasmon linewidth is given by the imaginary part of the bare structure factor at $`\omega =\omega _q`$,

$$\gamma _q=\frac{1}{2}V_q\omega _q|\mathrm{}\chi (\omega _q,q)|.$$
(40)
For the plasmon to be an undamped excitation its frequency should lie above the electron-hole continuum, $`\omega _q>qv`$, which requires $`\kappa >2q`$. As was already pointed out in Ref. MCE,  the plasmon acquires damping when $`\omega _q2\mathrm{\Delta }`$. Since $`q(q\kappa )^{1/2}m\alpha `$ at this range, the direct transitions (37) make the principal contribution to Eq. (40).

Impurity scattering. The equations presented in this paper assume ballistic electron motion. The absence of impurities allows one to write kinetic equation as a closed set of equations, Eq. (Transport equations for a two-dimensional electron gas with spin-orbit interaction), for the density matrix integrated over the energy variable $`ϵ`$, \[see Eq. (11)\], i.e. at coinciding times. In the presence of disorder the self-energy due to impurity scattering should be added to the right-hand side of Eq. (Transport equations for a two-dimensional electron gas with spin-orbit interaction). In general, since plain waves are no longer eigenstates of the system with impurities, the equations for the distribution function depending on the momentum $`𝐩`$ (and not on the energy $`ϵ`$) become not very convenient. More natural (though more complicated) equations would result from integration over $`\xi _p`$, similar to the usual spin-degenerate caseRS . Such equations are beyond the scope of the present paper.

Special case $`\alpha =\pm \beta `$. Recently, Schliemann et. al. SEL  proposed a spin field-effect transistor based on a particular tuning of the spin-orbit coupling constants such that $`\alpha =\beta `$ (or $`\alpha =\beta `$). This special system is expected to preserve spin coherence even in the presence of disorder. This is due to the fact that the spin eigenstates (5) are independent of the electron momenta, $`\chi _𝐩=\text{const}`$, therefore a scalar impurity potential does not result in the intersubband transitions. The same observation holds for the structure factor. Since the matrix elements of the density are identically zero for the transitions between different subbands the second line of Eq. (31) is absent in this case and the structure factor is intact by the presence of spin-orbit interaction (up to higher order corrections).

Conclusions. To summarize, we have derived transport equations for the distribution function of a two-dimensional electron gas with spin-orbit interaction of both the Bychkov-Rashba and the Dresselhaus mechanisms. The distribution function is a 2$`\times `$2-matrix in the spin space. General expressions for the particle and energy currents and densities are available in terms of the density $`f_𝐩`$ and spin $`𝐠_𝐩`$ distribution functions. The obtained equations are applied to the wave-packet propagation in a ballistic 2DEG and to the calculation of the density-density correlation function $`\chi (\omega ,q)`$. We observe that for $`q>m\alpha `$ the structure factor $`\chi (\omega ,q)`$ is almost not affected by the spin-orbit interaction, but it reveals new features when $`qm\alpha `$ due to the direct transitions between different spin-orbit subbands.

Fruitful discussions with A. Andreev, K. Flensberg, L. Glazman, and C. Marcus are gratefully appreciated. This material is based on work supported by the NSF under grant PHY-01-17795 and by the Defence Advanced Research Programs Agency (DARPA) under Award No. MDA972-01-1-0024.


text:
# Heralded Two-Photon Entanglement from Probabilistic Quantum Logic Operations on Multiple Parametric Down-Conversion Sources

## Abstract

An ideal controlled-NOT gate followed by projective measurements can be used to identify specific Bell states of its two input qubits. When the input qubits are each members of independent Bell states, these projective measurements can be used to swap the post-selected entanglement onto the remaining two qubits. Here we apply this strategy to produce heralded two-photon polarization entanglement using Bell states that originate from independent parametric down-conversion sources, and a particular probabilistic controlled-NOT gate that is constructed from linear optical elements. The resulting implementation is closely related to an earlier proposal by Sliwa and Banaszek \[quant-ph/0207117\], and can be intuitively understood in terms of familiar quantum information protocols. The possibility of producing a “pseudo-demand” source of two-photon entanglement by storing and releasing these heralded pairs from independent cyclical quantum memory devices is also discussed.

The entanglement of two (or more) particles remains one of the most fascinating aspects of quantum mechanics . In addition to its relevance in a number of fundamental issues , quantum entanglement has recently been identified as a valuable resource for a variety of quantum information processing tasks . One promising approach to the practical realization of many of these tasks relies on qubits that are encoded in the polarization states of single photons. Consequently, a reliable source of heralded or even “push-button” two-photon polarization entanglement is of great interest at the present time.

One recent suggestion for such a source involves photon-pair production from controlled biexciton emission of a single quantum dot . Alternatively, one can consider a source based on the photon-pairs produced in parametric down-conversion . However, the random nature of this spontaneous emission process prohibits the direct use of a single two-photon down-conversion source for the production of on-demand pairs, and schemes based on conventional entanglement swapping  between two down-conversion sources can lead to false heralding signals due to double pair emission from one of the sources . In fact, Kok and Braunstein have shown that the production of one heralded polarization-entangled photon pair using only conventional down-conversion sources, linear optical elements, and projective measurements requires at least three initial pairs .

In this brief paper, we describe a method for producing heralded two-photon entanglement along theses lines. We consider the use of our previously proposed probabilistic controlled-NOT gate , which consumes one entangled photon pair as a resource, to essentially perform an entanglement swapping procedure on two additional entangled photon pairs. The end result is a unique detection signature of four photons that heralds the presence of one remaining polarization-entangled pair, with a potentially negligible probability of false heralding due to undesired multiple pair emission from the three independent initial pair sources. As will be shown, the resulting set-up is closely related to an earlier proposal by Sliwa and Banaszek that relies on triple-pair emission from a single down-conversion source . Other related work includes proposals for generating heralded two-photon polarization entanglement from multiple sources of single photons , and methods for generating heralded photonic “path-entanglement” from single-photon sources  or down-converted photon pairs .

Our approach is based on the well-known fact that a controlled-NOT gate, followed by certain measurements, can be used to identify specific Bell states of its two input qubits . As shown in Figure 1, we therefore consider the interesting case when the two input qubits are each members of independent Bell states. Specifically, we envision $`A`$ and $`B`$ as independent down-conversion sources designed  to produce two-photon Bell states of the form $`\varphi _A^+=\frac{1}{\sqrt{2}}[0_10_2+1_11_2]`$, $`\varphi _B^+=\frac{1}{\sqrt{2}}[0_30_4+1_31_4]`$. Assuming that exactly one pair has been emitted by each source, it can be seen that the total four-photon state after the operation of the controlled-NOT gate is of the form :

$$\psi =\frac{1}{2}[0000+0011+1101+1110]$$
(1)
where the modes are sequentially labelled, 1,2,3,4. By considering the heralding signal to be the joint detection of a single photon with logical value 0 in mode $`4`$, and a single photon with logical value $`\frac{1}{\sqrt{2}}[0+1]`$ (or, say, a Hadamard followed by detection of value 0) in mode $`2`$, it can be seen from equation 1 that the heralded output in modes $`1`$ and $`3`$ is the Bell state $`\varphi _{out}^+=\frac{1}{\sqrt{2}}[0_10_3+1_11_3]`$.

This measurement procedure can be understood as follows: the detection of a single photon by each (ideal) detector assures that exactly one pair has been emitted from each random down-conversion source, while the requirement that these detected photons have the logical values described above serves to post-select the $`\varphi ^+`$ Bell state of the controlled-NOT inputs. Because these inputs are members of $`\varphi ^+`$ Bell states, this post-selected entanglement is swapped  onto the remaining photons in modes $`1`$ and $`3`$.

From a practical point of view, the key feature of this protocol is that a probabilistic controlled-NOT gate  (as opposed to a deterministic one) may be sufficient for the production of the heralded entanglement. We therefore consider the implementation of the protocol using our proposed probabilistic controlled-NOT gate, as shown in Figure 2. The details of this probabilistic quantum logic gate are described in reference , and some limited aspects of its operation have been experimentally demonstrated .

To briefly review, our controlled-NOT gate is constructed from two polarizing beam splitters (PBS) which transmit horizontally polarized photons (qubit value 0), and reflect vertically polarized photons (qubit value 1) when aligned in the computational basis. As seen in Figure 2, the upper PBS is aligned in the computational basis, while the lower PBS is rotated by $`45^o`$ and therefore transmits polarization states corresponding to $`|+=\frac{1}{\sqrt{2}}[0+1]`$, while reflecting polarization state corresponding to $`|=\frac{1}{\sqrt{2}}[01]`$. These two PBS’s are linked by a source of two entangled ancilla photons (labelled $`C`$) which, once again, we envision as a down-conversion source designed  to produce a two-photon Bell state of the form $`\varphi _C^+=\frac{1}{\sqrt{2}}[0_50_6+1_51_6]`$. Given the entangled ancilla state, this controlled-NOT gate signals its successful operation on two input qubits by the joint detection of a single photon in the $`+/`$ basis above the upper PBS, and a single photon in the computational basis below the lower PBS . The probability of successful operation is $`\frac{1}{4}`$ provided that feed-forward control  is used to perform various single qubit corrections in the outputs.

Within the framework of the protocol outlined in Figure 1, the heralding signal for the desired polarization entanglement in modes $`1`$ and $`3`$ would therefore consist of the detection of a single photon by each of the four ideal detectors. In the first step, the successful operation of the controlled-NOT gate is signalled by the detection of single photons in modes $`5`$ and $`6`$. In the second step, the detection of a single photon in each of the controlled-NOT gate outputs is used to perform the post-selection and swapping of the entanglement as described earlier.

One advantage of this particular arrangement is that it is nearly immune to false heralding signals due to undesired multiple pair emission by the three initial pair sources. For example, double (or even triple) pair emission from either source $`A`$ or $`B`$ simply cannot distribute four single photons among the four detectors. Given ideal photon-number resolving detectors, the same is true for double pair emission from any one of the three sources combined with single pair emission from either of the two remaining sources.

Furthermore, double pair emission from source $`C`$ alone cannot result in the correct distribution of photons as a consequence of the two PBS’s being aligned $`45^{}`$ with respect to one another. This latter interference effect can be understood by the fact that amplitudes (expressed in the computational basis) describing photon pairs in mode $`5`$ must be orthogonally polarized to be correctly distributed by the upper PBS (ie. one into each output port), while those in the mode $`6`$ must be parallel to have some probability of being correctly distributed by the lower PBS. In contrast to these requirements, double pair emission from the $`\varphi ^+`$ Bell state source $`C`$ is of the form $`\frac{1}{2}[0000+0101+1010+1111]`$, where the first two qubits correspond to photons in mode $`5`$ while the second two qubits correspond to photons in mode $`6`$. In other words, amplitudes that would be successfully distributed by the upper PBS are correlated with amplitudes that can not be correctly distributed by the lower PBS, and vice-versa.

Even with ideal detectors, however, a false heralding signal could be obtained from simultaneous double pair emission from both source $`A`$ and source $`B`$. Although the probability of this quadruple pair process can be made significantly smaller than the desired triple pair emission process in the limit of weak down-conversion, such contamination of the heralded Bell state may limit the usefulness of this method for applications in which these errors can not be tolerated.

As we mentioned earlier, the overall set-up which results from implementing the general protocol outlined in Figure 1 with our particular probabilistic controlled-NOT gate  is closely related to an earlier proposal by Sliwa and Banaszek . For comparative purposes, a loose interpretation of that proposal in terms of our probabilistic controlled-NOT gate is shown in Figure 3. Roughly speaking, additional (non-polarizing) beam splitters are inserted in modes 5 and 6, and the required three initial photon pairs are generated by triple-pair emission from the single source $`C`$. In this configuration, the heralding signal (of four single-photon detections) is analogous to that described above, and the production rate of heralded-pairs (in modes $`7`$ and $`8`$) was found to be optimized when the reflectivity of the additional beam splitters is $`\frac{2}{3}`$ .

Finally, it is instructive to consider the conversion of a heralded source of two-photon polarization entanglement, produced by any means, into a source of entanglement on “pseudo demand”. A pseudo-demand source approximates a “push-button” source in the sense that the heralding signal is used to activate a storage mechanism for the heralded resource, which can then be released on-command at some chosen later time. Pseudo-demand sources of this kind are particularly relevant in the optical realization of quantum information protocols that can operate on a fixed clock cycle such as that provided by the repetition rate of a master pulsed laser. For these applications the storage mechanisms can be based on the use of simple repetitive optical loops whose round-trip propagation times are equal to the period of the clock cycle.

We previously demonstrated the conversion of a heralded single-photon source into a source of single-photons on pseudo demand in this manner . In that experiment, the detection of one member of a down-conversion photon pair heralded the presence of its twin photon, and activated a switching mechanism to reroute the twin into a storage loop. By analogy, the conversion of a source of heralded two-photon polarization entanglement into a pseudo-demand source of entanglement can be accomplished by storing each member of the heralded pair in its own storage loop, as shown in Figure 4. The storage loops and associated switching mechanisms, however, must not measure or alter the value of their stored qubits in any way, and in this sense are equivalent to a cyclical (or periodic) quantum memory device. A proof-of-principle demonstration of one particular implementation of a suitable cyclical quantum memory device for single-photon qubits is described in reference .

Whereas loss in either of the quantum memory devices will cause an obvious reduction in the quality of the stored and released Bell state, the effects of phase shifts are not as detrimental. In fact, because one member of the heralded pair is contained in each cyclical quantum memory device of Figure 4, overall phase shifts between the two devices essentially factor out of the two-photon state in such a way that phase-locking of the two memory loops is not required to maintain the coherence of the stored Bell state. However, any relative phase shifts between qubit values 0 and 1 (ie. birefringence for the case of polarization encoded qubits) in the individual memory devices will result in a final state of the form $`\frac{1}{\sqrt{2}}[0_10_3+e^{i\varphi }1_11_3]`$, where $`\varphi `$ represents the total of these phase shifts. One of the tangible benefits of the particular cyclical quantum memory device described in reference  is a robustness against these types of relative phase-shift errors due to repeated bit-flipping of the stored qubits.

In summary, we have described a protocol for producing heralded two-photon polarization entangled states by performing probabilistic quantum logic operations on multiple parametric down-conversion sources. We then described a method for converting a source of heralded entanglement into a source of entanglement on pseudo demand by controlled storage of the entangled photon pairs in two independent cyclical quantum memory devices . Such a pseudo-demand source of entanglement is expected to be valuable for the optical realization of a variety of quantum information processing tasks, particularly within the framework of linear optics quantum computation .

For example, the primary motivation for this work was to provide a pseudo-demand source of two-photon polarization entanglement to use as the entangled ancilla pair required for a complete (as opposed to post-selected) demonstration of our previously proposed non-deterministic controlled-NOT gate . It is interesting to see that the controlled-NOT gate itself, run with a source of random entanglement, can be used to produce the deterministic ancilla entanglement needed for a complete implementation of the gate. Such “bootstrapping” of resources and lower efficiency logic gates to produce higher efficiency gates is a fundamental feature of linear optics quantum computation .

An experimental implementation of the protocol described in this paper presents a formidable challenge in several regards. By far the most significant difficulty appears to be the implementation of the high efficiency photon-number resolving detectors that are required to reject certain false heralding signals during the post-selection procedure. In addition, the suppression of a false heralding signal due to double-pair emission from source $`C`$ is based on an experimentally delicate destructive four-photon interference effect. Because an error of that kind would correspond to double-pair emission while the effects of interest correspond to triple-pair emission, the visibility of this destructive interference would have to be very near unity. Although the required mode-matching would be facilitated to some degree by the fact that source $`C`$ is a single down-conversion crystal , so far the highest reported visibilities for any type of multi-photon interference effects have only been in the range of roughly $`0.85`$ to $`0.95`$ (see, for example, ).

TBP, MMD, MJF, BCJ, and JDF were supported by Army Research Office, National Security Agency, Advanced Research and Development Activity, Office of Naval Research, and Independent Research & Development funding. JPD, PK, and HL would like to acknowledge that this work was performed in part by the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration, with additional support by the Office of Naval Research, the Defense Advanced Research Projects Agency, the National Security Agency, the Advanced Research and Development Activity, and the National Reconnaissance Office.


text:
# Density–functional study of Cu atoms, monolayers, and coadsorbates on polar ZnO surfaces

## I Introduction

The ternary composite system Cu/ZnO/Al2O3 is an important commercial catalysts for the methanol synthesis, the water–gas shift reaction and the methanol–steam reforming process.cat1 ; cat2  Although pure ZnO is a good high temperature catalyst for hydrogenation and dehydrogenation reactions, only the combination with Cu makes it an effective and selective catalyst for low temperature applications on an industrial scale. In particular the role of copper, including its electronic and atomistic structure on ZnO surfaces, is not clearly understood. There exist many controversial models to explain how copper enhances the activity of such ZnO–based catalysts.rodriguez ; campbell0  Different “active species” are proposed in the literature, like the existence of Cu+ impurities at the ZnO surfaces, but the most commonly accepted picture is that Cu forms small metallic clusters which are dispersed on defect–rich and heterogeneous ZnO surfaces.

The initial growth of vapor–deposited copper on different ZnO surfaces at ultra–high vacuum (UHV) conditions has been studied experimentally by various techniques such as x–ray photoelectron spectroscopy (XPS),campbell0 ; campbell1 ; campbell2 ; campbell3  low–energy ion–scattering spectroscopy (LEIS),campbell0 ; campbell1 ; campbell2 ; campbell3  surface x–ray diffraction (SXRD),jedrecy  and scanning tunneling microscopy (STM).diebold1  For the two polar surface terminations of ZnO a growth mode between a Volmer–Weber and Stranski–Krastonov behavior was found. Initially, the deposited Cu atoms form two–dimensional islands (i.e. clusters with a thickness of one atom) until a critical coverage is reached, at which further Cu atoms begin to add on top of the islands to form three–dimensional clusters, whereas the metal–free regions between the clusters fill only very slowly.campbell1 ; campbell2  This behavior was attributed to purely kinetic limitations at low temperatures since a thickening of the islands was observed upon annealing the surfaces. The critical coverages at which three–dimensional growth of the clusters begins, depends slightly on temperature and was estimated to be about 1/2 monolayer for the O–terminated (000$`\overline{1}`$)–O surfacecampbell1  and 1/3 of a monolayer for the Zn–terminated (0001)–Zn polar face of ZnO.campbell2  In contrast, for clean, adsorbate–free nonpolar (10$`\overline{1}`$0) surfaces exclusively three–dimensional islands were observed at all coverages.diebold1

The main focus of theoretical investigations of the bonding of metal films and clusters on oxide surfaces has been on the MgO (001) surfaceroesch1 ; musolino1 ; musolino2 ; musolino3 ; roesch2 ; roesch3  and to a lesser extent on the surfaces of more complex oxides such as TiO2,tio2  Al2O3,finnis1 ; finnis2  and MgAl2O4.koestlmeier  Almost all studies were restricted to nonpolar surface terminations. Polar surfaces, on the other hand, like the two basal planes of ZnO, may show a completely different behavior compared to nonpolar surfaces. The reason is that polar surfaces are characterized by an excess of one atomic speciesfinnis3  which leads to an interesting interplay between their electronic and atomic structure: In order to be charge–neutral the polar surfaces either have to adopt an electronic structure with partially occupied surface bands (sometimes referred to as “metallization” of polar surfaces), or they have to reconstruct (preferentially by the loss of suitable surface atoms) in order to avoid partial band filling, or charged species have to be adsorbed.nosker ; noguera ; zno37  Depending on which mechanism is realized in a particular system, very different properties of the interaction with metal films and coadsorbates may occur.

Studies of the interaction of copper with ZnO surfaces are very scarce. In an early semiempirical quantum–chemical studyrodriguez  the electronic structure and atomic charges were calculated for various ZnO clusters with adsorbed and substituted Cu atom, but no information on binding energies and relaxed atomic structures was obtained. In addition, clusters without embedding were used which can not catch the peculiarities of the polar surfaces described above. This was taken into account in a recent quantum–chemical ab–initio study of the adsorption of single Cu atoms on different cluster modelsfink  to which we will compare our periodic density–functional results.

In the present study we have calculated the adhesion energy and geometry of single Cu atoms and copper monolayers at various coverages on the two polar surface terminations of ZnO. The surfaces are described by periodically repeated slabs so that the different stabilization mechanisms for the polar surfaces, in particular the existence of partially occupied surface bands, can be handled without any restrictions. Additionally, atomic relaxations are fully taken into account. The main focus of this study will be on how the charge neutralization process of the polar surfaces and the filling of the surface bands influences the adsorption properties of copper spanning the range from single atoms to full monolayers. In addition to the very recent quantum–chemical finite cluster study, Ref. fink, , these are the first ab–initio results for the adhesion energies and geometries of copper on ZnO surfaces.

## II Computational Method

Self–consistent total–energy calculations based on density–functional theory (DFT)hks  were carried out to determine the adsorption energies and geometries of Cu atoms and monolayers on the polar ZnO surfaces. The exchange–correlation energy and potential were treated within the generalized–gradient approximation (GGA) using the functional of Perdew, Burke and Ernzerhof (PBE).pbe  Normconserving pseudopotentialsvan-pp  were employed together with a mixed–basis set consisting of plane waves and non–overlapping localized orbitals for the O–$`2p`$, the Zn–$`3d`$ and the Cu–$`3d`$ electrons.mb  A plane–wave cut–off energy of 20 Ry was sufficient to get well converged results.

The polar surfaces are represented by periodically repeated slabs. No mirror symmetry parallel to the polar surfaces is present which implies that all slabs are inevitably Zn–terminated on one side and O–terminated on the other. For the study of Cu monolayers (1$`\times `$1) surface unit cells and very thick slabs consisting of 8 Zn–O double-layers were used to reduce the residual internal electric field.zno34 . All atomic configurations were fully relaxed by minimizing the atomic forces. The calculations for single Cu atoms were performed on (2$`\times `$2) surface unit cells (corresponding to a Cu surface coverage of 1/4 monolayer), and the slab thickness was reduced to 4 ZnO double-layers. In a (2$`\times `$2) surface unit cell arrangement the Cu atoms are more than 6.5 Å apart. The interaction energy between the Cu atoms is calculated to be less than 0.02 eV so that this setup is sufficient to describe very well isolated Cu atoms.

Three different high–symmetry adsorption sites are present at the polar ZnO surfaces (see Fig. 1): an ‘on-top’ position which would be the regular lattice site for the next atomic layer, a ‘hcp-hollow site’ above atoms in the second surface layer and a ‘fcc-hollow site’ with no atoms beneath. All three adsorption sites are considered in the present study. Copper atoms are only adsorbed on one of the two slab terminations. The electrostatic potential is calculated as for a genuinely isolated slabbengtsson ; bm  so no artificial interactions between the periodically repeated images of the slabs is present. Monkhorst–Pack k–point meshesmp  with a density of at least (6$`\times `$6$`\times `$4) points in the Brillouin–zone of the primitive ZnO unit cell are chosen. For more details on convergence parameters, the construction of appropriate supercells as well as the calculated bulk and clean surface structures of ZnO we refer to Ref. zno34, , where the same computational settings as in the present study were used.

## III Results and Discussion

### III.1 Validation of the Method: Finite Cluster Model

Before we start our investigation of copper on the polar ZnO surfaces, we address the question how well our DFT method and in particular the PBE functional are able to describe the interaction of Cu and ZnO surfaces. For Cu on MgO(001) significant differences between various wave–function based quantum–chemical methods and DFT calculations using different GGA functionals were observed.roesch3  As a general trend it was found that DFT slightly overestimates the adhesion energy of Cu on MgO clusters compared to best coupled–cluster calculations, with noticeable variations between the GGA functionals.

To test the reliability of our DFT results we have calculated the binding energy and equilibrium distance of Cu on a small, isolated Zn4O4 cluster which we can compare with recently published results of an accurate coupled–cluster–type calculation (multi–configuration coupled electron–pair approximation MC–CEPA).fink  A wurtzite–like cluster geometry was chosen in order to mimic the stacking sequence of polar slabs (see Fig. 2). The Cu atom is adsorbed to the top O and bottom Zn atom of the cluster which are three–fold coordinated like at a ZnO surface.

The Zn4O4 cluster itself is a closed–shell system, but with the additional Cu atom it becomes open-shell. Therefore, in the first calculation for a single Cu atom on-top of the three–fold O atom of the cluster we used the spin–polarized version of the PBE functional to obtain the total energy of the adsorbate system as well as for the energy of the single Cu reference atom. The binding energy and the equilibrium distance are given in Table 1. No significant differences in the adsorption properties are found if the spin–polarization is completely neglected and both, the calculations for the adsorbate system and the separate Zn4O4 cluster and the isolated Cu atom are done with the nonpolarized PBE functional. Also using a Cu dimer instead of a single Cu atom gave a very similar result. Since the inclusion of spin polarization doubles the computational cost, we only used a Cu dimer for the study of the interaction of Cu with the three–fold coordinated Zn atom of the finite cluster.

Compared to MC–CEPA results we find a satisfying agreement for the binding energies and equilibrium distances. Like in the case of Cu on MgO, the DFT calculation (here with the PBE functional) slightly overestimates the bond strength. However, the adsorption energy is ‘equally’ to large by 0.15 to 0.25 eV, so that the relative stability of the O and the Zn adsorption site is quite well described. Overall we find that Cu preferentially binds to oxygen atoms with an adhesion energy of almost 1 eV. The same behavior was also found for MgOroesch3 , but with a lower binding energy of 0.4 eV (coupled–cluster) and 0.6 eV (DFT-GGA).

### III.2 Cu Atoms on Ideal Polar Surfaces

For a long time it was believed that the polar ZnO surfaces exist in an unreconstructed, truncated–bulk–like state.noguera ; wander  If no reconstruction by loss of surface atoms and no adsorption of charged species occurs, the polar surfaces are inevitably stabilized by partially filled surface states with a 3/4 and a 1/4 filled band for the O and the Zn termination, respectively.nosker ; noguera ; zno37 ; wander

Assuming such a situation in our DFT calculations we find a very strong binding of Cu atoms on both polar ZnO surfaces (see Table 2). The Cu binding energies are much larger than what is found for the finite Zn4O4 cluster (see Table 1) and, more importantly, for the nonpolar ZnO(10$`\overline{1}`$0) surfacezno39 . In addition, it is also larger than what is typically known from other oxide surfaces.roesch3  The reason for this behavior is that at low copper coverages, the $`4s`$ electrons of the Cu atoms fill the partially occupied surface bands and thereby stabilize the polar surfaces. The same mechanism was found for the adsorption of atomic hydrogen on the polar O–terminated surface (see Ref. zno37, ): at low H coverages the binding energy of an H atom relative to a reference of isolated H2 molecules is more than 2 eV. However, as soon as a H coverage of 1/2 monolayer is reached, the surface bands are filled and the adsorption energy for additional H atoms decreases rapidly. As a consequence, surface structures with more than 1/2 monolayer of H are not stable in thermodynamic equilibrium with a H2 gas phase at any temperature and partial pressure.zno37  Therefore, we can expect that the strong binding of copper on the polar ZnO surfaces is only possible as long as the surface bands are not fully occupied. If the bands are filled, either by higher copper coverages or by the presence of coadsorbates like hydrogen, the Cu adsorption properties may change drastically. However, for copper there is one significant difference compared to hydrogen: at higher Cu coverages strong Cu–Cu bonds will form so that Cu will also adsorb at coverages exceeding 1/2 monolayer.

### III.3 Cu Atoms on Adsorbate Pre–Covered Polar Surfaces

Several very recent studies have created considerable doubt that the polar ZnO surfaces really exist as ideal, truncated–bulk–like surfaces. With scanning tunneling microscopy (STM) nanoscaled island and pit structures were found at the Zn–terminated facediebold2  which lead to an overall loss of roughly 1/4 of the surface Zn atoms. Such a stabilization mechanism was confirmed in DFT calculationsdiebold2 ; kresse  which furthermore predict high concentrations of adsorbed OH–groups in hydrogen–rich environments. For the O–terminated surface it was shown by He atom scattering,woell  in a study of the CO adsorption energetics,zno35 ; zno36 ; staemmler  and also by DFT calculations zno37  that they are H covered for a wide range of temperatures and H partial pressures. Both mechanisms, loss of Zn atoms and adsorption of hydrogen, lead to a filling of the surface bands. The filling need not to be perfect, however, in order to study the consequences of such a stabilization of the polar surfaces on the Cu (co)adsorption properties and also in order to support our prediction from the last subsection, we have considered two idealized surface structures of the polar ZnO surfaces for which the surface bands are fully occupied.

For the O– and the Zn–terminated surface we assume that 1/2 monolayer of H atoms and OH groups, respectively, is present. We place the H–atoms/OH–groups on every second atomic row at the surface. At the adsorbate–free rows in between these adatoms the Cu atoms are placed and the Cu adsorption energies are computed as in Sec. III.2. The results are given in Table 3. We find indeed that the Cu binding energy is significantly reduced compared to the result of the ideal surface termination with partially occupied (“metallized”) surface bands. The adsorption energy of roughly 1 eV of Cu atoms on the H covered O–terminated surface is only slightly larger than what was found for the finite closed-shell Zn4O4 cluster and for the nonpolar ZnO surface (having fully occupied bands, i.e. nonmetallic character).zno39  In addition, a similarly small value was also found in a quantum–chemical embedded cluster study of the polar O–terminated surfacefink  in which the embedding scheme was chosen in such a way that it leads to fully occupied surface bands, similar to our hydrogen covered surface.

### III.4 Cu Monolayers on Polar Surfaces

We turn now to the study of Cu monolayers on the polar ZnO surfaces. ZnO and bulk Cu have very different lattice constants. The equilibrium distance of the surface atoms at the polar surfaces is computed to be 3.28 Å (experimentally: 3.25 Å), and for the Cu bulk lattice constant we found 3.65 Å which corresponds to a nearest–neighbor separation between the Cu atoms of 2.58 Å (experimentally: 3.615 Å and 2.56 Å, respectively). Epitaxial growth with a (111) orientation of copper would therefore lead to a lattice mismatch of 21 %. Thus, we can not expext that expitaxial growth occurs. Experimentally it was found that Cu on the polar surfaces is basically unstrained, even at very low coverages.jedrecy  Exclusively (111) Cu planes were observed with a fixed orientation given by the ZnO substrate.

The realistic treatment of such unstrained copper layers on the ZnO surfaces would be computationally very demanding since we would have to use prohibitively large surface unit cells for our slab. However, in order to get a first idea on how the Cu adsorption strength and the preferred adsorption sites change if a high concentration of copper is present on the polar surfaces we have nevertheless studied ideal expitaxial Cu (111) monolayers on the ZnO surfaces. For such hypothetical films we have calculated the equilibrium distance and the work of separation which is defined as the energy difference between the adsorbed monolayer and separate ZnO slabs and metal films. The results are given in Table 4. For the O–terminated surface we find that Cu preferentially binds on-top of the oxygen atoms. The equilibrium distance and the binding energy of 1.23 eV per Cu atom are very similar to the results obtained for single atoms on the finite Zn4O4 cluster and the surface with filled surface states. Also for the Zn–terminated surface we find a good binding of the copper monolayers with a work of separation of almost 0.9 eV per Cu atom. This is in agreement with experiment where a good adhesion of Cu on both polar surface terminations was found.campbell1 ; campbell2  In contrast to the O–terminated surface, the hollow positions for the Cu atoms are more stable than the on-top site. The two hollow positions are very similar in energy and may both be populated for non–epitaxial films. This may explain the experimental observation that the orientation of deposited Cu relative to the substrate is less strict for the Zn–terminated face than for the O–termination.campbell2 ; jedrecy  Importantly, we also find that the structural relaxation of the ZnO substrate can not be neglected. The relaxation energy contributes 20 % to 30 % to the binding energy of these copper films.

### III.5 Cu Films on Polar Surfaces

Finally, we briefly study the effect higher copper coverages than full monolayer on the bonding at the Cu/ZnO interface. This will give some insight into the nature of the chemical bond at the interface. The higher coverages were simulated by adding a second and a third monolayer of Cu on-top of the polar surface. For the two and three layer films different stacking sequences are possible. As shown in Table 5 and Table 6, the energy differences between different stacking sequences turn out to be very small so that the specific stacking does not play a role in the following arguments.

Overall we find that the same adsorption sites as for the monolayer coverage are the most stable ones for the thicker Cu films: the on-top position for the O–terminated surface and the hollow-sites for the Zn–terminated face. For the O–terminated surface the work of separation increases with the thickness of the copper film, whereas for the Zn–termination a small decrease in the binding energy is observed. This behavior as well as the preferred adsorption sites can be explained as follow: For the O–terminated surface the first copper layer adopts the role of the missing next Zn layer. Therefore the on-top position, which is the next regular lattice site of the wurtzite structure, is preferred. Charge is donated to the surface O ions to fill the partially occupied surface bands, which can more easily occur for thicker copper films, and a partially ionic/covalent bond like in bulk ZnO is formed.

On the other hand, for the Zn–terminated surface, the topmost Zn layer behaves more like the first layer of the adsorbed metal film. Therefore the hollow positions, which correspond to an extension of the metal film, are the low-energy adsorption sites. The $`4s`$ electrons of the partially filled Zn surface band contribute to the metallic state of the metal film so that the Zn–Cu bond is more metallic in character. However, the “integration” of the topmost Zn layer becomes weaker for thicker copper films.

## IV Summary and Conclusions

Using DFT calculations and periodically repeated slabs to represent the two polar ZnO surfaces, i.e. the (000$`\overline{1}`$)–O and (0001)–Zn faces, we find a pronounced difference in the adsorption energy of Cu atoms depending on the coverage. Crucial to the proposed understanding of this phenomenon is the fact that the ideal, unreconstructed O– and Zn–terminated surfaces feature partially filled surface states with a 3/4 and a 1/4 filled band, respectively (often referred to as “metallization” of the polar surfaces). Up to 1/2 monolayer coverage the adhesion of copper on such polar surfaces is found to be very strong since the $`4s`$ electrons of the Cu atoms can fill these partially occupied bands and thereby stabilize the polar surface itself. However, as soon as these bands are filled due to copper (or due to coadsorbates!) already present on the surface, the binding of additional Cu to the surface becomes weaker and less favorable in comparison to the cohesion energy of Cu in small Cu clusters and thin Cu films.

Considering the experimentally known coverage–dependent change of the growth mode of copper on the polar ZnO surfaces (i.e. formation of two–dimensional islands until critical coverages of about 50 % and 30 % for the O– and the Zn–terminated surface, respectively, is reached at which three–dimensional growth of Cu clusters sets in) it may very well be that this sudden change is not caused by kinetic limitations at low temperatures alone. Instead, it might also be driven thermodynamically by the significant change in the adsorption energetics of copper on such polar surfaces due to successive band filling, which in turn is directly coupled to increasing the surface coverage. The lower critical coverage found experimentally for the Zn–terminated surface may then be explained quite naturally by the lower adsorption energy of copper found for this particular surface. Additionally, the presence of coadsorbates during or prior to copper deposition may strongly influence both the critical coverage and the morphology of the grown Cu particles.

Thus, it becomes evident that the polar ZnO surfaces allow for a wealth of quite different adsorption scenarios depending on parameters such as coverage, defects, coadsorbates etc. Since these parameters typically change in complex ways during catalytic cycles one is tempted to speculate that the polar character might be at the very heart of the catalytic activity of ZnO.

## V Acknowledgments

We are grateful to Karin Fink for making the results of the cluster calculations available to us prior to publication, and we wish to thank her and Volker Staemmler for fruitful discussions. The work was supported by SFB 558 and FCI.