- Email: [email protected]
Theory of H− in SiO2
Journal of Non-Crystalline Solids 289 (2001) 42±52
www.elsevier.com/locate/jnoncrysol
Theory of H in SiO2 Arthur H. Edwards *, W.M. Shedd, R.D. Pugh Air Force Research Laboratory, Space Vehicles Directorate, 3550 Aberdeen Ave. SE, KAFB, NM 87117-5776, USA Received 11 August 2000; received in revised form 18 May 2001
Abstract We present an ab initio molecular orbital study of H in SiO2 with special emphasis on H . We have calculated equilibrium geometries, vibrational spectra, binding energies and electrical levels. For the electrical levels, we have included long-range polarization eects in three approximations. We compare our results to those of Yokozawa and Miyamoto and others. We ®nd that H0 is unstable in SiO2 . However, we ®nd that, rather than disproportionate, it prefers to dimerize. We also predict that H0 is a deep electron trap. Finally, we ®nd that long-range polarization eects are crucially important for obtaining even qualitatively correct values for electrical levels. Ó 2001 Published by Elsevier Science B.V.
1. Introduction Radiation-induced charge traps in metal oxide semiconductor (MOS) insulators, especially, in SiO2 , have been studied intensely over the last three decades (see, for example, [1±4]). It is known that in a radiation environment, MOS oxides become positively charged, leading to reductions in n-channel and increases in p-channel threshold voltages. As a result of this experimental fact, and of the common wisdom that electrons, being much more mobile, are much less easily trapped, most of the literature has focused on the nature of hole traps in SiO2 . In fact, several hole traps have been unambiguously identi®ed in SiO2 including the E0 centers, the non-bridging oxygen center, and the self-trapped hole [5±7]. * Corresponding author. Tel.: +1-505 853 6042; fax: +1-505 846 4558. E-mail addresses: [email protected] (A.H. Edwards), [email protected] (W.M. Shedd), robert.pugh@ kirtland.af.mil (R.D. Pugh).
Using thermally stimulated current measurements, Fleetwood and co-workers [8] recently have shown that there is a very large density of radiation-induced electron traps in both gate and buried oxides. These observations have been con®rmed in conjunction with Leray and co-workers [9]. Also, Stahlbush et al. [10] have observed shallow electron traps in buried oxides using standard C±V measurements. The large density of traps, along with a large range of depths in energy has important conceptual and technological consequences. The conceptual importance arises, at least in part, because the measured positive charge is net charge. Thus, there is a much larger density of hole traps than was originally thought. Technologically, the charging and discharging character of these traps may be important if we are to understand longterm radiation eects in gate oxides and in buried oxides. While the presence of large numbers of electron traps has been con®rmed, the understanding of the microscopic nature of these traps is limited. To date, no intrinsic electron traps have been
0022-3093/01/$ - see front matter Ó 2001 Published by Elsevier Science B.V. PII: S 0 0 2 2 - 3 0 9 3 ( 0 1 ) 0 0 6 4 9 - 4
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
experimentally identi®ed. Several theoretical studies in the literature predict that the oxygen vacancy [11,12] and the non-bridging oxygen (NBO) center [13] can trap electrons. In the former case, the predicted binding energy was small (0.5 eV). In the latter case the predicted binding energy is 1.2 eV. However, c-induced NBOs have only been observed in bulk, fused silica. Only after a-irradiation has it been observed in MOS oxides [14]. Yokozawa and Miyamoto [15] (hereafter YM) were the ®rst to show that H was stable in bcristobalite. The bonding structure is similar to that in Fig. 1. In fact, based on calculations using
43
density functional theory (DFT) applied to a fairly small super cell, this group predicted that hydrogen would exhibit a negative Ueff . This energy measures the stability of two neutral hydrogen atoms relative to one H and one H . When Ueff < 0, the charged species are favored over the neutral species. This is discussed in more detail in Section 2. Since then, other groups have obtained the same prediction using DFT applied to super cells [16,17]. This result is so interesting and so surprising that we decided to study it using a completely dierent formalism. By doing this we have deepened our understanding of the origin of the disproportionation. Furthermore, while there are almost certainly several dierent electron traps, we believe that this defect may be important enough, because of the persistent presence of hydrogen in MOS oxides, to study in depth. We present a cluster-based, Hartree±Fock (HF) and post Hartree±Fock, molecular orbital study of this defect using the clusters shown in Fig. 1. We have calculated a variety of properties, including equilibrium geometries, Ueff , the vibrational spectrum, and the electrical level structure. In the study of electrical levels, we used several approximate techniques for treating long-range polarization eects. We ®nd that polarization effects are fundamentally important; they change qualitatively our predictions about the stability of H0 . 2. Theoretical methods We have used GAMESS [18] and Gaussian-94 [19], ab initio, Hartree±Fock molecular orbital packages capable of treating electron correlation in several approximations. 2.1. Geometries and vibrational spectra
Fig. 1. Clusters used in present study. Clusters in (a) and (b) are referred to as I and II in the text. In (a) Oc indicates the central oxygen. In (a) and (b) H indicates added hydrogen atom. In (b) Oa and Oe denote axial and equatorial oxygen atoms, respectively.
For both clusters shown in Fig. 1, we used the 6-31G all-electron atomic orbital basis set. Here the 6 indicates that all core atomic orbitals are represented by six Gaussian functions, the 31 indicates that the valence orbitals are represented by two atomic functions, one, expanded in three
44
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
Gaussian functions, with a smaller spatial extent, and one, expanded in only one Gaussian, with a larger spatial extent. The two asterisks indicate that we used d-polarization functions on the heavy atoms and p-polarization functions on the hydrogen atoms. We allowed all atoms to relax freely in the search for energy minima. For the smaller cluster, we have used perturbative many-body theory (Mùller±Plesset theory [20]) during the geometry optimizations to study the eects of electron correlation. These calculations are labeled MP2. The charged states were calculated by adding or removing an electron and allowing the system to adjust self-consistently. For the neutral case, we calculated equilibrium geometries without the extra hydrogen atom, and added the total energies of this relaxed, neutral cluster to the energy of a free hydrogen atom. We chose this method for two reasons. First, when we performed neutral calculations with the added hydrogen, the equilibrium geometry was asymptotic; the hydrogen atom merely drifted away from the cluster. This was a practical diculty because the gradient would be above the tolerance value even if the main portion of the cluster was essentially converged. The second reason is that YM found that in b-cristobalite, the equilibrium position of the H0 atom was at the center of a large axial cavity, and that the total energy of the system was, to a very good approximation, the sum of the perfect crystal energy and the energy of an isolated hydrogen atom. This is not surprising, since the 1 H hyper®ne value in SiO2 is approximately unperturbed from the gas-phase value, indicating that the wave function was little perturbed by the solid [21,22]. We have calculated the vibrational spectrum for H two ways. We have calculated a one-dimensional potential surface from the Hartree± Fock and MP2 total energies without inclusion of polarization eects. Using our own cubic ®nite element program, we have calculated the ®rst few (10) vibrational levels and used these to calculate the standard vibrational parameters, xe , xe xe , and xe ye , representing the quadratic, cubic and quartic terms in the Taylor series expansion of the potential surface, by least-squares ®tting to Eq. (1)
2 1 1 x e xe n En hc xe n 2 2 3 ! 1 xe ye n : 2
1
We have also used Gaussian-94 to calculate the full vibrational spectrum, along with the IR intensities, in the harmonic approximation. 2.2. Electrical levels and Ueff We have calculated electrical levels approximately by using Etot n 1=n Etot n 1
Etot n v;
2
where v is the electron anity of SiO2 , and Etot n is the total energy for the system with n excess electrons, in the appropriate geometrical conformation. We have used the experimental value, 0.9 eV, for v. The ®rst two terms on the right-hand side merely give the energy required to move an electron to in®nity from a defect with n 1 excess electrons. The addition of v in Eq. (2) moves the reference energy from the vacuum level to the bottom of the conduction band. Experimentally, this de®nition gives the position in the forbidden gap where the defect will undergo a change in occupation from n to n 1 electrons at T 0 K. We note that on band diagrams the electrical levels are labeled by the actual charge states, rather than by the number of excess electrons. So, for example, the electrical level for the transition from neutral to negative would be Etot =0 Etot
Etot 0 v:
3
Of particular interest in this study is the calculation of Ueff , de®ned by Ueff Etot Etot
2Etot 0:
4
Eq. (4) shows clearly that when Ueff is negative, it is energetically favorable for two neutral hydrogen atoms spontaneously to disproportionate into a positive and a negative ion. As YM and others [23±25] have pointed out, there are intrinsic problems in using DFT with periodic boundary conditions to calculate electri-
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
cal levels of defects. (Here we assume that the defect in the solid is represented by a super cell and that the defect wave function is localized enough to be negligibly small outside this super cell.) For charged states, the relevant lattice sums diverge unless the unit cell is neutral. It is typical to introduce a uniform background charge of opposite sign. While this does yield the required neutrality, it also leads to spurious lowering of total energy because of the interaction between the background charge and the added (or removed) electron in the band states. Furthermore, it leaves higher-order multipole ®elds that are unphysical. Even in neutral defects the ®nite size of the central cell leads to unphysical electrostatic interactions between defects. There are several methods to overcome the intrinsic diculties with periodic boundary conditions [23±27]. YM used the transition state theory of Slater [26], while Bunson et al. [16] used Janak's theorem [27]. Both transition state theory and Janak's theorem use local density eigenenergies calculated with some fraction of an electron in the highest energy unoccupied state to extrapolate to the total energy for the system. It is well known that local density theory leads to band gaps that are systematically too small [28], so that use of these energies is suspect. Furthermore the fractional electron requires neutralization as well. The other techniques used to improve on the simple jellium background use sets of ®ctitious charges in the central unit cell to cancel the multipole ®eld generated either by the neutral defect, or by the charged defect and the background charge [23±25]. One of these [24] has been used recently by Bl ochl [17]. As mentioned in Section 1, all of these techniques give the same qualitative prediction for the sign of Ueff . However, the various calculations gave very large variations (between )7 eV [16] and )2.6 eV [15]) in the calculated Ueff . There are major diculties in calculating charge states in a cluster formalism as well. In general, one must be concerned that the defect wave function is not arti®cially con®ned by the cluster size. This is less crucial for an insulator, such as SiO2 , than for a semiconductor, as the defects are typically more strongly localized.
45
However, in both cases, charged defects lead to long range polarization eects that are absent in ®nite clusters. This problem has been recognized for many years and Mott and Littleton [29] had introduced approximate adjustments that have been used by many other groups (see, for example, [30±32]). In this paper, we have modeled the dielectric eects using a dielectric continuum approximation. We have used a simple non-self-consistent model due to Jost [33], as well as two self-consistent techniques. In Jost theory, we model the cluster as a sphere of the same dielectric as the bulk with a point charge at its center. In this case, the polarization ®eld can be calculated exactly. Furthermore, we can calculate the contribution to the total energy of the system due to the interaction of the polarization of the rest of the solid with the point charge in the sphere. This energy is EJost
1
1 0
Q2 ; 2RJost
5
where RJost is the radius of the cavity, which should approximate the radius of the cluster, Q is the net charge (+1 or )1 in our case), and 0 is the static dielectric constant in SiO2 . This energy is always negative and is added to the total energies for the charged species. We obtained our value for the for dielectric correction by assuming RJost 3 A, the smaller cluster and RJost 4:2 A, for the larger cluster. We should emphasize that the greater value of this model is heuristic. It shows immediately that the polarization is proportional to the square of the charge on the defect, and that the convergence with sphere size and, hence, cluster size, is rather mild 1=RJost . Finally, it gives the dependence on the dielectric constant. We have also calculated the polarization corrections using two self-consistent continuum methods. In both cases, the cluster is immersed in the dielectric during the self-consistent calculation. The simpler method, called the self-consistent reaction ®eld method (SCRF), also uses a spherical enclosure for the cluster. In the more sophisticated method, called the polarizable continuum model (PCM) [34] the cavity is formed by taking the union of the set of spheres centered at each of the
46
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
atomic positions in the cluster with radii slightly larger than the Van der Waals radii. For both PCM and SCRF the method of apparent surface charge is used, where the potential everywhere can be determined from the charge density inside the cavity and the apparent surface charge, given by r
1 rUin n; 4po
6
where r is the surface charge density, Uin is the potential inside the cavity, n is the normal vector to the surface at the point where r is de®ned, and o is the relative dielectric constant of the material. Details of the method can be found in [34]. We are required to perform PCM and SCRF calculations in all three charge states, as the charge density in the neutral case has multipole terms that will induce charges on the cavity surface. Finally, the PCM and SCRF methods as implemented in GAMESS do not currently allow calculations at the MP2 level. However, Gaussian-94 has that feature. We should note, however, that the implementation of PCM in Gaussian-94 is signi®cantly dierent and can give dramatically dierent, even unphysical, results if not used correctly. First, the sphere sizes for each atom type are signi®cantly dierent and, in our case, signi®cantly smaller. Second, GAMESS recognizes when the union of spheres has an unphysical shape and introduces other spheres to smooth the cavity. To remedy the short-comings in the Gaussian-94 implementation of PCM, one can use the isodensity surface PCM
(IPCM) model. IPCM gives results within a few tenths of an eV of the GAMESS PCM results.
3. Results 3.1. Geometries In Table 1 we present the relevant geometrical parameters associated with Fig. 1. We see that the cluster results are in good agreement with each other, with the local density results of YM, and, in more detail, with those of Bl och. We can understand the geometry simply by noting that Si is isoelectronic with phosphorus, so that the ®ve fold co-ordination is analogous to, for example, PCl5 [35]. That is, there are three equatorial ligands and two axial ligands, indicative of sp2 hybridization with bonding to the lone pair orbital. Fig. 1 illustrates the terms axial and equatorial. We should point out that Bl ochl predicts that in a-quartz, the hydrogen atom prefers to be equatorial, while both Bl ochl and YM predict that it prefers to be axial in cristoballite. In the small cluster, we have calculated equilibria for both the equatorial and the axial con®guration and the axial con®guration is favored by 0.15 eV. The dierence is probably that our clusters have no constraints. It is interesting that the distortions induced by the hydrogen ion are remarkably local. We say this, even though we are using relatively small clusters, because the introduction of the hydrogen ion distorts only the
Table 1 Equilibrium geometrical parameters for H for the two clusters and for the two approximations employed in this study RSi1 ±H (A) RSi1 ±Oc (A) RSi2 ±Oc (A) RSi1 ±Ob 1 (A) RSi2 ±Ob 2 (A) \H±Si1 ±Ob (deg.) \Ob ±Si1 ±Oc (deg.) \Ob ±Si2 ±Oc (deg.) \Si1 ±Oc ±Si2 (deg.)
I (RHF)
I (MP2)
II (RHF)
Ref. [15]
Ref. [17]
1.56 1.77 1.60 1.70 1.65 90.0 90.0 111.7 130.6
1.55 1.79 1.63 1.73 1.67 90.5 89.5 110.0 128.2
1.56 1.77 1.59 1.70 1.65 89.4 90.6 113.0 147.5
1.54
1.51 1.76 1.694
As given in Fig. 1, Si1 is the silicon atom to which the extra hydrogen ion is attached and Oc is the central oxygen atom for the small cluster. The set of Ob ±1 2 are the three other oxygen atoms attached to Si1 2 , respectively. For Cluster I we show only the parameters for axial hydrogen.
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
47
nearest neighbors of Si1 . Other geometrical parameters are very close to those found in a-quartz. 3.2. Vibrational spectra and potential surfaces We show vibrational spectra in Table 2, along with the the value obtained by YM. We see that the ®nite element solution is in good agreement with the calculation of the full vibrational spectrum for the cluster, and is in reasonable agreement with the local density result, but that there is signi®cant anharmonicity. The infrared activity of 2 for the RHF this mode is 578 Debye2 = amu A 2 for the calculation and 487 Debye2 = amu A MP2 calculation, both of which are very strong. Thus, we agree with YM that this mode should be observable, as it is far from any other calculated modes for this system, and suciently far 100±300 cm 1 from any observed lines [36]. The HF and MP2 potential surfaces are shown in Fig. 2. We include this to show that over a very large portion of the potential surface, the qualitative shape of the surface is unaltered with the perturbative inclusion of correlation eects. We cannot draw very many other conclusions from such surfaces. For instance, it would be tempting to calculate a binding energy. In fact, we would estimate a minimum for the binding energy of 4.9 eV. However, in an oxide, it is highly unlikely that a hydrogen ion would ®nd a cavity with a 5 A radius. A more realistic maximum distance would Also, at large distances we would expect be 2±3 A. that the lattice will relax leading to a much smaller binding energy. We have, in fact, calculated an equilibrium geometry and total energy for Not surprisingly, the rest of the RSi±H 2:5 A. cluster does relax to nearly the neutral con®guraTable 2 Vibrational parameters for the Si±H bond-stretch vibration FE (HF) FS (HF) FE (MP2) FS (MP2) Ref. [15]
xe cm 1
xe xe
xe y e
1818 1870 1887 1929 1763
34.32
0.194
32.65
0.51
HF Hartree±Fock, MP2 second-order perturbative con®guration interaction, FE ®nite element, FS full spectrum.
Fig. 2. Si±H stretching mode potential surfaces obtained from ab initio electronic structure calculations. Only the hydrogen atom is allowed to move. MP2 surface has been shifted rigidly so that the energy at equilibrium coincides with that of the RHF surface.
tion. Using that as an assumed mid-point in a cavity, the estimated binding energy is, including zero-point energy, is 0.9 eV. We have performed MP2 optimizations and vibrational analyses on the small cluster. We include the vibrational parameters from the MP2 calculations in Table 2. It is interesting that while the qualitative shapes are the same, there are signi®cant shifts to higher wave numbers in the vibrational frequencies. 3.3. Electrical levels 3.3.1. Smaller cluster results We calculate the electrical levels in four approximations. On the small cluster, we calculate total energies in both the standard Hartree±Fock theory and using second-order perturbative con®guration interaction. We also compare the results with and without two separate estimates of dielectric corrections, and with previous DFT calculations. The values from [15,16] were taken immediately from the stabilization diagrams in these papers by assuming an experimental band gap of 9.0 eV. Those from [17] were deduced from the thermodynamic level for the (+/)) transition
48
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
Table 3 Electrical levels and Ueff for hydrogen in SiO2 in eV, relative to the conduction band edge. All results are for the smaller cluster E =0 Ec E =0 Ec E =0 Ec E =0 Ec E 0= Ec E 0= Ec E 0= Ec E 0= Ec Ueff Ueff (Jost) Ueff (SCRF) Ueff (PCM)
(Jost) (SCRF) (PCM) (Jost) (SCRF) (PCM)
RHF
MP2
Ref. [16]
Ref. [15]
Ref. [17]
)3.55 )1.7 )2.02 )2.1 )0.82 )2.62 )2.43 )2.32 2.73 )0.92 )0.41 )0.22
)3.77 )1.9
0.0
)4.4
)3.3
)1.96 )3.76
)7.0
)7.0
)4.4
1.81 )1.86
)7.0
)2.6
)1.1
and from the corrected value of Ueff . The results are shown in Table 3. The DFT results should be compared to the Jost-MP2 cluster results and with each other. As stated in Section 1, there is wide variability within the DFT results. In our opinion, the results of [17] are probably the most reliable because they depend only on calculated total energies, they have been corrected more carefully for spurious electrostatic interactions, and they include the GGA. The cluster results give a signi®cantly larger magnitude for Ueff than Bl ochl's result, although the origins of these dierences could be a series of small eects. The positions in the gap seem to indicate that the larger discrepancy is for the (+/0) level. We will comment further below. There are several features of the cluster results that bear notice. First, inclusion of electron correlation leads to a much deeper ()/0) electrical level but leads to insigni®cant changes in the (0/+) electrical level. Second, inclusion of correlation in the PCM model leads to insigni®cant changes in the polarization contribution. This is to be expected because polarization is a long range eect that is not strongly eected by local variations in the charge density. There was also little change in the polarization contribution when we tested an augmented basis set (DZP). Third, the agreement between the Jost estimation and the self-consistent calculations is very good for both charge states. This is not so apparent in Table 3 because, as mentioned above, there is a correction for the neutral charge state in the PCM theory that is
absent in the Jost theory. The energy dierences for the positive and negative charge states are given in Table 4. Note especially that the adjustment for the positive charge state is nearly equal to that for the negative charge state. This comment deserves discussion. Fourth, and most important, the inclusion of dielectric eects leads to qualitatively dierent conclusions about Ueff . Without the corrections, both approximations lead to a positive Ueff , while their inclusion leads to a negative Ueff . In all cases, Ueff is much smaller than the experimental value, 12 eV, for H in vacuum. In hindsight, this is a rather arti®cial comparison, as the two charged states in the solid involve strong bonding to the lattice. It is interesting, however, that this added bond energy would not be enough to destabilize the neutral charge state. Rather, the polarization energy is required. The Jost theory is based on completely classical electrostatics. It yields an electrostatic potential that is added to the total energies of both of the charged states. It is intrinsically symmetrical beTable 4 Changes in total energy for the three charge states and for both clusters for Jost theory, for SCRF theory and for PCM theory. All energies are in eV Approximation
+
)
Jost SCRF PCM (GAMESS-RHF) IPCM (Gaussian-RHF) IPCM(Gaussian-MP2)
)1.78 )1.79 )1.88
)1.78 )2.18 )1.92 )1.99 )1.94
0 0.0 0.0 0.495 )0.91
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
cause we have assumed the same value in both cases for RJost , and we have assumed the same charge distribution (a point charge at the center). The PCM calculation mixes a classical calculation of the surface charge with a quantum mechanical calculation for the total energy of the cluster. The eects of the surface charge are included self-consistently. Note, however, that the sign of the surface charge should be, and is, dierent for the two charge states and this has dramatic eects on the eigenenergies. For a positively charged defect, the polarization charge on the cavity surface is negative and repels electrons, localizing them arti®cially and driving the eigenenergies higher. One might expect, then, that the quantum mechanical interactions that would tend to raise the total energy might cancel the purely electrostatic interactions, leading to a much smaller eect than one would predict using Jost theory. For a negatively charged defect, the polarization charge would be positive and would attract electrons and lower the eigenenergies. So we would expect that the quantum and the purely electrostatic interactions would superpose, leading to a much larger eect than predicted by Jost theory. Inspection of the eigenvalues con®rms at least part of this argument. For the negative charge state, the eigenenergies are pushed down by 3.0 eV in the presence of the dielectric continuum, while for the positive charge state the eigenenergies are pushed up by 2.7 eV. However, the total energy is not simply the sum of the occupied one-electron eigenenergies. Rather, there are sums of Coulomb and exchange integrals that void such a simple argument. 3.3.2. Larger cluster results It is of interest to compare the smaller cluster results to those for a larger cluster. In a larger cluster, more of the polarization energy should be absorbed in the quantum cluster. Thus, we would expect that the RHF estimation for the (+/0) (0/ )1) level would be closer to (further from) the conduction band. However, if the PCM representation were perfect, the position of the embedded levels would be unchanged. In Table 5 we show the results for cluster II. The agreement between the larger and smaller clusters is much better than we should expect.
49
Table 5 Electrical levels and Ueff for hydrogen in SiO2 in eV, relative to the conduction band edge. All results are for the larger cluster RHF E 0= Ec E 0= Ec E 0= Ec E =0 Ec E =0 Ec E =0 Ec Ueff Ueff (Jost) Ueff (PCM)
(Jost) (PCM) (Jost) (PCM)
)3.28 )2.0 )2.2 )0.78 )2.08 )1.96 2.5 )0.08 )0.24
First, a larger cluster should exhibit greater relaxation energy because more atoms are allowed to relax. Second, it is highly unlikely that the quantum cluster would exhibit the same dielectric constant as bulk SiO2 . Nevertheless, the near invariability of the embedded electrical level with cluster size and the changes in the bare cluster electrical levels indicate that embedding gives a reasonable representation of the solid. Furthermore, the essential qualitative and quantitative agreement between the several methods of approximation for long range polarization eects and between cluster sizes lends strength to the our estimates for the electrical levels, and for Ueff .
4. Discussion We are in qualitative agreement with the conclusion from several DFT calculations for two isolated hydrogen atoms. That is, when long range polarization eects are included, we agree that two isolated H0 atoms should disproportionate into H and H . We reiterate that these cluster studies, wherein we can isolate the various energetic contributions, have led to a deeper insight into the origin of the negative Ueff . Local rebonding alone will not lead to disproportionation. Long range polarization will always contribute signi®cantly to this mechanism because it tends to lower the energy of charged species independent of the sign of the charge. We should point out that this has been discussed by Stoneham and Sangster [37,38] and by Stoneham and Ramos [39] in highly ionic systems where no stable Ueff defects have been identi®ed.
50
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
While there seems to be complete consensus that neutral hydrogen atoms are thermodynamically unstable, we must now ask about the practical implications of this result. We know, for example, that at low temperatures isolated hydrogen atoms are easily observed in a-quartz, fused quartz and in thermal oxides [21,22]. On warming above 130 K the signal rapidly disappears [21]. Based on the results shown here, one might be tempted to argue that they simply disproportionate. However, it has long been thought that hydrogen actually dimerizes. Griscom [52] has incorporated this ansatz to interpret successfully the annealing of NBOHC's in a-SiO2 . When we consider the possibility that two hydrogen atoms could dimerize, we ®nd that it is far more favorable energetically to form H2 (DE 4:35 eV Jost-MP2) than to form H and H DE 1:86 eV, Jost-MP2 estimate). Furthermore, we contend that there will be strong kinetic factors that will favor dimerization. Ueff is a thermodynamic quantity. There will be an activation energy to form H and H . This energy should be signi®cant because to form the charged species a hydrogen atom has to lose its electron to the conduction band of SiO2 . On the other hand, it is well known that H0 diuses rapidly in SiO2 indicating that the barrier to motion is very small. Thus dimerization should involve very small activation barriers making the reaction even more advantageous than disproportionation. Recently, Bl ochl [17] has considered the transformation of H0 to H . Using high quality DFT-GGA calculations, he concluded that the neutral hydrogen atom would readily give up its electron. However, his energies were referenced to a Fermi level at the Si-mid gap point. 1 Hydrogen atoms within a tunneling distance of the Si±SiO2 interface 2 would be able to release their electron. However, if the
1 Of course, this choice is appropriate whenever there is a large density of states at that energy level and within a tunneling distance of a hydrogen atom. These could be oxide defects or states associated with a silicon substrate. 2 In fact, by the same argument, hydrogen atoms near the interface could pick up an electron from the Si Fermi level to become negatively charged, so that disproportionation is a clear possibility.
distance from a hydrogen atom to the interface is large enough so that tunneling times are long, this reference level is not appropriate. For an estimate of the activation energy for ionization far from the interface, we can use Bl ochl's results from Fig. 4 of [17]. If we assume that the hydrogen atom must lose its electron to the SiO2 conduction band, then, by shifting the positive charge-state potential upward by 3.7 eV (half the Si band-gap plus the Si±SiO2 conduction band oset), we estimate a barrier of 1.7 eV. We should comment on the individual charge states. The proton has been studied fairly extensively and there is strong evidence that it exists in thermal oxides on silicon [17,41±45]. While the negative charge state has now received some theoretical attention, there is no compelling evidence for its existence. We predict that H0 should be a deep electron trap. That is to say, if there is atomic hydrogen in the network and an electron in the conduction band, the system can gain about 3.8 eV (Jost-MP2) by forming an H . We would like to comment on the general implications of these embedding results on the validity of previously published cluster calculations on point defects in SiO2 . As discussed in Section 1, several groups, including us, have published estimates of electrical levels for the E0 centers [11,46,47], the super-peroxide radical [48] and the NBOHC [13]. All of these published results ignored polarization eects. The clusters are small enough that the results are qualitatively incorrect. When these results are included, all electron trapping would be much deeper. Thus, we would expect the various E0 centers to be a fairly deep electron trap, though not as deep as an isolated hydrogen atom. Also, the NBOHC should be a very deep electron trap. On the other hand, the recent calculations of optical energies by Pacchioni and Ierano [49,50] should be little changed by the inclusion of polarization eects because the charge state does not change. We have seen that changes in the charge distribution for a given charge state have led to small changes in the polarization energy. Thus, we would expect that the excitation energies calculated for localized ground and excited states should change by a few tenths of an eV.
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
We have couched our discussion of these results in the assumptions that the defect wave function and the lattice distortions were naturally localized within the cluster. While in many cases the former is probably true to a good approximation, the latter has been demonstrated to be a poor approximation in highly ionic systems, such as MgO [51]. This must clearly be tested in SiO2 by developing more sophisticated, atomistic representations of polarization and especially of relaxation for a larger portion of the solid. Nevertheless, the calculations presented here show clearly and simply the crucial role of polarization. 5. Conclusions We have presented a set of cluster calculations to study the physics of H in SiO2 . We have found that when dielectric eects are included H0 is a relatively deep electron trap. Its electrical level is predicted to be 3.76 eV below the conduction band edge. Finally, we have considered the Ueff for hydrogen in SiO2 and, when dielectric eects are included, we agree with the DFT results that hydrogen should be a negative Ueff system. Our value for Ueff , )1.8 eV, is to be compared with a range of )7 to )1.1 eV from the DFT calculations. However, we predict that the spontaneous disproportionation of two neutral hydrogens will compete with the energetically more favorable dimerization. The rates for these two reactions will be determined largely by the activation energies for diusion of H0 and for formation of H and H . In equilibrium, however, our calculations predict that the vast majority of hydrogen atoms will dimerize, in agreement with previous models for radiation damage in SiO2 [40,52]. Finally, we have shown unequivocally, the crucial importance of including polarization in obtaining the correct qualitative physics. Acknowledgements It is a pleasure to thank Professors W.B. Fowler, A.M. Stoneham, and A. Shluger and Drs. P. Bl ochl, D.M. Fleetwood and J.-L. Leray for
51
stimulating discussions during the course of this work. A.H.E wishes to thank AFOSR for their support during the summer 1998 under the Summer Faculty Research Program. We acknowledge generous grants of computer time from the Albuquerque High Performance Computing Center and the Maui High Performance Computing Center.
References [1] K.G. Aubuchon, IEEE Trans. Nucl. Sci. NS-18 (1971) 117. [2] D.L. Griscom, D.B. Brown, N.S. Saks, in: C.R. Helms, B.E. Deal (Eds.), The Physics and Chemistry of SiO2 and the Si±SiO2 Interface, Plenum, New York, 1988, p. 287. [3] A.J. Lelis, H.E. Boesch, T.R. Oldham, F.B. MacLean, IEEE Trans. Nucl. Sci. NS-35 (1988) 1186. [4] R.K. Freitag, D.B. Brown, C.M. Dozier, IEEE Trans. Nucl. Sci., 1993. [5] D.L. Griscom, E.J. Friebele, G.H. Sigel Jr., Solid State Commun. 15 (1974) 479. [6] M. Stapelbroek, D.L. Griscom, E.J. Friebele, J.G.H. Sigel, J. Non-Cryst. Solids 32 (1979) 313. [7] D.L. Griscom, Phys. Rev. B 40 (1989) 4224. [8] D.M. Fleetwood, P.S. Winokur, R.A. Reber, T.L. Meisenheimer, J.R. Schwank, M.R. Shaneyfelt, L.C. Riewe, J. Appl. Phys. 73 (1993) 5058. [9] D.M. Fleetwood, P.S. Winokur, O. Flament, J.L. Leray, Appl. Phys. Lett 74 (1999) 2969. [10] R.E. Stahlbush, G.J. Campisi, J.B. McKitterick, W.P. Maszara, P. Roitman, G.A. Brown, IEEE Trans. Nucl. Sci. NS-39 (1992) 2086. [11] J.K. Rudra, W.B. Fowler, Phys. Rev. B 35 (1987) 8223. [12] M. Boero, A. Pasquarello, R. Car, Phys. Rev. Lett. 78 (1997) 887. [13] A.H. Edwards, W.B. Fowler, in: G.E. Walrafen, A.G. Revesz (Eds.), Structure and Bonding in Non-Crystalline Solids, Plenum, New York, 1986, p. 139. [14] E.H. Poindexter, Z. Phys. Chem. NF 151 (1987) 165. [15] A. Yokozawa, Y. Miyamoto, Phys. Rev. B 55 (1997) 13783. [16] P.E. Bunson, M.D. Ventra, S.T. Pantelides, R.D. Schrimf, K.F. Galloway, IEEE Trans. Nucl. Sci. 46 (1999) 1568. [17] P.E. Bl ochl, Phys. Rev. B (in press). [18] M.W. Schmidt, K.K. Baldridge, J.A. Boatz, S.T. Elbert, M.S. Gordon, J.H. Jensen, S. Koseki, N. Natsunaga, K.A. Nguyen, S.J. Su, T.L. Windus, M. Dupuis et al., J. Comput. Chem. 14 (1993) 1347. [19] M.J. Frisch, G.W. Trucks, H.B. Schlegel, P.M.W. Gill, B.G. Johnson, M.A. Robb, J.R. Cheeseman, T. Keith, G.A. Petersson, J.A. Montogometry, K. Raghavachari, M.A. Al-Laham, et al., Gaussian 94, Revision B.3, Gaussian Inc., Pittsburgh, PA, 1995. [20] C. Mùller, M.S. Plesset, Phys. Rev. 46 (1934) 618.
52
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
[21] K.L. Brower, P.M. Lenahan, P.V. Dressendorfer, Appl. Phys. Lett. 41 (1982) 251. [22] J.A. Weil, Phys. Chem. Miner. 10 (1984) 149. [23] P.A. Schultz, Phys. Rev. Lett. 84 (1999) 1942. [24] P.E. Bl ochl, J. Chem. Phys. 103 (1995) 7422. [25] G. Makov, M.C. Payne, Phys. Rev. B 51 (1995) 4014. [26] J.C. Slater, The Self-consistent Field for Molecules and Solids, McGraw-Hill, New York, 1974. [27] J.F. Janak, Phys. Rev. B 18 (1978) 7165. [28] M.S. Hybertsen, S.G. Louie, Phys. Rev. B 34 (1986) 5390. [29] N.F. Mott, M.J. Littleton, Trans. Faraday Soc. 34 (1938) 485. [30] A.B. Lidiard, J. Chem. Soc. Faraday Trans. 2 (85) (1989) 341. [31] M. Leslie, J. Chem. Soc. Faraday Trans. 2 (1989) 407. [32] A. Harker, J. Chem. Soc. Faraday Trans. 2 85 (1989) 471. [33] W. Jost, J. Chem. Phys. 1 (1933) 466. [34] J. Tomasi, M. Persico, Chem. Rev. 94 (1994) 2027. [35] P.C.V.D. Voorn, R.S. Drago, J. Chem. Phys. 88 (1966). [36] A.G. Dias, L. Guimaraes, M. Brunel, in: R.A.B. Devine (Ed.), The Physics and Technology of Amorphous SiO2 , Plenum, New York, 1987, p. 359. [37] M.J.L. Sangster, A.M. Stoneham, Philos. Mag. B 43 (1981) 597. [38] A.M. Stoneham, M.J.L. Sangster, Radiat. E. 73 (1983) 267. [39] A.M. Stoneham, M.M.D. Ramos, J. Solid State Chem. 2± 12 (1993).
[40] D.L. Griscom, J. Electron. Mater. 21 (1992) 763. [41] A.H. Edwards, Phys. Rev. Lett. 71 (1993) 3190. [42] A.M. Ferreira, S.P. Karna, C.P. Brothers, R.D. Pugh, B.B. Singaraju, K.Vanheusden, W.L. Warren, R.A.B. Devine, in: W.L. Warren, R.A.B. Devine, M. Matsumura, S. Cristoloveanu, Y. Homma, J. Kanicki (Eds.), Amorphous and Crystalline Insulating Thin Films, Materials Research Society, Pittsburgh, PA, 1997, p. 187. [43] N.S. Saks, R.B. Klein, D.L. Griscom, IEEE Trans. Nucl. Sci. 35 (1988) 1234. [44] K. Vanheusden, A. Stesmans, Appl. Phys. Lett. 64 (1994) 2575. [45] R.E. Stahlbush, R.K. Lawrence, H.L. Hughes, IEEE Trans. Nucl. Sci. 45 (1999) 2398. [46] W.B. Fowler, J.K. Rudra, M.E. Zvanut, F.J. Feigl, Phys. Rev. B 41 (1990) 8313. [47] S.P. Karna, H.A. Kurtz, W.M. Shedd, R.D. Pugh, B.K. Singaraju, IEEE Trans. Nucl. Sci. 46 (1999) 1544. [48] A.H. Edwards, W.B. Fowler, Phys. Rev. B 26 (1982) 6649. [49] G. Pacchioni, G. Ierano, Phys. Rev. Lett., 1997. [50] G. Pacchioni, G. Ierano, Phys. Rev. B 56 (1997) 1. [51] P.V. Sushko, A.L. Shluger, C.R.A. Catlow, Surf. Sci. 450 (2000) 153. [52] D.L. Griscom, D.B. Brown, N.S. Saks, in: C.R. Helms, B.E. Deal (Eds.), The Physics and Chemistry of SiO2 and the Si±SiO2 Interface, Plenum, New York, 1988, p. 287.
www.elsevier.com/locate/jnoncrysol
Theory of H in SiO2 Arthur H. Edwards *, W.M. Shedd, R.D. Pugh Air Force Research Laboratory, Space Vehicles Directorate, 3550 Aberdeen Ave. SE, KAFB, NM 87117-5776, USA Received 11 August 2000; received in revised form 18 May 2001
Abstract We present an ab initio molecular orbital study of H in SiO2 with special emphasis on H . We have calculated equilibrium geometries, vibrational spectra, binding energies and electrical levels. For the electrical levels, we have included long-range polarization eects in three approximations. We compare our results to those of Yokozawa and Miyamoto and others. We ®nd that H0 is unstable in SiO2 . However, we ®nd that, rather than disproportionate, it prefers to dimerize. We also predict that H0 is a deep electron trap. Finally, we ®nd that long-range polarization eects are crucially important for obtaining even qualitatively correct values for electrical levels. Ó 2001 Published by Elsevier Science B.V.
1. Introduction Radiation-induced charge traps in metal oxide semiconductor (MOS) insulators, especially, in SiO2 , have been studied intensely over the last three decades (see, for example, [1±4]). It is known that in a radiation environment, MOS oxides become positively charged, leading to reductions in n-channel and increases in p-channel threshold voltages. As a result of this experimental fact, and of the common wisdom that electrons, being much more mobile, are much less easily trapped, most of the literature has focused on the nature of hole traps in SiO2 . In fact, several hole traps have been unambiguously identi®ed in SiO2 including the E0 centers, the non-bridging oxygen center, and the self-trapped hole [5±7]. * Corresponding author. Tel.: +1-505 853 6042; fax: +1-505 846 4558. E-mail addresses: [email protected] (A.H. Edwards), [email protected] (W.M. Shedd), robert.pugh@ kirtland.af.mil (R.D. Pugh).
Using thermally stimulated current measurements, Fleetwood and co-workers [8] recently have shown that there is a very large density of radiation-induced electron traps in both gate and buried oxides. These observations have been con®rmed in conjunction with Leray and co-workers [9]. Also, Stahlbush et al. [10] have observed shallow electron traps in buried oxides using standard C±V measurements. The large density of traps, along with a large range of depths in energy has important conceptual and technological consequences. The conceptual importance arises, at least in part, because the measured positive charge is net charge. Thus, there is a much larger density of hole traps than was originally thought. Technologically, the charging and discharging character of these traps may be important if we are to understand longterm radiation eects in gate oxides and in buried oxides. While the presence of large numbers of electron traps has been con®rmed, the understanding of the microscopic nature of these traps is limited. To date, no intrinsic electron traps have been
0022-3093/01/$ - see front matter Ó 2001 Published by Elsevier Science B.V. PII: S 0 0 2 2 - 3 0 9 3 ( 0 1 ) 0 0 6 4 9 - 4
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
experimentally identi®ed. Several theoretical studies in the literature predict that the oxygen vacancy [11,12] and the non-bridging oxygen (NBO) center [13] can trap electrons. In the former case, the predicted binding energy was small (0.5 eV). In the latter case the predicted binding energy is 1.2 eV. However, c-induced NBOs have only been observed in bulk, fused silica. Only after a-irradiation has it been observed in MOS oxides [14]. Yokozawa and Miyamoto [15] (hereafter YM) were the ®rst to show that H was stable in bcristobalite. The bonding structure is similar to that in Fig. 1. In fact, based on calculations using
43
density functional theory (DFT) applied to a fairly small super cell, this group predicted that hydrogen would exhibit a negative Ueff . This energy measures the stability of two neutral hydrogen atoms relative to one H and one H . When Ueff < 0, the charged species are favored over the neutral species. This is discussed in more detail in Section 2. Since then, other groups have obtained the same prediction using DFT applied to super cells [16,17]. This result is so interesting and so surprising that we decided to study it using a completely dierent formalism. By doing this we have deepened our understanding of the origin of the disproportionation. Furthermore, while there are almost certainly several dierent electron traps, we believe that this defect may be important enough, because of the persistent presence of hydrogen in MOS oxides, to study in depth. We present a cluster-based, Hartree±Fock (HF) and post Hartree±Fock, molecular orbital study of this defect using the clusters shown in Fig. 1. We have calculated a variety of properties, including equilibrium geometries, Ueff , the vibrational spectrum, and the electrical level structure. In the study of electrical levels, we used several approximate techniques for treating long-range polarization eects. We ®nd that polarization effects are fundamentally important; they change qualitatively our predictions about the stability of H0 . 2. Theoretical methods We have used GAMESS [18] and Gaussian-94 [19], ab initio, Hartree±Fock molecular orbital packages capable of treating electron correlation in several approximations. 2.1. Geometries and vibrational spectra
Fig. 1. Clusters used in present study. Clusters in (a) and (b) are referred to as I and II in the text. In (a) Oc indicates the central oxygen. In (a) and (b) H indicates added hydrogen atom. In (b) Oa and Oe denote axial and equatorial oxygen atoms, respectively.
For both clusters shown in Fig. 1, we used the 6-31G all-electron atomic orbital basis set. Here the 6 indicates that all core atomic orbitals are represented by six Gaussian functions, the 31 indicates that the valence orbitals are represented by two atomic functions, one, expanded in three
44
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
Gaussian functions, with a smaller spatial extent, and one, expanded in only one Gaussian, with a larger spatial extent. The two asterisks indicate that we used d-polarization functions on the heavy atoms and p-polarization functions on the hydrogen atoms. We allowed all atoms to relax freely in the search for energy minima. For the smaller cluster, we have used perturbative many-body theory (Mùller±Plesset theory [20]) during the geometry optimizations to study the eects of electron correlation. These calculations are labeled MP2. The charged states were calculated by adding or removing an electron and allowing the system to adjust self-consistently. For the neutral case, we calculated equilibrium geometries without the extra hydrogen atom, and added the total energies of this relaxed, neutral cluster to the energy of a free hydrogen atom. We chose this method for two reasons. First, when we performed neutral calculations with the added hydrogen, the equilibrium geometry was asymptotic; the hydrogen atom merely drifted away from the cluster. This was a practical diculty because the gradient would be above the tolerance value even if the main portion of the cluster was essentially converged. The second reason is that YM found that in b-cristobalite, the equilibrium position of the H0 atom was at the center of a large axial cavity, and that the total energy of the system was, to a very good approximation, the sum of the perfect crystal energy and the energy of an isolated hydrogen atom. This is not surprising, since the 1 H hyper®ne value in SiO2 is approximately unperturbed from the gas-phase value, indicating that the wave function was little perturbed by the solid [21,22]. We have calculated the vibrational spectrum for H two ways. We have calculated a one-dimensional potential surface from the Hartree± Fock and MP2 total energies without inclusion of polarization eects. Using our own cubic ®nite element program, we have calculated the ®rst few (10) vibrational levels and used these to calculate the standard vibrational parameters, xe , xe xe , and xe ye , representing the quadratic, cubic and quartic terms in the Taylor series expansion of the potential surface, by least-squares ®tting to Eq. (1)
2 1 1 x e xe n En hc xe n 2 2 3 ! 1 xe ye n : 2
1
We have also used Gaussian-94 to calculate the full vibrational spectrum, along with the IR intensities, in the harmonic approximation. 2.2. Electrical levels and Ueff We have calculated electrical levels approximately by using Etot n 1=n Etot n 1
Etot n v;
2
where v is the electron anity of SiO2 , and Etot n is the total energy for the system with n excess electrons, in the appropriate geometrical conformation. We have used the experimental value, 0.9 eV, for v. The ®rst two terms on the right-hand side merely give the energy required to move an electron to in®nity from a defect with n 1 excess electrons. The addition of v in Eq. (2) moves the reference energy from the vacuum level to the bottom of the conduction band. Experimentally, this de®nition gives the position in the forbidden gap where the defect will undergo a change in occupation from n to n 1 electrons at T 0 K. We note that on band diagrams the electrical levels are labeled by the actual charge states, rather than by the number of excess electrons. So, for example, the electrical level for the transition from neutral to negative would be Etot =0 Etot
Etot 0 v:
3
Of particular interest in this study is the calculation of Ueff , de®ned by Ueff Etot Etot
2Etot 0:
4
Eq. (4) shows clearly that when Ueff is negative, it is energetically favorable for two neutral hydrogen atoms spontaneously to disproportionate into a positive and a negative ion. As YM and others [23±25] have pointed out, there are intrinsic problems in using DFT with periodic boundary conditions to calculate electri-
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
cal levels of defects. (Here we assume that the defect in the solid is represented by a super cell and that the defect wave function is localized enough to be negligibly small outside this super cell.) For charged states, the relevant lattice sums diverge unless the unit cell is neutral. It is typical to introduce a uniform background charge of opposite sign. While this does yield the required neutrality, it also leads to spurious lowering of total energy because of the interaction between the background charge and the added (or removed) electron in the band states. Furthermore, it leaves higher-order multipole ®elds that are unphysical. Even in neutral defects the ®nite size of the central cell leads to unphysical electrostatic interactions between defects. There are several methods to overcome the intrinsic diculties with periodic boundary conditions [23±27]. YM used the transition state theory of Slater [26], while Bunson et al. [16] used Janak's theorem [27]. Both transition state theory and Janak's theorem use local density eigenenergies calculated with some fraction of an electron in the highest energy unoccupied state to extrapolate to the total energy for the system. It is well known that local density theory leads to band gaps that are systematically too small [28], so that use of these energies is suspect. Furthermore the fractional electron requires neutralization as well. The other techniques used to improve on the simple jellium background use sets of ®ctitious charges in the central unit cell to cancel the multipole ®eld generated either by the neutral defect, or by the charged defect and the background charge [23±25]. One of these [24] has been used recently by Bl ochl [17]. As mentioned in Section 1, all of these techniques give the same qualitative prediction for the sign of Ueff . However, the various calculations gave very large variations (between )7 eV [16] and )2.6 eV [15]) in the calculated Ueff . There are major diculties in calculating charge states in a cluster formalism as well. In general, one must be concerned that the defect wave function is not arti®cially con®ned by the cluster size. This is less crucial for an insulator, such as SiO2 , than for a semiconductor, as the defects are typically more strongly localized.
45
However, in both cases, charged defects lead to long range polarization eects that are absent in ®nite clusters. This problem has been recognized for many years and Mott and Littleton [29] had introduced approximate adjustments that have been used by many other groups (see, for example, [30±32]). In this paper, we have modeled the dielectric eects using a dielectric continuum approximation. We have used a simple non-self-consistent model due to Jost [33], as well as two self-consistent techniques. In Jost theory, we model the cluster as a sphere of the same dielectric as the bulk with a point charge at its center. In this case, the polarization ®eld can be calculated exactly. Furthermore, we can calculate the contribution to the total energy of the system due to the interaction of the polarization of the rest of the solid with the point charge in the sphere. This energy is EJost
1
1 0
Q2 ; 2RJost
5
where RJost is the radius of the cavity, which should approximate the radius of the cluster, Q is the net charge (+1 or )1 in our case), and 0 is the static dielectric constant in SiO2 . This energy is always negative and is added to the total energies for the charged species. We obtained our value for the for dielectric correction by assuming RJost 3 A, the smaller cluster and RJost 4:2 A, for the larger cluster. We should emphasize that the greater value of this model is heuristic. It shows immediately that the polarization is proportional to the square of the charge on the defect, and that the convergence with sphere size and, hence, cluster size, is rather mild 1=RJost . Finally, it gives the dependence on the dielectric constant. We have also calculated the polarization corrections using two self-consistent continuum methods. In both cases, the cluster is immersed in the dielectric during the self-consistent calculation. The simpler method, called the self-consistent reaction ®eld method (SCRF), also uses a spherical enclosure for the cluster. In the more sophisticated method, called the polarizable continuum model (PCM) [34] the cavity is formed by taking the union of the set of spheres centered at each of the
46
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
atomic positions in the cluster with radii slightly larger than the Van der Waals radii. For both PCM and SCRF the method of apparent surface charge is used, where the potential everywhere can be determined from the charge density inside the cavity and the apparent surface charge, given by r
1 rUin n; 4po
6
where r is the surface charge density, Uin is the potential inside the cavity, n is the normal vector to the surface at the point where r is de®ned, and o is the relative dielectric constant of the material. Details of the method can be found in [34]. We are required to perform PCM and SCRF calculations in all three charge states, as the charge density in the neutral case has multipole terms that will induce charges on the cavity surface. Finally, the PCM and SCRF methods as implemented in GAMESS do not currently allow calculations at the MP2 level. However, Gaussian-94 has that feature. We should note, however, that the implementation of PCM in Gaussian-94 is signi®cantly dierent and can give dramatically dierent, even unphysical, results if not used correctly. First, the sphere sizes for each atom type are signi®cantly dierent and, in our case, signi®cantly smaller. Second, GAMESS recognizes when the union of spheres has an unphysical shape and introduces other spheres to smooth the cavity. To remedy the short-comings in the Gaussian-94 implementation of PCM, one can use the isodensity surface PCM
(IPCM) model. IPCM gives results within a few tenths of an eV of the GAMESS PCM results.
3. Results 3.1. Geometries In Table 1 we present the relevant geometrical parameters associated with Fig. 1. We see that the cluster results are in good agreement with each other, with the local density results of YM, and, in more detail, with those of Bl och. We can understand the geometry simply by noting that Si is isoelectronic with phosphorus, so that the ®ve fold co-ordination is analogous to, for example, PCl5 [35]. That is, there are three equatorial ligands and two axial ligands, indicative of sp2 hybridization with bonding to the lone pair orbital. Fig. 1 illustrates the terms axial and equatorial. We should point out that Bl ochl predicts that in a-quartz, the hydrogen atom prefers to be equatorial, while both Bl ochl and YM predict that it prefers to be axial in cristoballite. In the small cluster, we have calculated equilibria for both the equatorial and the axial con®guration and the axial con®guration is favored by 0.15 eV. The dierence is probably that our clusters have no constraints. It is interesting that the distortions induced by the hydrogen ion are remarkably local. We say this, even though we are using relatively small clusters, because the introduction of the hydrogen ion distorts only the
Table 1 Equilibrium geometrical parameters for H for the two clusters and for the two approximations employed in this study RSi1 ±H (A) RSi1 ±Oc (A) RSi2 ±Oc (A) RSi1 ±Ob 1 (A) RSi2 ±Ob 2 (A) \H±Si1 ±Ob (deg.) \Ob ±Si1 ±Oc (deg.) \Ob ±Si2 ±Oc (deg.) \Si1 ±Oc ±Si2 (deg.)
I (RHF)
I (MP2)
II (RHF)
Ref. [15]
Ref. [17]
1.56 1.77 1.60 1.70 1.65 90.0 90.0 111.7 130.6
1.55 1.79 1.63 1.73 1.67 90.5 89.5 110.0 128.2
1.56 1.77 1.59 1.70 1.65 89.4 90.6 113.0 147.5
1.54
1.51 1.76 1.694
As given in Fig. 1, Si1 is the silicon atom to which the extra hydrogen ion is attached and Oc is the central oxygen atom for the small cluster. The set of Ob ±1 2 are the three other oxygen atoms attached to Si1 2 , respectively. For Cluster I we show only the parameters for axial hydrogen.
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
47
nearest neighbors of Si1 . Other geometrical parameters are very close to those found in a-quartz. 3.2. Vibrational spectra and potential surfaces We show vibrational spectra in Table 2, along with the the value obtained by YM. We see that the ®nite element solution is in good agreement with the calculation of the full vibrational spectrum for the cluster, and is in reasonable agreement with the local density result, but that there is signi®cant anharmonicity. The infrared activity of 2 for the RHF this mode is 578 Debye2 = amu A 2 for the calculation and 487 Debye2 = amu A MP2 calculation, both of which are very strong. Thus, we agree with YM that this mode should be observable, as it is far from any other calculated modes for this system, and suciently far 100±300 cm 1 from any observed lines [36]. The HF and MP2 potential surfaces are shown in Fig. 2. We include this to show that over a very large portion of the potential surface, the qualitative shape of the surface is unaltered with the perturbative inclusion of correlation eects. We cannot draw very many other conclusions from such surfaces. For instance, it would be tempting to calculate a binding energy. In fact, we would estimate a minimum for the binding energy of 4.9 eV. However, in an oxide, it is highly unlikely that a hydrogen ion would ®nd a cavity with a 5 A radius. A more realistic maximum distance would Also, at large distances we would expect be 2±3 A. that the lattice will relax leading to a much smaller binding energy. We have, in fact, calculated an equilibrium geometry and total energy for Not surprisingly, the rest of the RSi±H 2:5 A. cluster does relax to nearly the neutral con®guraTable 2 Vibrational parameters for the Si±H bond-stretch vibration FE (HF) FS (HF) FE (MP2) FS (MP2) Ref. [15]
xe cm 1
xe xe
xe y e
1818 1870 1887 1929 1763
34.32
0.194
32.65
0.51
HF Hartree±Fock, MP2 second-order perturbative con®guration interaction, FE ®nite element, FS full spectrum.
Fig. 2. Si±H stretching mode potential surfaces obtained from ab initio electronic structure calculations. Only the hydrogen atom is allowed to move. MP2 surface has been shifted rigidly so that the energy at equilibrium coincides with that of the RHF surface.
tion. Using that as an assumed mid-point in a cavity, the estimated binding energy is, including zero-point energy, is 0.9 eV. We have performed MP2 optimizations and vibrational analyses on the small cluster. We include the vibrational parameters from the MP2 calculations in Table 2. It is interesting that while the qualitative shapes are the same, there are signi®cant shifts to higher wave numbers in the vibrational frequencies. 3.3. Electrical levels 3.3.1. Smaller cluster results We calculate the electrical levels in four approximations. On the small cluster, we calculate total energies in both the standard Hartree±Fock theory and using second-order perturbative con®guration interaction. We also compare the results with and without two separate estimates of dielectric corrections, and with previous DFT calculations. The values from [15,16] were taken immediately from the stabilization diagrams in these papers by assuming an experimental band gap of 9.0 eV. Those from [17] were deduced from the thermodynamic level for the (+/)) transition
48
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
Table 3 Electrical levels and Ueff for hydrogen in SiO2 in eV, relative to the conduction band edge. All results are for the smaller cluster E =0 Ec E =0 Ec E =0 Ec E =0 Ec E 0= Ec E 0= Ec E 0= Ec E 0= Ec Ueff Ueff (Jost) Ueff (SCRF) Ueff (PCM)
(Jost) (SCRF) (PCM) (Jost) (SCRF) (PCM)
RHF
MP2
Ref. [16]
Ref. [15]
Ref. [17]
)3.55 )1.7 )2.02 )2.1 )0.82 )2.62 )2.43 )2.32 2.73 )0.92 )0.41 )0.22
)3.77 )1.9
0.0
)4.4
)3.3
)1.96 )3.76
)7.0
)7.0
)4.4
1.81 )1.86
)7.0
)2.6
)1.1
and from the corrected value of Ueff . The results are shown in Table 3. The DFT results should be compared to the Jost-MP2 cluster results and with each other. As stated in Section 1, there is wide variability within the DFT results. In our opinion, the results of [17] are probably the most reliable because they depend only on calculated total energies, they have been corrected more carefully for spurious electrostatic interactions, and they include the GGA. The cluster results give a signi®cantly larger magnitude for Ueff than Bl ochl's result, although the origins of these dierences could be a series of small eects. The positions in the gap seem to indicate that the larger discrepancy is for the (+/0) level. We will comment further below. There are several features of the cluster results that bear notice. First, inclusion of electron correlation leads to a much deeper ()/0) electrical level but leads to insigni®cant changes in the (0/+) electrical level. Second, inclusion of correlation in the PCM model leads to insigni®cant changes in the polarization contribution. This is to be expected because polarization is a long range eect that is not strongly eected by local variations in the charge density. There was also little change in the polarization contribution when we tested an augmented basis set (DZP). Third, the agreement between the Jost estimation and the self-consistent calculations is very good for both charge states. This is not so apparent in Table 3 because, as mentioned above, there is a correction for the neutral charge state in the PCM theory that is
absent in the Jost theory. The energy dierences for the positive and negative charge states are given in Table 4. Note especially that the adjustment for the positive charge state is nearly equal to that for the negative charge state. This comment deserves discussion. Fourth, and most important, the inclusion of dielectric eects leads to qualitatively dierent conclusions about Ueff . Without the corrections, both approximations lead to a positive Ueff , while their inclusion leads to a negative Ueff . In all cases, Ueff is much smaller than the experimental value, 12 eV, for H in vacuum. In hindsight, this is a rather arti®cial comparison, as the two charged states in the solid involve strong bonding to the lattice. It is interesting, however, that this added bond energy would not be enough to destabilize the neutral charge state. Rather, the polarization energy is required. The Jost theory is based on completely classical electrostatics. It yields an electrostatic potential that is added to the total energies of both of the charged states. It is intrinsically symmetrical beTable 4 Changes in total energy for the three charge states and for both clusters for Jost theory, for SCRF theory and for PCM theory. All energies are in eV Approximation
+
)
Jost SCRF PCM (GAMESS-RHF) IPCM (Gaussian-RHF) IPCM(Gaussian-MP2)
)1.78 )1.79 )1.88
)1.78 )2.18 )1.92 )1.99 )1.94
0 0.0 0.0 0.495 )0.91
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
cause we have assumed the same value in both cases for RJost , and we have assumed the same charge distribution (a point charge at the center). The PCM calculation mixes a classical calculation of the surface charge with a quantum mechanical calculation for the total energy of the cluster. The eects of the surface charge are included self-consistently. Note, however, that the sign of the surface charge should be, and is, dierent for the two charge states and this has dramatic eects on the eigenenergies. For a positively charged defect, the polarization charge on the cavity surface is negative and repels electrons, localizing them arti®cially and driving the eigenenergies higher. One might expect, then, that the quantum mechanical interactions that would tend to raise the total energy might cancel the purely electrostatic interactions, leading to a much smaller eect than one would predict using Jost theory. For a negatively charged defect, the polarization charge would be positive and would attract electrons and lower the eigenenergies. So we would expect that the quantum and the purely electrostatic interactions would superpose, leading to a much larger eect than predicted by Jost theory. Inspection of the eigenvalues con®rms at least part of this argument. For the negative charge state, the eigenenergies are pushed down by 3.0 eV in the presence of the dielectric continuum, while for the positive charge state the eigenenergies are pushed up by 2.7 eV. However, the total energy is not simply the sum of the occupied one-electron eigenenergies. Rather, there are sums of Coulomb and exchange integrals that void such a simple argument. 3.3.2. Larger cluster results It is of interest to compare the smaller cluster results to those for a larger cluster. In a larger cluster, more of the polarization energy should be absorbed in the quantum cluster. Thus, we would expect that the RHF estimation for the (+/0) (0/ )1) level would be closer to (further from) the conduction band. However, if the PCM representation were perfect, the position of the embedded levels would be unchanged. In Table 5 we show the results for cluster II. The agreement between the larger and smaller clusters is much better than we should expect.
49
Table 5 Electrical levels and Ueff for hydrogen in SiO2 in eV, relative to the conduction band edge. All results are for the larger cluster RHF E 0= Ec E 0= Ec E 0= Ec E =0 Ec E =0 Ec E =0 Ec Ueff Ueff (Jost) Ueff (PCM)
(Jost) (PCM) (Jost) (PCM)
)3.28 )2.0 )2.2 )0.78 )2.08 )1.96 2.5 )0.08 )0.24
First, a larger cluster should exhibit greater relaxation energy because more atoms are allowed to relax. Second, it is highly unlikely that the quantum cluster would exhibit the same dielectric constant as bulk SiO2 . Nevertheless, the near invariability of the embedded electrical level with cluster size and the changes in the bare cluster electrical levels indicate that embedding gives a reasonable representation of the solid. Furthermore, the essential qualitative and quantitative agreement between the several methods of approximation for long range polarization eects and between cluster sizes lends strength to the our estimates for the electrical levels, and for Ueff .
4. Discussion We are in qualitative agreement with the conclusion from several DFT calculations for two isolated hydrogen atoms. That is, when long range polarization eects are included, we agree that two isolated H0 atoms should disproportionate into H and H . We reiterate that these cluster studies, wherein we can isolate the various energetic contributions, have led to a deeper insight into the origin of the negative Ueff . Local rebonding alone will not lead to disproportionation. Long range polarization will always contribute signi®cantly to this mechanism because it tends to lower the energy of charged species independent of the sign of the charge. We should point out that this has been discussed by Stoneham and Sangster [37,38] and by Stoneham and Ramos [39] in highly ionic systems where no stable Ueff defects have been identi®ed.
50
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
While there seems to be complete consensus that neutral hydrogen atoms are thermodynamically unstable, we must now ask about the practical implications of this result. We know, for example, that at low temperatures isolated hydrogen atoms are easily observed in a-quartz, fused quartz and in thermal oxides [21,22]. On warming above 130 K the signal rapidly disappears [21]. Based on the results shown here, one might be tempted to argue that they simply disproportionate. However, it has long been thought that hydrogen actually dimerizes. Griscom [52] has incorporated this ansatz to interpret successfully the annealing of NBOHC's in a-SiO2 . When we consider the possibility that two hydrogen atoms could dimerize, we ®nd that it is far more favorable energetically to form H2 (DE 4:35 eV Jost-MP2) than to form H and H DE 1:86 eV, Jost-MP2 estimate). Furthermore, we contend that there will be strong kinetic factors that will favor dimerization. Ueff is a thermodynamic quantity. There will be an activation energy to form H and H . This energy should be signi®cant because to form the charged species a hydrogen atom has to lose its electron to the conduction band of SiO2 . On the other hand, it is well known that H0 diuses rapidly in SiO2 indicating that the barrier to motion is very small. Thus dimerization should involve very small activation barriers making the reaction even more advantageous than disproportionation. Recently, Bl ochl [17] has considered the transformation of H0 to H . Using high quality DFT-GGA calculations, he concluded that the neutral hydrogen atom would readily give up its electron. However, his energies were referenced to a Fermi level at the Si-mid gap point. 1 Hydrogen atoms within a tunneling distance of the Si±SiO2 interface 2 would be able to release their electron. However, if the
1 Of course, this choice is appropriate whenever there is a large density of states at that energy level and within a tunneling distance of a hydrogen atom. These could be oxide defects or states associated with a silicon substrate. 2 In fact, by the same argument, hydrogen atoms near the interface could pick up an electron from the Si Fermi level to become negatively charged, so that disproportionation is a clear possibility.
distance from a hydrogen atom to the interface is large enough so that tunneling times are long, this reference level is not appropriate. For an estimate of the activation energy for ionization far from the interface, we can use Bl ochl's results from Fig. 4 of [17]. If we assume that the hydrogen atom must lose its electron to the SiO2 conduction band, then, by shifting the positive charge-state potential upward by 3.7 eV (half the Si band-gap plus the Si±SiO2 conduction band oset), we estimate a barrier of 1.7 eV. We should comment on the individual charge states. The proton has been studied fairly extensively and there is strong evidence that it exists in thermal oxides on silicon [17,41±45]. While the negative charge state has now received some theoretical attention, there is no compelling evidence for its existence. We predict that H0 should be a deep electron trap. That is to say, if there is atomic hydrogen in the network and an electron in the conduction band, the system can gain about 3.8 eV (Jost-MP2) by forming an H . We would like to comment on the general implications of these embedding results on the validity of previously published cluster calculations on point defects in SiO2 . As discussed in Section 1, several groups, including us, have published estimates of electrical levels for the E0 centers [11,46,47], the super-peroxide radical [48] and the NBOHC [13]. All of these published results ignored polarization eects. The clusters are small enough that the results are qualitatively incorrect. When these results are included, all electron trapping would be much deeper. Thus, we would expect the various E0 centers to be a fairly deep electron trap, though not as deep as an isolated hydrogen atom. Also, the NBOHC should be a very deep electron trap. On the other hand, the recent calculations of optical energies by Pacchioni and Ierano [49,50] should be little changed by the inclusion of polarization eects because the charge state does not change. We have seen that changes in the charge distribution for a given charge state have led to small changes in the polarization energy. Thus, we would expect that the excitation energies calculated for localized ground and excited states should change by a few tenths of an eV.
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
We have couched our discussion of these results in the assumptions that the defect wave function and the lattice distortions were naturally localized within the cluster. While in many cases the former is probably true to a good approximation, the latter has been demonstrated to be a poor approximation in highly ionic systems, such as MgO [51]. This must clearly be tested in SiO2 by developing more sophisticated, atomistic representations of polarization and especially of relaxation for a larger portion of the solid. Nevertheless, the calculations presented here show clearly and simply the crucial role of polarization. 5. Conclusions We have presented a set of cluster calculations to study the physics of H in SiO2 . We have found that when dielectric eects are included H0 is a relatively deep electron trap. Its electrical level is predicted to be 3.76 eV below the conduction band edge. Finally, we have considered the Ueff for hydrogen in SiO2 and, when dielectric eects are included, we agree with the DFT results that hydrogen should be a negative Ueff system. Our value for Ueff , )1.8 eV, is to be compared with a range of )7 to )1.1 eV from the DFT calculations. However, we predict that the spontaneous disproportionation of two neutral hydrogens will compete with the energetically more favorable dimerization. The rates for these two reactions will be determined largely by the activation energies for diusion of H0 and for formation of H and H . In equilibrium, however, our calculations predict that the vast majority of hydrogen atoms will dimerize, in agreement with previous models for radiation damage in SiO2 [40,52]. Finally, we have shown unequivocally, the crucial importance of including polarization in obtaining the correct qualitative physics. Acknowledgements It is a pleasure to thank Professors W.B. Fowler, A.M. Stoneham, and A. Shluger and Drs. P. Bl ochl, D.M. Fleetwood and J.-L. Leray for
51
stimulating discussions during the course of this work. A.H.E wishes to thank AFOSR for their support during the summer 1998 under the Summer Faculty Research Program. We acknowledge generous grants of computer time from the Albuquerque High Performance Computing Center and the Maui High Performance Computing Center.
References [1] K.G. Aubuchon, IEEE Trans. Nucl. Sci. NS-18 (1971) 117. [2] D.L. Griscom, D.B. Brown, N.S. Saks, in: C.R. Helms, B.E. Deal (Eds.), The Physics and Chemistry of SiO2 and the Si±SiO2 Interface, Plenum, New York, 1988, p. 287. [3] A.J. Lelis, H.E. Boesch, T.R. Oldham, F.B. MacLean, IEEE Trans. Nucl. Sci. NS-35 (1988) 1186. [4] R.K. Freitag, D.B. Brown, C.M. Dozier, IEEE Trans. Nucl. Sci., 1993. [5] D.L. Griscom, E.J. Friebele, G.H. Sigel Jr., Solid State Commun. 15 (1974) 479. [6] M. Stapelbroek, D.L. Griscom, E.J. Friebele, J.G.H. Sigel, J. Non-Cryst. Solids 32 (1979) 313. [7] D.L. Griscom, Phys. Rev. B 40 (1989) 4224. [8] D.M. Fleetwood, P.S. Winokur, R.A. Reber, T.L. Meisenheimer, J.R. Schwank, M.R. Shaneyfelt, L.C. Riewe, J. Appl. Phys. 73 (1993) 5058. [9] D.M. Fleetwood, P.S. Winokur, O. Flament, J.L. Leray, Appl. Phys. Lett 74 (1999) 2969. [10] R.E. Stahlbush, G.J. Campisi, J.B. McKitterick, W.P. Maszara, P. Roitman, G.A. Brown, IEEE Trans. Nucl. Sci. NS-39 (1992) 2086. [11] J.K. Rudra, W.B. Fowler, Phys. Rev. B 35 (1987) 8223. [12] M. Boero, A. Pasquarello, R. Car, Phys. Rev. Lett. 78 (1997) 887. [13] A.H. Edwards, W.B. Fowler, in: G.E. Walrafen, A.G. Revesz (Eds.), Structure and Bonding in Non-Crystalline Solids, Plenum, New York, 1986, p. 139. [14] E.H. Poindexter, Z. Phys. Chem. NF 151 (1987) 165. [15] A. Yokozawa, Y. Miyamoto, Phys. Rev. B 55 (1997) 13783. [16] P.E. Bunson, M.D. Ventra, S.T. Pantelides, R.D. Schrimf, K.F. Galloway, IEEE Trans. Nucl. Sci. 46 (1999) 1568. [17] P.E. Bl ochl, Phys. Rev. B (in press). [18] M.W. Schmidt, K.K. Baldridge, J.A. Boatz, S.T. Elbert, M.S. Gordon, J.H. Jensen, S. Koseki, N. Natsunaga, K.A. Nguyen, S.J. Su, T.L. Windus, M. Dupuis et al., J. Comput. Chem. 14 (1993) 1347. [19] M.J. Frisch, G.W. Trucks, H.B. Schlegel, P.M.W. Gill, B.G. Johnson, M.A. Robb, J.R. Cheeseman, T. Keith, G.A. Petersson, J.A. Montogometry, K. Raghavachari, M.A. Al-Laham, et al., Gaussian 94, Revision B.3, Gaussian Inc., Pittsburgh, PA, 1995. [20] C. Mùller, M.S. Plesset, Phys. Rev. 46 (1934) 618.
52
A.H. Edwards et al. / Journal of Non-Crystalline Solids 289 (2001) 42±52
[21] K.L. Brower, P.M. Lenahan, P.V. Dressendorfer, Appl. Phys. Lett. 41 (1982) 251. [22] J.A. Weil, Phys. Chem. Miner. 10 (1984) 149. [23] P.A. Schultz, Phys. Rev. Lett. 84 (1999) 1942. [24] P.E. Bl ochl, J. Chem. Phys. 103 (1995) 7422. [25] G. Makov, M.C. Payne, Phys. Rev. B 51 (1995) 4014. [26] J.C. Slater, The Self-consistent Field for Molecules and Solids, McGraw-Hill, New York, 1974. [27] J.F. Janak, Phys. Rev. B 18 (1978) 7165. [28] M.S. Hybertsen, S.G. Louie, Phys. Rev. B 34 (1986) 5390. [29] N.F. Mott, M.J. Littleton, Trans. Faraday Soc. 34 (1938) 485. [30] A.B. Lidiard, J. Chem. Soc. Faraday Trans. 2 (85) (1989) 341. [31] M. Leslie, J. Chem. Soc. Faraday Trans. 2 (1989) 407. [32] A. Harker, J. Chem. Soc. Faraday Trans. 2 85 (1989) 471. [33] W. Jost, J. Chem. Phys. 1 (1933) 466. [34] J. Tomasi, M. Persico, Chem. Rev. 94 (1994) 2027. [35] P.C.V.D. Voorn, R.S. Drago, J. Chem. Phys. 88 (1966). [36] A.G. Dias, L. Guimaraes, M. Brunel, in: R.A.B. Devine (Ed.), The Physics and Technology of Amorphous SiO2 , Plenum, New York, 1987, p. 359. [37] M.J.L. Sangster, A.M. Stoneham, Philos. Mag. B 43 (1981) 597. [38] A.M. Stoneham, M.J.L. Sangster, Radiat. E. 73 (1983) 267. [39] A.M. Stoneham, M.M.D. Ramos, J. Solid State Chem. 2± 12 (1993).
[40] D.L. Griscom, J. Electron. Mater. 21 (1992) 763. [41] A.H. Edwards, Phys. Rev. Lett. 71 (1993) 3190. [42] A.M. Ferreira, S.P. Karna, C.P. Brothers, R.D. Pugh, B.B. Singaraju, K.Vanheusden, W.L. Warren, R.A.B. Devine, in: W.L. Warren, R.A.B. Devine, M. Matsumura, S. Cristoloveanu, Y. Homma, J. Kanicki (Eds.), Amorphous and Crystalline Insulating Thin Films, Materials Research Society, Pittsburgh, PA, 1997, p. 187. [43] N.S. Saks, R.B. Klein, D.L. Griscom, IEEE Trans. Nucl. Sci. 35 (1988) 1234. [44] K. Vanheusden, A. Stesmans, Appl. Phys. Lett. 64 (1994) 2575. [45] R.E. Stahlbush, R.K. Lawrence, H.L. Hughes, IEEE Trans. Nucl. Sci. 45 (1999) 2398. [46] W.B. Fowler, J.K. Rudra, M.E. Zvanut, F.J. Feigl, Phys. Rev. B 41 (1990) 8313. [47] S.P. Karna, H.A. Kurtz, W.M. Shedd, R.D. Pugh, B.K. Singaraju, IEEE Trans. Nucl. Sci. 46 (1999) 1544. [48] A.H. Edwards, W.B. Fowler, Phys. Rev. B 26 (1982) 6649. [49] G. Pacchioni, G. Ierano, Phys. Rev. Lett., 1997. [50] G. Pacchioni, G. Ierano, Phys. Rev. B 56 (1997) 1. [51] P.V. Sushko, A.L. Shluger, C.R.A. Catlow, Surf. Sci. 450 (2000) 153. [52] D.L. Griscom, D.B. Brown, N.S. Saks, in: C.R. Helms, B.E. Deal (Eds.), The Physics and Chemistry of SiO2 and the Si±SiO2 Interface, Plenum, New York, 1988, p. 287.