ARTICLE pubs.acs.org/JPCC
Lattice Stability of Si[100] Wires From First Principles Wai-Leung Yim*,† and Thorsten Kl€uner‡ † ‡
Institute of High Performance Computing, 1 Fusionopolis Way, #16-16 Connexis, Singapore 138632 Institut f€ur Reine und Angewandte Chemie, Theoretische Chemie, Carl von Ossietzky Universit€at Oldenburg, 26129 Oldenburg, Germany
bS Supporting Information ABSTRACT: Quantum many-body effects such as electronphonon coupling are generally enhanced in reduced dimensions. However, the influence of many-body effects on the structural stability is often underestimated. We used ab initio phonon and molecular dynamics calculations to illustrate the strong electron-phonon coupling in metallic Si[100] wires. We found that the smallest crystalline Si[100] wire is a sharpcorner Si25 wire. Electron-phonon coupling favors a tilted orientation of Si-Si bond pairs at the walls of the wires. Both sharp- and cut-corner Si[100] wires are thermally stable. A partial charge density near the Fermi energy was calculated to reveal the chemical reactivity of different Si sites. Surprisingly, the most reactive sites of cut-corner Si[100] were located at the subsurface atoms at the groove sites adjacent to the corners of the walls. Our work predicts the crystallinity of nanowire structures which will find further applications particularly for the systems sensitive to crystalline-to-amorphous transition.
’ INTRODUCTION Physical properties of silicon nanowires (SiNWs) have been studied intensively.1 Doped or binary nanowires of smaller diameters are of great interest because of their expectedly high transconductance when compared to the performance of traditional materials.2 The selective growth and purification processes of silicon wires and the band gap engineering are the keys to manipulate SiNW devices.3 Moreover, lattice stability/instability and crystallinity which govern the durability of SiNWs are the focus of this work. SiNWs were synthesized in various ways, for example, by using nanocluster-catalyzed reactions in the presence of SiH4 or Si2H6. SiNWs with different diameters and orientations were formed, in which the most probable orientation was found to be Si[110].4 Recently, mesoporous and branched silicon nanowires can be synthesized, and the crystallinity of the silicon networks was also examined.5,6 In addition to the well-ordered crystalline phase of SiNWs, Al-doped p-type SiNWs were obtained by catalyzing the growth on Si(111) substrates using Al as catalyst. In such a case, the growth was thought to involve liquid Al, which may involve crystalline and noncrystalline phases of Al-doped SiNWs.7 Crystallinity is also an important issue for tin nanowires which are amorphous at small diameters and for the durability of nanowires resisting structural distortions upon repeated mechanical loading.8 The understanding of crystallinity and lattice stability of one-dimensional Group IV nanostructures is lacking. The computational approach allows us to examine lattice crystallinity. The lattice instability mechanism is usually driven by electronic instability that cannot be examined by common experimental techniques. The absence of sharp-corner atoms of Si[100] wires was rationalized from first principles, but some problems are yet to be answered, r 2011 American Chemical Society
notably the chemistry of the corner atoms and the crystallinity of ultrathin SiNWs.9 Originally, a cross section of Si[100] wire exhibited a square shape, in which the four atoms at the corners of the square are connected to the other three Si atoms underneath. The angle between two Si-Si bonds was bent to a smaller value which was energetically disfavored. To compare the stability of sharp-corner (sc) and cutcorner (cc) SiNWs, their respective formation energies were computed. The results implied that the sharp-corner atoms were not favored in Si[100] wires of larger diameters. However, the computation of formation energies can only give the information concerning the chemical stability but not chemical kinetics nor the chemical reactivity. Therefore, we have examined the lattice stability by ab initio phonon calculations and performed chemical bonding analyses on the Si[100] lattice. Quantum many-body effects have to be taken into account when describing low dimensional systems such as surfaces and wires.10 Bare Si[100] wires are metallic, and the electronphonon interaction is particularly strong. In a commonly used approach in computational materials science, a system being studied is simulated as a supercell. However, sometimes the simulations may be unreliable, because the computational setup may induce some unrealistic constraints. A typical example is the surface reconstruction of Si(100), where the correct structure will be obtained only with a physically correct model in which Si-Si dimers are formed.11 In general, this might be a very challenging task especially when the system is new and the available experimental data are inadequate or not existing at all. In this case, it is mandatory to carry out ab initio phonon calculations to examine potential structural instabilities. Received: November 3, 2010 Revised: December 30, 2010 Published: February 04, 2011 3286
dx.doi.org/10.1021/jp110500d | J. Phys. Chem. C 2011, 115, 3286–3290
The Journal of Physical Chemistry C
ARTICLE
Figure 1. Cross sections of Si[100] wires. cc-SiNWs: Si21 and Si57. sc-SiNWs: Si25 and Si61.
Electron-phonon coupling cannot be obtained from a standard DFT approach. By performing phonon calculations the structural (in-)stability can be obtained by monitoring the phonon band structure. Usually the imaginary phonon structure will reveal a structural distortion leading to a structure of reduced symmetry. In our case, the Si[100] wire would reduce from C4 to C2 symmetry. The computation of electron-phonon coupling parameters was, however, discouraged by the very big size of the silicon wire which we were investigating. Quite often the electron-phonon calculations have been reported in the literature for a system involving several atoms in a primitive cell only. Vibrational frequency calculations have been used extensively to characterize the optimized molecular structures. The total energy of the system can be expressed as: E(r) = E0(r) þ V0 (r)Δr þ ΔrTH(r)Δr þ higher-order terms, where V0 is the energy gradient and H is the Hessian matrix. The Hessian matrix of a stationary structure is positive definite; i.e., all force constants (k) are positive. When there is one or more negative force constants, the structure is at a saddle point and is subject to instabilities which result in structural transformations. In analogy to the dispersion of electronic states in crystalline structures, quantized vibrational energy levels do also exhibit dispersion relations which lead to phonon dispersion structures. A structurally unstable compound exhibits an imaginary frequency branch in its phonon dispersion spectrum. In this work, we used an interdisciplinary approach—ab initio phonon and molecular dynamics simulations—to examine the structural stability/instability of bare Si[100] lattices.12 Such an approach has been used to study chemistry under high pressure,13 and it provides key quantities for studying ab initio thermodynamics.14 The electronic structure is also illustrated by means of electronic band structure calculations and partial charge density analysis. Our method is not limited to Si[100] wires but also applicable to other one-dimensional systems like sandwiched molecular wires and other Group IV nanowires.15 In this report, we will be able to answer: (1) What is the smallest size of a crystalline Si[100] nanowire that can exist by theory? (2) What is the source of potential structural instabilities? (3) What are the reactive sites of bare Si[100] nanowires (corner or wall atoms)?
’ COMPUTATIONAL DETAILS We used the Vienna Ab Initio Simulation Package (VASP) to perform calculations within density functional theory.12,16,17 Throughout the study, the PBE exchange-correlation functional was used.18 Projected augmented wave (PAW) pseudopotentials were adopted for Si atoms.19 The plane wave energy cutoff and the augmentation
charge cutoff were set to 307 and 322 eV, respectively. Geometry optimizations were carried out using the conjugate gradient minimization scheme in VASP. The convergence thresholds for electronic structure calculations and geometry optimizations were set to 1 10-8 eV and 1 10-4 eV/Å, respectively. We used a 1 1 12 Monkhorst-Pack (MP) scheme for k-point sampling in geometry optimizations. The axis of SiNW was aligned along the c-axis of the tetragonal unit cell. The separation between the closest images of SiNWs was larger than 10 Å to avoid spurious interactions. Lattice dynamics were carried out by using the direct force constant approach.20 A supercell containing three primitive unit cells was used, and a 1 1 4 MP mesh was selected. The atomic displacement was set to 0.015 Å for nonequivalent Si atoms. Furthermore, the structural stability/instability was tested by AIMD. Supercells containing three primitive units were chosen. Plane wave and charge density cutoffs were set to 184 and 322 eV, respectively. The AIMD was set up as a canonical ensemble which was controlled by a Nose-Hoover thermostat.21,22 The stability of SiNWs was monitored for 10 ps for smaller wires and for 4 ps for larger wires. The time interval was set to 1 fs.
’ RESULTS AND DISCUSSIONS Sharp- and cut-corner Si[100] wires (sc/cc) of different sizes have been taken into account (cf. Figure 1 and Supporting Information). Figure 2 shows the phonon dispersion relations of selected sc- and cc-wires. It was proposed by Cao et al. that sc-Si25 was preferred as compared to corresponding cc-wires.9 In this work we found that the cc-Si21 wire lattice is indeed unstable and therefore disordered. In analogy to the vibrational frequencies in the gas phase, an unstable structure would be identified by imaginary vibrational frequency modes. The phonon spectrum of an ordered cc-Si21 wire clearly contains imaginary frequency bands down to 200i cm-1 (Figure 2a). The structural instability is due to the presence of unsaturated surface atoms which were under high strain. Most of the outer atoms have a coordination number of 3. However, some Si-Si-Si linkages have been bent seriously, and the Si-SiSi angle was as large as 146.7°. Furthermore, our ab initio molecular dynamics (AIMD) calculations confirmed that cc-Si21 does not exhibit an ordered structure, and the structural skeleton became seriously distorted in the first 2 ps of the molecular dynamics simulation (cf. Figure 3a). Some of the strained subsurface Si atoms have changed their coordination number from 4 to 2, resulting in bicoordinated unsaturated Si atoms. Atomic coordinates of cc-Si21 wire are provided in the Supporting Information. 3287
dx.doi.org/10.1021/jp110500d |J. Phys. Chem. C 2011, 115, 3286–3290
The Journal of Physical Chemistry C
ARTICLE
Figure 2. Phonon dispersion relation of selected SiNWs. (a) cc-Si21; (b) sc-Si25; (c) cc-Si57; (d) sc-Si61.
Figure 3. Ab initio molecular dynamics simulations on selected SiNWs. (a) cc-Si21; (b) sc-Si25; (c) cc-Si57; (d) sc-Si61.
In sharp contrast, sc-Si25 was found to be stable from our ab initio phonon calculations. The phonon spectrum (cf. Figure 2b) of sc-Si25 was free of imaginary modes revealing a stable structure. The structural stability has further been demonstrated by AIMD simulations. The surfaces of sc-Si25 were stabilized by a “benign” reconstruction in which the number of unsaturated Si atoms was reduced. The Si atoms oscillate around their equilibrium positions and do not show permanent distortions in the molecular dynamics simulation. Remarkably, the smallest Si[100] wire which contains sharp corner atoms is sc-Si25. The diameter of sc-Si25 was 1.08 nm, which was comparable to the smallest size of H-terminated Si[120] wire with its diameter at 1.3 nm.23 The standard Si-H bond length is 1.47 Å. By subtracting the two Si-H bond lengths from the diameter of H-terminated Si[120], the diameter of bare Si[120] wire can be estimated roughly to be 1.06 nm. It is very
interesting that the smallest size of the silicon diameter is not sensitive to the orientation of the wire. Figure 2c,d shows the lattice dynamics results for cc-Si57 and sc-Si61, respectively. cc-Si57 exhibits real phonon dispersion relations, except for a small phonon softening near the Γ point implying a long-range structural instability. The cc-Si57 wire exhibits a C4 symmetry along the axis of the wire. However, the phonon dispersion relations of the sc-Si61 wire contain an imaginary frequency band around 100i cm-1 which has a subtle effect on its structural stability. Originally, sc-Si61 was optimized, and its structure could be classified by the C4 point group. The visualization of the imaginary frequency mode showed that this mode originates from the dangling motion of Si atoms at the wall. Very importantly, sharp corner atoms have no contribution to the atomic displacements of the imaginary modes. AIMD simulations showed that wall atoms 3288
dx.doi.org/10.1021/jp110500d |J. Phys. Chem. C 2011, 115, 3286–3290
The Journal of Physical Chemistry C
ARTICLE
Figure 4. Partial charge density of sc-Si25, cc-Si57, and sc-Si61. Partial density was integrated from (εf - 0.2) to εf, where εf is Fermi energy. Isosurface thresholds were set to 0.0025 and 0.016 for the left and right panels, respectively.
become tilted (similar to Si-Si dimers on the reconstructed Si(100) surface). Thus, the reactivity of the Si atoms is governed by the overall symmetry in sc-Si61, where the wire of C2 symmetry is more stable than that of C4 symmetry. Nevertheless, the skeleton of sc-Si61 wire remains intact for 10 ps. So, we conclude that both ccSi57 and sc-Si61 wires are thermally stable. In Figure 4, we displayed the partial charge density of SiNWs close to the Fermi level. In detail, the partial density was integrated from (εf - 0.2 eV) to εf, where εf is the Fermi energy of SiNWs. This illustrates the chemical reactivity of different SiNWs. We have used two isosurface values to plot the partial charge density, one at 0.0025 (left panels) and the other at 0.016 (right panels). In sc-Si25, ccSi57, and sc-Si61, the cores of the wires contribute to the partial density (left panels of each wire in Figure 4) because of quantum confinement effects existing in these wires and also because of their small diameters. For larger SiNWs, the main contribution to the partial charge density near the Fermi level was predominantly located at the outer surface atoms (see Supporting Information). Figure 4 highlights the most reactive sites on cc-Si57 and scSi61. Interestingly, among all cc-SiNWs, the most reactive sites were not the outermost corner atoms but atoms which are located at the groove sites adjacent to the corners. This result reveals interesting insight into the growth mode of the nanowires. As the cc-SiNWs grow laterally, the next incoming Si atom would react with these wider groove sites rather than the narrower sites at the corners. Moreover, there are eight wider groove sites versus four sites at the corners. Therefore, the next incoming Si atom would preferentially allocate to the wider groove sites of cc-SiNWs as facilitated by energetic and statistical reasons. The reactivity of sublayer Si atoms was reported in a recent theoretical study on the chemisorption of C2H4 and C2H6 on Si(100) surfaces.24 In contrast to sc-Si57, the situation for the sc-Si61 nanowires is altogether different. For the sc-Si61 SiNW, the C4 molecular
symmetry is reduced, and the partial density near the Fermi level has a substantial contribution from the wall atoms. The two adjacent atoms at the walls were spatially allowed to interact with each other. As a result, the Si-Si bond pair would be tilted. Among sc-SiNWs, the most reactive sites are the corner atoms. Larger nanowires such as cc-Si109, sc-Si113, cc-Si177, and scSi181 have also been considered (cf. Figure S1 in Supporting Information). Depending on the width of the walls, cc-Si109 and scSi113 belong to the C2 point group, while cc-Si177 and sc-Si181 exhibit a C4 symmetry. AIMD calculations were performed on cc-Si109 and scSi113 (cf. Figure S2). We observed that the structural skeleton keeps intact, and the wall Si-Si dimers prefer a tilted orientation.
’ CONCLUSIONS In summary, the structural instability of several Si[100] nanowires has been rationalized. The Si-Si bond pairs at the walls prefer a tilted orientation, which was attributed to a strong electron-phonon coupling particularly in these metallic Si[100] wires. Both sc- and cc-SiNWs are thermally stable. We analyzed the partial charge density near the Fermi level and identified the most reactive sites of the SiNWs. Among cc-SiNWs, the most reactive sites are located along the groove sites next to the corners of the wires with important implications on the growth properties of these wires. Apparently, incoming Si atoms will be attached to these reactive sites as favored by both energetic and statistical reasons. The reactivity of the subsurface Si atom is supported by the recent theoretical study on the chemisorption of small hydrocarbons on two-dimensional Si(100) surfaces. The smallest bare Si[100] wire is found to be sc-Si25, with its diameter at 1.08 nm. Very interestingly, the size is comparable to the experimentally reported smallest H-terminated Si[120] wire. It is appealing that the size of the smallest Si nanowire is not sensitive to the choice of the growth orientation. Here, ab initio phonon calculations have successfully 3289
dx.doi.org/10.1021/jp110500d |J. Phys. Chem. C 2011, 115, 3286–3290
The Journal of Physical Chemistry C been used to study the lattice stability of pure Si[100] wires. The method will find more applications on nanowires of heavier Group IV elements like germanium and tin. It is especially important to study tin wires because it has been postulated as a component for lithium storage in which Sn-Li alloys are formed. The crystallinity of tin will be an interesting issue which determines the diffusion of lithium and also the overall Sn-Li structures. Chemical bonding analysis beautifully identified the reactive sites of both cc- and sc-SiNWs, which provided great insight into the growth mechanism. It should be noted that the full description of the growth mechanism is out of the scope of this study and will also be very specific depending on experimental conditions such as the choice of catalyst and substrates and the ratio of SiH4:H2. Nevertheless, the real-time growth mode is attainable, for example, by metadynamics method based on a density functional tight-binding approach.25,26
ARTICLE
(15) Shen, L.; Jin, H. M.; Ligatchev, V.; Yang, S. W.; Sullivan, M. B.; Feng, Y. P. Phys. Chem. Chem. Phys. 2010, 12, 4555. (16) Kresse, G.; Furthm€uller, J. Phys. Rev. B 1996, 54, 11169. (17) Kresse, G.; Furthm€uller, J. Comput. Mater. Sci. 1996, 6, 15. (18) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1997, 78, 1396. (19) Kresse, G.; Joubert, D. Phys. Rev. B 1999, 59, 1758. (20) Parlinski, K.; Li, Z. Q.; Kawazoe, Y. Phys. Rev. Lett. 1997, 78, 4063. (21) Nose, S. J. Chem. Phys. 1984, 81, 511. (22) Nose, S. Mol. Phys. 1984, 52, 255. (23) Ma, D. D. D.; Lee, C. S.; Au, F. C. K.; Tong, S. Y.; Lee, S. T. Science 2003, 299, 5614. (24) Zhang, Q. J.; Fan, X. L.; Lau, W. M.; Liu, Z. F. Phys. Rev. B 2009, 79, 195303. (25) Martonak, R.; Laio, A.; Parrinello, M. Phys. Rev. Lett. 2003, 90. (26) Aradi, B.; Hourahine, B.; Frauenheim, T. J. Phys. Chem. A 2007, 111, 5678.
’ ASSOCIATED CONTENT
bS
Supporting Information. Dynamical properties, electronic band structures, and partial charge density of SiNWs and distorted cc-Si21 wire in the end of AIMD trajectory are provided. This material is available free of charge via the Internet at http:// pubs.acs.org.
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
’ ACKNOWLEDGMENT Computations have been performed on the computers allocated at the A*STAR Computational Resource Centre (A*CRC) and on the national supercomputer NEC SX-8 at the High Performance Computer Center Stuttgart (HLRS) under the grant WLYIM. ’ REFERENCES (1) Cui, Y.; Lieber, C. M. Science 2001, 291, 851. (2) Hu, Y. J.; Xiang, J.; Liang, G. C.; Yan, H.; Lieber, C. M. Nano Lett. 2008, 8, 925. (3) Lu, W.; Xie, P.; Lieber, C. M. IEEE Trans. Electron Devices 2008, 55, 2859. (4) Il Park, W.; Zheng, G. F.; Jiang, X. C.; Tian, B. Z.; Lieber, C. M. Nano Lett. 2008, 8, 3004. (5) Doerk, G. S.; Ferralis, N.; Carraro, C.; Maboudian, R. J. Mater. Chem. 2008, 18, 5376. (6) Hochbaum, A. I.; Gargas, D.; Hwang, Y. J.; Yang, P. D. Nano Lett. 2009, 9, 3550. (7) Ke, Y.; Weng, X. J.; Redwing, J. M.; Eichfeld, C. M.; Swisher, T. R.; Mohney, S. E.; Habib, Y. M. Nano Lett. 2009, 9, 4494. (8) Lugstein, A.; Steinmair, M.; Steiger, A.; Kosina, H.; Bertagnolli, E. Nano Lett. 2010, 10, 3204. (9) Cao, J. X.; Gong, X. G.; Zhong, J. X.; Wu, R. Q. Phys. Rev. Lett. 2006, 97, 136105. (10) Kroger, J. Rep. Prog. Phys. 2006, 69, 899. (11) Healy, S. B.; Filippi, C.; Kratzer, P.; Penev, E.; Scheffler, M. Phys. Rev. Lett. 2001, 87, 016105. (12) Kresse, G.; Hafner, J. Phys. Rev. B 1994, 49, 14251. (13) Tse, J. S.; Li, Z. Q.; Uehara, K.; Ma, Y. M.; Ahuja, R. Phys. Rev. B 2004, 69, 132101. (14) Alapati, S. V.; Johnson, J. K.; Sholl, D. S. Phys. Rev. B 2007, 76, 104108. 3290
dx.doi.org/10.1021/jp110500d |J. Phys. Chem. C 2011, 115, 3286–3290