Size, vacancy and temperature effects on Young’s modulus of silicene nanoribbons

Abstract

We report results on the Young’s modulus (YM) of defect-free and defective silicene nanoribbons (SNRs) as a function of length and temperature. In this study, we perform molecular dynamics simulations using the Environment-Dependent Interatomic Potential (EDIP) to describe the interaction of the Si atoms. We find that the YM of pristine and defective SNRs increases with the ribbon length in both chirality directions. It is shown that the YM of defective SNRs exhibits a complex dependence on the combinations of vacancies. With respect to temperature, we find that YM for SNRs with and without vacancy defects shows nonlinear behavior and it could be tailored for a given length and chirality.

---

1 Introduction

At the beginning of the 21st century a new class of materials, called two-dimensional materials, was introduced.^1,2 The first identified material in this new category was graphene,³ which is a single atomic layer of carbon with a honeycomb structure. Since its discovery, in 2004, several studies have revealed the extraordinary properties of graphene, making it one of the most promising materials for applications in electronics,⁴ optics,⁵ sensors,⁶ energy storage,⁷ water purification,⁸ and biodevices.⁹ Besides graphene, a large variety of these layered materials have been studied such as boron nitride,^1,10,11 dichalcogenides,^1,2,11–14 germanene,^11,15–20 silicene,^{11,15–19,21–24} stanene,^15,16,25 phosphorene,^16,26–28 borophene,^29,30 and most recently titanium trisulfide,³¹ arsenene and antimonene.³² From this growing family of 2D materials, silicene has emerged to replace not only graphene, but also bulk silicon in current electronics, because studies have shown that silicene has electronic properties similar to those of graphene.^17,18 Silicene, analogous to graphene, is made of Si atoms arranged in a honeycomb structure. Nevertheless, it differs from graphene in the degree of flatness. While graphene has a perfectly flat stable structure, silicene shows a buckled stable structure. It should be noted that silicene is a material whose structure has not been obtained yet, like graphene. So, the existence of freestanding silicene sheets are so far hypothetical.¹⁷ Despite this, research shows that silicene can be synthesized on ordered substrates like Ag,^21,33–37 ZrBr₂,³⁸ and Ir.³⁹

Nowadays, the application of these layered materials in nanoelectromechanical systems (NEMS) has attracted much attention, so it is essential to understand their mechanical properties in order to exploit their potential. Thus far, the outstanding electronic^40,41 and mechanical⁴² properties of graphene, make it one of the most studied materials for such applications.^43–48 In this context, there are quite a few studies of the most fundamental mechanical property, which is the Young’s modulus (YM) or elastic modulus defined as the ratio of the stress to the strain as a material is stretched. Experimentally the YM of graphene has been obtained for monolayer, bilayer and multiple layers. Reported values for a monolayer are 1.2 TPa,⁴² 0.89 TPa⁴⁹ and 2.4 TPa;⁵⁰ for a bilayer are 0.39 TPa⁴⁹ and 2.0 TPa,⁵⁰ and finally for a multilayer (less than 5 layers) is 0.5 TPa.⁵¹ From a theoretical point of view, not only density functional theory (DFT) and tight binding (TB) calculations have been made, but also molecular dynamics (MD) and molecular mechanics (MM). Employing different kinds of potentials, these methods have been used to determine the YM of graphene, with reported values of for example 1.05,⁵² 0.91 TPa,⁵³ 1.01 and 1.19 TPa, respectively.^53,54 These values are in good agreement with those reported by Lee et al.⁴² However, the TB calculation shows a slightly lower value than others. As NEMS and their components (nanostructures) have structural dimensions around or below 100 nm, graphene nanoribbons (GNRs) are good candidates for their development. Thus, research on ribbon mechanical behavior is of increasing interest. Results have shown that the YM depends on the ribbon width as well as their chirality.^53,55,56 Since silicene could be equivalent to graphene and therefore has the advantage of easy integration into present technology, the interest of using silicene nanoribbons (SNRs) has opened a wide field of study. Currently, several studies of the mechanical properties of pristine silicene^11,57–61 and SNRs^62–64 have been reported.

A very important mechanical property in this kind of material is the YM because many design problems in nanoscale mechanical devices, like NEMS, are related to it. Mechanical properties of nanomaterials, as is well known, are size dependent. YM could also be affected by the presence of defects in the nanomaterial. In honeycomb structures, such as graphene and silicene, vacancy defects are the main defects observed.⁶⁵ So, a good understanding of the behavior of YM in nanomaterials with and without vacancy defects is helpful for the design and fabrication of reliable nanoscale mechanical devices.

In this paper, we apply classical molecular dynamics simulations to investigate the size and temperature dependence of the YM of free-standing SNRs. YM of pristine and defective SNRs as a function of length and temperature are reported. We found that the mechanical properties for pristine and defective SNRs exhibit a size dependence in both chirality directions. YM shows not only a size dependence but also a meaningful dependence of temperature and combinations of vacancy defects. The rest of the paper is organized as follows: in Section 2, we briefly describe our calculation method; in Section 3, the results are analyzed; and finally, conclusions are given in Section 4.

2 Methodology

To study the YM of low-buckled SNRs,^21,64 we perform molecular dynamics simulations of SNRs under uniaxial tension using LAMMPS.⁶⁶ Si atoms interact via Environment-Dependent Interatomic Potential (EDIP)⁶⁷ as implemented in LAMMPS. All simulations are performed at room temperature for three square-shaped SNRs, with and without defects, whose lengths are shown in Table 1. As the chirality matches with the direction of deformation, each ribbon length is represented by its chirality as ach or zz, and number labels indicate the size of the ribbon. In what follows, this nomenclature will be used to specify chirality direction and size. An NVE ensemble is used with 30 Å vacuum space on each side along the z-direction to relax the ideal structure (see Fig. 1). Fig. S1 in the ESI† shows front and lateral views of a relaxed structure, where it can be seen that the structure presents variable buckling. Similar behavior is found in all structures considered in this study. After relaxing the structure, NVT simulations are performed with a time step of 1 fs. We apply periodic boundary conditions (PBCs) along the direction of tension, which is also the chirality direction, as mentioned before. PBCs are needed in order to avoid any boundary effect in the structure along the tensile direction. Perpendicular edges to the tensile direction are free, see Fig. 1. A constant strain rate of 0.005 ps⁻¹ is applied to the structure to simulate uniaxial tension every time step during 90 ps. For integration of the motion equations we use the velocity-Verlet algorithm with an integration time of 1 fs. For each simulation, we obtain a stress–strain curve, since the YM is generally defined as the corresponding slope of the linear regime (elastic region) curve. A linear interpolation is performed in order to find YM.⁵⁶

Fig. 1 Top: ideal defect-free silicene, along with the definitions of armchair and zig-zag directions, to illustrate the simulation procedure. The dashed box indicates the boundaries of the simulation box. Arrows indicate the armchair or zig-zag straining direction. Bottom: side view of silicene structure with buckling distance h = 0.75 Å.^21,64

Table 1 Armchair and zig-zag square-shaped SNRs. Chirality is represented by ach or zz and number labels indicate the size of the ribbon

Length/Å	Chirality
	Armchair	Zig-zag
32.48 × 32.15	ach1	zz1
60.32 × 60.28	ach2	zz2
81.20 × 80.36	ach3	zz3

---

3 Results

First, we perform simulations at room temperature. A bulk (sheet) silicene with 4050 Si atoms is considered. YM obtained for this sheet are 154.7 and 152.5 GPa for zz and ach directions, respectively. These values are in good agreement with 154.8 and 153.8 GPa, obtained by Jing et al.⁶³ In order to study the effect of different vacancies, simulations are performed for different lengths of NRs, with and without vacancies. Later on, temperature dependence is also analyzed (Section 3.2).

3.1 Size effects

3.1.1 Pristine nanoribbons. Fig. 2(a) shows the strain–stress relationship for the zz3 (red on line) and ach3 (blue on line) pristine SNRs under uniaxial tensile test. From the figure we can see the well defined elastic (linear) region at small strain, followed by yielding until the ultimate tensile strength is reached (maximum stress). Necking then continues until fracture (the latter one is not included here). Fig. S2 in the ESI† shows the atomic configurations of a relaxed pristine SNRs at different strains. YMs obtained for each defect-free ribbon are plotted in Fig. 2(b). YM for both directions increases with ribbon length. This trend corresponds with that reported in ref. 63. Furthermore, zz-SNRs show larger values than those obtained for ach. Therefore, YM exhibits size as well as chirality dependence. This behavior has also been found in GNRs.⁵³ The difference between the YM values with chirality occurs due to the orientation of the bonds along the tension direction. When a tensile load is applied along to the zz direction, i.e., strain is applied perpendicular to the ach direction, the stretched Si–Si bonds are inclined, (see Fig. 1). Similar to graphene, bonds of this kind carry a part of the strain by changing the bond angles as well as increasing the bond length. Hence, deformations in the zz direction can be more supported. On the other hand under ach deformation, strain is applied in such a way that bonds parallel to the strain direction control deformation and their lengths increase more easily than in the case of zz, leading to a smaller slope and, consequently, a smaller YM.⁶⁸ For the same reason, ach-SNRs show a stronger size effect. We can also observe that for both chiralities, YM approaches pristine sheet values around 100 Å.

Fig. 2 (a) Strain–stress relations of ach(zz)3 defect-free SNRs under tensile test. (b) Young’s modulus for pristine SNRs as a function of length.

3.1.2 Defective nanoribbons. Vacancy defects are created by removing atoms. We consider three kinds: monovacancies (mv) in which one atom is removed, and two types of bivacancies (bv) where two adjacent atoms are eliminated. A bv is created by removing two adjacent atoms along the ach direction (bvp or parallel) or the zz one (bva or angular). These defects are schematically illustrated in Fig. 3. YM may be affected by the position of the missing bonds since vacancies near to the edge are expected to make the ribbon more fragile. Three cases for the position are included. Vacancies at the center of the ribbon are labeled with “c” (mvc, bvac, and bvpc), while those close to the boundaries with “e” (mve, bvae, and bvpe). Next, combinations of one central and four externals are built. Notice that four externals are included (close to each corner of the NR) to have a symmetric ribbon, making the applied tensile stress along both chiral directions comparable. Since YM for pristine NRs depends on length, from now on results will be presented in all cases normalized to the free-defective value and we will refer to them as normalized values or simply YM.

Fig. 3 Schematic representation of the vacancy-defects in SNRs.

3.1.2.1 Mono- and bi-vacancy defects. Fig. 4 is a schematic representation of the central and external (close to the ribbon boundaries) mono- and bi-vacancy defects set into the studied ribbons. As explained before, external vacancies are situated symmetrically near the corners.

Fig. 4 Schematic representation of the vacancy defects corresponding to the graphs shown in Fig. 5. Broken bonds are indicated with thin red lines.

In Fig. 5 we observe the effect on YM of type and position of the vacancy, as a function of length. On the top of each plot the value for the pristine nanoribbon is shown. As mentioned before, for a given length results are normalized to this value. Left (right) graphs correspond to ach(zz) ribbons. Panels a, b and c correspond to bvp, bva and mv. Numbers in the legend indicate the configuration of vacancies in the NRs shown in Fig. 4. It can be seen that dependence on length is similar to that of pristine case, i.e., the longer the ribbon the larger the YM. We also observe that the change in YM due to vacancies will be negligible when the length approaches 100 Å. Le et al. have reported similar behavior of YM for a pristine sheet and one with a monovacancy in its center.⁶⁹

Fig. 5 Young’s modulus for defect-free and defective armchair and zig-zag SNRs vs. length. Defects correspond to central (c) and external (e) vacancies. The numerical values at the top correspond to the pristine SNR, in GPa.

We first analyze the ach case. Each panel shows the effect of the vacancy position (see Fig. 5). Same trend with length as in pristine NRs is observed for all vacancies. YM for central vacancies (blue triangles) is larger than for external ones, as for the latter the NRs have more missing bonds (5 for central vs. 20 for external bv, while 3 vs. 12 for mv, see Fig. 4 for the structure). More missing bonds induces a larger deterioration of the structure, allowing the bonds to carry less strain making the structure easier to deform. Regarding the type, for bivacancies (panels a and b in Fig. 5) YM does not show a significant difference since the number of missing bonds is the same for them (structures 3–6 in Fig. 4). The small difference could be explained in terms of the orientation of the missing bonds. While for the angular one (bvac) we have two missing bonds parallel to the deformation direction, for bvpc there is only one. Thus, bvpc has a larger relative value.

As for the zz direction (right graphs), the same behavior is observed in general. The main difference is a larger separation between external (down triangle, red on line) and central (up triangle, blue on line) bivacancies for the smallest ribbon (panels d and e in Fig. 5). Along this orientation, the externals have a larger number of missing bonds (Fig. 4, structures 3–4 and 5–6). Therefore, YM are smaller than those for central vacancies.

3.1.2.2 Combinations of vacancy defects. In Fig. 6 we illustrate the different combinations of vacancy defects placed into the ribbons.

Fig. 6 Schematic representation of the vacancy-defects corresponding to the graphs shown in Fig. 7.

In order to analyze the behavior of normalized YM in the presence of combinations of vacancy-defects, Fig. 7 shows different kinds of YM values. As before, the left (right) graphs correspond to ach(zz) direction. Each panel contains three columns of data for a given length. The first column shows YM for NR with non-combined vacancies (shown in Fig. 5). The second and third columns show, respectively, YM as determined by MD (filled squares) for the indicated combination, and the average of the corresponding non-combined vacancies (empty circles and crosses). Panels a, b, and c deal with mv, bvp, bva located at the center of the NR. The overall behavior of bv follows the non-combined NRs, i.e., the shorter the NR, the more sensitive the YM is to length, see Section 3.1.2. In all cases, values given by the MD calculation are smaller than the average. The difference between them decreases as length increases, except when the ach-NR has an angular bivacancy at the center (bvac + bvpe and bvac + mve, panel c) since both combinations have the same number of missing bonds along the deformation direction. This is overlooked for long NRs given their large number of bonds, thus both cases approach the pristine case. The same number of missing bonds in the mentioned combination also causes that the MD values to be close to those of the non-combined NRs, opposite to the general trend of the rest for which the MD values are smaller. As for the zz direction, we observe a similar trend in the type of vacancy for corresponding panels (a and d, as well as b and e). Thus, YM for a combination of vacancies in general follows the same trend with size as in the pristine case. However, YM is not a simple average of the values of the non-combined NRs, giving a complex relation between type and position of vacancies. Thus, tailoring YM can be possible by means of a combination of vacancies of different types, for a given length and chirality.

Fig. 7 Young’s modulus for defect-free and defective armchair and zig-zag SNRs vs. length. Triangles correspond to non-combined vacancies, square symbols represent values obtained from the MD simulation and empty circles and crosses give the average (avg) of the corresponding non-combined NRs. The numerical values at the top correspond to the pristine SNR, in GPa.

3.1.2.3 Remarks. In general: (a) Relaxed structures with mve, show that vacancies do not coalesce into a multivacancy (see ESI,† Fig. S3(a)). This behavior is preserved even in the presence of two nearby monovacancies, where vacancies do not move to form a bivacancy (see ESI,† Fig. S3(b)). These results are explained by the self-healing mechanism of vacancy-defects, which has already been studied in graphene and silicene. In this process, the atoms surrounding the vacancy move towards the centre of the vacancy and the bonds are reoriented to close it. In the case of a monovacancy, the three dangling bonds from the atoms around the vacancy-defect form bonds that are stable and have the same length. For the bivacancy case, the atoms around the defect have four dangling bonds, which close the vacancy forming two Si–Si bonds.⁶⁵ This reconstruction mechanism of the atoms avoids the coalescence of vacancies. (b) Critical strains for SNRs under uniaxial tension are 19.5(17.5), 17.0(12.0), and 14.5(11.5)% in the zz(ach) direction. As we can see, the critical strain decreases as length increases. These results show that the critical strain depends on the ribbon length. Similar behaviour has been observed in square-shaped graphene-NRs (GNRs). On the one hand, zz-GNRs critical strain decreases with the ribbon width, and on the other critical strain for ach-GNRs slightly depends on the ribbon width.⁷⁰ Our results are similar to those reported for graphene as we are also studying square-shaped structures, which by increasing its length also increases its width. (c) We test the case where the vacancy concentration is the same in two NRs. In one, length is taken as L (ach(zz)1) and the other has length 2L. The first one has one mvc and the second one has four mv (see ESI,† Fig. S4). YM values obtained for the NR with length of 2L are 121.27 and 137.96 GPa for ach and zz directions, respectively. While, YM for the NR with length of L are 99.18 GPa for the ach-direction and 121.40 GPa for the zz-direction. These results show that not only the density of vacancy-defects are important as length increases, but also the number and position of vacancies strongly influences the value of YM. (d) Fig. S5(a) in ESI† shows the average buckling behaviour of a pristine SNR as function of the time step. It can be seen that the buckling amplitude varies with the time step, this dispersion is due to the displacement of the atoms in the x-direction, which is compressed by exerting a tension in the y-direction. At around 40 ps the structure has reached the necking region, and the buckling dispersion decreases as a result of the strain as expected, because with strain the structure tends to be planar. In Fig. S5(b)† the buckling for each atom confirms a qualitative change after necking.

3.2 Temperature effects

3.2.1 Mono- and bi-vacancy defects. NRs of 60 × 60 Å represent an intermediate size regime for our nanoribbons. For this reason, we have chosen them to analyze the temperature effect. Left (right) graphs of Fig. 8 shows the normalized YM of defective ach(zz) 2-SNRs at temperatures between 100 and 1000 K. It is worth mentioning that one atom misses during the calculation at 1000 K. This might be interpreted as if one phase transition could take place. For this reason, this is the highest temperature we study. Panels a, b, and c correspond to bvp, bva and mv types, respectively. On top of the panel, the value of pristine NR is shown at the given temperature. As expected, the larger the temperature the smaller the YM because at high temperatures the bond lengths are larger due to the thermal motion of atoms, therefore the deformation resistance is lower. A similar behavior has been observed in pristine ach graphene, which even at 2400 K has 90% of its value at room temperature.⁷¹ For silicene instead, the value has diminished to 61(90)% of its value for ach(zz) at 1000 K. This is consistent with the known facts that graphene is stronger than silicene and pristine zz is in turn stronger than ach. However, while the pristine value follows a systematic decrease, defective NRs show an overall complex behavior. YM decreases up to 500 K as with the pristine case, then it increases (depending on type of vacancy) and finally it drops. Central vacancies (blue triangles) always have larger relative values than the external (red triangles). As mentioned before, this is a consequence of the missing bonds. Some particular characteristics for vacancy types are described now. At some temperatures, YM does not depend on the type of vacancy as seen at 600 and 800 for the ach-bva, whereas for ach-mv only at 800 K (panels b and c, respectively). As for the zz case, this occurs only at 1000 K for mv (panel f). On the other hand, YM values for central and externals are similar around 800 K, but at higher temperatures they tend to separate. However, an opposite behavior is observed for the zz-bv case (d). This shows that thermal motion causes atoms to arrange in a complex structure, giving rise to similar temperature dependence but with small differences determined by the vacancy type and chirality.

Fig. 8 Young’s modulus for defect-free and defective armchair and zig-zag SNRs vs. temperature. Defects correspond to central and close to the ribbon’s boundary vacancies (see Fig. 4). The numerical values at the top correspond to the pristine SNR, in GPa.

3.2.2 Combinations of vacancy defects. In Fig. 9, normalized YMs are plotted for different combinations of vacancy defects as a function of temperature. As before, left (right) graphs correspond to ach(zz)-NRs. The decrease of YM with increasing temperature remains, as in the previous case. However, the results do not show a strong dependence on the type of combined vacancies up to 600 K. At higher temperatures, differences with respect to the pristine case are significant only for combinations with mvc (panel c). In contrast, all zz-NRs (right panels) show a decrease and a clear difference for the combinations is observed at high temperatures. This is due to the stronger bonds along the zz direction such that atoms rearrange less easily compared to the ach. Therefore, chirality plays an important role in this case. A more detailed study of the structural behavior is in progress.

Fig. 9 Young’s modulus for defect-free and defective armchair and zig-zag SNRs vs. temperature. Defects correspond to different combinations of vacancies (see Fig. 6). The numerical values at the top correspond to the pristine SNR, in GPa.

4 Conclusions

We have performed MD simulations to study the YM of SNRs, with and without vacancy defects, as a function of length and temperature. At room temperature, the YM obtained for pristine nanoribbons nonlinearly increases with ribbon length. YM values for zz-SNRs are larger than those obtained for ach-SNRs. The difference arises due to the different type of bonds along the tension direction, with zz being the stronger as reported before. Three cases of vacancies (mono-, parallel and angular vacancies) were analyzed and, in all the studied cases, YM increases with length for the given number of vacancies. Thus, mechanical properties of defective SNRs exhibit a size dependence. As length increases, YM tends to the corresponding pristine value, which means that the vacancy defects in the SNRs considered in this work do not affect the YM as the size approaches 100 Å. Combinations of vacancies give rich and complex behavior. Values obtained by numerical calculations do not coincide with the average of the YMs of non-combined vacancies. As for the temperature, calculations were performed for an intermediate length, up to 1000 K. YM for the pristine NRs nonlinearly decrease with increasing temperature. When vacancies are included, it shows a similar nonlinear-dependence regardless of the type of vacancy, but different to the pristine case. However, some combinations show particular features which depend significantly on the temperature. Thus, we have shown that tailoring of YM can be done by means of a combination of vacancies of different types, for a given length and chirality. Our results can be helpful for the design and performance of devices based on these nanomaterials.

A main subject for further investigations, is to find the YM of silicene on a substrate. As we mentioned before, silicene could form different structures on metal surfaces.^72–74 Since Ag and Si atoms do not easily form an alloy, silver has been considered to be the appropriate substrate for the deposition of silicene.^34,75 However, to study the YM of silicene on a Ag-substrate we have to consider not only the strains that the substrate introduces to the material due to the lattice mismatch, but also the interactions of silicon atoms with the substrate, the effect of the substrate temperature and the hybridization between silicon and silver atoms.^35,76 All these factors can strongly influence the value of the YM, because they should modify the strain tensor components giving rise to a different value. Thus, YM values obtained here are to be considered as lower bounds.

Acknowledgements

Partially supported by CONACYT and VIEP-BUAP, México.

References

1. K. S. Novoselov, D. Jiang, F. Schedin, T. J. Booth, V. V. Khotkevich, S. V. Morozov and A. K. Geim, Proc. Natl. Acad. Sci. U. S. A., 2005, 102, 10451–10453.
2. P. Miro, M. Audiffred and T. Heine, Chem. Soc. Rev., 2014, 43, 6537–6554.
3. K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, Y. Zhang, S. V. Dubonos, I. V. Grigorieva and A. A. Firsov, Science, 2004, 666–669.
4. M. Woszczyna, A. Winter, M. Grothe, A. Willunat, S. Wundrack, R. Stosch, T. Weimann, F. Ahlers and A. Turchanin, Adv. Mater., 2014, 26, 4831–4837.
5. G. Jo, M. Choe, S. Lee, W. Park, Y. H. Kahn and T. Lee, Nanotechnology, 2012, 23, 112001.
6. Y. Wang, L. Wang, T. Yang, X. Li, X. Zang, M. Zhu, K. Wang, D. Wu and H. Zhu, Adv. Funct. Mater., 2014, 24, 4666–4670.
7. J. Zhu, D. Yang, Z. Yin, Q. Yan and H. Zhang, Small, 2014, 10, 3480–3498.
8. Y. Han, Z. Xu and C. Gao, Adv. Funct. Mater., 2013, 23, 3693–3700.
9. Y. Hu, F. Li, D. Han and L. Niu, in Biocompatible Graphene for Bioanalytical Applications, Springer Berlin Heidelberg, 2015, pp. 11–33.
10. L. Song, L. Ci, H. Lu, P. B. Sorokin, C. Jin, J. Ni, A. G. Kvashnin, D. G. Kvashnin, J. Lou, B. I. Yakobson and P. M. Ajayan, Nano Lett., 2010, 10, 3209–3215.
11. H. Şahin, S. Cahangirov, M. Topsakal, E. Bekaroglu, E. Akturk, R. T. Senger and S. Ciraci, Phys. Rev. B: Condens. Matter Mater. Phys., 2009, 80, 155453.
12. M. Chhowalla, H. S. Shin, G. Eda, L.-J. Li, K. P. Loh and H. Zhang, Nat. Chem., 2013, 5, 263–275.
13. Q. H. Wang, K. Kalantar-Zadeh, A. Kis, J. N. Coleman and M. S. Strano, Nat. Nanotechnol., 2012, 7, 699–712.
14. S. Z. Butler, S. M. Hollen, L. Cao, Y. Cui, J. A. Gupta, H. R. Gutiérrez, T. F. Heinz, S. S. Hong, J. Huang, A. F. Ismach, E. Johnston-Halperin, M. Kuno, V. V. Plashnitsa, R. D. Robinson, R. S. Ruoff, S. Salahuddin, J. Shan, L. Shi, M. G. Spencer, M. Terrones, W. Windl and J. E. Goldberger, ACS Nano, 2013, 7, 2898–2926.
15. J. C. Garcia, D. B. de Lima, L. V. C. Assali and J. F. Justo, J. Phys. Chem. C, 2011, 115, 13242–13246.
16. S. Balendhran, S. Walia, H. Nili, S. Sriram and M. Bhaskaran, Small, 2015, 11, 640–652.
17. S. Cahangirov, M. Topsakal, E. Aktürk, H. Şahin and S. Ciraci, Phys. Rev. Lett., 2009, 102, 236804.
18. S. Lebègue and O. Eriksson, Phys. Rev. B, 2009, 79, 115409.
19. C.-C. Liu, W. Feng and Y. Yao, Phys. Rev. Lett., 2011, 107, 076802.
20. M. E. Dávila, L. Xian, S. Cahangirov, A. Rubio and G. L. Lay, New J. Phys., 2014, 16, 095002.
21. P. Vogt, P. de Padova, C. Quaresima, J. Avila, E. Frantzeskakis, M. C. Asensio, A. Resta, B. Ealet and G. Le Lay, Phys. Rev. Lett., 2012, 108, 155501.
22. X. Lin and J. Ni, Phys. Rev. B, 2012, 86, 075440.
23. M. Houssa, A. Dimoulas and A. Molle, J. Phys.: Condens. Matter, 2015, 27, 253002.
24. D. Jose and A. Datta, Acc. Chem. Res., 2014, 47, 593–602.
25. Y. Xu, B. Yan, H.-J. Zhang, J. Wang, G. Xu, P. Tang, W. Duan and S.-C. Zhang, Phys. Rev. Lett., 2013, 111, 136804.
26. H. Liu, A. T. Neal, Z. Zhu, Z. Luo, X. Xu, D. Tománek and P. D. Ye, ACS Nano, 2014, 8, 4033–4041.
27. L. Li, Y. Yu, G. J. Ye, Q. Ge, X. Ou, H. Wu, D. Feng, X. H. Chen and Y. Zhang, Nat. Nanotechnol., 2014, 9, 372–377.
28. F. Xia, H. Wang and Y. Jia, Nat. Commun., 2014, 5, 5458.
29. Z. A. Piazza, H.-S. Hu, W.-L. Li, Y.-F. Zhao, J. Li and L.-S. Wang, Nat. Commun., 2014, 5, 4113.
30. H.-J. Zhai, Y.-F. Zhao, W.-L. Li, Q. Chen, H. Bai, H.-S. Hu, Z. A. Piazza, W.-J. Tian, H.-G. Lu, Y.-B. Wu, Y.-W. Mu, G.-F. Wei, Z.-P. Liu, J. Li, S.-D. Li and L.-S. Wang, Nat. Chem., 2014, 6, 727–731.
31. J. Dai and X. C. Zeng, Angew. Chem., Int. Ed., 2015, 54, 7572–7576.
32. S. Zhang, Z. Yan, Y. Li, Z. Chen and H. Zeng, Angew. Chem., Int. Ed., 2015, 54, 3112–3115.
33. C.-L. Lin, R. Arafune, K. Kawahara, N. Tsukahara, E. Minamitani, Y. Kim, N. Takagi and M. Kawai, Appl. Phys. Express, 2012, 5, 045802.
34. H. Enriquez, S. Vizzini, A. Kara, B. Lalmi and H. Oughaddou, J. Phys.: Condens. Matter, 2012, 24, 314211.
35. H. Jamgotchian, Y. Colignon, N. Hamzaoui, B. Ealet, J. Y. Hoarau, B. Aufray and J. P. Bibérian, J. Phys.: Condens. Matter, 2012, 24, 172001.
36. L. Chen, C.-C. Liu, B. Feng, X. He, P. Cheng, Z. Ding, S. Meng, Y. Yao and K. Wu, Phys. Rev. Lett., 2012, 109, 056804.
37. L. Chen, H. Li, B. Feng, Z. Ding, J. Qiu, P. Cheng, K. Wu and S. Meng, Phys. Rev. Lett., 2013, 110, 085504.
38. A. Fleurence, R. Friedlein, T. Ozaki, H. Kawai, Y. Wang and Y. Yamada-Takamura, Phys. Rev. Lett., 2012, 108, 245501.
39. L. Meng, Y. Wang, L. Zhang, S. Du, R. Wu, L. Li, Y. Zhang, G. Li, H. Zhou, W. A. Hofer and H.-J. Gao, Nano Lett., 2013, 13, 685–690.
40. K. S. Novoselov, S. V. Morozov, T. M. G. Mohinddin, L. A. Ponomarenko, D. C. Elias, R. Yang, I. I. Barbolina, P. Blake, T. J. Booth, D. Jiang, J. Giesbers, E. W. Hill and A. K. Geim, Phys. Status Solidi B, 2007, 244, 4106–4111.
41. A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov and A. K. Geim, Rev. Mod. Phys., 2009, 81, 109–162.
42. C. Lee, X. Wei, J. W. Kysar and J. Hone, Science, 2008, 321, 385–388.
43. B. Standley, W. Bao, H. Zhang, J. Bruck, C. N. Lau and M. Bockrath, Nano Lett., 2008, 8, 3345–3349.
44. D. Garcia-Sanchez, A. M. van der Zande, A. S. Paulo, B. Lassagne, P. L. McEuen and A. Bachtold, Nano Lett., 2008, 8, 1399–1403.
45. C. Chen, S. Rosenblatt, K. I. Bolotin, W. Kalb, P. Kim, I. Kymissis, H. L. Stormer, T. F. Heinz and J. Hone, Nat. Nanotechnol., 2009, 4, 861–867.
46. X. Liu, J. W. Suk, N. G. Boddeti, L. Cantley, L. Wang, J. M. Gray, H. J. Hall, V. M. Bright, C. T. Rogers, M. L. Dunn, R. S. Ruoff and J. S. Bunch, Adv. Mater., 2014, 26, 1571–1576.
47. X. Zang, Q. Zhou, J. Chang, Y. Liu and L. Lin, Microelectron. Eng., 2015, 132, 192–206.
48. X. Wang, H. Tian, W. Xie, Y. Shu, W.-T. Mi, M. Ali Mohammad, Q.-Y. Xie, Y. Yang, J.-B. Xu and T.-L. Ren, NPG Asia Mater., 2015, 7, e154.
49. Y. Zhang and C. Pan, Diamond Relat. Mater., 2012, 24, 1–5.
50. J.-U. Lee, D. Yoon and H. Cheong, Nano Lett., 2012, 12, 4444–4448.
51. I. W. Frank, D. M. Tanenbaum, A. M. van der Zande and P. L. McEuen, J. Vac. Sci. Technol., B: Microelectron. Nanometer Struct.--Process., Meas., Phenom., 2007, 25, 2558–2561.
52. B. Hajgató, S. Güryel, Y. Dauphin, J.-M. Blairon, H. E. Miltner, G. Van Lier, F. De Proft and P. Geerlings, J. Phys. Chem. C, 2012, 116, 22608–22618.
53. H. Zhao, K. Min and N. R. Aluru, Nano Lett., 2009, 9, 3012–3015.
54. L. Zhou, J. Xue, Y. Wang and G. Cao, Carbon, 2013, 63, 117–124.
55. Z. Xu, J. Comput. Theor. Nanosci., 2009, 6, 625–628(4).
56. Q. Lu, W. Gao and R. Huang, Modell. Simul. Mater. Sci. Eng., 2011, 19, 054006.
57. R. Qin, C.-H. Wang, W. Zhu and Y. Zhang, AIP Adv., 2012, 2, 022159.
58. H. Zhao, Phys. Lett. A, 2012, 376, 3546–3550.
59. Q. Peng, X. Wen and S. De, RSC Adv., 2013, 3, 13772–13781.
60. R. E. Roman and S. W. Cranford, Comput. Mater. Sci., 2014, 82, 50–55.
61. Q.-X. Pei, Z.-D. Sha, Y.-Y. Zhang and Y.-W. Zhang, J. Appl. Phys., 2014, 115, 023519.
62. M. Topsakal and S. Ciraci, Phys. Rev. B, 2010, 81, 024107.
63. Y. Jing, Y. Sun, H. Niu and J. Shen, Phys. Status Solidi B, 2013, 250, 1505–1509.
64. R. Ansari, S. Rouhi and S. Ajori, Superlattices Microstruct., 2014, 65, 64–70.
65. V. O. Özçelik, H. H. Gurel and S. Ciraci, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 88, 045440.
66. S. Plimpton, J. Comp. Physiol., 1995, 117, 1–19.
67. J. F. Justo, M. Z. Bazant, E. Kaxiras, V. V. Bulatov and S. Yip, Phys. Rev. B: Condens. Matter Mater. Phys., 1998, 58, 2539–2550.
68. M. A. N. Dewapriya, A. Srikantha Phani and R. K. N. D. Rajapakse, Modell. Simul. Mater. Sci. Eng., 2013, 21, 065017.
69. M.-Q. Le and D.-T. Nguyen, Appl. Phys. A, 2015, 118, 1437–1445.
70. Y. Chu, T. Ragab and C. Basaran, Comput. Mater. Sci., 2014, 81, 269–274.
71. H. Zhao and N. R. Aluru, J. Appl. Phys., 2010, 108, 064321.
72. D. Kaltsas, L. Tsetseris and A. Dimoulas, Appl. Surf. Sci., 2014, 291, 93–97.
73. H. Liu, J. Gao and J. Zhao, J. Phys.: Conf. Ser., 2014, 491, 012007.
74. R. Quhe, Y. Yuan, J. Zheng, Y. Wang, Z. Ni, J. Shi, D. Yu, J. Yang and J. Lu, Sci. Rep., 2014, 4, 1–8.
75. B. Feng, Z. Ding, S. Meng, Y. Yao, X. He, P. Cheng, L. Chen and K. Wu, Nano Lett., 2012, 12, 3507–3511.
76. Y. Yuan, R. Quhe, J. Zheng, Y. Wang, Z. Ni, J. Shi and J. Lu, Phys. E, 2014, 58, 38–42.

---

Footnote

† Electronic supplementary information (ESI) available: Snapshots of atomic configurations of silicene nanoribbons. See DOI: 10.1039/c5ra15312c

---