The Elastic Property of Bulk Silicon Nanomaterials through an Atomic Simulation Method

Abstract

This paper reports a systematic study on the elastic property of bulk silicon nanomaterials using the atomic finite element method. The Tersoff-Brenner potential is used to describe the interaction between silicon atoms, and the atomic finite element method is constructed in a computational scheme similar to the continuum finite element method. Young's modulus and Poisson ratio are calculated for [100], [110], and [111] silicon nanowires that are treated as three-dimensional structures. It is found that the nanowire possesses the lowest Young's modulus along the [100] direction, while the [110] nanowire has the highest value with the same radius. The bending deformation of [100] silicon nanowire is also modeled, and the bending stiffness is calculated.

Academic Editor:Bohayra Mortazavi

Department of Civil Engineering and Architecture, Zhongyuan University of Technology, Zhengzhou, China

Received 19 July 2016; Accepted 5 September 2016

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Over the past decade, one-dimensional nanostructures such as nanowires have been extensively studied through both theoretical and experimental methods. The most important one is the silicon nanowire (SiNW) that has been successfully applied in the nanoelectromechanical devices such as field effect transistors (FETs) [1-3] due to their unique optoelectronic and mechanical properties [4, 5]. Researchers have made significant contributions in the study of the elastic property of SiNWs, in which they are considered with different diameters and different axial or surface orientations. It is generally accepted that Young's modulus tends to increase with the increasing diameter [6, 7].

The molecular dynamics and density functional theory [6-9] are the common atomic scale simulating methods to study the elastic property of nanostructures. However, their computational cost is very huge and they are valid only for the small size structures. In order to overcome this shortcoming, Liu et al. [10, 11] proposed an order-N atomic finite element method (AFEM) which can achieve the same accuracy with molecular mechanics but is much faster than the latter one. Sun and Tao [12] adopted this method to investigative the elastic properties of carbon, boron nitride, and silicon carbide nanotubes, and Young's moduli are accurately obtained and the buckling behavior is better displayed. The present work extends its application to three-dimensional nanostructures in order to study the fine elastic property of SiNWs. An appropriate type of Tersoff-Brenner potential is employed to describe the atomic interaction. An atomic element contains 17 atoms for the bulk silicon crystal, and the global stiffness matrix and nonequilibrium force vector are assembled similar to the continuum finite element method (FEM). The nonlinear iteration is implemented to determine the equilibrium.

2. Atomic Scale Modeling Method for Bulk Silicon

Bulk silicon has the diamond lattice structure as illustrated in Figure 1, and its elastic property has a great dependence on their growth directions. The present work focuses mainly on SiNWs grown along the [100], [110], and [111] crystallographic orientations, and Figure 2 shows the representative cross section and side view of the SiNWs with different growth directions.

Figure 1: The diamond lattice structure of bulk silicon.

[figure omitted; refer to PDF]

Figure 2: Cross sections and side view for [100], [110], and [111] SiNWs.

(a) [ 100]

[figure omitted; refer to PDF]

(b) [ 110]

[figure omitted; refer to PDF]

(c) [ 111]

[figure omitted; refer to PDF]

In order to describe the Si-Si interaction, the Tersoff-Brenner potential is adopted with the parameter set refined by Erhart and Albe [13]. Tersoff-Brenner potential is a multibody potential which can be used to describe the interaction between atoms including C, B, H, Si, and N: [figure omitted; refer to PDF] where _VR and _VA are the repulsive pair potential and attractive pair potential, respectively; _rIJ is the distance from atom I to atom J; _B¯IJ is the bond order function. All above functions have analytic forms [figure omitted; refer to PDF] where _fc (r) is the cutoff function [figure omitted; refer to PDF] The angular function g(θ) is given by [figure omitted; refer to PDF] The potential parameters are listed in Table 1.

Table 1: Parameter sets for silicon.

D 0 (eV)	r 0 (Å)	S	β (Å-1 )	γ	c	d	h	2 μ	R (Å)	D (Å)
3.24	2.222	1.57	1.4760	0.09253	1.13681	0.63397	0.335	0.0	2.90	0.15

The basic idea of AFEM is to treat atoms as the continuum FEM nodes, and each element is characterized by a set of discrete atoms [10, 11]. For a system of N atoms, the energy stored in the atomic bonds can be evaluated using Tersoff-Brenner potential, and the atomic positions are determined by minimizing the system energy [figure omitted; refer to PDF] in which x=(_x1 ,_x2 ,...,_xn^)T , _xi is the position of atom i, and _Fi is the external force exerted on atom i.

Giving Taylor expansion of _Et (x) and substituting it into [figure omitted; refer to PDF] we have [figure omitted; refer to PDF] where Δu is displacement increment, K and P are, respectively, the stiffness matrix and nonequilibrium force vector given by [figure omitted; refer to PDF]

For the present nonlinear system, (7) is solved iteratively until P reaches zero. Here, Newton method is used to solve the problem.

In AFEM, the choice of element depends on the atomic structure and nature of atomistic interaction. In the studies for carbon nanotubes by Liu et al. [10, 11], an AFEM element containing 10 atoms is developed in order to capture the multibody interaction. For the current diamond lattice structure and multibody interaction, an AFEM element should contain 17 atoms that include 4 nearest-neighbor atoms (red balls) and 12 second nearest-neighbor atoms (blue balls) as shown in Figure 3. Therefore, the element stiffness matrix and the nonequilibrium force vector are given by [figure omitted; refer to PDF]

Figure 3: The AFEM element for silicon.

[figure omitted; refer to PDF]

In the present study, all work, including assembling the stiffness matrix and force vector and solving the equation system, is performed with our Fortran codes. The following steps are used. Firstly, the initial configuration of SiNWs with a uniform bond length is constructed, and the obtained coordinates and bond information are stored in an array. Using the obtained array, the information about the first and second neighbor atoms for each atom can be found and is stored in another array. Secondly, boundary conditions are applied to SiNWs to get the initial equilibrium coordinates of the system. Here, we restrain one side of SiNWs and make the other side free until the system returns to the equilibrium configuration. Table 2 gives the bond lengths of different types of SiNWs when the potential parameters in [13] are employed. These results are similar to the work in [14]. In the third step, displacement field is applied to the initial equilibrium coordinates in order to simulate the axial tension or compression. It is achieved by completely fixing one side of SiNWs and incrementally imposing an axial movement at the other end. The length of the nanowire is changed by 0.01 nm per loading step. In the final step, the total energy of SiNWs is obtained, and it is treated as a function of strain. With the first- and second-order derivatives, Newton iteration method is thus applied to obtain the equilibrium state, and 3-5 iterative steps can generally achieve a good convergence. Therefore, the present method is far faster than molecular dynamic method and is very efficient for the nanostructures with a larger number of atoms.

Table 2: Bond lengths of [100], [110], and [111] SiNWs.

Growth directions	100	110	111
Diameter (nm)	1.21	1.83	2.43	1.21	1.47	1.77	0.83	1.20	1.59
Bond length (nm)	0.233	0.234	0.234	0.233	0.233	0.234	0.232	0.233	0.233

3. Results and Discussions

In continuum mechanics, the material can be represented by two independent constants, namely, Young's modulus Y and Poisson's ratio ν. For a material undergoing a uniaxial deformation, Y is defined as [11, 15] [figure omitted; refer to PDF] where Ω is the volume corresponding to the initial equilibrium configurations of SiNWs; _Vt is the total strain energy under tensile or compressive deformation, and [straight epsilon] is strain.

The Poisson ratio ν is calculated as [figure omitted; refer to PDF] where _D0 and D are the mean diameters of SiNWs corresponding to the initial and deformed equilibrium configurations, respectively. They can be derived from S=π^D2 /4, where S is the cross-sectional area.

Young's modulus for the bulk silicon crystal under tension along [100] directions is first calculated in this paper, and it has the bulk structure with the side length 2.16 nm. The total potential energy is plotted versus strain in Figure 4. Using the polynomial curve fitting, the obtained Y is 134.03 GPa. For [110] and [111] directions, the present computational results are, respectively, 190.2 GPa and 178.3 GPa, which drop in the literature reports. Ma et al. [6] used all-electron density functional theory to calculate the binding energy, heat of formation, and Young's modulus of hydrogen-passivated SiNWs with various diameters and orientations, and Young's moduli of bulk silicon for the [100], [110], and [111] directions are 133.24 GPa, 157.66 GPa, and 163.58 GPa, respectively. Lee and Rudd [15] found that the bulk value of Young's modulus is 122.53 GPa for the [100] direction by the density functional theory.

Figure 4: Potential versus the length change of [100] bulk Si crystal for axial tensile/compressive deformation.

[figure omitted; refer to PDF]

Next, the numerical computation is carried out for Young's modulus and Poisson ratio of SiNWs with different dimensions. Table 3 shows the Poisson ratio for [100], [110], and [111] SiNWs, which have a small dependence on diameter. The calculated Young's moduli are plotted in Figure 5. Leu et al. [5] have given a systematic study of the mechanical property of strained small diameter silicon nanowires using ab initio density functional theory, and the values of Young's modulus are 139-153 GPa and 43-131 GPa for [110] and [111] SiNWs, respectively. In the researches by Ma et al. [6] and Lee and Rudd [15], Young's modulus increases while the wire mean diameter increases, and at certain critical size it approaches the bulk value. The present results agree well with these previous studies. It can also be seen from Figure 5 that, with the same mean diameter, nanowires with different grown directions have distinctly different Young's moduli, and the wires along [110] direction have the highest value while [100] wires have the lowest value.

Table 3: Poisson ratio of [100], [110], and [111] SiNWs.

Growth directions	100	110	111
Wire number	1	2	1	2	1	2
Poisson ratio	0.14	0.14	0.19	0.22	0.22	0.25

Figure 5: Young's modulus versus the mean diameter.

[figure omitted; refer to PDF]

Besides the axial strain, the simulation has also been carried out for the bending of SiNWs. The bending stiffness of a SiNW can be calculated as [figure omitted; refer to PDF] where L is the length and κ is the curvature of the bent SiNWs and is related to the bending angle θ as κ=θ/L [16].

The bending deformation can be achieved by incrementally rotating the two end planes in the opposite directions with respect to a line that is perpendicular to the axis of the undeformed SiNW and passing through its center. During the simulating process, the two end planes are reversely rotated 0.1 rad per loading step. Two [100] SiNWs with different diameters are studied in this paper. The change of the total energy is plotted as a function of the bending angle in Figure 6. Using polynomial curve fitting, we can obtain [figure omitted; refer to PDF] The bending stiffness can be calculated by replacing θ with κL. The bending stiffness is, respectively, 5432.7 and 1964.2 ev·nm for SiNWs with the diameters 2.43 nm and 1.83 nm. Menon et al. [17] investigated the bending stiffness of T and C type nanowires, and their computed results are 12.2 × 10^-6 ev·m and 8.5 × 10^-6 ev·m, respectively.

Figure 6: The change of the total energy as a function of bending angle for [100] SiNWs.

[figure omitted; refer to PDF]

4. Conclusions

The AFEM element is developed for the silicon nanowire, and a basic element contains 17 atoms while it is treated as the bulk structure. Young's modulus and Poisson ratio are calculated for [100], [110], and [111] SiNWs, and the bending stiffness is also obtained. AFEM is proved to be an accurate atomic simulation method for the bulk nanomaterials, and it is found that the [110] nanowire has the largest Young's modulus while the [100] nanowire has the smallest value. Newton iteration method is applied to obtain the equilibrium state, in which the first- and second-order derivatives are used and 3-5 iterative steps can achieve a good convergence, and the developed method is much faster than the molecular dynamic. In order to take advantage of the order-N nature of AFEM, further studies can be made to combine AFEM and the classical FEM software to solve problem with a large number of degrees, for example, the nanoindentation of bulk monocrystalline silicon.

Acknowledgments

The work in this research was supported by the Natural Science Foundation of China (Grant no. 11472316), Innovation Scientists and Technicians Troop Construction Projects of Henan Province (Project no. 164200510020), and Natural Science Foundation from Department of Education of Henan Province (Project no. 13B560304).

References

[1] K. Eom, H. S. Park, D. S. Yoon, T. Kwon, "Nanomechanical resonators and their applications in biological/chemical detection: nanomechanics principles," Physics Reports , vol. 503, no. 4-5, pp. 115-163, 2011.

[2] S.-W. Chung, J.-Y. Yu, J. R. Heath, "Silicon nanowire devices," Applied Physics Letters , vol. 76, no. 15, pp. 2068-2070, 2000.

[3] Y. Huang, X. Duan, Y. Cui, L. J. Lauhon, K.-H. Kim, C. M. Lieber, "Logic gates and computation from assembled nanowire building blocks," Science , vol. 294, no. 5545, pp. 1313-1317, 2001.

[4] T. Vo, A. J. Williamson, G. Galli, "First principles simulations of the structural and electronic properties of silicon nanowires," Physical Review B , vol. 74, no. 4, 2006.

[5] P. W. Leu, A. Svizhenko, K. Cho, "Ab initio calculations of the mechanical and electronic properties of strained Si nanowires," Physical Review B , vol. 77, no. 23, 2008.

[6] L. Ma, J. Wang, J. Zhao, G. Wang, "Anisotropy in stability and Young's modulus of hydrogenated silicon nanowires," Chemical Physics Letters , vol. 452, no. 1-3, pp. 183-187, 2008.

[7] G. Yun, H. S. Park, "Bridging the gap between experimental measurements and atomistic predictions of the elastic properties of silicon nanowires using multiscale modeling," Finite Elements in Analysis & Design , vol. 49, no. 1, pp. 3-12, 2012.

[8] G. H. Cocoletzi, N. Takeuchi, "First principles calculations of the structural and electronic properties of zinc sulfide nanowires," Quantum Matter , vol. 2, no. 5, pp. 382-387, 2013.

[9] K. Kang, W. Cai, "Brittle and ductile fracture of semiconductor nanowires--molecular dynamics simulations," Philosophical Magazine , vol. 87, no. 14-15, pp. 2169-2189, 2007.

[10] B. Liu, Y. Huang, H. Jiang, S. Qu, K. C. Hwang, "The atomic-scale finite element method," Computer Methods in Applied Mechanics & Engineering , vol. 193, no. 17-20, pp. 1849-1864, 2004.

[11] B. Liu H, Y. Huang, S. Qu, M. F. Yu, "Atomic-scale finite element method in multiscale computation with applications to carbon nanotubes," Physical Review B: Condensed Matter , vol. 72, no. 3, 2005.

[12] Y. Z. Sun, J. Tao, "The mechanical properties of tubular nanostructures through an atomic simulation method," Nanoscience & Nanotechnology Letters , vol. 7, no. 8, pp. 648-654, 2015.

[13] P. Erhart, K. Albe, "Analytical potential for atomistic simulations of silicon, carbon, and silicon carbide," Physical Review B , vol. 71, no. 3, 2005.

[14] A. K. Singh, V. Kumar, R. Note, Y. Kawazoe, "Effects of morphology and doping on the electronic and structural properties of hydrogenated silicon nanowires," Nano Letters , vol. 6, no. 5, pp. 920-925, 2006.

[15] B. Lee, R. E. Rudd, "First-principles study of the Young's modulus of Si (001) nanowires," Physical Review B , vol. 75, no. 4, 2007.

[16] D. Srivastava, C. Wei, K. Cho, "Nanomechanics of carbon nanotubes and composites," Applied Mechanics Reviews , vol. 56, no. 2, pp. 215-229, 2003.

[17] M. Menon, D. Srivastava, I. Ponomareva, L. A. Chernozatonskii, "Nanomechanics of silicon nanowires," Physical Review B , vol. 70, no. 12, 2004.