Roger Haydock and Ronald L. Te

Department of Physics and Materials Science Institute

University of Oregon, Eugene OR 97403-1274, USA

**Abstract**

Recent observations of transitions in the non-linear polarizability of submonolayers of metal atoms on silicon, and between insulating and metallic behavior for electrons in silicon inversion layers, indicate the need for further study of the electronic states in such systems. In this work, the projected density of states and normalizations for individual states are used to characterize the localization of independent electrons in a model for a layer of metal atoms randomly deposited at low density on a crystalline surface of silicon. We find a localization transition, whose sharpness in energy is limited only by the sample size of over 80,000 sites for a model of the Si(001)2x1 surface. This is consistent with previous analytic and numerical calculations, by similar methods, of the localization properties of other models for electronic states in weakly disordered layers, and with the observed transition in the non-linear polarizability of submonolayers of metal atoms on silicon.

The nature of electronic states at disordered interfaces is a long standing problem. In the 1970s it was thought that the independent electronic spectrum of a disordered layer had mobility edges which separated energies where the states were exponentially localized from energies where the states were plane-wave-like. Many theoretical and numerical attempts were made to determine the location of these edges for the Anderson [1] and related models, but the mobility edges were elusive. During the 1980's, many opinions changed to the view that all states were exponentially localized and that the localization length varied smoothly with energy, becoming very large at energies where the states were previously thought to have been delocalized. Recent observations by Kravchenko and others [2] of an apparent metal-insulator transition in silicon inversion layers has raised doubts about this view.

A variety of theoretical and numerical approaches have been applied to the behavior of electrons at disordered interfaces. The scaling theory of Abrahams and others [3] persuaded many that the electronic states in such systems must all be exponentially localized, and this has been supported by field theoretic approaches[4] (for a review see Lee and Ramakrishna[5]) as well as the extensive numerical work of McKinnon, Kramer, and collaborators, see Ref. [6] for a review. On the other hand, Haydock[7] presented analytical arguments, and Godin and Haydock (GH) [8] produced numerical evidence that there is a sharp transition between exponential and power-law localized states in disordered two-dimensional systems. Starting from high magnetic fields where individual Landau levels are known to carry currents, Azbel [9] argued that there is a metal-insulator transition at zero field. In addition to the work on independent electrons, there have been efforts to take electron-electron interactions into account, and it has been suggested [10] that interactions can delocalize electrons in circumstances where independent electrons would be localized.

Most of the experimental information about electron localization at interfaces has come from silicon inversion layers, devices in which electrons are trapped by an electric field at the (001) interface between p-type silicon and a layer of silicon oxide. Uren and others [11] found two distinct kinds of localization, and Davies and Pepper [12] extended this work to show evidence for a transition between exponential and power-law localization. Kravchenko and others [2] saw evidence similar to that of Ref. [11] of two kinds of localization in low mobility devices, and found a much more dramatic transition in high mobility devices. While electron-electron interactions clearly play an important role in the Kravchenko devices, as long as the electronic excitations close to the Fermi level can be expanded in localized excitations, no matter how complicated, the electronic Hamiltonian expanded in these excitations is of the same form as the one used below.

In a different experimental approach Arakat, Kevan, and Richmond (AKR) [13] deposited alkali metal atoms on a 2x1 reconstructed Si(001) surface, and measured the second harmonic produced by an optical, in-plane electric field. The advantages of this approach are that the distribution of metal atoms, the disorder, can be observed directly, and that an electric field is applied directly to the electrons rather than a potential difference being applied between contacts to a tenuous electron gas. The disadvantage of this approach is that the measurements are at optical frequencies rather than DC. The second harmonic showed a sharp threshold at about 1/6th of a monolayer for each of three different alkali metals. Although the AKR suggested that this threshold is due to a metal-insulator transition induced by electron-electron interactions, we were intrigued to see if a localization transition could be excluded.

In the first of the following four Secs., we present an independent particle model for the electrons of metal atoms deposited on the 2x1 reconstructed surface of (001) silicon. Section 3 contains a description of the calculations for different coverages, of normalizations for states at the Fermi level, and of the projected densities of states at various atoms. In the subsequent Sec. the results of these calculation are compared with those obtained by GH using the Anderson Model to simulate an inversion layers complete with contacts. In the last Sec we discuss the relation between the results of these calculations and the metallization of overlayers.

Our model for the metallic electrons is similar to that of Debney[14] with a basis of a single s-orbital centered on each metal atom. We assume that at low coverages only the most binding site in each 2x1 surface cell is occupied. We take the energies of these orbitals to be independent of whether nearby surface cells are occupied, and allow the electrons to hop from one metal atom to another with matrix elements decreasing exponentially with the distance between cells. In our calculations, metal atoms occupy each surface cell randomly, independent of the occupation of any other cells.

The bulk terminated Si(001) surface has two dangling silicon bonds for each surface atom. For each surface silicon atom, one of the dangling bonds dimerizes with that of a neighboring atom to produce a 2x1 reconstruction which has two dangling bonds per surface cell . The remaining dangling bonds further hybridize and transfer electrons to produce one band of doubly occupied surface states, and one band of unoccupied surface states separated in energy by about 1 eV. See Ref. [13] for further discussion of the Si(001) surface.

The valence s-orbital of the metal adatom hybridizes most strongly with bands of surface states to produce a singly occupied orbital whose energy must consequently lie in the gap between the occupied and unoccupied surface bands of the silicon. The energy of the orbital is determined by the chemical bond formed between the metal and silicon, and this in turn fixes the exponent with which the orbital decays with distance outside the cell in which it is centered. It is this hybridized orbital of the isolated metal adatom which forms our basis for the electronic states of the metal overlayer, and because the spatial decay of the orbital outside its central cell is determined by the chemistry of the metal-silicon bond, we expect the different alkalis to behave similarly.

In the electronic Hamiltonian for this model, the energies of each orbital are the same for the reasons described in the previous paragraph. Outside the central cell, each orbital decays exponentially, and we take this to be independent of direction so that the matrix element of the Hamiltonian between metal orbitals in cells i and j takes the Hückel form:

, (1)

where *E**0* is the on-site
energy of the equivalent orbitals, is the separation between the ith and the
jth orbitals, and *a* describes the spatial extent of the orbital wave
function.

For the calculations which follow, we took
a two-dimensional grid of sites with a 2:1 ratio of its primitive
vectors corresponding to the 2x1 reconstruction of the Si(001)
surface. The samples consisted of square regions of over 200 of
the large primitive vectors by over 400 of the small primitive
vectors, containing over 80,000 adsorption sites which were
populated randomly for the different coverages studied. The
resulting Hamiltonians were represented as matrices of dimension
over 80,000 by 80,000, and matrix elements of the Hamiltonian of
less than 10-6*E0* were neglected.

Normalizations and projected densities of states were calculated for the model using the recursion method[15]. The electronic Hamiltonian was made tridiagonal by constructing a new basis starting with a single orbital, usually near the center of the model lattice. From a tridiagonal matrix it is easy, to generate eigenvectors and calculate normalization integrals, or to construct the electronic Green function whose imaginary part is the projected density of states.

Starting with an atomic orbital, say *u**0*, the
tridiagonalization proceeds by constructing *u**1*, *u**2*, *u**3*, ... *u**n*, ...
using the recurrence,

L (2)

where .

The new basis elements *u**1*, *u**2*, *u**3*, ... *u**n*, ...
consists of linear combinations of the metal orbitals. The second
element of the new basis, *u**1* , is a linear
combination of orbitals with non-zero coupling to *u**0*, *u**2* is a
linear combination of orbitals coupled to *u**1*, and so
forth. In exact arithmetic, the new basis would be orthonormal,
but in finite precision it is not. However, this loss of
orthonormality during the computation does not affect the
accuracy results[16].

An un-normalized eigenstate of the
Hamiltonian with energy *E* can then be expressed as a
linear combination of the tridiagonal basis {*un*}, i.e., , where the coefficients
{*P**n**(E)*} form a set of polynomials in *E*,
orthogonal with respect to the projected density of states for *u**0*, and
satisfying the same recurrence relation as the basis: .

Since the metal atoms have one electron per
orbital, the Fermi level of the adlayer must lie in the band of
states derived from these orbitals. For a calculation which
begins with the orbital *u**0*, we approximate
the Fermi level by requiring that this orbital contain a single
electron. Although this approximate Fermi level varies a little
with the choice of *u**0* [17], this ensures that the state at the Fermi
level has a large amplitude on *u**0* and so is
centered close to *u**0*.

From the expansion of the eigenstates in
the tridiagonal basis, we can relate the normalization of the
eigenstate at the Fermi level, , to the sum of the squares of the orthogonal
polynomials *P**n**(E**f** )* [15,18], i.e.,

. (3)

The numerical loss of orthogonality of the tridiagonal basis has a negligible effect on the normalization[16].

Since strongly localized states have a
normalization that is independent of the size of the system while
weakly localized states have normalizations which scale with the
size of the system, we expect Eq. 3 to give: (i) , for strongly or
exponentially localized states, i.e., the slope of the
normalization as a function of recursion level, N, goes to zero ,
and (ii) , for
weakly localized states, where *p* >0; a positive slope
for the normalization as a function of N. This serves as a
criterion for distinguishing states with different degrees of
localization.

Individual Hamiltonian configurations were
generated for each of the different coverages (*r* = 1, 1/4, 1/5,
1/6, 1/7, 1/8, 1/9, 1/10, 1/20, and 1/100). Taking a localized
atomic orbital near the center of the lattice as our starting
state, *u**0*, the normalization of the eigenstate at the Fermi
level (Eq. 3) was calculated. Results are shown in Fig. 1. (For
clarity, only the cases where *r* =1/8, 1/9, 1/10, and 1/20 are given here. See Ref.
[18] for more data.) Observe that for low adatom densities (*r* = 1/10), the
normalization of the wave function is consistent with strong electron
localization - its value remains constant as a function of
recursion level. For high adatom coverage (*r* >1/9), the
states show qualitatively weaker localization with having a positive slope
for all N. The transition occurs at a density of 1/9, and is
resolvable to as small a change as 1 part in 90 of the adsorbate
density. We have examined between five and ten different
realizations for each adsorbate coverage and all are consistent
with a sharp transition for the infinite system, to within the
fluctuations expected for these finite samples,.

The changes in the projected density of
states correspond to the above changes in the normalizations of
states. From Fig. 2, it is clear that for *r * 1/10, the projected densities of states are
composed of discrete bound states, whereas for *r* > 1/9,
projected densities of states become smooth near the Fermi energy
while states at the band edges remain discrete. This smoothing of
the projected densities of states indicates the formation of a
metallic subband corresponding to the onset of the adlayer's
metallization. The transition from strong to weak localization of
states at *E**f* with increasing adatom density suggests that the
Fermi level crossed a localization edge.

Note that the normalizations in Fig.1 and projected densities of states in Fig. 2 are typical values rather than averages in the sense that they were calculated for single realizations of the random system rather than averaged over many realizations. The typical values used in this work and in the work of Kramer and MacKinnon [6] avoid the problem that averages of many properties of disordered systems do not converge with the number of realizations sampled, see further Draeger and Bunde[19].

The critical coverage *r**c*, at which the states at the Fermi level become
weakly localized, depends mainly on *a*, the exponent with which the adatom orbitals
decrease at large distances, although there could be weak
dependence on the relationships between lattice parameters. Since
the coverage scales as an inverse area for a two-dimensional
system, and *a* scales as an inverse distance, we expect,

. (4)

This is consistent with data compiled from
our numerical samples and plotted in Fig. 3. Note that the line
corresponding to a constant of proportionality of 0.12 in Eq. 4
separates the strongly localized and weakly localized states for
all values of *a.*

The functional relationship between the
critical density and the atomic parameter (Eq. 4) can be used to
extract a value of *a* from experiment. Using a value of *a*=1 in our
calculations, gives a critical density of 1/9 for the
metallization transition (Fig. 1). From the second harmonic
experiments [13], the onset for the metallization, as indicated
by the dramatic increase in the in-plane second harmonic signal,
is around 1/6 monolayer (ML). Since there are two silicon atoms
per unit cell (7.68Å¥3.84Å), 1/6 ML is equivalent to a value of *r*c =1/3 in our model. This gives a value of *a*=v3 or 0.451Å-1.
The corresponding orbital radius of an adatom with this value of *a* is 2.217Å,
which compares well with the metallic radii of the alkali atoms
(e.g. Na=1.83Å, K=2.26Å, Cs=2.62Å).

The most surprising aspect of this calculation is how clearly the transition from weak to strong localization can be seen in the normalizations. In contrast to this, the localization edge is a relatively subtle feature in the results of the calculations by GH. The main differences between the two calculations are the models for disorder, and in the method used to distinguish strongly and weakly localized states.

As a check on the method Wolfram Arnold[20] kindly calculated the normalizations of states for the Anderson model on square lattices. The interesting result of these calculations is that it is impossible to see any transition in the Anderson model by this method, even for samples of over one million sites. In contrast to the results presented here, the normalizations for the Anderson model fluctuate enormously both with energy and from one realization to another.

This qualitative difference in the numerical behavior of the two models can perhaps be understood in terms of the range of the disordered part of the Hamiltonians. In the Anderson model the disorder is in the diagonal matrix elements of the Hamiltonian, a disorder in the energies of the individual electronic orbitals, while in this work the disorder is in the off-diagonal matrix elements of the Hamiltonian, and is topological in the sense that it is in the coupling between orbitals which can be several lattice spacings apart. The range of the topological disorder is distinct from the range which might be introduced into the Anderson model by correlating the energies of nearby orbitals leaving the way orbitals are connected to one another, the topology, unchanged. The off-diagonal matrix elements of the Hamiltonian determine the way the tridiagonal basis {un} spreads out over the layer, and the randomness in these seems to suppress fluctuations relative to randomness in the diagonal matrix elements.

In GH the energies of localization edges were determined from calculations of the energy dependence of the conductance between probes placed at opposite corners of the sample. The localization edges appears in these calculations as a cusp where the conductivity goes from a fluctuating but roughly constant value on the weakly localized side to an exponential decrease with energy on the strongly localized side. Calculating the conductance is much more difficult than calculating normalizations, but we suppose that such calculations for the current model would produce results similar to those of GH.

Another computational approach is that of Kramer and MacKinnon(KM)[6] which when applied to the two-dimensional Anderson model, shows exponential localization of the state at the center of the band for all disorders. In this method the increase in the localization length with width of long strips of sites is used to estimate the localization length of a layer. Although the Anderson model is different from the model used in this paper, it is hard to see how this difference could eliminate the transition. So, the disagreement between KM and GH remains a puzzle.

In their paper, AKR suggest that the transition they saw in the optical second harmonic was due to electron-electron repulsion rather than disorder. Their picture is that at low coverage the extra Coulombic energy of doubly occupied orbitals splits the metal-derived band in two: the lower Hubbard band being fully occupied with one electron per adatom and an antiferromagnetic ordering of spins, and the upper Hubbard band whose states span the second spin orbital of each adatom being empty. As the coverage increases so does the electron density and the kinetic energy of the electrons until it becomes comparable with their Coulomb repulsion so that the gap between the two bands goes to zero, and the adlayer becomes metallic.

Hubbard [21] showed that for ordered
systems, the upper and lower bands meet when the ratio of the
hopping bandwidth to the energy of repulsion between two electron
in the same orbital increases to about 1.15. The bandwidth due to
hopping is easily obtained from the calculations presented above,
and the energy of repulsion can be estimated for the metal
orbitals using their simple exponential form and dielectric
constant for silicon. The result is that the Hubbard criterion is
satisfied at a coverage of about 1/25 when *a* is unity, a much
lower coverage than 1/9 where the independent electron model
shows a transition for the same value of *a*.

The observed transition occurs at a
coverage of about 1/3 of the surface cells, so we must scale the
results of the calculation for *a* = 1. As explained above, the experimental result
corresponds to *a*= if the
effects of electron repulsion are completely neglected, and
taking repulsion to promote localization, serves as an upper bound on the value of *a*. The Coulomb
energy of repulsion scales with *a* and so is greater for *a*=. The
hopping matrix elements scale as e^{-}^{a}^{R}
according to Eq. 1, so the hopping bandwidth scales similarly,
and the average R is taken to be 1/. Consequently for *a*= the coverage at which the
upper and lower Hubbard bands meet should be about 1/7, still
much smaller than the coverage of 1/3 where the transition is
seen.

There are, however, two effects present in
the experiments which, we argue, lower the effective value of *a* from . The first of these is
the same electron repulsion as was discussed in the previous
paragraph. Since electron repulsion is surely present in the
experimental system and since it promotes localization, it pushes
the transition to a higher coverage than if the electrons were
truly independent, hence decreasing the value of *a* for which the
model is consistent with the experiment.

The second effect in the experimental
system is the tendency of alkali metal atoms to form chain
structures on the silicon surface[22] This decreases the Coulomb
repulsion between electrons by allowing them to spread over
several atoms, and it increases the effective disorder by
increasing the width of the distribution of distances between
metal atoms over that resulting from truly independent deposition
of each atom. Since some atoms are closer due to the formation of
chains, others must be further away to keep the same coverage.
This increase in the effective disorder again reduces the value
of *a*
estimated from experiment, and the reduction in repulsive energy
between electrons reduces the estimate of the coverage where the
upper and lower Hubbard bands meet and weakens any interaction
effects on the localization.

Another observation by AKR is that the metallization transition seems to be associated with ordering of the metal atoms into a 2x3 structure. It is difficult to think of forces which could produce such an ordering, other than those arising from the electronic states of the metal electrons. For the metal atoms to order over large regions, the forces must be long-ranged which implies that the electronic states must be more delocalized than at lower coverages. This is also consistent with a localization transition of the kind we propose. Speculation by AKR that the 2x3 structure is insulating seems to be inconsistent with our estimate that the Hubbard bands meet at much lower coverages.

Our interpretation of the results of AKR as a transition in the ground state of the metal electrons must be qualified by the observation that the experiments were carried out at optical frequencies and therefore could involve excitation of electrons out of the ground state, for example into a band of delocalized states an eV or so higher in energy. This possibility cannot be excluded without further experiments in which the frequency of the electric field is varied, however the variation of the phase shift in the second harmonic with coverage, reported by AKR, is inconsistent with resonant excitation to a higher band. As the coverage increases, the phase of the second harmonic advances, whereas it would be expected to lag if the transition were due to excitation of electrons to a band such as the upper Hubbard band whose energy decreases with increasing coverage.

The advancing phase of the second harmonic
itself poses a challenge to any theory of these experiments. If
we associate the increasing second harmonic intensity with
delocalized electronic states on the surface, then in the absence
of scattering, the second harmonic currents should lag the
applied field. It seems that the only way the phase could advance
with increasing coverage would be if the scattering rate of the
delocalized electrons were to also increase. This is possible
because the increase in the number of delocalized electronic
states both above and below the Fermi level would increase the
phase space available to scattering and hence the scattering
rate; the matrix elements between delocalized and localized
states being suppressed by their relatively small overlap. It
might even be possible to explain the advance in phase beyond 180^{o}
as evidence that at high coverages the second harmonic currents
in the surface saturate before the electric field reaches its
peak in each cycle.

In summary we have compared the experiments of AKR with a model of independent electrons hopping between metal atoms deposited randomly on an insulating surface. The model shows a localization transition at a coverage not dissimilar from the coverage at which a transition was seen experimentally whereas an estimate of the coverage for a transition due to interactions is much lower. Making allowances for interactions and other effects present in the experiment suggests that disorder may dominate the observed transition and that this is consistent with suggestions of ordering at the transition and possibly even the advance of the second harmonic phase.

This work was supported by the National Science Foundation Materials Theory Program under grant no. DMR- 9319246, and by the Office of Science and Technology Infrastructure under grant no. STI-9413532.

[1] P.W. Anderson, Phys. Rev.109, 1492(1958).

[2] S.V.Kravchenko, W.E. Mason, G.E. Bowker, J.E. Furneaux, V.M. Pudalov, and M. D'Iorio, Phys. Rev. B51, 7038(1995), and S.V.Kravchenko, G.V.Kravchenko, J.E. Furneaux, V.M. Pudalov, and M. D'Iorio, Phys. Rev. B50, 8039(1994).

[3] E. Abrahams, P.W. Anderson, D.C. Licciardello, and T.V. Ramakrishnan, Phys. Rev. Lett.42, 673(1979).

[4] D. Vollhardt and P. Woelfle, Phys. Rev. Lett.45, 482(1980); F. Wegner, Z. Phys. B35, 207(1979).

[5] P.A. Lee and T.V. Ramakrishnan, Rev. Mod. Phys.57, 287(1985).

[6] B. Kramer and A. MacKinnon, Rep. on Prog. in Phys.56, 1469(1993).

[7] R. Haydock, Philos. Mag. B43, 203(1981); and53,554(1986).

[8] T.J. Godin and R. Haydock, Phys. Rev. B,461528(1992).

[9] M.Ya. Azbel, Phys. Rev. B45, 4208(1992).

[10] Ph. Jacquod and D.L. Shepelyansky, Phys. Rev. Lett.78, 4986(1997).

[11] M.J. Uren, R.A. Davies, M. Kaveh, and M. Pepper, J. Phys. C14, 5737(1981).

[12] R.A. Davies and M. Pepper, J. Phys. C15, L371(1982).

[13] S. Arekat S., S.D. Kevan, and G.L. Richmond G. L., Europhys. Lett.22, 377(1993); and references therein.

[14] B.T. Debney. J. Phys. C10, 4719(1977).

[15] R. Haydock, inSolid State Physics,edited by H. Ehrenreich, F. Seitz and D. Turnbull, Vol. (Academic, New York, NY, 1980), vol. 35 p. 215.

[16] R. Haydock and R.L. Te, Phys. Rev. B49, 10845(1994).

[17] R. Haerle, R. Haydock, and R.L. Te, Comp. Phys. Commun.90, 81(1995).

[18] R.L. Te, Ph. D. Thesis, University of Oregon, Eugene (1994).

[19] J. Draeger and A. Bunde, Phys. Rev. E,54, 4597(1996).

[20] W. Arnold, private communication.

[21] J. Hubbard, Proc. Roy. Soc. A281, 401(1964).

[22] Y. Hasegawa, I. Kamiya, T. Sakurai, H. Tochihara, M. Kubota, and Y. Murata, Phys. Rev. B41, 9688(1990).

**Figures**

Figure 1. Normalization of
the wave functions at E*f* for different adlayer density.

Figure 2. Projected densities of states for various adatom coverages.

Figure 3. The dependence of localization on coverage and the atomic parameter, compared with Eq. 4.

ÿ