Excitation process simulation for atoms leaving a metal surface

Excitation process simulation for atoms leaving a metal surface

Vacuum 66 (2002) 149–155 Excitation process simulation for atoms leaving a metal surface V.G. Drobnich*, S.Yu. Medvedev, I.S. Sharodi Uzhgorod Nation...

127KB Sizes 0 Downloads 18 Views

Vacuum 66 (2002) 149–155

Excitation process simulation for atoms leaving a metal surface V.G. Drobnich*, S.Yu. Medvedev, I.S. Sharodi Uzhgorod National University, Voloshin str. 54, 88000 Uzhgorod, Ukraine Received 1 September 2001; accepted 4 February 2002

Abstract Recently, experimental data about the dependences of secondary atom excitation probability on atom velocity were obtained. In the present work an attempt is made to give a theoretical explanation. We used a model with the Anderson-type Hamiltonian, where all possible resonant and non-resonant single-electron processes, taking part in ‘‘metal–emitted multi-level particle’’ system, are taken into account. Calculations which are carried out in the frame of the proposed model are found to be in good quantitative agreement with the experimental data and explain facts such as: (1) the high efficiency of the electron-exchange mechanism in secondary atom excitation, (2) the relatively slow dependence of the probability of excitation of secondary atoms on their velocity, (3) the non-trivial characteristic type of such dependence for the secondary atom excitation and ionization probabilities. r 2002 Elsevier Science Ltd. All rights reserved.

1. Introduction During recent years, in connection with the development of quantitative optical ion beam surface diagnostics methods, different apparatus for investigation of the most informative differential characteristics of excited particle emission (EPE) were developed and reliable experimental data on the regularities and mechanisms of this phenomenon were obtained (see, for example, Refs. [1–4]). A subsequent necessary stage for application of these methods is connected with the development of the EPE theory, which would allow quantitative estimation of excited state formation probabilities for the particles which leave the surface. Secondary emission of excited atoms from a metal surface may be considered as *Corresponding author.

the top priority of the theoretical investigations, because it is well investigated experimentally, and the excited atom data offer a basic contribution to EPE [1]. The available experimental data on the integral and differential characteristics of EPE testify that secondary atom excited states are formed in an electron-exchange (EE) mechanism mode [1–4]. This mechanism is known to be responsible also for the emission of positive single charged ground state secondary ions. Recently, we developed the quantum-mechanical model of excitation and ionization of secondary atoms [5] in the EE framework. It is based on the more precise multielectron Anderson-type Hamiltonian than those used earlier, where all possible resonant and nonresonant single-electron processes, taking part in the ‘‘metal-emitted particle’’ system and also the Coulomb repulsion of atomic electrons are taken

0042-207X/02/$ - see front matter r 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 4 2 - 2 0 7 X ( 0 2 ) 0 0 1 7 7 - X

V.G. Drobnich et al. / Vacuum 66 (2002) 149–155

150

into account. The main purpose of this work is to test whether it is enough to consider the mentioned processes to explain non-trivial type of experimentally obtained dependences of secondary atom excitation and ionization probability on their velocity [3,4,6,7]. The model [5] is developed in two directions in this work: (a) the surface periodic potential is introduced in a Hamiltonian model; (b) the procedure for the numerical solution of basic model equations is outlined.

2. Periodic potential Let us designate the potential of aggregate electric field of surface ions as V ðrÞ; where r ¼ ðq; zÞ are the coordinates of the electron in the frame of reference related to one of the surface atoms (OZ is the surface normal). Our aim is to calculate the matrix element V01 ðtÞ ¼ /c0 ðtÞjV jc1 ðtÞS; where jc0 S and jc1 S are wavefunctions of the ground and excited states of the moving atom (in the indicated frame of reference). To simplify the task, we shall be dealing with a rectangular lattice; moreover, the crystallographic axes with periods dx and dy will be oriented along the coordinate axes x and y; respectively. As the surface ionic cores form a two-dimensional periodic structure, i.e. V ðq þ d; zÞ ¼ V ðq; zÞ; where d is the spatial period of the surface lattice, so V ðrÞ may be expanded as the Fourier series V ðrÞ ¼

X

Vg ðzÞ expðigqÞ;

ð1Þ

g

where g is the inverse lattice vector. Let us note that only the item with g ¼ 0 from series (1) were taken into consideration in the work [5], i.e. only the potential V0 ðzÞ was used instead of V ðrÞ: The form of Fourier components Vg ðzÞ is determined with the electric potential wðr  dn Þ; which is created by every nth surface ion. We approximate the potential wðr  dn Þ by means of a Thomas–Fermi–Firsov potential in a Molie" re approximation [8]. Using this approximation and turning to the reference system connected with the moving atom, we obtain

V01 ðtÞ ¼

3 X

Bj =A

j¼1

P



Cng ðza ðtÞÞ ¼

g

Cgj ðZa ðtÞÞ exp½bgj Za ðtÞ ; bgj exp½igvr t

Z

ð2Þ

d3 r /c1 jrS/rjc0 S

zþza Xx

ð3Þ exp½bng z  igr ; qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi where bng ¼ b2n þ g2 ; bn and Bn are known parameters of the Molie" re potential, A is the surface cell square of the crystal, za ðtÞ ¼ z0 þ vz t is the distance between the atom and the surface at the time moment t; x is the parameter determined from the condition of seaming the complete surface potential with the bottom of the metal conductivity zone. Expression (2) illustrates the substantial difference in the roles played by the normal, vz ; and parallel to the surface, vr ; velocity components of the moving atom velocity: vr determines the time perturbation frequency gvr ; and vz determines the rapidity of decrease of this periodical perturbation while the atom moves away from the surface. The main difficulty of computation of the matrix elements V10 ðtÞ was connected with computation of integral (3). In the case of Coulomb wave functions, we have obtained analytical expressions for Cng ðza Þ suiting various transitions in the hydrogen and hydrogen-like atoms. For multi-electron atoms, we have developed the procedure of numerical computation of Cng ðza Þ with the use of Hartree–Fock wave functions for multi-electron atoms. In particular, for the transitions np03n0 s Z za x Cng ðzaÞ ¼ 31=2 =2 dl l exp½bng l Fg ðlÞ; ð4Þ N

Fg ðlÞ ¼

Z

N

dr Rn0 s ðrÞRnp ðrÞJ0 ðgðr2 þ l2 Þ1=2 Þ;

ð5Þ

jlj

where Rn0 s ðrÞ and Rnp ðrÞ are radial wave functions of jn0 sS and jnp0S states, and J0 ðxÞ is a zero-order Bessel function. The testing of the software performance of the procedure was carried out by

V.G. Drobnich et al. / Vacuum 66 (2002) 149–155

collating the Cng ðza Þ dependences calculated according to the above-mentioned analytical expressions, and according to formulae (4) and (5) where Coulomb wave functions Rn0 s ðrÞ and Rnp ðrÞ were used. The finally obtained matrix elements V10 ðtÞ were used for numerical solution of the system of differential equations of our model.

3. Basic model equations According to Ref. [5] we chose the following Hamiltonian for the ‘‘metal–moving atom’’ electron system HðtÞ ¼ HM þ HA ðtÞ X Xaa ðtÞna ðtÞ þ ½X01 ðtÞCþ þ 0 ðtÞC1 ðtÞ þ h:c: a¼0;1

þ

X

½Yka ðtÞCþ k Ca ðtÞ þ h:c: þ Jn0 ðtÞn1 ðtÞ;

k;a

ð6Þ where HM is the Hamiltonian of the metal which does not interact with the atom; index k denotes the stationary single-electron states jkS of this system; HA ðtÞ is the Hamiltonian of a moving free atom; a ¼ 0 or 1 denotes the ground and excited states of this particle (|0S and |1S; respectively), þ Cþ a ðtÞ; Ca ðtÞ and Ck ; Ca are creation and annihilation operators for electrons of states jaðtÞS and jkS; respectively; na ¼ Cþ a Ca ; Xaa0 and Yka are the matrix elements of those operators which are responsible for the transitions between the jaS; ja0 S and jkS; jaS basic states, respectively: J is the Coulomb repulsion energy of electrons captured by the atom into the states j0S and j1S simultaneously. Differences of Hamiltonian (6) from known Hamiltonians, which were used in models of ion emission and emission of excited atoms, are considered in detail in work [5]. Here, we shall note only that Hamiltonian (6) allows to simulate evolution of an electron shell of secondary particle, as a result of set of resonant and non-resonant one-electron transitions jkS2jaðtÞS and jaðtÞS2ja0 ðtÞS: Elementary two-electron transitions are presented in (6) by one process. This is

151

Coulomb interaction of atomic electrons which is presented by the last item of Hamiltonian. The given interaction results in partial blocking of oneelectron transition jkS-jaS in a case when other atomic state is already occupied. Interaction of atomic electrons can be taken into account only approximately. For big Coulomb repulsion the approximation of full blocking is valid. Otherwise, they apply another calculation scheme [9–14]. Because of complexity, Hamiltonian (6) we will use the most simple of such schemes, namely timedependent Hartree–Fock approximation (TDHF) [9]. Use of this approximation allows to reasonably limit efficiency of jkS-jaS transitions and to obtain the correct order of value of population of the excited level of emitted atom [5]. In the given approximation n0 ðtÞn1 ðtÞ ¼ /n0 ðtÞSn1 ðtÞ þ n0 ðtÞ/n1 ðtÞS  /n0 ðtÞS/n1 ðtÞS; where /na ðtÞS is the value of the operator na averaged over the states of multi-electron ‘‘metalmoving atom’’ system at the moment of time t: The state of this N-electron system within the TDHF approximation is a Slater determinant composed of N single-electron states jCðtÞ >¼ jj1 ðtÞ; j2 ðtÞ ; y; jN ðtÞS; . where each state jjj ðtÞS satisfies the Schrodinger equations with Hamiltonian (6). Therefore the calculation of the N-electron state jCðtÞS is reduced to the solving of N single-electron Shr. odinger equations for jjj ðtÞS-states under the known initial conditions. Further, we will consider the typical case when free atom levels E0 and E1 are located below and above the Fermi-level EF of the metal. Besides this, being limited with the case of secondary emission of excited atoms (transition to the case of scattering does not cause difficulties [15]), we can formulate our problem and initial conditions in the following way. When t ¼ t0 ) N the atom being initially located far away from the surface in the ground state j0S begins to approach the metal surface with velocity v0 ) 0: The electron system of the metal is in the ground state too when t ¼ t0 ; i.e. only the states jkS with energies ek pEF are occupied. When t ¼ 0 the distance between atom and surface

V.G. Drobnich et al. / Vacuum 66 (2002) 149–155

152

is the shortest. When t > 0 the atom moves from the surface with velocity v: The mentioned electron processes form the electron state of a twolevel particle flying off the surface. In the observation act at t ) N this state is reduced to one of the following alternative states: a ground state (with probability P0 ), an excited atomic state (P1 ), a positive ion state (Pþ ) and negative ion state (P ). In accordance with such formulation of our problem and initial conditions, it is reasonable to indicate one of the functions jjj ðtÞS in the Slater determinant as jj0 ðtÞS and the rest of the N  1 functions as jjk ðtÞS; taking into consideration that jj0 ðt0 ÞS ¼ j0S and jjk ðt0 ÞS ¼ jkS (ek pEF ). Then all the mentioned probabilities can be expressed through the scalar products /ajj0 S and /ajjk S: In particular P1 ¼ /n1 S  /n0 n1 S; Pþ ¼ 1  /n0 S  /n1 S þ /n0 n1 S;

ð7Þ

where /na S ¼ j/ajj0 Sj2 þ

X

j/ajjk Sj2 ;

k

/n0 n1 S ¼ /n0 S/n1 S  j

X

fore, we designate X01 ¼ X10 X ), we have # dj0 ðsÞ=ds ¼ iAðsÞj 0 ðsÞ;

ð9aÞ

# djk ðsÞ=ds ¼ iAðsÞj k ðsÞ þ f k ðsÞ:

ð9bÞ

Here, the vector columns uj and f k are ! /0jjj S jj ¼ ; /1jjj S ! Y0k f k ¼ i expðiek sÞ Y1k # is and the matrix AðsÞ # ¼ AðsÞ

e0 ðsÞ  i D0 ðsÞ

X ðsÞ  i DðsÞ

X ðsÞ  i DðsÞ

e1 ðsÞ  i D1 ðsÞ

! ;

where ea ¼ Ea þ Xaa þ J/na S; Da play the role of half-width of the a-levels, D ¼ ðD0 D1 Þ1=2 : In the present work, the standard assumption that Da ðsÞ ¼ Da ð0Þ expðsÞ is used. Here, the dimensionless variable s is equal gvz t with typical values ( 1. All energy characterof parameter g ¼ 223 A istics are already divided by _gvz ; i.e. they are dimensionless.

/0jjj S/jj j1Sj2 :ð8Þ

j¼0;k

Thereby, the task is reduced to the determination of the values /ajj0 S and /ajjk S: Starting . from the Shrodinger equation with Hamiltonian (6) we can find that the dependence of these values on time is described with the system of integrodifferential equations. To our regret, because of the complexity of Hamiltonian (6) and correspondingly, complexity of the obtained equations, their direct resolving faces the significant calculation obstacles. Therefore, we will transform them into system of usual differential equations using the well-known broad band approximation [16]. It is necessary to note, that adequacy of this approximation in the case of simulation of ion emission from the surface of simple metals is demonstrated in numerous works (see for example [17–24]). Taking into account that in cases of practical interest the matrix elements Xaa0 are real (there-

4. The solution of the main equations Let us discuss the procedure of solution of Eqs. (9a) and (9b) for the practically important case of sputtering. In this case u0 ðsÞ 0 for sX0 [5]. Therefore, now we should compute uk ðsÞ from (9b) for sX0: In accordance with the theory of ordinary differential equations, the solution at s-N is as follows:  1 # jk ðsÞ ¼ UðsÞ U# ð0Þjk ð0Þ Z s 1 # dt U ðtÞf k ðtÞ ; ð10Þ þ 0

# is the fundamental solution matrix for where UðsÞ the homogenous Eq. (9a). Thus, the problem reduces to determination of uk ð0Þ and matrix # UðsÞ: At first, let us deal with the latter. The # analysis shows (see below) that to compute UðsÞ; it is advisable to solve not Eq. (9a) itself, but one

V.G. Drobnich et al. / Vacuum 66 (2002) 149–155

equivalent to it for the vectors zðsÞ ¼ cðsÞu0 ðsÞ; where Z s

e1 ðtÞ þ e0 ðtÞ  i DðtÞ : cðsÞ ¼ exp i dt 2 0 The equation for the vectors zðsÞ is as follows: ! il iw # zðsÞ; OðsÞ # ¼ dzðsÞ=ds ¼ OðsÞ ; iw il ð11Þ

153

vectors uk ðs-Þ according to formula (10). Before that, it would be necessary to determine the initial conditions, i.e. uk ð0Þ: This would not be difficult provided that the equations under review are solved at t0 oto0 (s0 oso0) and vz ) 0 [5]. It is clear that, in this case, dependence of the matrix O# elements upon s may be assigned conventionally, provided they acquire correct values at the point s ¼ 0: In particular, they may be considered constant. Then, according to the general theory # # ¼ exp½Oð0Þðs # 0 Þ: uðsÞ  s0 Þ uðs

ð13Þ

where l ¼ ½e1 ðsÞ  e0 ðsÞ =2; wðsÞ ¼ X ðsÞ  i DðsÞ: If # is the fundamental solution matrix of Eq. (11) uðsÞ # ¼ c1 ðsÞuðsÞ: # then UðsÞ From Eq. (11) one may draw the following conclusion that the determi# does not depend on s: It can be chosen to nant uðsÞ be det u# ¼ 1 and, when s ) N ! eils 0 # ¼ uðsÞ : ð12Þ 0 eils

The proper values P0 and P1 of matrix O# differ from one another at the point s ¼ 0; therefore, the exponent in formula (13) may be entered according to Silvester’s theorem [25]. As a result

# are As can be seen, the elements of matrix uðsÞ rapidly oscillating functions. The oscillation frequency is approximately the same for all elements, and it approximately equals l: By means of transformations, one may achieve a more ‘‘smooth’’ behaviour of one or two elements, although always at the cost of oscillation frequency of the others. Hence, Eq. (11) is the most suitable one for numerical computations as it is the only way to solve the basic model equations without the use of simplifying assumptions. One may use the asymptotic form (12) of the matrix # uðsÞ fixing a certain rather high s ¼ sass ; as the initial conditions to determine the solution in the region sosass : To check the correctness of the procedure of numerical computation, one may also take advantage of the analytical solution obtained above [5], and also the condition # ¼ 1: det uðsÞ

where em ð0Þ ¼ ðe1 ð0Þ þ e0 ð0ÞÞ=2: So we have vectors uk ð0Þ: The integration was carried out numerically within the interval ½0; sass and analytically within the region ½sass ; NÞ during the estimation of the uk ðNÞ vectors. Numerical integration was carried out simultaneously for all ek with the help of the fast Fourier transformation algorithm with multiple precision developed in this study [5]. By substituting uk ðNÞ into (7) and (8), we would compute /n0 S; /n1 S; /n0 n1 S and, further, P1 and Pþ : Summation against k was substituted with integration against ek from the bottom of the conductivity zone eb up to EF :

5. Solution of non-homogeneous Eq. (9b) # in the interval Upon computing the matrix uðsÞ ½0; sass with the chosen values of horizontal and vertical components of secondary atom velocity vr and vz ; one may proceed to the computation of

/ajjk ð0ÞS ¼ Yk ð0Þ 

ea ð0Þ  i Dð0Þ þ wð0Þ  ek ; l2 ð0Þ þ w2 ð0Þ  ½em ð0Þ  i Dð0Þ  ek 2

6. Results and discussion The calculations have been carried out for the experimentally studied case of the secondary emission of excited Al atoms from an Al surface. We treated the ground 3p and excited 4s states of the Al atom with Hartree–Fock wave functions j3p0i and j4si; respectively. In these calculations, we took into account secondary atom electron interaction with: (1) the surface atoms (Thomas–Fermi approximation),

V.G. Drobnich et al. / Vacuum 66 (2002) 149–155

154

(2) the constant potential Eb inside the metal, (3) the self-electric image in the metal, and (4) electric image of the ionic core of the secondary atom. The sum of potentials (1) and (3) that form the surface potential barrier, would seam with constant potential (2) (the value Eb ¼ 15:9 eV at EF ¼ 4:25 eV, determined in free electron approximation, was used here). For the energy of Coulumb repulsion of the electrons captured onto Hartree–Fock orbitals j3p0i and j4si; we used the value 3.79 eV found earlier in [5]. For the half-width Dð0Þ and ( 1, respectively, parameter g values 3 eV and 2.5 A are used. Chosen values of Dð0Þ and g are considered to be typical [1]. Let us note, that they closely correspond to results of precision calculation of the half-width of 3s-level of Na

atom flying off the Al surface [26]. In the present work, we have computed the oscillating part of the function X ðsÞ; i.e. the matrix element V01 from Eq. (2). Fig. 1 presents the calculation results of the probability of formation of the excited 4s-state of secondary Al atoms P4s ; as a function of vz at vr ¼ 0: To demonstrate the possibilities of the model, the probabilities Pþ of the formation of secondary Al ions are obtainable simultaneously with P4s. As seen from the figure, at a characteristic velocity vz ¼ 106 cm/s, the probability P4s ¼ 6:2  106 ; this being well-harmonized with the experimental data we possess [3,4]. From the figure, it is also clear that, in the region 0.6  106 cm/ sovz o1:7  106 cm/s, the graphs are well approximated by a linear dependence, i.e. P4s ðvz ; vr ¼ 0ÞBexpðC=vz Þ; where we obtain C ¼ 8:2 106 cm/s. As shown by the calculations, at

-1

10

-2

10

-3

10

-4

Probabilities P4s , P+

10

-5

10

-6

10

-7

10

2

-8

10

1 -9

10

-10

10

0.5

1.0

1.5

2.0

2.5

-6

1/vz , 10 s/cm Fig. 1. Calculated probabilities of excitation P4s and single positive ionization Pþ of the secondary Al atoms (curves 1 and 2, respectively). Dashed curve refers approximation of P4s ðvz Þ dependence with exp½C4s =vz function where the coefficient C4s ¼ 8:2  106 cm/s.

V.G. Drobnich et al. / Vacuum 66 (2002) 149–155

vz > 1:7  106 cm/s, the growth in P4s accelerates with the increase of vz while, at vz o0:6  106 cm/s, the slump in P4s slows down with the decrease in vz : As it could be seen from Fig. 1, the dependence of probability of ionization Pþ on vz in the region vz o1.7  106 cm/s is similar to P4s probability, but for high velocities these probabilities are anticorrelate. This fact is in agreement with experimental data on the P4s ðvz Þ [27] and Pþ ðvz Þ [6] functions. We obtained diagrams similar to those shown in Fig. 1, in our earlier work [5] where an appropriate # expression was used for uðsÞ: However, the earlier diagrams have significantly greater slope, i.e. they correspond with a greater dependence of P4s ðvz Þ upon vz : The plot shown in Fig. 1, matches with experimental data on the dependence of the probability of the excitation of Al atoms versus their velocity much better. As shown by our calculations, the slope depends substantially upon the choice of energy of the bottom of conductivity zone eb : with a decrease of eb ; both the slope and the values P4s also decrease. In particular, if eb ¼ 9:9 eV, then the slope would decrease down to 0 (P4s is independent of vz ). In this case, the value of excitation probability would become equal to 1.2  106. We have also computed the dependence of P4s upon vr at fixed vz ¼ 1:5  106 cm/s. It has appeared that with the increase of vr from 0 up to 107 cm/s, the excitation probability would decrease and fall to approximately half its initial value. These and another results obtained in this work allow us to make the following conclusion: the consideration of the processes mentioned above for the formation of excited and charge states of secondary atoms is quite enough to allow an explanation of such essential facts as: (a) the high efficiency of EE, (b) the relatively slow dependence of the probability of excitation of secondary atoms versus their velocity and (c) the non-trivial characteristic type of such dependence for the secondary atom excitation and ionization probabilities.

155

References [1] Pop SS, Belykh SF, Drobnich VG, Ferleger VKh. Ionnophotonnaja emissija metallov. Tashkent: FAN, 1989. 200pp [in Russian]. [2] Lediankin DV, Urazgildin IF, Yurasova VE. ZETP 1988;94:90. [3] Drobnich VG, Mastyugin VA, Pop SS. Nucl Instrum Meth B 1991;58:443. [4] Drobnich VG, Medvedev SYu, Nady CsB, Sharodi IS. Physico-chimicheskije,structurnyje i emissionnyje svojstva tonkikh plenok i poverkhnosti tverdogo tela. Kiev: UMK VO, 1992. p. 221 [in Russian]. [5] Drobnich VG, Medvedev SYU. Nucl Instrum Meth B 1993;78:148. " Makarenko BN, Shergin AP. Izv Akad Nauk [6] Popov AB, SSSR Ser Fiz 1994;58:47 [in Russian]. [7] Popov AB, Makarenko BN, Shergin AP. Nucl Instrum Meth B 1995;100:267. [8] Parilis ES, Turajev NYu, Umarov FF, Nizhnaya SL. Teoriya rasseyaniya atomov srednikh energij poverchnosty tverdogo tela. Tashkent: FAN, 1987. 210pp [in Russian]. [9] Grimley TB, Juothi Bhasu VC, Sebastian KL. Surf Sci 1983;124:305. [10] Sebastian KL. Phys Lett A 1983;98:39. [11] Sebastian KL. Phys Rev B 1985;31:6976. [12] Zanardi P, Calandra C. Surf Sci 1993;290:141. [13] Nordlander P, Shao H, Langreth DC. Nucl Instrum Meth B 1993;78:11. [14] Dunning FB, Nordlander P. Nucl Instrum Meth B 1995;100:245. [15] Medvedev SYu, Drobnich VG. Ukr Phys J 1995;40:966 [in Ukrainian]. [16] Blandin A, Nourtier A, Hone DW. J Phys 1976;37:369. [17] Brako R, Makoshi K, Newns DM. Proceedings of the IXth IVC-V ICSS, Madrid, Spain, 1983. p. 108. [18] Kato M, Mix W, Snowdon KJ. Surf Sci 1993;294:429. [19] Pan L, Wang Ya, Huang F, Fang D, Tang J, Yang F. Surf Sci 1994;311:281. [20] Shao H, Nordlander P, Langreth DC. Nucl Instrum Meth B 1995;100:260. [21] Klushin DV, Gusev MYu, Urazgil’din IF. Nucl Instrum Meth B 1995;100:316. [22] Mii T, Makoshi K. Surf Sci 1996;363:145. [23] Kato M, O’Connor DJ, Yamamoto K, Souda R. Surf Sci 1996;363:150. [24] Onufriev AV, Marston JB. Phys Rev B 1996;53:13340. [25] Korn G, Korn T. Spravochnik po matematike. Moscow: Nauka, 1984. 831pp [in Russian]. [26] Borisov AG, Teilet-Billy D, Gauyacq JP, Winter H, Dierkes G. Phys Rev B 1996;54:17166. [27] Drobnich VG. Izv Akad Nauk SSSR Ser Fiz 1998;62(7):1382–96 [in Russian].