text:
# First-principles study of the polar O–terminated ZnO surface in thermodynamic equilibrium with oxygen and hydrogen

## I Introduction

To understand the structure, composition and stability of polar surfaces on a solid theoretical basis is one of the challenging problems in surface science.noguera  The most interesting polar surfaces are so called “Tasker-type(3)” surfacestasker  which are formed by alternating layers of oppositely charged ions. Assuming a purely ionic modelcomment1  in which all ions are in their formal bulk oxidation state, such a stacking sequence creates a dipole moment perpendicular to the surfaces which diverges with slab thickness, and with simple electrostatic arguments it can be shown that the surface energy will diverge with sample size.tasker  To quench the dipole moment and to make the polar surfaces stable, a redistribution of charges in the surface layers has to take place.nosker  For various polar surfaces different mechanisms to accomplish the charge compensation have been observed,fks  however, in many cases the underlying stabilization mechanism is very controversially discussed in the literature.

One of the most widely investigated examples of Tasker-type(3) polar surfaces are the two basal planes of ZnO: the O–terminated (000$`\overline{1}`$)–O and the Zn–terminated (0001)–Zn surface. The two surfaces are the terminating planes of a stacking sequence of hexagonal Zn and O layers along the crystallographic $`c`$-axis with alternating distances of $`R_1`$=0.61 Å and $`R_2`$=1.99 Å. In this case, as can be easily shownnoguera ; zno34 , the polar ZnO surfaces are only stable if the O–terminated face is less negative and the Zn–terminated surface layer less positively charged compared to the formal bulk oxidation state by a factor of $`R_1/(R_1+R_2)`$$``$1/4. In principle, three different scenarios are conceivable to accomplish this charge redistribution: The ionic charge of the surface ions may be reduced from $`\pm `$2 to $`\pm `$3/2, which may be regarded as an “electron transfer” from the O– to the Zn–terminated surface (Ia). As a result, partially occupied surface bands will appear with a 3/4 filled O–$`2p`$ band at the (000$`\overline{1}`$)–O and a 1/4 filled Zn–$`4s`$ band at the (0001)–Zn surface. This is often referred to as “intrinsic surface state compensation”nosker  or as “metalization of the polar surfaces”.wander  However, whether a true metallic state is present will depend on the dispersion of the partially occupied bands. Additionally, in a second step, the surface may reconstruct and undergo a distortion in which, for example for the O–terminated surface, four surface atoms combine in such a way that an unoccupied $`2p`$–band splits from the other eleven occupied $`2p`$–bands and the surface becomes insulating again (Ib). Secondly, the charge reduction of the surface layers may take place by removing 1/4 of the surface ions (II). The so created vacancies may be ordered and form a reconstruction or may be randomly distributed. Finally, charged species may be adsorbed to reduce the formal oxidation state of the surface ions (III). For example, water may dissociate and protons (H+) and hydroxyl groups (OH-) adsorb on every second O and Zn surface ion, respectively.comment2  All three scenarios represent idealizations of the charge compensation process. In general, any combinations of the three mechanisms are conceivable, like the simultaneous formation of vacancies and partially filled bands, as long as the charge compensation rule is obeyed. The surface structure finally realized for a specific polar surface will then be the one with the lowest surface energy.

For ZnO it was believed for a long time that both polar surfaces exist in an unreconstructed, truncated-bulk-like state. After standard preparation procedures both surfaces show regular (1$`\times `$1) pattern in low-energy electron diffraction (LEED)duke  and other diffraction experiments.wander ; jedrecy ; overbury  Some evidence for missing Zn ions on the (0001)–Zn surface was found in gracing incidence X-ray diffraction (GIXD)jedrecy , however, for the (000$`\overline{1}`$)–O surface no evidence for substantial amounts of surface oxygen vacancies was detected in GIXDjedrecy  and low-energy alkali-ion scattering (LEIS).overbury

For ideal, truncated-bulk-like surface terminations only mechanism (Ia) can explain the stability of the polar ZnO surfaces. Consequently, in most theoretical first-principles studies of the polar ZnO surfaces ideal surface terminations together with partially filled surface bands were assumed.wander ; carlsson ; zno34  Studies exploring the other two stabilizations mechanisms are very scarce. In a pioneering ab-initio study Wander and Harrisonwander\_h  investigated whether the polar surfaces may be stabilized by the dissociation of water and the adsorption of H+ and OH- groups according to mechanism (III). They found this energetically unfavorable compared to situation (Ia). However, instead of 1/2 monolayers, which would be the ideal configurations for band filling, only full monolayers of H+ and OH- were considered, thereby overcompensating the needed charge transfer between the polar surfaces. In addition, only one specific adsorption site for the H+ and OH- groups was probed.

Meanwhile two recent experiments have created considerable doubt that the polar ZnO surfaces really exist in a clean, unreconstructed state. With scanning tunneling microscopy (STM) it was showndiebold1 ; diebold2  that the Zn–terminated surface is characterized by the presence of nanoscaled, triangular islands with a height of one ZnO double-layer. The shape of the islands is size-dependent and typical pit structures appear for larger islands. Since the step edges are O–terminated, the high step concentration leads to a significant decrease of Zn ions in the surface. A rough analysis of the island and pit size distribution yielded that approximately 1/4 of the Zn ions is missing, in agreement with mechanism (II). With detailed density-functional theory (DFT) calculationsdiebold1 ; kresse  it was confirmed that a crystal termination with triangular shaped islands and pits is indeed lower in energy than the perfect bulk-truncated surface for a wide range of oxygen and hydrogen chemical potentials. Under H rich conditions, structures with up to 1/2 monolayer of hydroxyl groups were even more stable, indicating that the actual surface morphology will sensitively depend on the chemical environment.

On the other hand, for the O-terminated polar surface no such island and pit structure was found with STM.diebold2  However, with He scattering (HAS) it was discoveredwoell  that at ultrahigh vacuum (UHV) conditions O–terminated surfaces with (1$`\times `$1) LEED and HAS diffraction pattern are always hydrogen covered, whereas after a careful removal of the hydrogen a (1$`\times `$3) structure is found. The (1$`\times `$3) spots are best visible in HAS, but under certain conditions can also be observed in LEED.woell  The H–free surface is very reactive and dissociates molecular hydrogen and water and therefore exits only for a limited time even at UHV conditions. In a subsequent studyzno35 ; zno36 ; staemmler ; fink  CO was used as a probe molecule to distinguish between clean and hydrogen saturated surfaces. By comparing calculated CO adsorption energies for different surface structures with experimental results it was confirmed that the polar O–terminated surface is usually hydrogen covered whereas a clean (1$`\times `$1) (000$`\overline{1}`$)–O surface is very unlikely to exist.

In previous studies the H coverage of the polar O–terminated surface was not observed most likely because LEED and X–rays are not sensitive to hydrogen. However, it is not clear how the structures found in the HAS study, Ref. woell, , lead to a stabilization of the polar O–terminated surface. A full hydrogen monolayer would overcompensate the charge transfer so that again partially occupied bands have to be present, and the nature of the H–free (1$`\times `$3) structure is still unknown. There is some evidence from X–ray photoelectron spectroscopy (XPS)woell  that oxygen vacancies are involved, but 1/3 of the oxygens missing is also not the expected vacancy concentration to stabilize the surface.

Motivated by these experimental findings we explore in the present paper the competition between the three stabilization mechanisms in a very general way. The main focus will be on the O–terminated (000$`\overline{1}`$)–O surface, and we take an approach very similar in spirit to the investigation of the Zn–terminated surface in Ref. diebold1,  and kresse, . For a series of surface models we determine the total energies and the fully relaxed atomic structures using a first-principles DFT approach. Surface structures with various oxygen vacancy concentrations and different amounts of adsorbed hydrogen are considered, including structures corresponding to the three ideal stabilization scenarios and structures compatible with the HAS observations.

Static total-energy DFT calculations only give results for zero temperature, zero pressure and for surfaces in contact with vacuum. However, the actual lowest-energy structure of the (000$`\overline{1}`$)–O surface will depend on the environment and can change with temperature $`T`$, pressure $`p`$ and exposure to O2 and H2 gas phases. Therefore, to determine the equilibrium structure and composition of the surface at finite temperature and oxygen and hydrogen partial pressures, we combine our DFT results with a thermodynamic description of the surfaces. To take deviations in surface composition and the presence of gas phases into account, we introduce appropriate chemical potentialschadi  and calculate an approximation of the Gibbs free surface energy.finnis1  Depending on the chemical potentials we then determine the surface structure with the lowest free energy which allows us to construct a phase diagram for the surface. If we assume that the surface is in thermodynamic equilibrium with the gas phases, we can relate the chemical potentials to a given temperature $`T`$ and pressure $`p`$. In this way we are able to extend our zero temperature and zero pressure DFT results to experimentally relevant environments, thereby bridging the gap between UHV-like conditions and temperatures and gas phase pressures that are typically applied, for example, in catalytic processes like the methanol synthesis.catalysis

## II Computational Approach

### II.1 Thermodynamics

In this section we will give a brief description of the thermodynamic formalism which we have used to determine the most stable structures of the polar O–terminated ZnO surface. The formalism has been successfully applied in several previous surface studiespadilla ; wang ; batyrev ; pojani ; reuter ; diebold1  and is described in more detail in Refs. finnis1,  and reuter, .

The general expression for the free energy of a surface in equilibrium with particle reservoirs at the temperature $`T`$ and pressure $`p`$ is given bycahn

$$\gamma (T,p)=\frac{1}{A}\left(G(T,p,\{N_i\})\underset{i}{}N_i\mu _i(T,p)\right),$$
(1)
where $`G(T,p,\{N_i\})`$ is the Gibbs free energy of the solid with the surface of interest, $`A`$ is the surface area, and $`\mu _i`$, $`N_i`$ are the chemical potentials and particle numbers of the various species. In contrast to the usual convention in macroscopic thermodynamics we define here the chemical potentials per atom rather than per mole. For the study of the polar O–terminated ZnO surface in contact with an oxygen and a hydrogen gas phase we have to consider the three chemical species $`i`$ = Zn, O and H.

For simplicity we have assumed two independent reservoirs for O2 and H2 with a common pressure $`p`$. Experimentally, it is more likely that a mixture of O2 and H2 is present. In this case, the pressure $`p`$ in Eq. (1) has to be replaced by appropriate partial pressures $`p_{\mathrm{O}_2}`$ and $`p_{\mathrm{H}_2}`$. However, in the present study we will keep the restriction of separate reservoirs in the sense that we do not allow O2 and H2 to react to H2O, which would be the case in full thermodynamic equilibrium. This is justified by arguing that the reaction barrier for the formation of H2O is high enough that the reaction plays no role on the time scales of interest.

In thermodynamic equilibrium the chemical potentials would be determined uniquely by the temperature $`T`$, the pressure $`p`$ and the total particle numbers of the solid and the gas phases. The surface structure, here represented by the particle numbers $`N_i`$, would then be determined by a unconstrained minimization of the surface free energy, Eq. (1). However, this is not very practical to do. Therefore we will take a different approach: We calculate the surface free energy of a series of model surfaces with different structures and compositions as a function of the chemical potentials. For given chemical potentials we predict which surface structure is the most stable one by searching for the surface model with the lowest surface free energy. In a second step, the chemical potentials are then related to actual temperature and pressure conditions by assuming that the surface is in thermodynamic equilibrium with the gas phases.

In our calculations all surfaces are represented by periodically repeated slabs so that the Gibbs free energy $`G(T,p,\{N_i\})`$ refers to the content of one supercell. Since ZnO is not centrosymmetric, slabs representing the polar ZnO surfaces are inevitably O–terminated on one side and Zn–terminated on the other side. It is therefore not possible to assign unique surface energies to the two polar surface terminations. Only the sum of the surface energies and thereby the cleavage energy are well defined quantities. However, in the present study we are only interested in the relative stability of O–terminated surfaces with different structures and compositions. The Zn-face of the slabs is unchanged in all calculations. Therefore we relate all energies relative to a reference state which we have taken to be the ideal, truncated-bulk termination. We define the change of the cleavage energy $`\mathrm{\Delta }\gamma `$ as the difference of Eq. (1) for a slab with a chosen surface structure and a slab with ideal surface terminations:

$`\mathrm{\Delta }\gamma (T,p)`$ $`=`$ $`{\displaystyle \frac{1}{A}}(G_{\mathrm{slab}}^{\mathrm{surf}}(T,p,N_\mathrm{V},N_\mathrm{H})G_{\mathrm{slab}}^{\mathrm{ref}}(T,p)`$ (2)
$`+N_\mathrm{V}\mu _\mathrm{O}(T,p)N_\mathrm{H}\mu _\mathrm{H}(T,p)).`$
Here, $`G_{\mathrm{slab}}^{\mathrm{surf}}`$ and $`G_{\mathrm{slab}}^{\mathrm{ref}}`$ are the Gibbs free energies of the supercells with the model surface and the reference configuration, respectively, and $`A`$ is now the area of the surface unit cell. Since we only consider structures of the O–terminated surface with O–vacancies and adsorbed H atoms, only the number of O–vacancies $`N_\mathrm{V}`$ (= difference of the number of O–atoms in the slab with the model surface and the reference state) and the number of adsorbed H atoms $`N_\mathrm{H}`$ appears in Eq. (2), i.e. the chemical potential of Zn drops out. The difference $`\mathrm{\Delta }\gamma `$ is negative if a model surface is more stable than the ideal, truncated-bulk-like surface termination and positive otherwise.

In principle we now have to calculate the Gibbs free energy of all slabs representing our surface models, including the contributions coming from changes in volume and in entropy. The entropy term may be calculated, for example, by evaluating the vibrational spectra in a quasiharmonic approximationfrank , but in practice this is computationally very demanding. As is apparent from Eq. (2), only the difference of the Gibbs free energy of two surface structures enters the expression for $`\mathrm{\Delta }\gamma `$. In Ref. reuter,  it was shown that the vibrational contributions to the entropy usually cancel to a large extent and that the influence of volume changes are even smaller. Therefore we will neglect all entropy and volume effects. The Gibbs free energies then reduce to the internal energies of the slabs and we can replace $`G_{\mathrm{slab}}^{\mathrm{surf}}`$ and $`G_{\mathrm{slab}}^{\mathrm{ref}}`$ in Eq. (2) by the energies as directly obtained from total–energy (e.g. DFT) calculations.

Finally we have to determine meaningful ranges in which we can vary the chemical potentials. First, there are upper bounds for all three chemical potentials $`\mu _\mathrm{O}`$, $`\mu _\mathrm{H}`$ and $`\mu _{\mathrm{Zn}}`$, beyond which molecular oxygen and molecular hydrogen would condensate and metallic Zn would crystallize at the surface. These bounds are given by the total energy of the isolated molecules $`E_{\mathrm{O}_2}`$, $`E_{\mathrm{H}_2}`$ and of bulk Zn $`E_{\mathrm{Zn}}^{\mathrm{bulk}}`$ (neglecting volume and entropy effects):

$$\mu _\mathrm{O}\frac{1}{2}E_{\mathrm{O}_2},\mu _\mathrm{H}\frac{1}{2}E_{\mathrm{H}_2},\mu _{\mathrm{Zn}}E_{\mathrm{Zn}}^{\mathrm{bulk}}.$$
(3)
In the following we will use these upper bounds as zero point of energy and relate the chemical potentials relative to the total energies of the isolated molecules:

$$\mathrm{\Delta }\mu _\mathrm{O}=\mu _\mathrm{O}\frac{1}{2}E_{\mathrm{O}_2},\mathrm{\Delta }\mu _\mathrm{H}=\mu _\mathrm{H}\frac{1}{2}E_{\mathrm{H}_2}.$$
(4)
Furthermore we impose that the surface is always in equilibrium with the ZnO bulk phase. Then the sum of $`\mu _\mathrm{O}`$ and $`\mu _{\mathrm{Zn}}`$ has to be equal to the total energy $`E_{\mathrm{ZnO}}`$ of bulk ZnO. Thus only one of the two chemical potentials $`\mu _\mathrm{O}`$ and $`\mu _{\mathrm{Zn}}`$ is independent, and together with Eq. (3) we introduce lower bounds for the chemical potentials. Using the energy of formation $`E_\mathrm{f}`$ of ZnO:

$$E_\mathrm{f}=E_{\mathrm{ZnO}}E_{\mathrm{Zn}}\frac{1}{2}E_{\mathrm{O}_2}$$
(5)
the allowed range for the chemical potential $`\mathrm{\Delta }\mu _\mathrm{O}`$ is given by:

$$E_\mathrm{f}\mathrm{\Delta }\mu _\mathrm{O}\mathrm{\hspace{0.33em}0}.$$
(6)
If we assume that the surfaces are in thermodynamic equilibrium with the gas phases we can relate the chemical potentials $`\mathrm{\Delta }\mu _\mathrm{O}`$ and $`\mathrm{\Delta }\mu _\mathrm{H}`$ to a given temperature $`T`$ and partial pressures $`p_{\mathrm{O}_2}`$ and $`p_{\mathrm{H}_2}`$. For ideal gases we can use the well-known thermodynamic expressionsreuter

$$\mathrm{\Delta }\mu _\mathrm{O}(T,p_{\mathrm{O}_2})=\frac{1}{2}\left(\stackrel{~}{\mu }_{\mathrm{O}_2}(T,p^0)+k_\mathrm{B}T\mathrm{ln}(p_{\mathrm{O}_2}/p^0)\right)$$
(7)
and

$$\mathrm{\Delta }\mu _\mathrm{H}(T,p_{\mathrm{H}_2})=\frac{1}{2}\left(\stackrel{~}{\mu }_{\mathrm{H}_2}(T,p^0)+k_\mathrm{B}T\mathrm{ln}(p_{\mathrm{H}_2}/p^0)\right)$$
(8)
in which $`p^0`$ is the pressure of a reference state and the temperature dependence of the chemical potentials $`\stackrel{~}{\mu }_{\mathrm{O}_2}(T,p^0)`$ and $`\stackrel{~}{\mu }_{\mathrm{H}_2}(T,p^0)`$ is tabulated in thermochemical reference tables.janaf  However, it should be noted, as pointed out by Finnisfinnis2 , that the equilibrium with the gas phase need not to be perfect. It is sufficient if the surface is in equilibrium with the bulk phase. In this case, the chemical potential is related to defect concentrations of the bulk. For example, if oxygen vacancies are the dominant defects we havemayer1 ; hagen ; mayer2

$$\mathrm{\Delta }\mu _\mathrm{O}(T,p)=\stackrel{~}{\mu }_0k_\mathrm{B}T\mathrm{ln}(c_\mathrm{V}/c_0)$$
(9)
with the vacancy concentration $`c_\mathrm{V}`$ and the oxygen occupancy of the O–lattice site of $`c_0`$.

### II.2 Ab-Initio Calculations

The total energies of slabs representing different model surfaces as well as the bulk and molecular reference energies were calculated within the framework of density-functional theory (DFT).hks  Exchange and correlation effects were included in the generalized–gradient approximation (GGA) using the functional of Perdew, Becke 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`$ and the Zn–$`3d`$ electrons.mb  A plane-wave cut-off energy of 20 Ry was sufficient to get well converged results. Monkhorst-Pack k-point meshesmp  with a density of at least (6$`\times `$6$`\times `$6) points in the primitive ZnO unit cell were chosen. A dipole correctionbengtsson ; bm  to the electrostatic potential was included in the calculations to eliminate all artificial interactions between the periodically repeated supercells due to the dipole moment of the slabs. For more details on convergence parameters, the construction of appropriate supercell 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.

All surfaces were modeled by periodically repeated slabs. Very thick slabs consisting of 8 Zn-O double-layers were used to reduce the residual internal electric field.zno34 . To represent different surface structures (1$`\times `$2), (1$`\times `$3) and (2$`\times `$2) surface unit cells with different combinations of O vacancies and H adatoms were considered. All atomic configurations were fully relaxed by minimizing the atomic forces.

In Table 1 we compare the computed binding energies of different bulk and molecular reference structures with experimental heat of formations. While the calculated binding energies for isolated $`H_2`$ molecules and bulk Zn agree quite well with experiment, there is a noticeable error of 1.2 eV in the binding energy of the free O2 molecule. This is a well known deficiency of DFT.batyrev ; reuter  The overestimation of the O2 binding energy is also reflected in a formation energy for ZnO which is to low by 0.6 eV. Such deviations would seriously influence our analysis of the surface energies, Eq. (2), and would alter the allowed range for the oxygen chemical potential, Eq. (6). Therefore, to circumvent errors introduced by a poor description of the O2 molecule, we have applied the following procedure: we take the experimental value for the formation energy $`E_\mathrm{f}`$ of ZnO from Table 1 and we use Eq. (5) to replace the total energy of O2 by $`E_\mathrm{f}`$ and the total energies of Zn and ZnO. $`E_{\mathrm{Zn}}`$ and $`E_{\mathrm{ZnO}}`$ are both bulk quantities which are typically more accurately described in DFT than molecular energies.

## III Results and Discussion

### III.1 Surface Distortions

First we explore the possibility whether the ideal, truncated-bulk-like (000$`\overline{1}`$)–O surface may lower its energy by breaking the symmetry of the surface layers, thereby adopting a distorted surface structure according to mechanism (Ib). A tendency for symmetry breaking reconstructions is often indicated by a strong nesting of the Fermi surface. In Fig. 1 we have plotted the two-dimensional Fermi surface formed by the partially occupied O–$`2p`$ bands. The plot represents a cut through the Brillouin zone of our supercell including only k-vectors with a zero component perpendicular to the surface. Figure 1 reveals that actually two surface bands cross the Fermi level. This is well known and in full agreement with band structure plots presented in Refs. wander,  and carlsson, . The corresponding wave functions are strongly localized at the oxygen atoms of the terminating surface layers and are mainly formed by the two O–$`2p`$ orbitals parallel to the surface. By integrating the occupied and unoccupied areas of the Brillouin-zone we find that indeed roughly 1/2 electron per surface atom is missing to fill the two surface bands. However, both Fermi contours are almost spherical and only a weak nesting is present.

As a second test, we did several calculations in which we randomly displaced the surface atoms in the top atomic layer of our slabs and started an atomic relaxation. Different slabs with (1$`\times `$2), (1$`\times `$3) and (2$`\times `$2) surface unit cells were used but in all cases the surfaces relaxed back toward a symmetric structure with a 3-fold symmetry. The surface energy was always higher than in the fully symmetric state, so that no hint for a tendency toward symmetry breaking reconstructions was found.

### III.2 Oxygen Vacancies

In the next step we investigate whether the (000$`\overline{1}`$)–O surface is stabilized by creating oxygen vacancies. From slabs with (1$`\times `$1), (1$`\times `$2), (1$`\times `$3) and (2$`\times `$2) surface unit cells we have removed one O–atom, thereby creating vacancy concentrations $`c_\mathrm{V}`$ of 1, 1/2, 1/3 and 1/4. In Fig. 2 we have plotted the change of the surface energy $`\mathrm{\Delta }\gamma `$ of the four defect structures relative to the defect-free surface as a function of the oxygen chemical potential $`\mathrm{\Delta }\mu _\mathrm{O}`$. As to be expected from mechanism (II) the defective surface with 1/4 of the O atoms missing is the most stable surface structure over a wide range of chemical potentials. Translating the chemical potential into temperature and pressure conditions (assuming that the surface is in equilibrium with an O2 gas phase, see upper $`x`$–axis in Fig. 2) we see that this will be the most stable surface at typical UHV-conditions of surface science experiments. However, at oxygen rich conditions (high pressure and low temperature), exceeding a chemical potential of $``$1.58 eV, the ideal, defect-free surface becomes the most stable structure. The other surfaces with 1, 1/2 and 1/3 vacancy concentrations are higher in energy for all chemical potentials and will therefore not be present in thermodynamic equilibrium. In particular, it is very unlikely that the (1$`\times `$3) structure observed in the HAS experimentwoell  corresponds to a simple missing-row structure with every third O atom removed from the surface.

At this point we should emphasize that plots like Fig. 2 strictly only allow to rule out surface structures which are higher in energy than other surface models. Since we can only perform calculations for a limited set of surface models, it is always possible that a not yet considered structure with a lower energy exists. For example, since we use periodically repeated surface unit cells of a specific size in our DFT approach, all our defect structures are perfectly ordered. It is however very well possible, that disordered or even incommensurate structures might lead to a lower energy. Additionally, there are hints that island and pit structures like the ones observed for the Zn–terminated surfacediebold1  may also for the O–terminated surface be lower in energy than the ideal surface and the surface with 1/4 vacancy concentration considered in our study.kresse2

In Table 2 we have listed O vacancy formation energies $`E_\mathrm{V}`$ which we have defined in the present context as the energy difference between the ideal surface and a vacancy surface structure of concentration $`c_\mathrm{V}`$ plus $`1/2E_{\mathrm{O}_2}`$, i.e. $`E_\mathrm{V}`$ is defined with respect to the total energy $`E_{\mathrm{O}_2}`$ of free oxygen molecules, and not, as usually done, with respect to bulk ZnO. $`E_\mathrm{V}`$ depends strongly on the vacancy concentration, indicating a strong interaction between the vacancies. Up to $`c_\mathrm{V}`$=1/4 the energy cost to remove O–atoms is modest. This is not surprising since up to a vacancy concentration of 1/4 the removal of O-atoms supports the charge compensation of the O–terminated surface and will result in a better filling of the partially occupied O–$`2p`$ band. For higher vacancy concentrations, however, we start to overcompensate the charge transfer which stabilizes the surface. The O–$`2p`$ band is full now, and we have to start to fill Zn–$`4s`$–states in the conduction band. Therefore the energy cost to remove more O–atoms increases rapidly.

### III.3 Hydrogen Adsorption

We turn now to a situation in which the (000$`\overline{1}`$)–O surface is in equilibrium with a H2 gas phase. H2 molecules can dissociate, and hydrogen atoms may adsorb at the surface thereby forming OH-groups with the surface oxygen ions. Before we start to calculate the surface energy for different surface models with H adatoms, we consider the possibility that the preferred adsorption site of these OH-groups is no longer the regular lattice site of the O ions. Three different high-symmetry adsorption sites are conceivable above the underlying Zn layer: an ‘on-top’ position, a ‘hcp-hollow site’, which is the regular lattice position for the O ions in the wurtzite structure, and a ‘fcc-hollow site’.

First we consider the clean surface without adsorbed hydrogen. Then we see from Table 3 that indeed the ‘hcp-hollow’ sites are the most stable positions for the O surface layer. However, moving the whole layer to ‘fcc-hollow’ sites costs only a small amount of energy. Turning to the hydroxylated surface with a full monolayer coverage of hydrogen we find that the OH-groups now prefer the ‘fcc-hollow site’. So by gradually adding hydrogen, the regular lattice site of the surface O ions becomes unstable relative to the ‘fcc-hollow site’. However, the energy difference between ‘fcc-’ and ‘hcp-hollow sites’ is so small that we have neglected this effect in all further calculations and have only considered ‘hcp-hollow sites’ for surface oxygen atoms and OH–groups.

In the next step we construct different surface models of hydrogen covered (000$`\overline{1}`$)–O surfaces using a similar procedure as in Sec. III.2. We take slabs with (1$`\times `$1), (1$`\times `$2), (1$`\times `$3) and (2$`\times `$2) surface unit cells, and we add different amounts of hydrogen to create hydrogen coverages of $`c_\mathrm{H}`$ of 1/4, 1/3, 1/2, 2/3, 3/4 and 1 monolayer. The calculated surface energy changes $`\mathrm{\Delta }\gamma `$ relative to the clean (000$`\overline{1}`$)–O surface are plotted in Fig. 3 as a function of the hydrogen chemical potential $`\mathrm{\Delta }\mu _\mathrm{H}`$. A very similar behavior as in Section III.2 arises: At H–rich conditions the structure with a 1/2 monolayer hydrogen coverage is the most stable surface, in agreement with mechanism (III). On the other hand, at H–poor conditions the clean surface without hydrogen becomes the most stable structure. As a new feature we find that between the two limiting cases the surfaces with 1/3 and 1/4 monolayer coverage are slightly more stable for small intervals of the chemical potentials.comment3  Thus, by lowering the chemical potential $`\mathrm{\Delta }\mu _\mathrm{H}`$ from H–rich to H–poor conditions it is possible to gradually reduce the hydrogen coverage from 1/2 monolayer to 1/3, 1/4 and zero coverage. Translating the chemical potential into temperature and pressure conditions we find that we will start to remove H at UHV-pressures at a temperature of roughly 750 K, which is in reasonable agreement with the experimental observation.woell  The other surface structures with hydrogen coverages larger than 1/2 are always higher in energy and will be unstable for all temperatures and hydrogen partial pressures. In particular a surface with a full monolayer of hydrogen as postulated from the results of the HAS experimentwoell  is not likely to exist in thermodynamic equilibrium. However, from the intensity of the He–specular peak it was deduced that only about 0.1 % of the (000$`\overline{1}`$)–O surface consists of flat terraces with diameters exceeding 50 Å which contribute to the (1$`\times `$1) HAS signal.woell  Therefore, the H covered surface with a (1$`\times `$1) HAS diffraction pattern may well be a minority phase which is formed under suitable kinetic conditions.

In Table 4 we have summarized the H binding energies $`E_\mathrm{b}`$ per atom which we have calculated as energy difference between the ideal (000$`\overline{1}`$)–O surface plus $`1/2E_{\mathrm{H}_2}`$ and surface structures with H coverages of $`c_\mathrm{H}`$. We see that it becomes rapidly unfavorable to adsorb more hydrogen as soon as the concentration for ideal charge compensation of the polar surface of 1/2 monolayer is reached. The reason for this behavior is the same as in the case of the oxygen vacancies: Up to 1/2 monolayer of hydrogen we fill the partially occupied O–$`2p`$–band, beyond 1/2 monolayer the O–$`2p`$–band is completely filled and we have to populate the conduction band. The decrease in $`E_\mathrm{b}`$ from 1/2 to 1/4 monolayer coverage indicates that also at low coverages a weak repulsive interaction between the hydrogen atoms remains. This is the reason why also the low-coverage structures appear in the surface phase diagram. If we extrapolate $`E_\mathrm{b}`$ toward zero coverage of the surface we can estimate an initial binding energy of about 2.3 eV per H atom for the dissociative adsorption of H2.

### III.4 Water Dissociation

As evident from Fig. 3 hydrogen adsorption plays a major role for the stabilization of the polar (000$`\overline{1}`$)–O surface and for almost all conceivable experimental conditions hydrogen will be present on the surface. But until now we have only considered molecular hydrogen as the reservoir for the surface hydrogen. However, in many chemical reactions and catalytic processes also water is present. Therefore we will briefly explore if also water can act as a source for surface hydrogen in the presence of the (0001)–Zn face.

In a DFT study using the hybrid B3LYP functional Wander and Harrisonwander\_h  found that dissociating water and forming a full H and OH monolayer on the O– and Zn–terminated surface, respectively, is energetically unfavorable by 0.1 eV. As adsorption sites for the OH–groups they assumed the Zn ‘on-top’ positions which would be the next lattice sites for O ions if the crystal is extended. However, on the polar surfaces two more high-symmetry adsorption sites exist: the ‘hcp-hollow site’ position above atoms in the second surface layer and a ‘fcc-hollow’ site with no atoms beneath.

Using the same adsorption geometry as Wander and Harrison we also find that dissociating water is energetically unfavorable with a slightly larger energy cost of 0.3 eV. However, as shown in Table 5, the configuration with the OH groups adsorbed at the ‘fcc-hollow sites’ instead of the ‘on-top’ positions is much lower in energy. Considering the correct adsorption positions for the OH groups we now find that even for the thermodynamically unstable monolayer coverages the dissociation of water is energetically preferred by about 0.4 eV. Taking only 1/2 monolayer coverages into account, this energy gain will be significantly larger.

### III.5 Surface Phase Diagram

Finally we combine the results of the previous subsections and assume that the polar (000$`\overline{1}`$)–O surface is now simultaneously in equilibrium with an O2 and a H2 gas phase. In addition to the surface models described in Sec. III.2 and Sec. III.3 we have furthermore considered various mixed structures of O vacancies and adsorbed H atoms in the (1$`\times `$2), (1$`\times `$3) and (2$`\times `$2) surface unit cells.

The surface free energy now depends on two chemical potentials $`\mathrm{\Delta }\mu _\mathrm{O}`$ and $`\mathrm{\Delta }\mu _\mathrm{H}`$. The graphs of Fig. 2 and 3 therefore have to be extended to a 3-dimensional plot. Such a diagram would be rather complex and hard to follow. The most important information contained in the plot of the surface free energies is which of the surface models has the lowest surface energy for a given combination of chemical potentials $`\mathrm{\Delta }\mu _\mathrm{O}`$ and $`\mathrm{\Delta }\mu _\mathrm{H}`$. This information is better visualized if we project the 3-dimensional diagram onto the ($`\mathrm{\Delta }\mu _\mathrm{O}`$,$`\mathrm{\Delta }\mu _\mathrm{H}`$) plane and only mark the regions for which a certain surface structure is the most stable one. The result is a phase diagram of the (000$`\overline{1}`$)–O surface which is shown in Fig. 4.

The surface phase diagram in Fig. 4 summarizes in a condensed fashion the essential results of our study. This phase diagram is dominated by two structures: a surface with 1/2 monolayer of adsorbed H and a hydrogen-free surface with 1/4 of the oxygen atoms removed. These are the two scenarios denoted as mechanism (II) and (III) in the Introduction, indicating that filling the O–$`2p`$–bands is a very important mechanism to stabilize the (000$`\overline{1}`$)–O surface. Next to these two phases we find two structures in which H is gradually removed from the surface and only at very H–poor conditions and when plenty of oxygen is available, the ideal O–terminated surface stabilized by mechanism (I) becomes the most stable structure.

None of the additionally considered mixed structures which simultaneously contain O vacancies and H adatoms appears in the phase diagram. This is not very surprising since for the given sizes of the surface unit cells no combination of O vacancies and H adatoms exists that leads to fully occupied O–$`2p`$–bands. However, it is very well conceivable that for larger surface unit cells mixed structures with, for example, an O vacancy concentration of 1/8 and a H coverage of 1/4, become more stable which would then appear as new phases between the H covered and the O vacancy structures.

Relating the chemical potentials to temperature conditions and partial pressures of the gas phase shows, see Fig. 4, that for almost all realistic experimental conditions hydrogen will be present at the (000$`\overline{1}`$)–O surface. Even at UHV-conditions with a low hydrogen partial pressure one has to go to rather high temperatures to fully remove all hydrogen. In this case, a surface structure with O vacancies will become the most stable one. In order to stabilize the ideal O–terminated surface an oxygen atmosphere with an extremely low content of hydrogen (and also water vapor) is necessary, which is basically not achievable in experiment.

In Table 6 we have summarized the surface relaxations for the three most important surface structures appearing in the surface phase diagram, Fig. 4. For the extended surface structures with H adatoms and O vacancies we have averaged in each atomic plane parallel to the surface the atomic displacements, and we define the interlayer distances $`d`$ as the separation of the averaged atomic positions. Depending on the surface structure and the charge compensation process, very different surface relaxations occur. The largest relaxations are found for the clean, defect-free surface termination with a compression of the first double-layer distance of almost 50 % and also a significant contraction of the second and subsequent double-layer spacings. This is in agreement with other previous ab-initio studies.noguera ; zno34 ; wander ; carlsson  After filling the partially occupied surface bands by adsorbing 1/2 monolayer of H or by removing 1/4 of the O ions, however, the surface layers relax back to almost truncated-bulk-like positions. Thus, for surfaces with a lower H adatom or O vacancy concentration, surface relaxations between the two extremes are conceivable. This may explain why experimentally very different results for the surface relaxations were found. GIXD measurementswander ; jedrecy  have predicted an inward relaxation of the upper O–layer with a contraction of the first double-layer distance of 40 % and 20 %, respectively, whereas from LEEDduke  and LEISoverbury  experiments it was concluded that the first double-layer spacing is close to its bulk value.

## IV Summary and Conclusions

By combining first-principles density-functional calculations with a thermodynamic formalism we have determined lowest-energy structures of the polar O–terminated (000$`\overline{1}`$)–O surface of ZnO in thermal equilibrium with an O2 and H2 gas phase. This scheme allows us to extend our zero-temperature and zero-pressure DFT results to more realistic temperature and pressure conditions which are usually applied in surface science experiments or in catalytic processes like the methanol synthesis, and thus to bridge computationally the ‘pressure gap’.

The essential result of this approach is a phase diagram of the (000$`\overline{1}`$)–O surface which is plotted in Fig. 4. From this surface phase diagram we predict that hydrogen is adsorbed at the (000$`\overline{1}`$)–O surface for a wide range of temperatures and H2 partial pressures, including UHV–conditions. This is in agreement with the recent observations of a HAS experimentwoell  and was also confirmed in a study of the CO adsorption on the polar ZnO surfaces.zno35

We find a H binding energy of roughly 2.3 eV per atom if molecular hydrogen dissociates and adsorbs at the clean O–terminated surface. Furthermore we predict that in situations where both polar surface terminations are present (for example for powder samples) also the dissociative adsorption of water with H and OH–groups being adsorbed at the O– and Zn–terminated surface, respectively, is energetically preferable. However, as soon as a coverage of 1/2 monolayer of hydrogen is reached, the energy gain of adsorbing more hydrogen on the (000$`\overline{1}`$)–O surface drops very rapidly with increasing hydrogen coverage. Therefore no stable phases with more than 1/2 monolayer H coverage appear in the surface phase diagram, Fig. 4. In particular, a structure with a full monolayer of H as predicted in Ref. woell,  is not very likely to exist globally in thermodynamic equilibrium (which does not exclude a kinetic or local stabilization).

Going to low hydrogen partial pressures and higher temperatures it is possible to gradually remove the hydrogen from the surface and to form stable phases with less than 1/2 monolayer coverage of hydrogen. If all hydrogen is removed, oxygen vacancies will be created as was speculated in Ref. woell, . However, we find that a surface with a vacancy concentration of 1/4 is much more stable than a missing-row structure where 1/3 of the oxygens has been removed. Therefore, based on the limited set of surface structures taken into consideration in our study, we currently do not understand the (1$`\times `$3) structure observed in Ref. woell, .

Finally, at higher oxygen partial pressures the O vacancies will be filled and the clean, defect-free O–terminated surface becomes the most stable structure. However, the hydrogen partial pressure has to be very low so that we consider it very unlikely that a clean, defect-free (000$`\overline{1}`$)–O surface can be observed in experiment.

## V Acknowledgments

The author would like to thank Dominik Marx, Christof Wöll, Georg Kresse, and Ulrike Diebold for fruitful discussions. The work was supported by SFB 558 and FCI.


format: text

text:
# Muon transfer from deuterium to helium

## I Introduction

Muon transfer from hydrogen to helium is a loss channel in muon catalyzed fusion ($`\mu `$CF), the muon induced fusion of hydrogen isotope nuclei Breunlich et al. (1989). In the $`\mu `$CF cycle, where in favorable cases a negative muon can catalyze up to 200 fusions, muon transfer to helium limits the fusion yield. Muon transfer from hydrogen to helium can happen during the cascade in muonic hydrogen (excited state transfer) or from the muonic hydrogen 1s ground state through the formation of an excited, metastable hydrogen–helium molecule $`(\mathrm{h}\mu \mathrm{He})^{}`$ (h = proton p, deuteron d, or triton t, and He $`={}_{}{}^{3}\mathrm{He}`$ or $`{}_{}{}^{4}\mathrm{He}`$), a reaction first proposed by Aristov Aristov et al. (1981). These molecules decay from the excited state to the unbound ground state mostly by x–ray emission ($`\mathrm{E}_x6.8`$ keV). Auger electron emission and $`\mathrm{h}\mu \mathrm{He}`$ breakup are also possible. The scheme of the principal transfer and decay processes is presented in Fig. 1. The muon entering a deuterium–helium mixture may be captured either by deuterium (with probability $`W_\mathrm{d}`$) or by helium (with probability $`W_{\mathrm{He}}`$) via direct capture. The two vertical arrows indicate the cascade of the muon to the 1s ground state. The $`\mathrm{q}_{1s}^{\mathrm{He}}`$ represents the probability for the $`\mu \mathrm{d}^{}`$ to reach the ground state in the presence of helium. Excited state transfer is shown by the upper horizontal arrow. Ground state transfer is shown with a rate $`\lambda _{\mathrm{d}\mu \mathrm{He}}`$ via the $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecule. The $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecule decay channels are shown with rates $`\lambda _e`$ for the Auger decay, $`\lambda _\gamma `$ for the x–ray channel, and $`\lambda _p`$ for the break up channel, respectively.

The energies and widths of the molecular states have been characterized by measuring the x–ray energy spectra. The most precise experiment on $`(\mathrm{p}\mu \mathrm{He})^{}`$, with its intrinsically low x–ray yield, was carried out by our collaboration Tresch et al. (1998). The $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecules were also studied by our collaboration Gartner et al. (2000) and recently an experiment was performed on $`(\mathrm{t}\mu \mathrm{He})^{}`$ Matsuzaki et al. (2002). In those publications, earlier less precise experiments were referenced and discussed in detail. Our precision of $`0.2`$% for the energies and $`1.2`$% for the widths of the $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecular deexcitations make detailed comparisons with calculations possible. Precise results on the excited state transfer probabilities were also obtained. The combined use of results obtained with standard detectors and CCD techniques allowed us to determine the radiative branching ratio of the $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecules, a value which has been of considerable theoretical interest in recent years due to its direct and unique connection to the wave function overlap in the muonic molecule Kino and Kamimura (1993). The experimental challenges in obtaining the results were overcome thanks to months of beam time, and the use of large area Charge Coupled Device (CCD) x–ray detectors, as well as germanium detectors Gartner et al. (2000). A more detailed description of the present work can be found in Augsburger’s thesis Augsburger (2001).

## II Experiment

The experiment was performed at the $`\mu `$E4 channel at the Paul Scherrer Institute (PSI), Villigen, Switzerland. Setup, Ge and Si(Li) detector, gas handling and target conditions can be taken from Tresch Tresch et al. (1998) and Gartner Gartner et al. (2000). Figure 2 shows the target setup with the detectors. Detailed information about the large area CCD x–ray detectors is found in Tresch Tresch et al. (1998).

### II.1 Experimental conditions

The experimental setup consisted of a gas target, Ge and Si(Li) detectors, scintillators, and CCD detectors, as shown in Fig. 2. Both Tresch Tresch et al. (1998) and Gartner Gartner et al. (2000) measured the molecular formations rates $`\lambda _{\mathrm{h}\mu \mathrm{He}}`$ in protium and deuterium. The CCDs were used by Tresch Tresch et al. (1998) for the protium measurement. The present work shows our results for the deuterium measurement using CCDs.

The deuterium and related measurements, as well as the gas handling and mixture analysis are described in great detail in Gartner Gartner et al. (2000). From this reference, we summarized in Table 1 the gas mixtures used for our analysis. The choice of helium concentration was dictated by the goal of Gartner Gartner et al. (2000) measurement, namely the molecular formation rates. Due to the different theoretical values for the rates in $`\mathrm{D}_2+{}_{}{}^{3}\mathrm{He}`$ and $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$, the relative concentrations of $`{}_{}{}^{3}\mathrm{He}`$ and $`{}_{}{}^{4}\mathrm{He}`$ are very different.

### II.2 CCDs as low energy x–ray detectors

CCDs are excellent x–ray detectors in the energy range from 1 keV to 15 keV (details can be found in Egger (1999)). In most cases, an x ray produces charge in only a single pixel, whereas charged particles produce cluster events or tracks with more than one adjacent pixel hit. The usual way to distinguish x ray event pixels from charged particle, neutron, and higher energy gamma–ray background is to require that none of the eight surrounding pixels have a charge that is considerably above the noise level. The CCD type used for this experiment was a silicon based MOS type, model CCD–05–20 by E2V.111E2V, Technologies Inc, Waterhouse Lane, Chelmsford, Essex, CM1 2QU, England (previously EEV and Marconi). Each CCD chip has a size of 4.5 cm2 ($`770\times 1152`$ pixels of area $`22.5\times 22.5\mu `$m2) and a depletion depth of $`30\mu `$m. In this experiment the energy resolution of the muonic deuterium $`\mu \mathrm{D}(21)`$ line was 130 eV FWHM and the muonic helium $`\mu \mathrm{He}(21)`$ transition had a resolution of 215 eV FWHM.

Unfortunately, our CCDs cannot be triggered so there is no timing information. The CCD data were read out approximately every 3 minutes by a data-acquisition system which operated independently from the data acquisition of the other detectors. Therefore, we cannot normalize the collected data to the incoming muon rate. The results of CCD’s measurements can only be analyzed using absolute numbers.

## III Analysis and results

### III.1 Analysis of the x–ray energy spectra

We present in this section the different spectra obtained and explain the fitting procedures. Two CCDs were used. Since each CCD half was read out separately, we have 4 sets of measurements. At first, the data from each half CCD were analysed separately in order to detect any possible malfunction and to perform the energy calibration and background reduction by single-pixel analysis. After checking that the separate treatment of each half CCD gave consistent results, the calibrated energy spectra were added and the fits performed on the summed spectra.

#### III.1.1 Pure element spectra

The x–ray spectra from single element targets were studied in detail to find as much information as possible about the detector response function and target related backgrounds. The final results required that the entire energy range be fit at once and this was accomplished in several steps.

The final result for the muonic deuterium x–ray spectrum is presented in Figs. 3 and 4, namely the $`\mu \mathrm{d}`$ x–ray transition to the 1s ground state and a series of contaminant peaks essentially at higher energies. Also visible are the electronic K$`\alpha `$ and K$`\beta `$ transitions of Si (CCD), Cr, Mn, Fe, Ni, and Cu (target), although only the positions of the K$`\alpha `$ peaks are indicated in Fig. 4 (except for Fe, where both lines are clearly visible). In addition, muonic aluminium at 10.66 keV from the $`\mu \mathrm{Al}(54)`$ transition is clearly visible. Other muonic aluminium transitions, $`\mu \mathrm{Al}(65)`$ at 5.79 keV and $`\mu \mathrm{Al}(76)`$ at 3.49 keV, are also present but much weaker than the $`\mu \mathrm{Al}(54)`$ line.

The first fits were made using Gaussian peak shapes and a standard CCD background with the goal of locating all lines and characterizing the continuous background. The CCD background in a high–noise environment was studied in detail in Ref. Augsburger (2001). The large hill starting at 7 keV seen in Fig. 4 is due essentially to energy deposited by electrons crossing the CCD. The relative importance of each contaminant was estimated with the K$`\beta `$/K$`\alpha `$ intensity ratio held fixed according to values given in Ref. Salem et al. (1974). Once the relative intensities of the contaminant peaks were known, the first fits to the full spectra were carried out. The variation of the contaminant intensities for the different spectra were found to be small, and hence could be well parametrized.

Figure 5 presents the spectrum for pure $`{}_{}{}^{4}\mathrm{He}`$. The Lyman series between 8 and 11 keV and the Balmer series around 2 keV are clearly seen and the contaminant peaks are the same as in the muonic deuterium spectrum.

To further refine the fit, the muonic helium $`\mu \mathrm{He}(21)`$ transition was examined in detail to fully understand the true CCD response function for the line shape. This peak was chosen, even though it contained a small copper contamination (less than 1%), since it has high statistics and is well separated from the other peaks. The small copper contribution was subtracted. Since we could not use an analytical function for fitting the remaining asymmetric peaks, we interpolated the asymmetric peak shape of the $`\mu \mathrm{He}(21)`$ line, shown in Fig. 6, to fit the data correctly. As one can see, it looks like a Gaussian with an asymmetry on the left side. This same peak shape was then used to represent all other lines, replacing the Gaussian shape, and the spectra refit. In particular, the FWHM of each peak was obtained by using a scaling factor from the FWHM given in Fig. 6. In addition, all peak positions were defined by the center of gravity (not by the maximum). The values for the resolution of the $`\mu \mathrm{He}(21)`$ transition as well as the electronic lines are shown in Fig. 7. The fitting procedure was repeated for both the muonic deuterium and muonic helium spectra until a minimum $`\chi ^2`$ was obtained. Replacing the Gaussian line shape with the final fit function including background parametrization and asymmetric line shapes reduced the $`\chi _{dof}^2`$ from 5 to about 1.4 for both spectra.

#### III.1.2 Spectra of the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ and $`\mathrm{D}_2+{}_{}{}^{3}\mathrm{He}`$ mixtures

Figure 8 presents the energy spectrum of the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ mixture. In addition to the peaks from muonic deuterium, muonic helium, and the contaminants, a large x–ray peak from the decay of the $`(\mathrm{d}\mu {}_{}{}^{4}\mathrm{He})^{}`$ molecule appears around 6.8 keV. Again, the fitting procedure outlined above was used for all peaks except the molecular peak, for which theoretical curve has been calculated and is given in Ref. Belyaev et al. (1997). The difference in shape corresponding to decays of the $`J=0`$ and $`J=1`$ state, respectively, is negligible in our case. Hence, the calculated shape of Ref. Belyaev et al. (1997) was taken for the shape of this molecular peak. The positions of the maximum is a free parameter. We used two scaling factors to determine the amplitude and FWHM relative to the theoretical shape. Figure 9 shows the fit of the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ mixture in the region of the molecular peak (the fit was carried out over the whole energy region, 1.6 to 11.25 keV). The results are given in Table 2.

Figure 10 presents the spectrum of the $`\mathrm{D}_2+{}_{}{}^{3}\mathrm{He}`$ mixture and Fig. 11 shows the same spectrum in the region of the molecular peak. The analysis of that spectrum was identical to the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ analysis with results also given in Table 2.

### III.2 Relative intensities of the K series transitions in muonic $`{}_{}{}^{4}\mathrm{He}`$

The relative intensities of the K series transitions in pure muonic $`{}_{}{}^{4}\mathrm{He}`$ are given in Table 3. The errors include a statistical part (fit) and a systematic part (CCD detection efficiency). The main error comes from the CCD efficiency which is not surprising since the fit parameters for “CCD depletion depth” and “CCD window thickness” converge in a range of ($`28\pm 2`$$`\mu `$m and ($`35\pm 1`$$`\mu `$m respectively. The results are compared to Tresch et al. (1998) where no isotopic effect ($`{}_{}{}^{3}\mathrm{He}`$ or $`{}_{}{}^{4}\mathrm{He}`$) was seen (last column of Table 3). The agreement is excellent for the K$`\alpha `$ transition, and the significant discrepancies of the other values are understood, since the measurement of Tresch Tresch et al. (1998) was carried out at a lower density, which explains the K$`\beta `$ decrease and the K$`\gamma `$ increase. In addition, our accumulated experience with CCD background and detection efficiency resulted in a better fit in this work but we realize that the errors given in Tresch et al. (1998) were underestimated with respect to the CCD efficiency correction.

### III.3 Excited state transfer and the $`\mathrm{q}_{1s}^{\mathrm{He}}`$ probability

The $`\mathrm{q}_{1s}`$ value represents the probability for a newly formed light muonic atom to fully deexcite to the 1s state when the muon also has the possibility of transferring directly from an excited state to a heavier nucleus (cf. Fig. 1). In binary mixtures, the notation often includes the identity of the heavier nucleus, $`\mathrm{q}_{1s}^{{}_{}{}^{3}\mathrm{He}}`$ for example.

We begin our analysis with the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ mixture. The number of events in the K series transitions in $`\mu \mathrm{d}`$ and $`\mu {}_{}{}^{4}\mathrm{He}`$ in the gaseous mixture of $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ is given in Table 4. The sum of $`\mu \mathrm{d}`$ events represents the total number of $`\mu \mathrm{d}`$ which reach the ground state and is called $`\mathrm{N}_{\mu \mathrm{d}}^{{}_{}{}^{4}\mathrm{He}}(1s)`$.

Part of the $`{}_{}{}^{4}\mathrm{He}`$ events come from the direct capture of the muon by helium, the other part by excited state transfer from muonic deuterium. It was shown in Tresch et al. (1998) that in gaseous $`\mathrm{H}_2+{}_{}{}^{4}\mathrm{He}`$ mixtures, excited–state transfer proceeds only to the levels $`n=3`$ and $`n=2`$ of $`\mu ^4\mathrm{He}`$. A detailed comparison of Figs. 8 and 5 also shows that in the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ mixture there is a large enhancement of the K$`\alpha `$ and K$`\beta `$ He lines over the higher transitions. In addition, some L transitions can also be seen in Fig. 5. The fact that the L$`\gamma `$ transition (in Fig. 8) contains approximately $`\mathrm{200\hspace{0.25em}000}`$ (efficiency corrected) events versus only 2000 for the L$`\delta `$ line further confirms the above hypothesis. Therefore the sum of the events from the $`\mu \mathrm{He}(41)`$, $`\mu \mathrm{He}(51)`$, and $`\mu \mathrm{He}(61)`$ transitions in $`\mu {}_{}{}^{4}\mathrm{He}`$ is due to direct capture. Taking this sum from Table 4, one gets a measured number of direct capture, $`\mathrm{N}_{\mathrm{dc}}^\mathrm{m}=(65.5\pm 6.2)\times 10^3`$, where we added the errors quadratically.

In the spectrum from pure $`\mu {}_{}{}^{4}\mathrm{He}`$ (see results in Table 3) the percentage sum of the Ki/Ktotal fractions for $`i4`$ is $`25.20\pm 1.81`$%. The $`\mathrm{N}_{\mathrm{dc}}^\mathrm{m}`$ therefore corresponds to 25.20% of the total number of K–series x rays in $`\mu {}_{}{}^{4}\mathrm{He}`$ ($`\mathrm{N}_{\mathrm{dc}}^{\mathrm{tot}}`$). Thus, we deduced the total number of direct capture events $`\mathrm{N}_{\mathrm{dc}}^{\mathrm{tot}}=(260\pm 43)\times 10^3`$. This number will allow us to differentiate between direct capture and excited state transfer events in the K$`\alpha `$ and K$`\beta `$ intensities, measured in the mixture of deuterium and $`{}_{}{}^{4}\mathrm{He}`$.

The numbers of K$`\alpha `$ and K$`\beta `$ events occurring with the direct capture, $`\mathrm{N}_{\mathrm{dc}}^{\mathrm{K}_\alpha }`$ and $`\mathrm{N}_{\mathrm{dc}}^{\mathrm{K}_\beta }`$, were obtained using the intensity ratios Ki/Ktotal determined in the pure $`\mu {}_{}{}^{4}\mathrm{He}`$ spectrum and $`\mathrm{N}_{\mathrm{dc}}^{\mathrm{tot}}`$. The total number of $`\mu {}_{}{}^{4}\mathrm{He}`$ K$`\alpha `$ and K$`\beta `$ events in the mixture is given in Table 4. The differences are due to excited state transfer from $`\mu \mathrm{d}^{}`$. The number of K$`\alpha `$ and K$`\beta `$ events coming from excited state transfer, $`\mathrm{N}_{\mathrm{exc}}^{\mathrm{K}_\alpha }`$ and $`\mathrm{N}_{\mathrm{exc}}^{\mathrm{K}_\beta }`$, are the difference between the first two lines of Table 4 and the previously determined $`\mathrm{N}_{\mathrm{dc}}^{\mathrm{K}_\alpha }`$ and $`\mathrm{N}_{\mathrm{dc}}^{\mathrm{K}_\beta }`$. Therefore the sum of events from excited state transfer is $`\mathrm{N}_{\mathrm{exc}}^{\mathrm{tot}}=(398\pm 47)\times 10^3`$. Now $`\mathrm{q}_{1s}^{{}_{}{}^{4}\mathrm{He}}`$ can be determined by

$$\mathrm{q}_{1s}^{{}_{}{}^{4}\mathrm{He}}=\frac{\mathrm{N}_{\mu \mathrm{d}}^{{}_{}{}^{4}\mathrm{He}}(1s)}{\mathrm{N}_{\mu \mathrm{d}}^{{}_{}{}^{4}\mathrm{He}}(1s)+\mathrm{N}_{\mathrm{exc}}^{\mathrm{tot}}}=90.1\pm 1.5\%$$
(1)
where $`\mathrm{N}_{\mu \mathrm{d}}^{{}_{}{}^{4}\mathrm{He}}(1s)`$ is the total number of $`\mu \mathrm{d}`$ Lyman x rays in the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ mixture (see Table 4).

The analysis carried out in the case of the $`\mathrm{D}_2+{}_{}{}^{3}\mathrm{He}`$ mixture was the same as in the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ mixture with the additional hypothesis that the muonic cascade was the same in both $`\mu {}_{}{}^{3}\mathrm{He}`$ and $`\mu {}_{}{}^{4}\mathrm{He}`$. Pure $`\mu {}_{}{}^{3}\mathrm{He}`$ was not measured (only $`\mu {}_{}{}^{4}\mathrm{He}`$) for this work. However, Tresch Tresch et al. (1998) has shown no isotopic effects between the two gases.

Therefore, the $`\mathrm{q}_{1s}^{{}_{}{}^{3}\mathrm{He}}`$ is determined as

$$\mathrm{q}_{1s}^{{}_{}{}^{3}\mathrm{He}}=\frac{\mathrm{N}_{\mu \mathrm{d}}^{{}_{}{}^{3}\mathrm{He}}(1s)}{\mathrm{N}_{\mu \mathrm{d}}^{{}_{}{}^{3}\mathrm{He}}(1s)+\mathrm{N}_{\mathrm{exc}}^{\mathrm{tot}}}=68.9\pm 2.7\%$$
(2)
where $`\mathrm{N}_{\mu \mathrm{d}}^{{}_{}{}^{3}\mathrm{He}}(1s)`$ is the total number of $`\mu \mathrm{d}`$ Lyman x rays in the $`\mathrm{D}_2+{}_{}{}^{3}\mathrm{He}`$ mixture (Table 4).

### III.4 Radiative branching ratio $`\kappa _{\mathrm{d}\mu \mathrm{He}}`$ of the $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecule

The radiative branching ratio $`\kappa _{\mathrm{d}\mu \mathrm{He}}`$ for the $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecular decay can be determined the same way as in Tresch et al. (1998) for the $`(\mathrm{p}\mu \mathrm{He})^{}`$ molecule. $`\kappa _{\mathrm{d}\mu \mathrm{He}}`$ is given by

$$\kappa _{\mathrm{d}\mu \mathrm{He}}W=\frac{\mathrm{N}_{\mathrm{d}\mu \mathrm{He}}}{\mathrm{N}_{\mu \mathrm{d}}^{\mathrm{He}}(1s)},$$
(3)
where $`\mathrm{N}_{\mathrm{d}\mu \mathrm{He}}`$ is the number of events in the molecular peak (see Table 2) and $`\mathrm{N}_{\mu \mathrm{d}}^{\mathrm{He}}(1s)`$ the total number of the $`\mu \mathrm{d}`$ Lyman series x rays in the mixture. $`W`$ is the probability of a $`\mu \mathrm{d}_{1\mathrm{s}}`$ forming a $`(\mathrm{d}\mu \mathrm{He})^{}`$ molecule and is given by the equation:

$$W=\frac{\varphi \mathrm{c}_{\mathrm{He}}\lambda _{\mathrm{d}\mu \mathrm{He}}}{\mathrm{\Lambda }_{\mu \mathrm{d}_{1s}}},$$
(4)
where $`\varphi `$ is the atomic density of the mixture, normalized to LHD, $`\mathrm{c}_{\mathrm{He}}`$ is the helium atom proportion, $`\lambda _{\mathrm{d}\mu \mathrm{He}}`$ the ground state transfer rate from $`\mu \mathrm{d}`$ to He, and $`\mathrm{\Lambda }_{\mu \mathrm{d}_{1s}}`$ the disappearance rate of muons from the $`(\mu \mathrm{d})_{1s}`$ level.

The so determined values for $`\kappa _{\mathrm{d}\mu \mathrm{He}}`$ and $`W`$ are given in Table 5 for both helium isotopes. The values for $`\lambda _{\mathrm{d}\mu \mathrm{He}}`$ and $`\mathrm{\Lambda }_{\mu \mathrm{d}_{1s}}`$ (in Table 5) were taken from Gartner Gartner et al. (2000). While the errors for $`\mathrm{N}_{\mathrm{d}\mu \mathrm{He}}`$ and $`\mathrm{N}_{\mu \mathrm{d}}^{\mathrm{He}}(1s)`$ include both the statistical and systematic uncertainties, the errors in $`\varphi `$ and $`\mathrm{c}_{\mathrm{He}}`$ given in Table 1 are purely systematic. The errors on $`W`$ and $`\kappa _{\mathrm{d}\mu \mathrm{He}}`$ were calculated by normal error propagation without specifying the type of error.

## IV Discussions

### IV.1 General features of the x–ray energy spectra

The spectra presented in Figs. 3, 4, 5, 8, and 10 deserve three general comments. First, the relative intensities of the muonic deuterium K$`\alpha `$ and K$`\beta `$ transitions are density dependent Lauss et al. (1998). Since the density changed between mixtures (see Table 1), the muonic deuterium K$`\alpha `$ peak is slightly enhanced over K$`\beta `$ in the $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ mixture.

Second, the x–ray count rate for the $`(\mathrm{d}\mu {}_{}{}^{3}\mathrm{He})^{}`$ molecule is smaller than for the $`(\mathrm{d}\mu {}_{}{}^{4}\mathrm{He})^{}`$ molecule since in the $`(\mathrm{d}\mu {}_{}{}^{3}\mathrm{He})^{}`$ case the two–particle breakup channel is more prominent. Third, the relative helium/deuterium line intensities depend on the helium concentration (see Table 1).

### IV.2 $`(\mathrm{d}\mu {}_{}{}^{3}\mathrm{He})^{}`$ and $`(\mathrm{d}\mu {}_{}{}^{4}\mathrm{He})^{}`$ molecules

In Table 6 our CCD results for the position of the maximum of the two molecular peaks are compared to theoretical predictions Belyaev et al. (1997); Czaplinski et al. (1997) and to results Gartner et al. (2000) obtained with Ge and Si(Li) detectors. The Ge and Si(Li) detector results seem to favor the J=1 state but the CCD results imply a preference for transitions from $`J=0`$. Even if the CCD results are more precise and the CCD statistics are significantly higher, it is difficult to decide for $`J=0`$ or $`J=1`$ since the results of both detector types are effectively compatible considering the errors. What can be said unambiguously is that the CCD results are in excellent agreement with both theoretical predictions for decay from the $`J=0`$ state. It should be stressed that both the Ge and Si(Li) and the CCD data were taken simultaneously during the experiment.

In Table 7 the CCD experimental FWHM widths for the molecular peaks are again compared with the theoretical predictions Belyaev et al. (1997); Czaplinski et al. (1997) and with the Ge and Si(Li) data Gartner et al. (2000). The CCD results are in very good agreement with theoretical predictions for both $`J=0`$ and $`J=1`$ states, but distinguishing between the two states is not possible due to the almost identical theoretical values. On the other hand, the results with the “classic” (Ge or Si(Li)) x–ray detectors are between 1.5 to 2.5 $`\sigma `$ away from theory. The somewhat smaller width of the ($`\mathrm{d}\mu {}_{}{}^{4}\mathrm{He}`$) molecule predicted by theory is also hinted at by our CCD data.

### IV.3 Radiative branching ratio $`\kappa `$

Table 8 presents the different theoretical values for the radiative branching ratio. The calculations are those of Kino Kino and Kamimura (1993), Gershtein Gershtein and Gusev (1993), Kravtsov Kravtsov et al. (1993), and Belyaev Belyaev et al. (1995, 1996). Except for Kravtsov Kravtsov et al. (1993) who includes the Auger decay channel $`\lambda _e`$ (see Fig. 1), only $`\lambda _p`$ (break up) and $`\lambda _\gamma `$ (x ray) are calculated. The x ray channel relates to $`\kappa _{\mathrm{d}\mu \mathrm{He}}`$ via the ratio

$$\kappa _{\mathrm{d}\mu \mathrm{He}}=\frac{\lambda _\gamma }{\lambda _p+\lambda _\gamma +\lambda _e}.$$
(5)
The different theoretical $`\kappa _{\mathrm{d}\mu \mathrm{He}}`$ values are compared with our experiment in Table 8. The large isotopic effect predicted by theory is seen by our experiment which is in contradiction to the $`(\mathrm{p}\mu \mathrm{He})`$ case Tresch et al. (1998). In general the agreement between theory and experiment is good, however, the experimental errors are sizable.

### IV.4 Ground state formation probabilities $`\mathrm{q}_{1s}^{\mathrm{He}}`$

The meaning of $`\mathrm{q}_{1s}^{\mathrm{He}}`$ has been described in Section III.3. Our results for the deuterium–helium mixtures, Eqs. (1) and (2), are listed together with those for the hydrogen–helium mixtures Tresch et al. (1998) in Table 9. It is interesting to note that in both cases (hydrogen or deuterium) $`\mathrm{q}_{1s}^{\mathrm{He}}`$ is smaller for $`{}_{}{}^{3}\mathrm{He}`$, and therefore, the excited state transfer is more probable for the lighter of the two helium isotopes. In the case of hydrogen–deuterium mixtures Lauss et al. (1996), $`\mathrm{q}_{1s}^{\mathrm{He}}`$ has been shown to depend on the concentration of the components of the mixture. The large difference seen in our case between $`\mathrm{D}_2+{}_{}{}^{3}\mathrm{He}`$ and $`\mathrm{D}_2+{}_{}{}^{4}\mathrm{He}`$ is therefore partially due to the differing helium concentrations (see Table 1). The second observation is that the $`\mathrm{q}_{1s}^{\mathrm{He}}`$ is significantly larger for deuterium, a result of consequence in the case of muon catalyzed fusion in deuterium–helium mixtures Knowles et al. (2001).

## V Conclusions

The use of CCDs for low energy x–ray detection allowed for a complete energy measurement of muonic deuterium, helium, and molecular $`(\mathrm{d}\mu \mathrm{He})`$ x rays with excellent energy resolution and low background. The large CCD surface resulted in an increased solid angle and therefore in better statistics when compared to traditional Ge or Si(Li) detectors. Of course, results like transfer rates still need the usual x–ray detectors since the CCDs give no timing information. The simultaneous use of CCDs and other x–ray detectors allows for systematic error checks of the experiment since the CCD electronics is completely independent. In conclusion, the addition of CCDs permitted a characterization of all transfer parameters and some high precision results. Only the use of CCD detectors allowed the determination of the radiative branching ratio, a result which was long awaited by theorists.

###### Acknowledgements.

Financial support by the Austrian Academy of Sciences, the Austrian Science Foundation, the Swiss Academy of Sciences, the Swiss National Science Foundation and the Beschleunigerlaboratorium of the University and the Technical University Munich is gratefully acknowledged.