invited papers\(\def\hfill{\hskip 5em}\def\hfil{\hskip 3em}\def\eqno#1{\hfil {#1}}\)

Journal logoJOURNAL OF
SYNCHROTRON
RADIATION
ISSN: 1600-5775

Path-integral approach to anharmonic vibration of solids and solid interfaces

CROSSMARK_Color_square_no_text.svg

aDepartment of Chemistry, Graduate School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan
*Correspondence e-mail: toshi@chem.s.u-tokyo.ac.jp

(Received 26 July 2000; accepted 28 November 2000)

The temperature dependence of EXAFS (extended X-ray absorption fine structure) cumulants was investigated for bulk Cu and a thin film of Cu by means of the path-integral effective classical potential method. By using the semi-empirical embedded-atom method as a potential, agreement between the experiments and calculations is found to be excellent for bulk Cu. In the thin Cu(111) film, anisotropic anharmonic vibration was clearly observed; the out-of-plane vibration is much more enhanced and more anharmonic than the lateral vibration. The results are semiquantitatively consistent with the previous experimental data on Cu(111)/graphite. Such a vibrational enhancement should be the driving force for the roughening transition and/or the surface pre-melting at higher temperature. The practical usefulness of the path-integral effective classical potential method combined with the embedded-atom method is demonstrated.

1. Introduction

EXAFS spectroscopy has been applied widely to the determination of local structures of crystalline and noncrystalline solids, amorphous materials and liquids, solutions, catalysts, biological and environmental materials, surface and interfaces, and so forth. EXAFS is commonly recognized to provide information on coordination numbers and interatomic distances. However, since EXAFS contains information on radial distribution functions around X-ray-absorbing atoms, one can obtain higher-order information on thermal fluctuation, including anharmonicity. Stimulated by interesting experimental observations by means of temperature-dependent EXAFS, existing theories have been exploited in recent years. The quantum-mechanical statistical perturbation theory gives analytical formulae of the EXAFS cumulants that are experimentally obtainable (Fujikawa & Miyanaga, 1993[Fujikawa, T. & Miyanaga, T. (1993). J. Phys. Soc. Jpn, 62, 4108-4122.]) and the perturbation formulae are useful to describe molecules and clusters (Yokoyama et al., 1996[Yokoyama, T., Yonamoto, Y., Ohta, T. & Ugawa, A. (1996). Phys. Rev. B, 54, 6921-6928.], 1997[Yokoyama, T., Ohta, T. & Sato, H. (1997). Phys. Rev. B, 55, 11320-11329.]; Yokoyama, 1999[Yokoyama, T. (1999). J. Synchrotron Rad. 6, 323-325.]). It is, however, rather complicated to apply the perturbation theory to large-dimensional systems, such as solids and solid surfaces, especially for higher-order cumulants that describe anharmonicity.

Another theoretical approach is the Feynman path-integral theory, which treats the motion of atoms in real space like a classical image for particles. The real-space representation should be more suitable to describe the radial distribution function, especially of solids observed by EXAFS, when compared with the perturbation theory, which is based on the eigenstate (or reciprocal) representation. Moreover, the path-integral theory can handle greater anharmonicity than the perturbation theory. Although the path-integral Monte Carlo method is established, the calculations usually require heavy computational loads. A more practical method is to employ the path-integral effective classical potential (PIECP) theory (Cuccoli et al., 1992[Cuccoli, A., Macchi, M., Neumann, M., Tognetti, V. & Vaia, R. (1992). Phys. Rev. B, 45, 2088-2096.], 1995[Cuccoli, A., Giachetti, R., Tognetti, V., Vaia, R. & Verrucchi, P. (1995). J. Phys. Condens. Matter, 7, 7891-7938.]; Kleinert, 1995[Kleinert, H. (1995). Path Integrals in Quantum Mechanics, Statistics and Polymer Physics. Singapore: World Scientific.]). The first application of the PIECP theory to EXAFS cumulants was conducted for a two-body potential (Fujikawa et al., 1997[Fujikawa, T., Miyanaga, T. & Suzuki, T. (1997). J. Phys. Soc. Jpn, 66, 2897-2906.]) and comparisons with experimental data of a diatomic molecule of Br2 and of solids such as face-centred cubic (f.c.c.) Kr and Ni, were subsequently performed (Yokoyama, 1998[Yokoyama, T. (1998). Phys. Rev. B, 57, 3423-3432.]).

In the former part of this article, a brief review of the PIECP theory is at first presented, including a practical usage of the potential function. The PIECP theory is based on the two-body potential, which is sufficiently accurate to describe rare-gas solids but is insufficient for metals or semiconductors. It is, however, found that one can apply the embedded-atom method (EAM) (Foiles, 1985[Foiles, S. M. (1985). Phys. Rev. B, 32, 3409-3415, 7685-7693.]) to the PIECP theory without modification of the formulae. The results for bulk Cu calculated using the PIECP theory combined with EAM are subsequently given.

In the latter part of this paper, a surface anisotropic and anharmonic vibration is discussed. An enhancement of vibrational amplitudes and anharmonicity of surface atoms is believed to be the driving force of surface roughening transitions, subsequent surface pre-melting and consequent bulk melting, the transition mechanisms of which have been a long-standing problem from a microscopic point of view. EXAFS provides valuable information on local disorder which cannot be obtained by diffraction techniques. In order to observe a surface anisotropic and anharmonic vibration, we have investigated thin Ni(111) and Cu(111) films grown epitaxically on graphite and have found that the out-of-plane surface vibration is significantly enhanced and becomes more anharmonic compared with the in-plane one (Kiguchi et al., 2000[Kiguchi, M., Yokoyama, T., Matsumura, D., Endo, O., Kondoh, H. & Ohta, T. (2000). Phys. Rev. B, 61, 14020-14027.]). In this work, the PIECP calculations of thin Cu(111) films are compared with the experimental observations to demonstrate the practical usefulness of the PIECP theory for temperature-dependent EXAFS studies of real systems.

2. Theory

2.1. Path integral effective classical potential method

In this section, the PIECP theory is reviewed briefly. Details of the theory can be found in the literature (Cuccoli et al., 1992[Cuccoli, A., Macchi, M., Neumann, M., Tognetti, V. & Vaia, R. (1992). Phys. Rev. B, 45, 2088-2096.], 1995[Cuccoli, A., Giachetti, R., Tognetti, V., Vaia, R. & Verrucchi, P. (1995). J. Phys. Condens. Matter, 7, 7891-7938.]; Kleinert, 1995[Kleinert, H. (1995). Path Integrals in Quantum Mechanics, Statistics and Polymer Physics. Singapore: World Scientific.]). According to the Feynman theory, the density function ρ(X) (X is the 3N-dimensional Cartesian coordinate) is expressed as a functional form:

[\eqalignno{\rho ({\bf X}) = \hskip.2em& (1/Z) \langle {\bf X} | \exp(-\beta {\scr H}) | {\bf X} \rangle \cr = \hskip.2em& (1/Z) \textstyle\int\limits _{({\bf X},0) \Rightarrow ({\bf X},\,\hbar\, \beta)} {\scr D}{\bf X} \exp\{-{\scr A}[{\bf X}(u)]/\,\hbar\}\, ,&\hfill\llap{(1)}}]

where β = (kBT)−1, [\hbar] is the Planck constant divided by 2π, [{\scr H}] is the Hamiltonian of the system, Z is the partition function and [{\scr A}][X(u)] is the Euclidean action,

[{\scr A}[{\bf X}(u)] = \textstyle\int\limits _{0}^{\hbar \,\beta} \,{\rm d}u \, \left [{\textstyle {1 \over 2}} \,{}^t {\dot{\bf X}}(u){\bf M}\dot{\bf X}(u) + V({\bf X})\right ] \, . \eqno(2)]

In the PIECP approximation, the Euclidean action [{\scr A}] is solved variationally. For a trial function of [{\scr A}][X(u)], a solution of a harmonic oscillator [{\scr A}_{0}] should be the best candidate since thermal properties of solids are being investigated:

[{\scr A}_{0}[{\bf X}(u)] = \textstyle\int\limits _{0}^{\hbar\, \beta} \,{\rm d}u \, \left [{\textstyle {1 \over 2}} \,{}^t{\dot{\bf X}}{\bf M}\dot{\bf X} + w(\bar{\bf X}) + {\textstyle {1 \over 2}} \,{}^t{({\bf X}-\bar{\bf X})} {\bf F} ({\bf X}-\bar{\bf X}) \right] \, ,\eqno(3) ]

where [\bar{\bf X}] is the average path,

[\bar{\bf X} = (1/ \,\hbar\, \beta) \textstyle\int\limits _{0}^{\hbar \,\beta} \,{\rm d}u\, {\bf X}(u) \, , \eqno(4)]

and the force constant matrix F and the scalar potential w are variational parameters.

The Cartesian coordinate X is transformed into a normal coordinate Q through Q = tUM1/2(X[\bar{\bf X}]), where U is the eigenvector matrix of M−1/2FM1/2. One can thus rewrite

[{\scr A}_{0}[{\bf X}(u)] = \textstyle\int\limits _{0}^{\hbar \,\beta} \,{\rm d}u \,\left [{\textstyle {1 \over 2}} \,{}^t{\dot{\bf Q}}\dot{\bf Q} + {\textstyle {1 \over 2}} \,{}^t{{\bf Q}}{\boldomega}^{2} {\bf Q} + w(\bar{\bf X}) \right] \, , \eqno(5) ]

where [\boldomega^{2}] is the diagonal eigenvalue matrix. The density matrix ρ0 of the harmonic oscillator is obtained as

[\eqalignno{\rho _{0} (\bar{\bf X}) = \hskip.2em& \exp[ -\beta w ({\bf X})] /{\rm det}{\bf M}^{-1/2} \cr &\!\times\textstyle\prod\limits_{\bf k} ( 2\pi \,\hbar\, ^{2} \beta)^{-1/2} \,(f_{\bf k}/\sinh f_{\bf k})\, (2\pi \alpha _{\bf k})^{-1/2}\cr &\!\times \textstyle\int\limits _{-\infty} ^{\infty} \,{\rm d}Q_{\bf k}\, \exp[- (Q_{\bf k}-\bar{Q}_{\bf k})^{2} / \alpha _{\bf k}] \, , &\hfill\llap{(6)}}]

where

[\alpha _{\bf k} = (\,\hbar\, / 2\omega _{\bf k}) (\coth f_{\bf k} - 1/f_{\bf k}) \; \; {\rm and} \; \; f_{\bf k} = \beta \,\hbar \,\omega _{\bf k}/2 \, . \eqno(7)]

αk is the difference between the quantum-mechanical and classical fluctuations of the phonon mode k. The expected value of a certain physical quantity [ \langle{\scr O} \rangle_{0}] is given as

[\eqalignno{\langle {\scr O} \rangle _{0} =\hskip.2em& (1/Z_{0})\,(1/{\rm det}{\bf M}^{-1/2})\, [1/(2 \pi \,\hbar{}^{\, 2} \beta )^{3N/2}] \cr &\!\times \textstyle\int\limits \,{\rm d} \bar{\bf X}\, \exp[-\beta V_{\rm eff}(\bar{\bf X}) ]\, \langle \langle {\scr O} (\bar {\bf X} + {\bf M} ^{1/2} {\bf UQ} ) \rangle \rangle \, , &\hfill\llap{(8)}}]

where [V_{\rm eff}(\bar{\bf X})] is the effective classical potential,

[V_{\rm eff}(\bar{\bf X}) = w(\bar{\bf X}) + (1/\beta)\textstyle\sum\limits _{\bf k} \ln [(\sinh f_{\bf k}) / f_{\bf k}] \, .\eqno(9)]

In order to optimize the variational parameters of w and F, the Jensen–Feynman inequality can be used:

[F \leq F_{0} + (1/ \,\hbar\, \beta) \langle {\scr A} - {\scr A}_{0} \rangle _{0} \, , \eqno(10)]

where F and F0 are respectively the true and trial free energies. The optimized results are as follows:

[\langle \langle V (\bar{\bf X} + {\bf M} ^{1/2} {\bf UQ} ) \rangle \rangle = w(\bar{\bf X}) + {\textstyle {1 \over 2}} \textstyle\sum\limits _{\bf k} \omega _{\bf k} ^{2} (\bar{\bf X})\, \alpha _{\bf k} (\bar{\bf X}) \eqno(11) ]

and

[\langle \langle \nabla \nabla V (\bar{\bf X} + {\bf M} ^{1/2} {\bf UQ} ) \rangle \rangle = {\bf F} \, , \eqno(12) ]

where [\nabla \nabla V] denotes the dyadic, the ij component of which is defined as

[[\nabla \nabla V({\bf X})]_{ij} = \partial \, ^{2}\,V({\bf X}) / \partial X_{i} \partial X_{j} \, . \eqno(13)]

Using equation (8)[link], the EXAFS cumulants can thus be evaluated.

2.2. Low-coupling approximation in two-body potential solids

The 3N-dimensional integrals in equations (6)[link] and (8)[link] cannot, however, be numerically calculated in the case of a large-dimensional system. We will here employ a further approximation: the so-called low-coupling approximation. This assumes that w and F are independent of [\bar{\bf X}]. Moreover, we will for simplicity confine ourselves in a monatomic Bravais lattice with an atomic mass m in which only two-body forces are active. The dynamical matrix D is defined as

[{\bf D} = \textstyle\sum\limits _{j} {\bf F} _{oj} \exp (i {\bf k} \cdot {\bf R}_{oj} ) \, ,\eqno(14)]

where Foj is the 3 × 3 matrix of F for atoms o and j, and Roj is the location vector of atom j with respect to atom o. When the eigenvalues and eigenvectors of the matrix D are written as [m \omega _{{\bf k} \mu}^{2}] and [{\bf e}_{{\bf k} \mu}^{2}] (μ = 1, 2 and 3 are the phonon branches), the effective classical potential [V_{\rm eff} (\bar{\bf X})] is obtained as

[\eqalignno{V_{\rm eff}(\bar{\bf X}) = \hskip.2em& V(\bar{\bf X}) + \textstyle\sum\limits _{i \neq j} \, [u ^{\prime \prime} (R_{ij}) - u ^{\prime \prime} (R_{ij}^{0}) ] \sigma _{ij}^{(2)L} \cr &\! + \textstyle\sum\limits _{i \neq j} \, [ u ^{\prime} (R_{ij}) / R_{ij} - u ^{\prime} (R_{ij}^{0}) / R_{ij}^{0} ) ] \sigma _{ij}^{(2)T}\cr &\! +(1 / \beta) \textstyle\sum\limits _{\bf k} \ln [(\sinh f_{\bf k}) / f_{\bf k}]\, , &\hfill\llap{(15)}}]

where [V(\bar{\bf X})] is the bare classical potential, u(Rij) is the two-body potential between atoms i and j with the distance Rij, and Rij0 is the equilibrium distance. [\sigma _{ij}^{(2)L}] and [\sigma _{ij}^{(2)T}] are respectively expressed as

[\sigma _{ij}^{(2)L} = (2 / Nm) \textstyle\sum\limits _{{\bf k},\mu} \, (1-\cos {\bf k} \cdot {\bf R}_{ij}^{0}) \,(\hat{{\bf R}}_{ij}^{0} \cdot {\bf e}_{{\bf k},\mu})^{2} \,\alpha _{{\bf k},\mu} \eqno(16)]

and

[\sigma _{ij}^{(2)T} = ( 2 / Nm) \textstyle\sum\limits _{{\bf k},\mu} \, (1-\cos {\bf k} \cdot {\bf R}_{ij}^{0}) \,[1 - (\hat{{\bf R}}_{ij}^{0} \cdot {\bf e}_{{\bf k},\mu})^{2}]\, \alpha _{{\bf k},\mu} \, , \eqno(17)]

where [\hat{\bf R}_{ij}^{0}] is the unit vector of [{\bf R}_{ij}^{0}]. [\sigma _{ij}^{(2)L}] is the difference of the EXAFS Debye–Waller factors between the quantum-mechanical and the classical schemes. The formula can be adopted for a monatomic Bravais lattice with a two-body potential, such as rare-gas solids.

2.3. Application of the embedded-atom method to PIECP

Metallic bonds are regarded as a many-body force between valence electron densities and ion cores, and two-body interaction potentials are inadequate to describe thermal properties of metals. The EAM has extensively been employed in classical simulations of metals, which assumes semi-empirical parameters based on the density-functional theory. The EAM potential V is written as

[V = \textstyle\sum\limits _{i} \Big [{\scr F}_{i}(\rho _{h,i}) + {\textstyle {1 \over 2}} \textstyle\sum\limits _{j \neq i} \varphi _{ij} (R_{ij}) \Big] \, , \eqno(18)]

where ρh,i is the host electron density,

[\rho _{h,i} = \textstyle\sum\limits _{j \neq i} \rho _{j}^{a} (R_{ij}) \, , \eqno(19)]

and [\rho _{j}^{a} (R_{ij})] is the electron density of atom j at the site of atom i. [{\scr F}_{i}(\rho _{h,i})] is the embedding energy of atom i into the host and corresponds to an attractive many-body force between valence electrons and ion cores. [{\scr F}_{i}(\rho _{h,i})] is in principle a functional but is replaced by a normal function Fi under the local-density approximation. [\varphi _{ij} (R_{ij})] is the repulsive potential between ion cores i and j. [\rho _{j}^{a} (R_{ij})] can be calculated by using atomic wavefunctions, and the function forms of Fi and φij are empirically parameterized.

The EAM potential seems not to be applicable to the PIECP formula (15)[link] since it contains the many-body force of Fi. When the EAM potential is Taylor expanded up to the second order, however, we obtain

[V \cong N [F(\bar{\rho}) - \bar{\rho} F ^{\,\prime} (\bar{\rho}) ] + {\textstyle {1 \over 2}} \textstyle\sum\limits _{i \neq j} u (R_{ij}) \, , \eqno(20) ]

where [\bar{\rho}] is given in the equilibrium geometry and

[u(R) = \varphi (R) + 2 F ^{\,\prime} (\bar{\rho}) \rho _{j}^{a}(R) + F ^{\,\prime \prime} [\rho ^{a} (R) ]^{2} \, .\eqno(21)]

The expansion includes a very important consequence: the EAM potential is two body within the harmonic approximation. We can thus employ equation (15[link]) of a two-body potential system within the EAM. This is ascribed to the absence of angular-dependent terms in the EAM potential, although it is inherently a many-body potential. We can replace the bare classical potential [V(\bar{\bf X})] in equation (15[link]) with the EAM potential V in equation (18[link]) for PIECP calculations of metals.

3. Computational details

For the PIECP calculations of bulk Cu, the above formalism is available without further approximation. The employed EAM potential was obtained from the literature (Foiles et al., 1986[Foiles, S. M., Basks, M. I. & Daw, M. S. (1986). Phys. Rev. B, 33, 7983-7991.]). The normal vibrational analysis was initially performed using a cubic Brillouine zone (−2π/a0, 2π/a0) (a0 is the f.c.c. lattice constant) by sampling about 106 phonons. This yielded [\sigma _{ij}^{(2)T}] and [\sigma _{ij}^{(2)L}] in Veff. The classical-like NPT (constant pressure and temperature in a closed system) Monte Carlo (MC) calculations were subsequently preformed using Veff instead of the bare potential for 256 Cu atoms (43 f.c.c. unit cells). Three-dimensional periodicity was imposed. 20000 MC steps were calculated to reach the equilibrium and a further 10000 MC steps were performed to obtain the EXAFS cumulants; each MC step contained 256 movements of atoms and one variation of the lattice constant.

For the thin Cu film, there exists no three-dimensional periodicity, thus leading to difficulties in applying the above formula. In a practical sense, however, for the estimation of the EXAFS cumulants for comparison with the experimental data, the quantum-mechanical correction in Veff can be replaced by the bulk values, although rather crudely from a theoretical point of view. This would give some overestimation of the vibrational amplitude of surface phonons at low temperature since the eigenfrequencies should be smaller than the bulk ones. Within this approximation, the PIECP calculations can similarly be performed. A six-layer Cu(111) film was assumed, each layer containing 48 atoms in a rectangular lattice; the lowest layer was assumed to be vibrationally fixed. Similar NPT MC simulations were performed to obtain the EXAFS cumulants of intra- and interlayer atom pairs.

4. Results and discussion

The calculated EXAFS cumulants for the first-nearest-neighbour Cu–Cu shell of bulk Cu are shown in Fig. 1[link], together with the experimental observations. Although the calculated values of the distance R are shifted by ∼0.01 Å compared with the EXAFS and X-ray diffraction experimental data, the temperature dependence of R (slope of the curves) indicates that the thermal expansions are comparable. The deviations of the absolute values are partly derived from the fact that the EAM parameters were determined in a classical scheme. Other temperature-dependent quantities, such as the mean-square relative displacement C2 and the mean-cubic relative displacement C3, are also in good agreement with the experimental EXAFS values. In the C2 plot, the Debye fit was evaluated by using the correlated Debye model (Beni & Platzman, 1976[Beni, G. & Platzman, P. M. (1976). Phys. Rev. B, 14, 1514-1518.]) with the Debye temperature as a fitting parameter. It can be concluded here that the PIECP calculations combined with the EAM is sufficiently adequate to describe the thermal properties of bulk Cu.

[Figure 1]
Figure 1
Temperature dependence of the interatomic distance R, the mean-square relative displacement C2 and the mean-cubic relative displacement C3 for the first-nearest-neighbour Cu–Cu shell in bulk Cu. The PIECP results are given as crosses and solid lines, while the experimental EXAFS data are shown as filled squares. X-ray diffraction data for the distance (dotted line) and the Debye fitting results for C2 (dashed line) are also shown.

Let us now turn to the results for the Cu(111) film. The temperature dependences of R, C2 and C3 for the first-nearest-neighbour Cu–Cu shells are depicted in Figs. 2[link], 3[link] and 4[link], respectively. In Figs. 2[link], 3[link] and 4[link], the obtained quantities for each layer are given separately. For instance, the curves denoted 12 correspond to the first-and-second interlayer quantities, while those denoted 11 correspond to the first-and-first intralayer quantities. Note that the first-and-second layer (12) distance in Fig. 2[link], as well as all the intralayer distances 11, 22, 33 and 44, are contracted compared with the bulk one, while the second-and-third (23) and third-and-fourth (34) layer distances are elongated. Thermal expansion (Fig. 2[link]), C2 (Fig. 3[link]) and C3 (Fig. 4[link]) are found to be significantly enhanced for the first-and-second layer.

[Figure 2]
Figure 2
Temperature dependence of the interatomic distance R for several first-nearest-neighbour Cu–Cu shells in a 6 ML Cu(111) film, together with the bulk values. 11 implies the first-and-first intralayer distance, 12 indicates the first-and-second interlayer distance, etc. (11: open squares; 22: open circles; 33: upward open triangles; 44: downward open triangles; 12: filled squares; 23: filled circles; 34: filled upward triangles; intralayer: dashed line; interlayer: solid line; bulk: crosses and dot-dashed line).
[Figure 3]
Figure 3
Temperature dependence of the mean-square relative displacement C2 for several first-nearest-neighbour Cu–Cu shells in a 6 ML Cu(111) film, together with the bulk values (see the legend of Fig. 2[link] for notations).
[Figure 4]
Figure 4
Temperature dependence of the mean-cubic relative displacement C3 for several first-nearest-neighbour Cu–Cu shells in a 6 ML Cu(111) film, together with the bulk values (see the legend of Fig. 2[link] for notations).

We will next compare the present theoretical results with experimental data (Kiguchi et al., 2000[Kiguchi, M., Yokoyama, T., Matsumura, D., Endo, O., Kondoh, H. & Ohta, T. (2000). Phys. Rev. B, 61, 14020-14027.]). In the previous work, polarization-dependent Cu (and also Ni) K-edge EXAFS spectra were recorded for ultra-thin Cu(111) [or Ni(111)] films grown epitaxically on HOPG (highly oriented pyrolitic graphite) by varying the thickness. By using the results for the thickness dependence of anisotropic C2 and C3, the surface out-of-plane and in-plane components can be extracted. Moreover, by assuming (roughly) that all the contributions, other than the in-plane first-and-first (11) and the out-of-plane first-and-second (12) contributions, are equal to the bulk contribution, three components of 11, 12 and bulk were successfully determined. For detailed procedures of the data analysis, see the work of Kiguchi et al. (2000[Kiguchi, M., Yokoyama, T., Matsumura, D., Endo, O., Kondoh, H. & Ohta, T. (2000). Phys. Rev. B, 61, 14020-14027.]). The results are summarized in Table 1[link]. In the experimental work, the parameters were obtained as differences between 100 and 300 K. For C2, the ΔC2 values are replaced by the effective Debye temperature θD, which is similarly evaluated by using the correlated Debye model (Beni & Platzman, 1976[Beni, G. & Platzman, P. M. (1976). Phys. Rev. B, 14, 1514-1518.]).

Table 1
Effective Debye temperature θD and difference of the third-order EXAFS cumulants ΔC3 (difference between 100 and 300 K) for the surface out-of-plane, surface in-plane and bulk vibrations

Both experimental (Kiguchi et al., 2000[Kiguchi, M., Yokoyama, T., Matsumura, D., Endo, O., Kondoh, H. & Ohta, T. (2000). Phys. Rev. B, 61, 14020-14027.]) and theoretical results (this work) are given.

  Out of plane In plane Bulk
  θD ΔC3 θD ΔC3 θD ΔC3
Calculation 272 3.65 290 2.47 313 1.48
Experiment 262 (25) 3.8 (8) 322 (30) 3.1 (6) 338 1.62

Although the difference between the out-of-plane and the in-plane vibrations is underestimated in the calculations, the calculated and experimental results agree well semiquantitatively. It is concluded, both experimentally and theoretically, that the surface out-of-plane vibration is significantly softer and more anharmonic than the inner, lateral or bulk vibration. As a reason for the underestimation in the calculations, one can suppose that the present theoretical calculations are valid for perfect films, while the actual films should contain many defects and some roughening might already have occurred at room temperature. This is indicated by the experimental findings of smaller coordination number and larger C2 (Kiguchi et al., 2000[Kiguchi, M., Yokoyama, T., Matsumura, D., Endo, O., Kondoh, H. & Ohta, T. (2000). Phys. Rev. B, 61, 14020-14027.]). Since the surface area becomes wider in the presence of defects, the surface Debye temperature would effectively be lowered. One should also note that the surface Debye temperature of Cu(100) determined by LEED (low energy electron diffraction) is 235 K (Müller et al., 1995[Müller, S., Kinne, A., Kottcke, M., Metzler, R., Bayer, P., Hammer, L. & Heinz, K. (1995). Phys. Rev. Lett. 75, 2859-2862.]), which is still lower than the EXAFS result of 262 K.

It is essential to recall intrinsic differences between EXAFS and diffraction. Diffraction techniques such as LEED provide information on absolute displacements with respect to the lattice, while EXAFS gives relative displacements between X-ray-absorbing and neighbouring atoms. Although within the correlated Debye model the value of the Debye temperature should be equal to that given by diffraction, C2 and C3 are obtained uniquely by EXAFS. It is natural that the out-of-plane vibration is enhanced at the surface as long as the absolute displacement is discussed. This is ascribed simply to a lack of atoms above the surface. On the contrary, from a local point of view, the out-of-plane motion is not always enhanced, as can be seen in the literature: in the p4g(2 × 2)N/Ni(001) surface, the out-of-plane N–Ni bond has been found to be much stiffer and less anharmonic than the lateral bonds (Wenzel et al., 1990[Wenzel, L., Arvanitis, D., Rabus, H., Lederer, T., Baberschke, K. & Comelli, G. (1990). Phys. Rev. Lett. 64, 1765-1768.]); similarly the out-of-plane S–Ni bond in c(2 × 2)S/Ni(110) has been found to be stiffer and less anharmonic than the lateral bonds (Yokoyama et al., 1994[Yokoyama, T., Hamamatsu, H., Kitajima, Y., Takata, Y., Yagi, S. & Ohta, T. (1994). Surf. Sci. 313, 197-208.]). In the present case, the surface out-of-plane bond is found to be weaker and more anharmonic than the in-plane bond.

The aim of this work is to obtain a hint of surface melting, which is an initial stage of bulk melting. The enhancements of the vibrational amplitude and anharmonicity along the surface normal give rise to an important conclusion, since a roughening transition should occur through the `hopping' of surface atoms, for which not only the amplitude but the anharmonicity of the vibration are essential.

5. Conclusions

In this article, the usefulness of the PIECP theory for the evaluation of the EXAFS cumulants or, in general, the radial distribution function is demonstrated. In the PIECP theory, the physical quantities converge to be harmonic at low temperature and to be classical at high temperature. Therefore, anharmonicity at low temperature cannot be estimated. In order to obtain such information, one has to perform quantum-statistical perturbation calculations or more sophisticated path-integral Monte Carlo simulations. Low-temperature anharmonicity is usually not very important for EXAFS researchers and the PIECP theory is thus practically very useful to simulate the EXAFS cumulants. In summary, the PIECP theory can handle the quantum effect, anharmonicity, many degrees of vibrational freedom, three-dimensional periodicity, higher-nearest-neighbour interactions and many-body interactions. It is essentially difficult or impossible to include all these contributions using other theories, at least in a practical sense of numerical calculations. It is worth noting again that the EAM matches the PIECP theory excellently, this allowing one to investigate the vibrational properties of metals quantum statistically.

In the simulations of thin Cu(111) film, thermal vibration and local thermal expansion are larger than in the bulk metal. The relative motions, focussing on the local bonds of the surface, are found to be enhanced in the surface normal direction. These findings are consistent with the EXAFS experiments. The enhanced vibrational amplitude and anharmonicity in the surface normal direction could be the trigger of a roughening transition, surface melting and consequent bulk melting.

Acknowledgements

The author is grateful for the financial support of the Grant-in-Aid for Scientific Research (No. 09640598).

References

First citationBeni, G. & Platzman, P. M. (1976). Phys. Rev. B, 14, 1514–1518. CrossRef CAS Web of Science
First citationCuccoli, A., Giachetti, R., Tognetti, V., Vaia, R. & Verrucchi, P. (1995). J. Phys. Condens. Matter, 7, 7891–7938. CrossRef CAS Web of Science
First citationCuccoli, A., Macchi, M., Neumann, M., Tognetti, V. & Vaia, R. (1992). Phys. Rev. B, 45, 2088–2096. CrossRef Web of Science
First citationFoiles, S. M. (1985). Phys. Rev. B, 32, 3409–3415, 7685–7693.
First citationFoiles, S. M., Basks, M. I. & Daw, M. S. (1986). Phys. Rev. B, 33, 7983–7991. CrossRef CAS Web of Science
First citationFujikawa, T. & Miyanaga, T. (1993). J. Phys. Soc. Jpn, 62, 4108–4122. CrossRef CAS Web of Science
First citationFujikawa, T., Miyanaga, T. & Suzuki, T. (1997). J. Phys. Soc. Jpn, 66, 2897–2906. CrossRef CAS Web of Science
First citationKiguchi, M., Yokoyama, T., Matsumura, D., Endo, O., Kondoh, H. & Ohta, T. (2000). Phys. Rev. B, 61, 14020–14027. Web of Science CrossRef CAS
First citationKleinert, H. (1995). Path Integrals in Quantum Mechanics, Statistics and Polymer Physics. Singapore: World Scientific.
First citationMüller, S., Kinne, A., Kottcke, M., Metzler, R., Bayer, P., Hammer, L. & Heinz, K. (1995). Phys. Rev. Lett. 75, 2859–2862. CrossRef PubMed Web of Science
First citationWenzel, L., Arvanitis, D., Rabus, H., Lederer, T., Baberschke, K. & Comelli, G. (1990). Phys. Rev. Lett. 64, 1765–1768. CrossRef PubMed CAS Web of Science
First citationYokoyama, T. (1998). Phys. Rev. B, 57, 3423–3432. Web of Science CrossRef CAS
First citationYokoyama, T. (1999). J. Synchrotron Rad. 6, 323–325. Web of Science CrossRef CAS IUCr Journals
First citationYokoyama, T., Hamamatsu, H., Kitajima, Y., Takata, Y., Yagi, S. & Ohta, T. (1994). Surf. Sci. 313, 197–208. CrossRef CAS Web of Science
First citationYokoyama, T., Ohta, T. & Sato, H. (1997). Phys. Rev. B, 55, 11320–11329. CrossRef CAS Web of Science
First citationYokoyama, T., Yonamoto, Y., Ohta, T. & Ugawa, A. (1996). Phys. Rev. B, 54, 6921–6928. CrossRef CAS Web of Science

© International Union of Crystallography. Prior permission is not required to reproduce short quotations, tables and figures from this article, provided the original authors and source are cited. For more information, click here.

Journal logoJOURNAL OF
SYNCHROTRON
RADIATION
ISSN: 1600-5775
Follow J. Synchrotron Rad.
Sign up for e-alerts
Follow J. Synchrotron Rad. on Twitter
Follow us on facebook
Sign up for RSS feeds