- Email: [email protected]
Molecular dynamics study on surface structure and surface energy of rutile TiO2 (1 1 0)
Applied Surface Science 255 (2009) 5702–5708
Contents lists available at ScienceDirect
Applied Surface Science journal homepage: www.elsevier.com/locate/apsusc
Molecular dynamics study on surface structure and surface energy of rutile TiO2 (1 1 0) Dai-Ping Song *, Ying-Chun Liang, Ming-Jun Chen, Qing-Shun Bai Precision Engineering Research Institute, Harbin Institute of Technology, Harbin 150001, China
A R T I C L E I N F O
A B S T R A C T
Article history: Received 28 August 2008 Received in revised form 1 December 2008 Accepted 23 December 2008 Available online 31 December 2008
The formula for surface energy was modified in accordance with the slab model of molecular dynamics (MDs) simulations, and MD simulations were performed to investigate the relaxed structure and surface energy of perfect and pit rutile TiO2(1 1 0). Simulation results indicate that the slab with a surface more than four layers away from the fixed layer expresses well the surface characteristics of rutile TiO2 (1 1 0) surface; and the surface energy of perfect rutile TiO2 (1 1 0) surface converges to 1:801 0:001 J m2. The study on perfect and pit slab models proves the effectiveness of the modified formula for surface energy. Moreover, the surface energy of pit surface is higher than that of perfect surface and exhibits an upper-concave parabolic increase and a step-like increase with increasing the number of units deleted along [0 0 1] and [1 1 0], respectively. Therefore, in order to obtain a higher surface energy, the direction along which atoms are cut out should be chosen in accordance with the pit sizes: [1¯ 1 0] direction for a small pit size and [0 0 1] direction for a big pit size; or alternatively the odd units of atoms along [1 1 0] direction are removed. ß 2008 Elsevier B.V. All rights reserved.
PACS: 31.15.xv 68.35.Md 68.35.bd 61.72.Ff Keywords: Molecular dynamics TiO2 Surface energy Surface structure Pit
1. Introduction Ti-materials have extensive biomedical applications [1], especially in medical implant and prosthesis [2], for their good biocompatibility, biological responses, osseointegration, excellent mechanical properties and corrosion resistance. In fact, there is a layer of TiO2on Ti and Ti-alloy surfaces, which has a strong direct effect on the success of implant [3,4]. The surface topographic characters of an implant at the micro/ nano-scale are related to the biological response of a host. The roughness of more than 10 mm of micro-rough surface will has an influence on the biocompatibility, mechanical characteristics and mechanical interlocking effect between implant and tissue. Due to the same size as a cell or a biomacromolecule, the micro-rough surface with a roughness ranging from 10 nm to 10 mm has remarkable influences on the biocompatibility and weakens the mechanical characteristics of the interface [5]. The micro-roughness of less than 10 nm of the surface has a significant influence on the interface structures of the implant. The reason is that the defects in a crystal structure, such as vacancy, grain boundary, step
* Corresponding author. Tel.: +86 451 86413840; fax: +86 451 86415244. E-mail address: [email protected] (D.-P. Song). 0169-4332/$ – see front matter ß 2008 Elsevier B.V. All rights reserved. doi:10.1016/j.apsusc.2008.12.062
and pit, are in this size range, which are active region for adsorption and influence further the integration of biomolecular and implant surface [6]. Moreover, the surface energy of biomaterials significantly influences the adhesion, spreading and growth of a cell. The surface with a higher surface energy promotes adhesion and spreading of a cell [7,8]. For example, good spreading only occurs when surface energy is higher than approximately 0.057 J m2[9]. Baier et al. [10,11] also showed the relationship between critical surface energy and cell adhesion: materials with a low surface energy show a low cell attachment. The rutile (1 1 0) surface is a very stable crystal face. Many scholars [12–17] have researched the surface structure of rutile (1 1 0) by ab initio calculations, such as density functional theory (DFT) local density approximation (LDA) [12,13], generalized gradient approximation (GGA) [14], full-potential linear augmented plane wave (FP-LAPW) [15,16] and linear combination of atomic orbitals (LCAO) [17]. On the other hand, Matsui and Akaogi [18] and Kim et al. [19] have developed the molecular dynamics (MDs) potential energy function and its parameters for TiO2 since 1991. A few scholars [20,21] studied the relaxed structure of the TiO2 (1 1 0) surface by MDs. However, there was a discrepancy between their theoretical [12–17,20–23] and experimental [24] results.
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5703
The surface energy of the rutile TiO2 (1 1 0) surface has been calculated recently by first principles [12–14,25,26], but there is few systematic calculations by MD possibly because the calculation accuracy of MD simulation is lower than that of first principles for small-scale models. However, the larger amount of computational work and the smaller atom amount of first principles make it very difficult to calculate the large-scale surface structures and evaluate the changes in surface energies resulting from the difference of the large-scale complex surface structures. Therefore, we modify the formula for surface energy in accordance with the slab model of MD simulations. Then, we study the relaxed structures of the rutile TiO2(1 1 0) surface by MD and investigate the relationship between surface energy and slab thickness. Finally, we also cut out some atoms to create pits, calculate the surface energies and relaxed structures of pit surfaces and study the influence of pit size on the surface energy. Fig. 1. Rutile TiO2(1 1 0) surface (1 1): O, big red ball; Ti, small grey ball. Atom labels follow as bridging oxygen O1, six-coordinate titanium Ti2, five-coordinate titanium Ti3, three-coordinate in-plane oxygen O4, O5, out-of-plane oxygen O6, and reference atom O7.
2. Methods and surface modeling 2.1. Molecular dynamics simulation method A survey of the literature indicated that several force fields have been published for titanium dioxide and titanium oxide [18,19,27– 32]. A detailed analysis of the different available force fields has been published by Smith and his co-worker [31]. They concluded that the force field of Matsui and Akaogi [18] is the most suitable of the available force fields for use in classical molecular dynamics simulation. Matsui and Akaogi [18] developed the following force fields for TiO2 in 1991. Vij ¼
qi q j C i C j Ai þ A j r i j 6 þ f Bi þ B j exp ri j Bi þ B j ri j
(1)
where the terms represent Coulomb, dispersion, and repulsion interactions, respectively. Here r i j is the interatomic distance between atoms i and j, and qi , Ai , Bi , and C i are the effective charge, repulsive radius, softness parameter, and van der Waals coefficient of atom i, respectively. Ai and Bi are relative to atom size and compressibility. The quantity f is a standard force of 4.184 kJ A˚1 mol1. For many oxides, Coulomb and repulsion interactions are necessary. Consequently, in this paper we used the reduced form for Eq. (1), i.e. Buckingham potential together with Coulombic interactions from partial charges ! r i j C i j qi q j V i j ¼ Ai j exp (2) 6 þ ri j ri j ri j where V i j is the interaction energy between atoms i and j, and the parameters Ai j , ri j , C i j , qTi and qO are shown in Table 1[18]. MD simulations were carried out at 310 K for further study on the interaction between TiO2 surface and biomolecules. Simulations were carried out in the canonical ensemble (i.e. NVT ensemble), and the thermostat derived by Hoover [34] was employed to maintain a constant temperature. The velocity Verlet algorithm was used to calculate the atomic motions and the particle–particle particle–mesh (PPPM) [35] solver was applied for the calculation of electrostatic interactions. A 10.0 A˚ cutoff distance was used for van der Waals interactions. Table 1 Interaction parameters for TiO2.
2.2. Perfect surface models The most stable (1 1 0) surface of rutile TiO2 was investigated in this p ffiffiffi paper. The (1 1 0) surface unit cell (1 1) has a dimension of 2a c which corresponds to 6.497 and 2.959 A˚ along [1¯ 1 0] and [0 0 1] directions, respectively [36,37]. The surface contains fiveand six-coordinate Ti atoms and two types of O atoms (in-plane and bridging). The six-coordinate Ti atoms are covered by the outermost bridging oxygen on the surface, while the fivecoordinate Ti atoms are coordinated to in-plane oxygen atoms on the surface. One Ti–O layer of rutile TiO2(1 1 0) surface unit cell (1 1) is built up of symmetrical three-plane O–Ti2O2–O unit along [1 1 0] direction, containing 6 atoms (2 Ti and 4 O), so that the sequence of Ti–O layers is O–Ti2O2–O–O–Ti2O2–O , as shown in Fig. 1. We constructed (10 18) surface slabs with the thickness of n Ti–O layers based on the rutile (1 1 0) surface unit cell (1 1), as shown in Table 2. For example, when n ¼ 13, this corresponds to a ˚ containing 14040 dimension of 64:969 A˚ 53:262 A˚ 42:229 A, atoms. Periodic boundary conditions were applied along [1¯ 1 0] and [0 0 1] directions of the surface, with a periodic vacuum gap of threefold slab thickness (3 3:2484n A˚) in the direction normal to the (1 1 0) surface. The simulations were carried out for 50 ps with time step of 0.5 fs, and the simulation results of last 30 ps were analyzed. 2.3. Pit surface models Based on the perfect TiO2(1 1 0) surface model, we cut out some atoms from the relaxed surface to create pit, as shown in Table 3. Fig. 2 shows the unrelaxed pit surface model (PitL16B2S1x (5 6 1)) in which the deleted atoms range over (5 6 1). For the pit surface models, due to the longer time to equilibrate the Table 2 Perfect surface slab models. Model set ID
Number of total layers n
Number of fixed layers nfixed 1 for n ¼ 2; 2 for n ¼3–6 1 for n ¼2–5; floor ðn=2Þ for n ¼6–16 1
Interation
Ai j (kcal mol1)
ri j (A˚)
C i j (kcal mol1 A˚6)
Perfect I Perfect II
2–16 2–16
Ti–Ti Ti–O O–O
717653.9571 391052.7442 271718.8311
0.154 0.194 0.234
120.9967 290.3920 696.9407
Perfect III
2–17
Notes: Atomic charges, qTi ¼ þ2:196, qO ¼ 1:098.
Notes: Type of fixed layers: Perfect I and II are one side fixed, and Perfect III are the central ceil(n=2)th layer fixed.
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5704 Table 3 Pit surface slab models. Model set ID a
PitL16B2S1x PitL16B2S1ya PitL16B2S1za PitL17C1S1xb PitL17C1S1yb PitL17C1S1zb PitL17C1S2xc PitL17C1S2yc PitL17C1S2zc
Pit size (½1¯ 1 0 ½0 0 1 ½1 1 0)
N delunit
N delunit 6 1 5 N delunit 1 5 6 N delunit N delunit 6 1 5 N delunit 1 5 6 N delunit N delunit 6 1 5 N delunit 1 5 6 N delunit
1–7 2–14 1–6 1–7 2–14 1–4 1–7 2–14 1–4
Notes:a is based on the 16-layer slab with the bottom 2 layers fixed and pit locates on the top surface.b,c is based on the 17-layer slab with the central 9th layer fixed.bPit locates on the top surface.cPits locate on the top and bottom surfaces.
surface atoms, the simulations were carried out for 200 ps, and the results of last 50 ps were analyzed. The other conditions of MD simulations for the pit surfaces are the same as those for the perfect surfaces given above. 2.4. Surface energy calculations According to the Gibbs definition [38], the surface energy density (s ) of a solid, which corresponds to the energy variation (per unit area) due to the creation of a surface, is given by
s¼
Enslab nEbulk Aslab
(3)
where Aslab is the total area of the surface considered, Enslab is the total energy of a n-layer slab and Ebulk is the bulk energy per layer of an infinite solid. However, the surface energy calculated in Eq. (3) is very sensitive to the accuracy in determination of the bulk energy term [39,40]. A small error in that term can make the calculated surface energies diverge linearly with increasing slab thickness. To avoid the divergence problem, we have employed the method [26,41,42] that makes use of Eq. (3) rewritten in the following form: Enslab ¼ Aslab s þ nEbulk
(4)
which implies that the bulk energy can be extracted from the slope of a linear fit of the slab_s total energy plotted versus n. This value is subsequently used in Eq. (3). For molecular dynamics simulations, a slab model, especially one with defects, usually needs to contain fixed layers, otherwise the slab structure will be distorted. To calculate the effect of the fixed layers on Enslab and surface energy, we simulated the following
perfect surface slab models which will not be distorted: (i) a sufficiently large number of layers are unconstrained and other layers are fixed, totally m layers; (ii) all m layers are unconstrained. The effect (DEfixed ) of fixed layers can be estimated from the difference between the total energies of the former and the latter. DEfixed is relative to the type and number of fixed layers (Tables 2 and 3). The effect of fixed layers on Enslab and surface energy can then be avoided by subtracting DEfixed from Enslab . Therefore, the surface energy of slab model containing fixed layers can be given by the following modified equation:
s¼
Enslab nEbulk DEfixed Aslab
(5)
For a large-scale complex surface structure, Eq. (5) can be rewritten as
s¼
unit EN slab NEbulk DEfixed Aslab
(6)
where N is the total number of units in slab. For our slab models, the numerators in Eqs. (3), (5) and (6) are the energy variation due to the creation of top and bottom surfaces, denoted as DE. Because the structure and area of the top surface is the same as those of the bottom surface in the perfect surface slab, the contribution of each perfect surface to DE is 50% of the energy variation (DEPerfM ) of the slab model with two-sided perfect surfaces. However, for the slab model with a perfect surface on one side and a pit surface on the other side, the contribution of the pit surface to DE is not 50% of the energy variation (DEHybridM ) of this slab model and the surface energy of the pit surface cannot be directly obtained using Eqs. (5) and (6). To calculate the surface energy of pit surface and compare DE and s with those of other slab models, we removed the contribution of the perfect surface to DE by transforming the slab model with one-sided pit (PitL16B2S1x, PitL17C1S1x, PitL16B2S1y, PitL17C1S1y, PitL16B2S1z and PitL17C1S1z) into the slab model with twosided pits (PitL16B2S1TS2x, PitL17C1S1TS2x, PitL16B2S1TS2y, PitL17C1S1TS2y, PitL16B2S1TS2z and PitL17C1S1TS2z, respectively). 2ðDEHybridM DEPerfM =2Þ and twice the surface area of the pit surface are used as the energy variation DE and the surface area Aslab of the transformed slab model, respectively. So, the surface energy of the pit surface of slab model with one-sided pit can be calculated using Eq. (6). The surface energy of the pit surface is the surface energy of the slab model with two-sided pits or the transformed slab model with two-sided pits. 3. Results and discussion 3.1. Analysis on perfect TiO2(1 1 0) surface structure
Fig. 2. Unrelaxed pit surface model PitL16B2S1x (5 6 1): O, big red ball; Ti, small grey ball. The box with green lines represents the pit.
The perfect surface unit cell (1 1) of the rutile TiO2 (1 1 0) contains bridging oxygen O1, six-coordinate titanium Ti2, fivecoordinate titanium Ti3, three-coordinate in-plane oxygen O4 and O5, and out-of-plane oxygen O6, as shown in Fig. 1. Supposing that O7 are reference atoms, Fig. 3 shows the relationship between the average displacements of (1 1 0) surface unit cells (10 18) along [1 1 0] and the slab thickness, which is a function of slab thickness. Although the average displacements of O1, Ti2 and O6 change considerably, the displacements of all atoms reach certain values if n 6 for Perfect I and n 9 for Perfect III. We compared further the average displacements of surface atoms in slab models Perfect I containing 6–16 layers with those of other scholars’ works in Table 4. The average displacements of all surface atoms are qualitatively in agreement with the experimental finding of Charlton [24] and the previous calculations [12– 15,17,20–23]. Moreover, our calculations coincide quantitatively
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5705
Moreover the average displacements of the atoms of (1 1 0) surface are less than 6e 4 A˚ along [0 0 1]. Thus it may seem that all atoms oscillate slightly in bulk positions along [1¯ 1 0] and [0 0 1] except the displacements of in-plane oxygen atoms along [1¯ 1 0]. 3.2. Surface energy of perfect TiO2(1 1 0) surfaces
Fig. 3. Calculated displacements of perfect surfaces.
with the experimental work [24] for the surface atoms except O1 and Ti3. On the outermost surface layer, the bridging oxygen atoms O1 and the five-coordinate titanium atoms Ti3 relax inward, while the six-coordinate titanium atoms Ti2 and the neighboring three-coordinate oxygen atoms O4 and O5 move outward, which causes a rumpled appearance of the surface. Although the average displacements of O1 and Ti3 are different considerably from those of experimental work [24], the average displacements of O6, O4, O5, and Ti2 along [1 1 0] are in good agreement with those of experimental work [24]. Langer et al. [43] reported that the bridging oxygen displaces laterally into an asymmetric position and the Ti–O distances are 1.75 and 2.1 A˚. However, in the experimental work [24], the bridging oxygen locates in the center of the neighboring six-coordinate titanium along [0 0 1], and the Ti–O distance is reduced to 1:71 0:07 A˚ mainly by relaxed inward largely (0:27 A˚). In our calculations, the O1 also locates laterally in a symmetric position and the Ti2–O1 average distance is 1.845 A˚, which is close to 1.82 and 1.83 A˚ for a five-layer slab and a three-layer slab, respectively by GGA calculations [37] and also in good agreement with 1.85 A˚ by MD simulation [18]. In addition, the average displacements of the three-coordinate in-plane oxygen atoms O4 and O5 along [1¯ 1 0] are very close to those of some works [12,14,15,22,23], and the in-plane oxygen atoms displace toward the neighboring five-coordinate titanium by 0.044 A˚, while the rest atoms displace by less than 1e 4 A˚.
Using 1828:6568 kcal mol1for Eunit bulk , as determined from a linear fit (Eq. (4)) to the slab energies (the energy of a thinnest fourlayer slab), yields very well converged DE and s (Figs. 4(a) and 5(a), respectively). Thus, the above values of Eunit bulk is adopted in all further calculations of s . Our calculations of energy variation DE and surface energies s versus slab thickness n are shown in Figs. 4(a) and 5(a), respectively. For the one-sided relaxed perfect surface models with different fixed layers and the two-sided relaxed perfect surface models (Table 2), their surface energies are converged to a value of 1:801 0:001 J m2as the slab thickness increases, which is very close to 1.78 J m2 from the MD simulations of Oliver et al. [44]. Although the surface energies obtained by MD simulations [33,44] are more than those by ab initio calculations, such as 0.89 J m2 from the average of the 5- and 6-layer LDA calculations of Ramamoorthy et al. [12], 0.73 J m2from the 7-layer GGA calculations of Bates et al. [14], 0.64 J m2from the average of the 9- and 10-layer PW–GGA calculations of Bredow et al. [45], 0.57 J m2by GGA–PW91 and 0.47 J m2by GGA–PBE from the average of the 10- and 11-layer of Kiejna et al. [26], large-scale surface slabs with complex morphologic surfaces can be evaluated more easily by MD simulations than by ab initio calculations. In addition, the calculation temperature results partly in the difference between our calculations and first-principles calculations [12,14,26,45]. The theoretical results of our reference literature of first-principles calculations are strictly valid only at zero temperature, while our MD simulations were carried out at 310 K. The surface energies of one-sided relaxed models are much the same as those of two-sided relaxed models by MD simulations, which is different from the result obtained by ab initio calculations from Kiejna et al. [26] that the surface energies of one-sided relaxed models are much more than those of two-sided relaxed models. Moreover, the slab model with one-sided relaxation is usually used as the adsorption substrate [46,47]. Therefore we will continue our study on the surface structure and surface energy of one-sided relaxed slab model with surface micro/nano-morphology.
Table 4 Atomic displacements for the perfect rutile TiO2 (1 1 0) surface (A˚). Atoma
9–18b
Expc[24]
4d[20]
6e[21]
7f[15]
5g[12]
5h[14]
3i[22]
4j[23]
O1 Ti2 Ti3 O4, O5 O4, O5 (1¯ 1 0) O6 O7
0.032 0.134 0.085 0.055 0:044 0.038 0
0:27 0:08 0:12 0:05 0:16 0:05 0:05 0:05 0:16 0:08 0:03 0:08 0 0:08
0.03 0.22 0.19 0.10 – 0.00 0.04
0.11 0.02 0.26 0.0 – 0.0 0.0
0.02 0.23 0.17 0.03 0:05 0.02 0
0.06 0.13 0.17 0.12 0:04 0.07 0.02
0.02 0.23 0.11 0.18 0:05 0.03 0.03
0.09 0.09 0.12 0.11 0:05 0.05 –
0.04 0.19 0.14 0.15 0:05 – –
a b c d e f g h i j
Atom labels refer to Fig. 1. This paper: MD, average displacements of 6–16 layers in Perfect I, S:D: 0:001 A˚. Charlton: SXRD experimental results [24]. Yin: MD, four layers [20]. San Miguel et al.: MD, six layers [21]. Harrison et al.: LCAO, seven layers [15]. Ramamoorthy et al.: PW–PP–LDA, five layers [12]. Bates et al.: PW–GGA, five layers [14]. Lindan et al.: PW–PP–GGA, three layers [22]. Maria et al.: GGA–PBE, four layers [23].
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5706
Fig. 4. Energy variations DE of rutile TiO2 surfaces: (a) perfect (1 1 0) surface slabs, (b) surface slabs with pit by cutting out atoms along [1¯ 1 0], (c) surface slabs with pit by cutting out atoms along [0 0 1] and (d) surface slabs with pit by cutting out atoms along [1 1 0].
To determine the contribution of each atomic layer to the surface energy, we calculated the energetic contribution of any atomic layer to the surface energy using the formula in Refs. [48– 50]. We found that, for the perfect TiO2(1 1 0) surface, s is stored mainly in the first four layers, especially the two outer layers, with 86.4, 9.9, 1.9 and 1.0% of the energy, respectively. As shown in Figs. 3 and 5(a), the changes in the atomic displacements govern the changes in the surface energy as a function of slab thickness, which is in agreement with the results of Bates et al. [14]. Therefore, it can be concluded that the slab with a surface more than four layers away from the fixed layer can express well the surface characteristics of rutile TiO2(1 1 0) surface. 3.3. Surface energy of pit TiO2(1 1 0) surfaces In order to describe the atomic displacements in the pit surface, according to the formula for slip vector [51], we defined the atomic displacement vector as sa ¼
k 1X ðxab Xab Þ k b 6¼ a
(7)
where k is the number of the nearest neighbors to atom a, and xab and Xab are the vector differences of atoms a and b at current and reference positions, respectively. The reference configuration is the arrangement of atomic positions with unrelaxation. The unstable atoms with the dangling bonds around the pit result in the maximum displacement vectors. Fig. 6 shows the
displacement vector of the relaxed pit surface model PitL16B2S1x (5 6 1) with atoms colored by jsa j, if jsa j 0:2. As shown in Figs. 4 and 5, due to the effect of the perfect surface on energy variation DE and surface energies s , DE and s of the slab models with one-sided pit are lower than those of the slab models with two-sided pits. For the pit surfaces, we also calculated the energetic contribution of any atomic layer to the surface energy. We found that s is stored mainly in the first two layers. For example, for the pit surface PitL16B2S1x (5 6 1), the energetic contribution of the first two layers is 97.9 and 1.8% of the energy, respectively. It can be seen from Figs. 4 and 5 that the energy variations and surface energies of the pit surfaces are higher than those of the perfect surfaces. The reason is that the pits created by removing some atoms from the perfect surfaces weaken the interaction of the atoms on the interface of the pits and many unstable atoms remain, and after relaxation these atoms offset their initial positions as shown in Fig. 6, which lead to a significant increase in the surface energy. Furthermore, according to the theory of solid surface physics [52], the closer the arrangement of atoms is, the stronger the interaction of atoms is and the lower the surface energy is, which is in agreement with our calculations. Figs. 4(b–d) and 5(b–d) show the relationship between DE and s of pit surfaces and the number (N delunit ) of units deleted along [1¯ 1 0], [0 0 1] and [1 1 0]. The energy variations and surface energies of the surfaces with pit by cutting out atoms along [1¯ 1 0] change slightly along a middle-convex parabola, and the changes in their surface energies are less than about 0.1 J m2. When the number of units deleted along [0 0 1] is 6 and that along [1 1 0] is 1,
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5707
Fig. 5. Surface energies s of rutile TiO2 surfaces: (a) perfect (1 1 0) surface slabs, (b) surface slabs with pit by cutting out atoms along [1¯ 1 0], (c) surface slabs with pit by cutting out atoms along [0 0 1] and (d) surface slabs with pit by cutting out atoms along [1 1 0].
the quadratic regressions of N delunit values along [1¯ 1 0] with DE and s values of the pit surfaces are shown as the solid lines in Figs. 4(b) and 5(b), respectively. The energy variations and surface energies of the surfaces with pit by cutting out atoms along [0 0 1] exhibit an upper-concave parabolic increase with the increasing number of units deleted. When the number of units deleted along [1¯ 1 0] is 5 and that along [1 1 0] is 1, the quadratic regressions of N delunit values along [0 0 1] with DE and s values of the pit surfaces are shown as the solid lines in Figs. 4(c) and 5(c), respectively. The energy variations DE of the surfaces with pit by cutting out atoms along [1 1 0] is linearly relative to the number of units deleted. When the number of units deleted along [1¯ 1 0] is 5 and that along [0 0 1] is 6, the linear regression of Ndelunit values along
[1 1 0] with DE values is shown as the solid line in Fig. 4(d). However, as their surface areas increases significantly, their surface energies s show a step-like increase with the increasing number of units deleted. When N delunit is odd, s increases significantly. The means of s at Ndelunit are shown as the solid line in Fig. 5(d). If the number of deleted atoms is more than 180 atoms as shown in Fig. 5, the surface energy of the pit surface from which the atoms are removed along [0 0 1] is higher than that along [1¯ 1 0]. In order to obtain a higher surface energy, the direction along which atoms are cut out should be chosen in accordance with the pit sizes: [1¯ 1 0] direction for a small pit size and [0 0 1] direction for a big pit size. We can also remove the odd units of atoms along [1 1 0] direction. Because the pit surface has higher surface energy and more adsorption sites than the perfect surface, one can design the surface structure of implant using Figs. 4 and 5, and the regression equations therein, and investigate further the effect of surface structure and surface energy on biomolecular adsorption on Ti-materials.
4. Conclusions We modified the formula used for calculation of surface energy in accordance with the slab model of MD simulation and studied the surface energy and relaxed structure of perfect and pit rutile TiO2(1 1 0) through MD simulation. The following are our findings:
Fig. 6. Displacement vector of pit surface PitL16B2S1x (5 6 1) with atoms colored by jsa j.
1. The outermost surface bridging oxygens and five-coordinate titaniums relax inward, whereas the six-coordinate titaniums and the in-plane oxygens move outward. The slab with a surface
5708
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
more than four Ti–O layers away from the fixed layer expresses well the surface characteristics of rutile TiO2 (1 1 0) surface. To reduce the amount of the computational work and improve computational efficiency, it is feasible to use the slab model with one-sided relaxation as the surface model. The unstable atoms with the dangling bonds around the pit in the pit surface have the maximum displacement vectors. Larger pit size generally results in an increase in the region of atoms with large displacement vectors. 2. The study on perfect and pit slab models proves the effectiveness of the modified formula for surface energy. The Surface energy of the perfect rutile TiO2 (1 1 0) surface converges to 1:801 0:001 J m2. 3. The surface energy of the pit surface is higher than that of the perfect surface and exhibits an upper-concave parabolic increase and a step-like increase with increasing the number of units deleted along [0 0 1] and [1 1 0], respectively. The surface energy of the surface with pit by cutting out atoms along [0 0 1] or [1 1 0] changes more significantly than that along [1¯ 1 0]. Moreover, if the number of deleted atoms is over a certain value (180 atoms in this paper), the surface energy of the pit surface from which the atoms are cut out along [0 0 1] is higher than that along [1¯ 1 0]. Thus, in order to obtain a higher surface energy, the direction along which atoms are cut out should be chosen in accordance with the pit sizes: [1¯ 1 0] direction for a small pit size and [0 0 1] direction for a big pit size; or alternatively the odd units of atoms along [1 1 0] direction are removed. It can therefore be concluded that the modified formula for surface energy and the displacement vector can be used to calculate the surface energy of the large-scale complex structures and evaluate these surface structures by MD, respectively, which will be aid to analyze further the effect of surface micro/nanomorphology of the TiO2 covering titanium implant surfaces on cell adhesion and protein adsorption. Acknowledgments This work is funded by the National Natural Science Foundation of China (No. 50675050) and the Multidiscipline Scientific Research Foundation of Harbin Institute of Technology (No. HIT. MD 2003. 10). References [1] C. Leyens, M. Peters, Translated by Z.H. Chen et al. Titanium and Titanium Alloys, Chemical Industry Press, Beijing (2005). [2] U. Diebold, Surf. Sci. Rep. 48 (2003) 53. [3] A. Hunter, C.W. Archer, P.S. Walker, G.W. Blunn, Biomaterials 16 (1995) 287. [4] B.D. Boyan, T.W. Hummert, D.D. Den, Z. Schwartz, Biomaterials 17 (1996) 137. [5] B. Kasemo, J. Prosth. Dent. 49 (1983) 832. [6] B. Kasemo, J. Lausmaa, Material selection surface characteristics and chemical processes at implant surface, in: P.I. Branemark, G.A. Zarb, T. Albrektsson (Eds.),
[7] [8] [9] [10]
[11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52]
Tissue Integrated Prostheses: Osseointergtation in Clinical Dentistry, Quintessence Publishing Co. Inc., Chicago, 1985, p. 99. T.G. van Kooten, J.M. Schakenraad, H.C. van der Mei, H.J. Busscher, Biomaterials 13 (1992) 897. S.A. Redey, S. Razzouk, C. Rey, D. Bernache-Assollant, G. Leroy, M. Nardin, G. Cournout, J. Biomed. Mater. Res. 45 (1999) 140. J.M. Schakenraad, H.J. Busscher, C.R.H. Wildevuur, J. Arends, J. Biomed. Mater. Res. 20 (1986) 773. R.E. Baier, Surface properties influencing biological adhesion, in: R.S. Mansly (Ed.), Adhesion in Biological Systems, Academic Press, New York, 1970 , pp. 15–48. R.E. Baier, A.E. Meyer, J.R. Natiella, R.R. Natiella, J.M. Carter, J. Biomed. Mater. Res. 18 (1984) 327. M. Ramamoorthy, D. Vanderbilt, R.D. King-Smith, Phys. Rev. B 49 (1994) 16721. M. Ramamoorthy, R.D. King-Smith, D. Vanderbilt, Phys. Rev. B 49 (1994) 7709. S.P. Bates, G. Kresse, M.J. Gillan, Surf. Sci. 385 (1997) 386. N.M. Harrison, X.G. Wang, J. Muscat, M. Schefflet, Faraday Discuss. Chem. Soc. 144 (1999) 305. D. Vogtenhuber, R. Podlouchy, A. Nechel, S.G. Steinemann, A.J. Freeman, Phys. Rev. B 49 (1994) 2099. P. Reinhardt, B.A. Heb, Phys. Rev. B 50 (1994) 12015. M. Matsui, M. Akaogi, Mol. Simul. 6 (1991) 239. D.W. Kim, N. Enomoto, Z. Nakagawa, K. Kawamura, J. Am. Ceram. Soc. 79 (1996) 1095. X. Yin, R. Miura, A. Endou, I. Gunji, R. Yamauchi, M. Kubo, A. Stirling, A. Fahmi, A. Miyamoto, Appl. Surf. Sci. 119 (1997) 199. M.A. San Miguel, C.J. Calzado, J.F. Sanz, Surf. Sci. 409 (1998) 92. P.J.D. Lindan, N.M. Harrison, M.J. Gillan, J.A. While, Phys. Rev. B 55 (1997) 15919. L.S. Maria, Y.G. Andrey, L.S. Alexander, J. Phys. Chem. B 110 (2006) 4853. G. Charlton, Phys. Rev. Lett. 78 (1997) 495. J. Goniakowski, J.M. Holender, L.N. Kantorovich, M.J. Gillan, Phys. Rev. Lett. 53 (1996) 957. A. Kiejna, T. Pabisiak, S.W. Gao, J. Phys.: Condens. Matter 18 (2006) 4207. M. Mostoller, J.C. Wang, Phys. Rev. B 32 (1985) 6773. V. Swamy, J.D. Gale, Phys. Rev. B 62 (2000) 5406. S. Ogata, H. Iyetomi, K. Tsuruta, F. Shimojo, R.K. Kalia, A. Nakano, P. Vashishta, J. Appl. Phys. 86 (1999) 3036. H. le Roux, L. Glasser, J. Mater. Chem. 7 (1997) 843. D.R. Collins, W. Smith, Council for the Central Laboratory of Research Councils (1996). F.H. Streitz, J.W. Mintmire, J. Adhes. Sci. Technol. 8 (1994) 853. P.K. Naicker, P.T. Cummings, H. Zhang, J.F. Banfield, Phys. Phys. Chem. B 109 (2005) 15243. W.G. Hoover, Phys. Rev. A 31 (1985) 1695. H. Eastwood, Computer Simulation Using Particles, Adam Hilger, NY, 1989. S.C. Abrahams, J.L. Bernstein, J. Chem. Phys. 55 (1971) 3206. A.V. Bandura, J.D. Kubicki, J. Phys. Chem. B 107 (2003) 11072. J.W. Gibbs, Scientific Papers, Thermodynamics, vol. 1, Dover Publictions, New York, 1926, pp. 55. J.C. Boettger, Phys. Rev. B 49 (1994) 16798. V. Fiorentini, M. Methfessel, J. Phys.: Condens. Matter 8 (1996) 6525. J.C. Boettger, Phys. Rev. B 53 (1996) 13133. J.G. Gay, J.R. Smith, R. Richter, F.J. Arlinghaus, R.H. Wagoner, J. Vac. Sci. Technol. A 2 (1984) 931. W. Langel, Surf. Sci. 496 (2002) 141. P.M. Oliver, G.W. Watson, E.T. Kelsey, S.C. Parker, J. Mater. Chem. 7 (1997) 563. T. Bredow, L. Giordano, F. Cinquini, G. Pacchioni, Phys. Rev. B 70 (2004) 035419. J.R.B. Gomes, J.P. Prates Ramalho, Phys. Rev. B 71 (2005) 235421. M. Prˇedota, A.V. Bandura, P.T. Cummings, J.D. Kubicki, D.J. Wesolowski, A.A. Chialvo, M.L. Machesky, J. Phys. Chem. B 108 (2004) 12049. O. Politano, S. Garruchet, J.M. Salazar, Mat. Sci. Eng. A 387–389 (2004) 749. S. Garruchet, O. Politano, J.M. Salazar, T. Montesin, Surf. Sci. 586 (2005) 15. S. Garruchet, O. Politano, J.M. Salazar, A. Hasnaoui, T. Montesin, Appl. Surf. Sci. 252 (2006) 5384. J.A. Zimmerman, C.L. Kelchner, P.A. Klein, J.C. Hamilton, S.M. Foiles, Phys. Rev. Lett. 87 (2001) 165507. L.-S. Wen, Physics Basis of Solid Materials Interface Research, Science Press, Beijing, 1991, p. 21.
Contents lists available at ScienceDirect
Applied Surface Science journal homepage: www.elsevier.com/locate/apsusc
Molecular dynamics study on surface structure and surface energy of rutile TiO2 (1 1 0) Dai-Ping Song *, Ying-Chun Liang, Ming-Jun Chen, Qing-Shun Bai Precision Engineering Research Institute, Harbin Institute of Technology, Harbin 150001, China
A R T I C L E I N F O
A B S T R A C T
Article history: Received 28 August 2008 Received in revised form 1 December 2008 Accepted 23 December 2008 Available online 31 December 2008
The formula for surface energy was modified in accordance with the slab model of molecular dynamics (MDs) simulations, and MD simulations were performed to investigate the relaxed structure and surface energy of perfect and pit rutile TiO2(1 1 0). Simulation results indicate that the slab with a surface more than four layers away from the fixed layer expresses well the surface characteristics of rutile TiO2 (1 1 0) surface; and the surface energy of perfect rutile TiO2 (1 1 0) surface converges to 1:801 0:001 J m2. The study on perfect and pit slab models proves the effectiveness of the modified formula for surface energy. Moreover, the surface energy of pit surface is higher than that of perfect surface and exhibits an upper-concave parabolic increase and a step-like increase with increasing the number of units deleted along [0 0 1] and [1 1 0], respectively. Therefore, in order to obtain a higher surface energy, the direction along which atoms are cut out should be chosen in accordance with the pit sizes: [1¯ 1 0] direction for a small pit size and [0 0 1] direction for a big pit size; or alternatively the odd units of atoms along [1 1 0] direction are removed. ß 2008 Elsevier B.V. All rights reserved.
PACS: 31.15.xv 68.35.Md 68.35.bd 61.72.Ff Keywords: Molecular dynamics TiO2 Surface energy Surface structure Pit
1. Introduction Ti-materials have extensive biomedical applications [1], especially in medical implant and prosthesis [2], for their good biocompatibility, biological responses, osseointegration, excellent mechanical properties and corrosion resistance. In fact, there is a layer of TiO2on Ti and Ti-alloy surfaces, which has a strong direct effect on the success of implant [3,4]. The surface topographic characters of an implant at the micro/ nano-scale are related to the biological response of a host. The roughness of more than 10 mm of micro-rough surface will has an influence on the biocompatibility, mechanical characteristics and mechanical interlocking effect between implant and tissue. Due to the same size as a cell or a biomacromolecule, the micro-rough surface with a roughness ranging from 10 nm to 10 mm has remarkable influences on the biocompatibility and weakens the mechanical characteristics of the interface [5]. The micro-roughness of less than 10 nm of the surface has a significant influence on the interface structures of the implant. The reason is that the defects in a crystal structure, such as vacancy, grain boundary, step
* Corresponding author. Tel.: +86 451 86413840; fax: +86 451 86415244. E-mail address: [email protected] (D.-P. Song). 0169-4332/$ – see front matter ß 2008 Elsevier B.V. All rights reserved. doi:10.1016/j.apsusc.2008.12.062
and pit, are in this size range, which are active region for adsorption and influence further the integration of biomolecular and implant surface [6]. Moreover, the surface energy of biomaterials significantly influences the adhesion, spreading and growth of a cell. The surface with a higher surface energy promotes adhesion and spreading of a cell [7,8]. For example, good spreading only occurs when surface energy is higher than approximately 0.057 J m2[9]. Baier et al. [10,11] also showed the relationship between critical surface energy and cell adhesion: materials with a low surface energy show a low cell attachment. The rutile (1 1 0) surface is a very stable crystal face. Many scholars [12–17] have researched the surface structure of rutile (1 1 0) by ab initio calculations, such as density functional theory (DFT) local density approximation (LDA) [12,13], generalized gradient approximation (GGA) [14], full-potential linear augmented plane wave (FP-LAPW) [15,16] and linear combination of atomic orbitals (LCAO) [17]. On the other hand, Matsui and Akaogi [18] and Kim et al. [19] have developed the molecular dynamics (MDs) potential energy function and its parameters for TiO2 since 1991. A few scholars [20,21] studied the relaxed structure of the TiO2 (1 1 0) surface by MDs. However, there was a discrepancy between their theoretical [12–17,20–23] and experimental [24] results.
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5703
The surface energy of the rutile TiO2 (1 1 0) surface has been calculated recently by first principles [12–14,25,26], but there is few systematic calculations by MD possibly because the calculation accuracy of MD simulation is lower than that of first principles for small-scale models. However, the larger amount of computational work and the smaller atom amount of first principles make it very difficult to calculate the large-scale surface structures and evaluate the changes in surface energies resulting from the difference of the large-scale complex surface structures. Therefore, we modify the formula for surface energy in accordance with the slab model of MD simulations. Then, we study the relaxed structures of the rutile TiO2(1 1 0) surface by MD and investigate the relationship between surface energy and slab thickness. Finally, we also cut out some atoms to create pits, calculate the surface energies and relaxed structures of pit surfaces and study the influence of pit size on the surface energy. Fig. 1. Rutile TiO2(1 1 0) surface (1 1): O, big red ball; Ti, small grey ball. Atom labels follow as bridging oxygen O1, six-coordinate titanium Ti2, five-coordinate titanium Ti3, three-coordinate in-plane oxygen O4, O5, out-of-plane oxygen O6, and reference atom O7.
2. Methods and surface modeling 2.1. Molecular dynamics simulation method A survey of the literature indicated that several force fields have been published for titanium dioxide and titanium oxide [18,19,27– 32]. A detailed analysis of the different available force fields has been published by Smith and his co-worker [31]. They concluded that the force field of Matsui and Akaogi [18] is the most suitable of the available force fields for use in classical molecular dynamics simulation. Matsui and Akaogi [18] developed the following force fields for TiO2 in 1991. Vij ¼
qi q j C i C j Ai þ A j r i j 6 þ f Bi þ B j exp ri j Bi þ B j ri j
(1)
where the terms represent Coulomb, dispersion, and repulsion interactions, respectively. Here r i j is the interatomic distance between atoms i and j, and qi , Ai , Bi , and C i are the effective charge, repulsive radius, softness parameter, and van der Waals coefficient of atom i, respectively. Ai and Bi are relative to atom size and compressibility. The quantity f is a standard force of 4.184 kJ A˚1 mol1. For many oxides, Coulomb and repulsion interactions are necessary. Consequently, in this paper we used the reduced form for Eq. (1), i.e. Buckingham potential together with Coulombic interactions from partial charges ! r i j C i j qi q j V i j ¼ Ai j exp (2) 6 þ ri j ri j ri j where V i j is the interaction energy between atoms i and j, and the parameters Ai j , ri j , C i j , qTi and qO are shown in Table 1[18]. MD simulations were carried out at 310 K for further study on the interaction between TiO2 surface and biomolecules. Simulations were carried out in the canonical ensemble (i.e. NVT ensemble), and the thermostat derived by Hoover [34] was employed to maintain a constant temperature. The velocity Verlet algorithm was used to calculate the atomic motions and the particle–particle particle–mesh (PPPM) [35] solver was applied for the calculation of electrostatic interactions. A 10.0 A˚ cutoff distance was used for van der Waals interactions. Table 1 Interaction parameters for TiO2.
2.2. Perfect surface models The most stable (1 1 0) surface of rutile TiO2 was investigated in this p ffiffiffi paper. The (1 1 0) surface unit cell (1 1) has a dimension of 2a c which corresponds to 6.497 and 2.959 A˚ along [1¯ 1 0] and [0 0 1] directions, respectively [36,37]. The surface contains fiveand six-coordinate Ti atoms and two types of O atoms (in-plane and bridging). The six-coordinate Ti atoms are covered by the outermost bridging oxygen on the surface, while the fivecoordinate Ti atoms are coordinated to in-plane oxygen atoms on the surface. One Ti–O layer of rutile TiO2(1 1 0) surface unit cell (1 1) is built up of symmetrical three-plane O–Ti2O2–O unit along [1 1 0] direction, containing 6 atoms (2 Ti and 4 O), so that the sequence of Ti–O layers is O–Ti2O2–O–O–Ti2O2–O , as shown in Fig. 1. We constructed (10 18) surface slabs with the thickness of n Ti–O layers based on the rutile (1 1 0) surface unit cell (1 1), as shown in Table 2. For example, when n ¼ 13, this corresponds to a ˚ containing 14040 dimension of 64:969 A˚ 53:262 A˚ 42:229 A, atoms. Periodic boundary conditions were applied along [1¯ 1 0] and [0 0 1] directions of the surface, with a periodic vacuum gap of threefold slab thickness (3 3:2484n A˚) in the direction normal to the (1 1 0) surface. The simulations were carried out for 50 ps with time step of 0.5 fs, and the simulation results of last 30 ps were analyzed. 2.3. Pit surface models Based on the perfect TiO2(1 1 0) surface model, we cut out some atoms from the relaxed surface to create pit, as shown in Table 3. Fig. 2 shows the unrelaxed pit surface model (PitL16B2S1x (5 6 1)) in which the deleted atoms range over (5 6 1). For the pit surface models, due to the longer time to equilibrate the Table 2 Perfect surface slab models. Model set ID
Number of total layers n
Number of fixed layers nfixed 1 for n ¼ 2; 2 for n ¼3–6 1 for n ¼2–5; floor ðn=2Þ for n ¼6–16 1
Interation
Ai j (kcal mol1)
ri j (A˚)
C i j (kcal mol1 A˚6)
Perfect I Perfect II
2–16 2–16
Ti–Ti Ti–O O–O
717653.9571 391052.7442 271718.8311
0.154 0.194 0.234
120.9967 290.3920 696.9407
Perfect III
2–17
Notes: Atomic charges, qTi ¼ þ2:196, qO ¼ 1:098.
Notes: Type of fixed layers: Perfect I and II are one side fixed, and Perfect III are the central ceil(n=2)th layer fixed.
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5704 Table 3 Pit surface slab models. Model set ID a
PitL16B2S1x PitL16B2S1ya PitL16B2S1za PitL17C1S1xb PitL17C1S1yb PitL17C1S1zb PitL17C1S2xc PitL17C1S2yc PitL17C1S2zc
Pit size (½1¯ 1 0 ½0 0 1 ½1 1 0)
N delunit
N delunit 6 1 5 N delunit 1 5 6 N delunit N delunit 6 1 5 N delunit 1 5 6 N delunit N delunit 6 1 5 N delunit 1 5 6 N delunit
1–7 2–14 1–6 1–7 2–14 1–4 1–7 2–14 1–4
Notes:a is based on the 16-layer slab with the bottom 2 layers fixed and pit locates on the top surface.b,c is based on the 17-layer slab with the central 9th layer fixed.bPit locates on the top surface.cPits locate on the top and bottom surfaces.
surface atoms, the simulations were carried out for 200 ps, and the results of last 50 ps were analyzed. The other conditions of MD simulations for the pit surfaces are the same as those for the perfect surfaces given above. 2.4. Surface energy calculations According to the Gibbs definition [38], the surface energy density (s ) of a solid, which corresponds to the energy variation (per unit area) due to the creation of a surface, is given by
s¼
Enslab nEbulk Aslab
(3)
where Aslab is the total area of the surface considered, Enslab is the total energy of a n-layer slab and Ebulk is the bulk energy per layer of an infinite solid. However, the surface energy calculated in Eq. (3) is very sensitive to the accuracy in determination of the bulk energy term [39,40]. A small error in that term can make the calculated surface energies diverge linearly with increasing slab thickness. To avoid the divergence problem, we have employed the method [26,41,42] that makes use of Eq. (3) rewritten in the following form: Enslab ¼ Aslab s þ nEbulk
(4)
which implies that the bulk energy can be extracted from the slope of a linear fit of the slab_s total energy plotted versus n. This value is subsequently used in Eq. (3). For molecular dynamics simulations, a slab model, especially one with defects, usually needs to contain fixed layers, otherwise the slab structure will be distorted. To calculate the effect of the fixed layers on Enslab and surface energy, we simulated the following
perfect surface slab models which will not be distorted: (i) a sufficiently large number of layers are unconstrained and other layers are fixed, totally m layers; (ii) all m layers are unconstrained. The effect (DEfixed ) of fixed layers can be estimated from the difference between the total energies of the former and the latter. DEfixed is relative to the type and number of fixed layers (Tables 2 and 3). The effect of fixed layers on Enslab and surface energy can then be avoided by subtracting DEfixed from Enslab . Therefore, the surface energy of slab model containing fixed layers can be given by the following modified equation:
s¼
Enslab nEbulk DEfixed Aslab
(5)
For a large-scale complex surface structure, Eq. (5) can be rewritten as
s¼
unit EN slab NEbulk DEfixed Aslab
(6)
where N is the total number of units in slab. For our slab models, the numerators in Eqs. (3), (5) and (6) are the energy variation due to the creation of top and bottom surfaces, denoted as DE. Because the structure and area of the top surface is the same as those of the bottom surface in the perfect surface slab, the contribution of each perfect surface to DE is 50% of the energy variation (DEPerfM ) of the slab model with two-sided perfect surfaces. However, for the slab model with a perfect surface on one side and a pit surface on the other side, the contribution of the pit surface to DE is not 50% of the energy variation (DEHybridM ) of this slab model and the surface energy of the pit surface cannot be directly obtained using Eqs. (5) and (6). To calculate the surface energy of pit surface and compare DE and s with those of other slab models, we removed the contribution of the perfect surface to DE by transforming the slab model with one-sided pit (PitL16B2S1x, PitL17C1S1x, PitL16B2S1y, PitL17C1S1y, PitL16B2S1z and PitL17C1S1z) into the slab model with twosided pits (PitL16B2S1TS2x, PitL17C1S1TS2x, PitL16B2S1TS2y, PitL17C1S1TS2y, PitL16B2S1TS2z and PitL17C1S1TS2z, respectively). 2ðDEHybridM DEPerfM =2Þ and twice the surface area of the pit surface are used as the energy variation DE and the surface area Aslab of the transformed slab model, respectively. So, the surface energy of the pit surface of slab model with one-sided pit can be calculated using Eq. (6). The surface energy of the pit surface is the surface energy of the slab model with two-sided pits or the transformed slab model with two-sided pits. 3. Results and discussion 3.1. Analysis on perfect TiO2(1 1 0) surface structure
Fig. 2. Unrelaxed pit surface model PitL16B2S1x (5 6 1): O, big red ball; Ti, small grey ball. The box with green lines represents the pit.
The perfect surface unit cell (1 1) of the rutile TiO2 (1 1 0) contains bridging oxygen O1, six-coordinate titanium Ti2, fivecoordinate titanium Ti3, three-coordinate in-plane oxygen O4 and O5, and out-of-plane oxygen O6, as shown in Fig. 1. Supposing that O7 are reference atoms, Fig. 3 shows the relationship between the average displacements of (1 1 0) surface unit cells (10 18) along [1 1 0] and the slab thickness, which is a function of slab thickness. Although the average displacements of O1, Ti2 and O6 change considerably, the displacements of all atoms reach certain values if n 6 for Perfect I and n 9 for Perfect III. We compared further the average displacements of surface atoms in slab models Perfect I containing 6–16 layers with those of other scholars’ works in Table 4. The average displacements of all surface atoms are qualitatively in agreement with the experimental finding of Charlton [24] and the previous calculations [12– 15,17,20–23]. Moreover, our calculations coincide quantitatively
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5705
Moreover the average displacements of the atoms of (1 1 0) surface are less than 6e 4 A˚ along [0 0 1]. Thus it may seem that all atoms oscillate slightly in bulk positions along [1¯ 1 0] and [0 0 1] except the displacements of in-plane oxygen atoms along [1¯ 1 0]. 3.2. Surface energy of perfect TiO2(1 1 0) surfaces
Fig. 3. Calculated displacements of perfect surfaces.
with the experimental work [24] for the surface atoms except O1 and Ti3. On the outermost surface layer, the bridging oxygen atoms O1 and the five-coordinate titanium atoms Ti3 relax inward, while the six-coordinate titanium atoms Ti2 and the neighboring three-coordinate oxygen atoms O4 and O5 move outward, which causes a rumpled appearance of the surface. Although the average displacements of O1 and Ti3 are different considerably from those of experimental work [24], the average displacements of O6, O4, O5, and Ti2 along [1 1 0] are in good agreement with those of experimental work [24]. Langer et al. [43] reported that the bridging oxygen displaces laterally into an asymmetric position and the Ti–O distances are 1.75 and 2.1 A˚. However, in the experimental work [24], the bridging oxygen locates in the center of the neighboring six-coordinate titanium along [0 0 1], and the Ti–O distance is reduced to 1:71 0:07 A˚ mainly by relaxed inward largely (0:27 A˚). In our calculations, the O1 also locates laterally in a symmetric position and the Ti2–O1 average distance is 1.845 A˚, which is close to 1.82 and 1.83 A˚ for a five-layer slab and a three-layer slab, respectively by GGA calculations [37] and also in good agreement with 1.85 A˚ by MD simulation [18]. In addition, the average displacements of the three-coordinate in-plane oxygen atoms O4 and O5 along [1¯ 1 0] are very close to those of some works [12,14,15,22,23], and the in-plane oxygen atoms displace toward the neighboring five-coordinate titanium by 0.044 A˚, while the rest atoms displace by less than 1e 4 A˚.
Using 1828:6568 kcal mol1for Eunit bulk , as determined from a linear fit (Eq. (4)) to the slab energies (the energy of a thinnest fourlayer slab), yields very well converged DE and s (Figs. 4(a) and 5(a), respectively). Thus, the above values of Eunit bulk is adopted in all further calculations of s . Our calculations of energy variation DE and surface energies s versus slab thickness n are shown in Figs. 4(a) and 5(a), respectively. For the one-sided relaxed perfect surface models with different fixed layers and the two-sided relaxed perfect surface models (Table 2), their surface energies are converged to a value of 1:801 0:001 J m2as the slab thickness increases, which is very close to 1.78 J m2 from the MD simulations of Oliver et al. [44]. Although the surface energies obtained by MD simulations [33,44] are more than those by ab initio calculations, such as 0.89 J m2 from the average of the 5- and 6-layer LDA calculations of Ramamoorthy et al. [12], 0.73 J m2from the 7-layer GGA calculations of Bates et al. [14], 0.64 J m2from the average of the 9- and 10-layer PW–GGA calculations of Bredow et al. [45], 0.57 J m2by GGA–PW91 and 0.47 J m2by GGA–PBE from the average of the 10- and 11-layer of Kiejna et al. [26], large-scale surface slabs with complex morphologic surfaces can be evaluated more easily by MD simulations than by ab initio calculations. In addition, the calculation temperature results partly in the difference between our calculations and first-principles calculations [12,14,26,45]. The theoretical results of our reference literature of first-principles calculations are strictly valid only at zero temperature, while our MD simulations were carried out at 310 K. The surface energies of one-sided relaxed models are much the same as those of two-sided relaxed models by MD simulations, which is different from the result obtained by ab initio calculations from Kiejna et al. [26] that the surface energies of one-sided relaxed models are much more than those of two-sided relaxed models. Moreover, the slab model with one-sided relaxation is usually used as the adsorption substrate [46,47]. Therefore we will continue our study on the surface structure and surface energy of one-sided relaxed slab model with surface micro/nano-morphology.
Table 4 Atomic displacements for the perfect rutile TiO2 (1 1 0) surface (A˚). Atoma
9–18b
Expc[24]
4d[20]
6e[21]
7f[15]
5g[12]
5h[14]
3i[22]
4j[23]
O1 Ti2 Ti3 O4, O5 O4, O5 (1¯ 1 0) O6 O7
0.032 0.134 0.085 0.055 0:044 0.038 0
0:27 0:08 0:12 0:05 0:16 0:05 0:05 0:05 0:16 0:08 0:03 0:08 0 0:08
0.03 0.22 0.19 0.10 – 0.00 0.04
0.11 0.02 0.26 0.0 – 0.0 0.0
0.02 0.23 0.17 0.03 0:05 0.02 0
0.06 0.13 0.17 0.12 0:04 0.07 0.02
0.02 0.23 0.11 0.18 0:05 0.03 0.03
0.09 0.09 0.12 0.11 0:05 0.05 –
0.04 0.19 0.14 0.15 0:05 – –
a b c d e f g h i j
Atom labels refer to Fig. 1. This paper: MD, average displacements of 6–16 layers in Perfect I, S:D: 0:001 A˚. Charlton: SXRD experimental results [24]. Yin: MD, four layers [20]. San Miguel et al.: MD, six layers [21]. Harrison et al.: LCAO, seven layers [15]. Ramamoorthy et al.: PW–PP–LDA, five layers [12]. Bates et al.: PW–GGA, five layers [14]. Lindan et al.: PW–PP–GGA, three layers [22]. Maria et al.: GGA–PBE, four layers [23].
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5706
Fig. 4. Energy variations DE of rutile TiO2 surfaces: (a) perfect (1 1 0) surface slabs, (b) surface slabs with pit by cutting out atoms along [1¯ 1 0], (c) surface slabs with pit by cutting out atoms along [0 0 1] and (d) surface slabs with pit by cutting out atoms along [1 1 0].
To determine the contribution of each atomic layer to the surface energy, we calculated the energetic contribution of any atomic layer to the surface energy using the formula in Refs. [48– 50]. We found that, for the perfect TiO2(1 1 0) surface, s is stored mainly in the first four layers, especially the two outer layers, with 86.4, 9.9, 1.9 and 1.0% of the energy, respectively. As shown in Figs. 3 and 5(a), the changes in the atomic displacements govern the changes in the surface energy as a function of slab thickness, which is in agreement with the results of Bates et al. [14]. Therefore, it can be concluded that the slab with a surface more than four layers away from the fixed layer can express well the surface characteristics of rutile TiO2(1 1 0) surface. 3.3. Surface energy of pit TiO2(1 1 0) surfaces In order to describe the atomic displacements in the pit surface, according to the formula for slip vector [51], we defined the atomic displacement vector as sa ¼
k 1X ðxab Xab Þ k b 6¼ a
(7)
where k is the number of the nearest neighbors to atom a, and xab and Xab are the vector differences of atoms a and b at current and reference positions, respectively. The reference configuration is the arrangement of atomic positions with unrelaxation. The unstable atoms with the dangling bonds around the pit result in the maximum displacement vectors. Fig. 6 shows the
displacement vector of the relaxed pit surface model PitL16B2S1x (5 6 1) with atoms colored by jsa j, if jsa j 0:2. As shown in Figs. 4 and 5, due to the effect of the perfect surface on energy variation DE and surface energies s , DE and s of the slab models with one-sided pit are lower than those of the slab models with two-sided pits. For the pit surfaces, we also calculated the energetic contribution of any atomic layer to the surface energy. We found that s is stored mainly in the first two layers. For example, for the pit surface PitL16B2S1x (5 6 1), the energetic contribution of the first two layers is 97.9 and 1.8% of the energy, respectively. It can be seen from Figs. 4 and 5 that the energy variations and surface energies of the pit surfaces are higher than those of the perfect surfaces. The reason is that the pits created by removing some atoms from the perfect surfaces weaken the interaction of the atoms on the interface of the pits and many unstable atoms remain, and after relaxation these atoms offset their initial positions as shown in Fig. 6, which lead to a significant increase in the surface energy. Furthermore, according to the theory of solid surface physics [52], the closer the arrangement of atoms is, the stronger the interaction of atoms is and the lower the surface energy is, which is in agreement with our calculations. Figs. 4(b–d) and 5(b–d) show the relationship between DE and s of pit surfaces and the number (N delunit ) of units deleted along [1¯ 1 0], [0 0 1] and [1 1 0]. The energy variations and surface energies of the surfaces with pit by cutting out atoms along [1¯ 1 0] change slightly along a middle-convex parabola, and the changes in their surface energies are less than about 0.1 J m2. When the number of units deleted along [0 0 1] is 6 and that along [1 1 0] is 1,
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
5707
Fig. 5. Surface energies s of rutile TiO2 surfaces: (a) perfect (1 1 0) surface slabs, (b) surface slabs with pit by cutting out atoms along [1¯ 1 0], (c) surface slabs with pit by cutting out atoms along [0 0 1] and (d) surface slabs with pit by cutting out atoms along [1 1 0].
the quadratic regressions of N delunit values along [1¯ 1 0] with DE and s values of the pit surfaces are shown as the solid lines in Figs. 4(b) and 5(b), respectively. The energy variations and surface energies of the surfaces with pit by cutting out atoms along [0 0 1] exhibit an upper-concave parabolic increase with the increasing number of units deleted. When the number of units deleted along [1¯ 1 0] is 5 and that along [1 1 0] is 1, the quadratic regressions of N delunit values along [0 0 1] with DE and s values of the pit surfaces are shown as the solid lines in Figs. 4(c) and 5(c), respectively. The energy variations DE of the surfaces with pit by cutting out atoms along [1 1 0] is linearly relative to the number of units deleted. When the number of units deleted along [1¯ 1 0] is 5 and that along [0 0 1] is 6, the linear regression of Ndelunit values along
[1 1 0] with DE values is shown as the solid line in Fig. 4(d). However, as their surface areas increases significantly, their surface energies s show a step-like increase with the increasing number of units deleted. When N delunit is odd, s increases significantly. The means of s at Ndelunit are shown as the solid line in Fig. 5(d). If the number of deleted atoms is more than 180 atoms as shown in Fig. 5, the surface energy of the pit surface from which the atoms are removed along [0 0 1] is higher than that along [1¯ 1 0]. In order to obtain a higher surface energy, the direction along which atoms are cut out should be chosen in accordance with the pit sizes: [1¯ 1 0] direction for a small pit size and [0 0 1] direction for a big pit size. We can also remove the odd units of atoms along [1 1 0] direction. Because the pit surface has higher surface energy and more adsorption sites than the perfect surface, one can design the surface structure of implant using Figs. 4 and 5, and the regression equations therein, and investigate further the effect of surface structure and surface energy on biomolecular adsorption on Ti-materials.
4. Conclusions We modified the formula used for calculation of surface energy in accordance with the slab model of MD simulation and studied the surface energy and relaxed structure of perfect and pit rutile TiO2(1 1 0) through MD simulation. The following are our findings:
Fig. 6. Displacement vector of pit surface PitL16B2S1x (5 6 1) with atoms colored by jsa j.
1. The outermost surface bridging oxygens and five-coordinate titaniums relax inward, whereas the six-coordinate titaniums and the in-plane oxygens move outward. The slab with a surface
5708
D.-P. Song et al. / Applied Surface Science 255 (2009) 5702–5708
more than four Ti–O layers away from the fixed layer expresses well the surface characteristics of rutile TiO2 (1 1 0) surface. To reduce the amount of the computational work and improve computational efficiency, it is feasible to use the slab model with one-sided relaxation as the surface model. The unstable atoms with the dangling bonds around the pit in the pit surface have the maximum displacement vectors. Larger pit size generally results in an increase in the region of atoms with large displacement vectors. 2. The study on perfect and pit slab models proves the effectiveness of the modified formula for surface energy. The Surface energy of the perfect rutile TiO2 (1 1 0) surface converges to 1:801 0:001 J m2. 3. The surface energy of the pit surface is higher than that of the perfect surface and exhibits an upper-concave parabolic increase and a step-like increase with increasing the number of units deleted along [0 0 1] and [1 1 0], respectively. The surface energy of the surface with pit by cutting out atoms along [0 0 1] or [1 1 0] changes more significantly than that along [1¯ 1 0]. Moreover, if the number of deleted atoms is over a certain value (180 atoms in this paper), the surface energy of the pit surface from which the atoms are cut out along [0 0 1] is higher than that along [1¯ 1 0]. Thus, in order to obtain a higher surface energy, the direction along which atoms are cut out should be chosen in accordance with the pit sizes: [1¯ 1 0] direction for a small pit size and [0 0 1] direction for a big pit size; or alternatively the odd units of atoms along [1 1 0] direction are removed. It can therefore be concluded that the modified formula for surface energy and the displacement vector can be used to calculate the surface energy of the large-scale complex structures and evaluate these surface structures by MD, respectively, which will be aid to analyze further the effect of surface micro/nanomorphology of the TiO2 covering titanium implant surfaces on cell adhesion and protein adsorption. Acknowledgments This work is funded by the National Natural Science Foundation of China (No. 50675050) and the Multidiscipline Scientific Research Foundation of Harbin Institute of Technology (No. HIT. MD 2003. 10). References [1] C. Leyens, M. Peters, Translated by Z.H. Chen et al. Titanium and Titanium Alloys, Chemical Industry Press, Beijing (2005). [2] U. Diebold, Surf. Sci. Rep. 48 (2003) 53. [3] A. Hunter, C.W. Archer, P.S. Walker, G.W. Blunn, Biomaterials 16 (1995) 287. [4] B.D. Boyan, T.W. Hummert, D.D. Den, Z. Schwartz, Biomaterials 17 (1996) 137. [5] B. Kasemo, J. Prosth. Dent. 49 (1983) 832. [6] B. Kasemo, J. Lausmaa, Material selection surface characteristics and chemical processes at implant surface, in: P.I. Branemark, G.A. Zarb, T. Albrektsson (Eds.),
[7] [8] [9] [10]
[11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52]
Tissue Integrated Prostheses: Osseointergtation in Clinical Dentistry, Quintessence Publishing Co. Inc., Chicago, 1985, p. 99. T.G. van Kooten, J.M. Schakenraad, H.C. van der Mei, H.J. Busscher, Biomaterials 13 (1992) 897. S.A. Redey, S. Razzouk, C. Rey, D. Bernache-Assollant, G. Leroy, M. Nardin, G. Cournout, J. Biomed. Mater. Res. 45 (1999) 140. J.M. Schakenraad, H.J. Busscher, C.R.H. Wildevuur, J. Arends, J. Biomed. Mater. Res. 20 (1986) 773. R.E. Baier, Surface properties influencing biological adhesion, in: R.S. Mansly (Ed.), Adhesion in Biological Systems, Academic Press, New York, 1970 , pp. 15–48. R.E. Baier, A.E. Meyer, J.R. Natiella, R.R. Natiella, J.M. Carter, J. Biomed. Mater. Res. 18 (1984) 327. M. Ramamoorthy, D. Vanderbilt, R.D. King-Smith, Phys. Rev. B 49 (1994) 16721. M. Ramamoorthy, R.D. King-Smith, D. Vanderbilt, Phys. Rev. B 49 (1994) 7709. S.P. Bates, G. Kresse, M.J. Gillan, Surf. Sci. 385 (1997) 386. N.M. Harrison, X.G. Wang, J. Muscat, M. Schefflet, Faraday Discuss. Chem. Soc. 144 (1999) 305. D. Vogtenhuber, R. Podlouchy, A. Nechel, S.G. Steinemann, A.J. Freeman, Phys. Rev. B 49 (1994) 2099. P. Reinhardt, B.A. Heb, Phys. Rev. B 50 (1994) 12015. M. Matsui, M. Akaogi, Mol. Simul. 6 (1991) 239. D.W. Kim, N. Enomoto, Z. Nakagawa, K. Kawamura, J. Am. Ceram. Soc. 79 (1996) 1095. X. Yin, R. Miura, A. Endou, I. Gunji, R. Yamauchi, M. Kubo, A. Stirling, A. Fahmi, A. Miyamoto, Appl. Surf. Sci. 119 (1997) 199. M.A. San Miguel, C.J. Calzado, J.F. Sanz, Surf. Sci. 409 (1998) 92. P.J.D. Lindan, N.M. Harrison, M.J. Gillan, J.A. While, Phys. Rev. B 55 (1997) 15919. L.S. Maria, Y.G. Andrey, L.S. Alexander, J. Phys. Chem. B 110 (2006) 4853. G. Charlton, Phys. Rev. Lett. 78 (1997) 495. J. Goniakowski, J.M. Holender, L.N. Kantorovich, M.J. Gillan, Phys. Rev. Lett. 53 (1996) 957. A. Kiejna, T. Pabisiak, S.W. Gao, J. Phys.: Condens. Matter 18 (2006) 4207. M. Mostoller, J.C. Wang, Phys. Rev. B 32 (1985) 6773. V. Swamy, J.D. Gale, Phys. Rev. B 62 (2000) 5406. S. Ogata, H. Iyetomi, K. Tsuruta, F. Shimojo, R.K. Kalia, A. Nakano, P. Vashishta, J. Appl. Phys. 86 (1999) 3036. H. le Roux, L. Glasser, J. Mater. Chem. 7 (1997) 843. D.R. Collins, W. Smith, Council for the Central Laboratory of Research Councils (1996). F.H. Streitz, J.W. Mintmire, J. Adhes. Sci. Technol. 8 (1994) 853. P.K. Naicker, P.T. Cummings, H. Zhang, J.F. Banfield, Phys. Phys. Chem. B 109 (2005) 15243. W.G. Hoover, Phys. Rev. A 31 (1985) 1695. H. Eastwood, Computer Simulation Using Particles, Adam Hilger, NY, 1989. S.C. Abrahams, J.L. Bernstein, J. Chem. Phys. 55 (1971) 3206. A.V. Bandura, J.D. Kubicki, J. Phys. Chem. B 107 (2003) 11072. J.W. Gibbs, Scientific Papers, Thermodynamics, vol. 1, Dover Publictions, New York, 1926, pp. 55. J.C. Boettger, Phys. Rev. B 49 (1994) 16798. V. Fiorentini, M. Methfessel, J. Phys.: Condens. Matter 8 (1996) 6525. J.C. Boettger, Phys. Rev. B 53 (1996) 13133. J.G. Gay, J.R. Smith, R. Richter, F.J. Arlinghaus, R.H. Wagoner, J. Vac. Sci. Technol. A 2 (1984) 931. W. Langel, Surf. Sci. 496 (2002) 141. P.M. Oliver, G.W. Watson, E.T. Kelsey, S.C. Parker, J. Mater. Chem. 7 (1997) 563. T. Bredow, L. Giordano, F. Cinquini, G. Pacchioni, Phys. Rev. B 70 (2004) 035419. J.R.B. Gomes, J.P. Prates Ramalho, Phys. Rev. B 71 (2005) 235421. M. Prˇedota, A.V. Bandura, P.T. Cummings, J.D. Kubicki, D.J. Wesolowski, A.A. Chialvo, M.L. Machesky, J. Phys. Chem. B 108 (2004) 12049. O. Politano, S. Garruchet, J.M. Salazar, Mat. Sci. Eng. A 387–389 (2004) 749. S. Garruchet, O. Politano, J.M. Salazar, T. Montesin, Surf. Sci. 586 (2005) 15. S. Garruchet, O. Politano, J.M. Salazar, A. Hasnaoui, T. Montesin, Appl. Surf. Sci. 252 (2006) 5384. J.A. Zimmerman, C.L. Kelchner, P.A. Klein, J.C. Hamilton, S.M. Foiles, Phys. Rev. Lett. 87 (2001) 165507. L.-S. Wen, Physics Basis of Solid Materials Interface Research, Science Press, Beijing, 1991, p. 21.