Defections induced hydrogenation of silicene: a density functional theory calculation study

Abstract

Based on density functional theory calculations, the dissociative adsorption of a H₂ molecule on silicene with a single atom vacancy defect (SV) or diatom vacancy defect (DV) has been investigated. The results show that the energy barrier for the dissociative adsorption of H₂ molecules on silicene can be significantly depressed with these two kinds of defects studied. The dissociation energy barrier for H₂ molecules is decreased from 2.23 eV on pristine silicene to 0.71 eV on silicene with single atom vacancy defects SV1(55|66), which is much smaller than the critical energy barrier of 0.91 eV, indicating that the reaction proceed at ambient temperature. The data in brackets indicate the number of silicon atoms for each member ring in the defective region. In addition, the dissociation energy barriers are 0.99 and 1.14 eV on silicene with divacancy defects DV1(5|8|5) and DV2(555|777), respectively, which are very close to the value of 0.91 eV. Therefore, we propose a promising method to tune the reaction activity of silicene and facilitate the hydrogenation of silicene, which can open the band gap of the zero band gap silicene, thus increasing the on/off ratio, which is essential for the potential applications of silicene in electronic devices.

---

1. Introduction

Stimulated by the excellent performance of graphene, silicene has attracted significant interest from both industrial and academic fields due to the outstanding electronic properties, such as the linear dispersion near the Fermi level and the quantum spin Hall effect.¹ Experimentally, silicene has been fabricated on many substrates, i.e. Ag(111),^2,3 ZrB₂(0001),⁴ Ir(111)⁵ and Au(110)⁶ surfaces. Unlike the flat honeycomb lattice of graphene, silicene has low buckling of about 0.44 Å,⁷ which allows more flexibility to functionalize their electronic and magnetic properties.^8–12 Moreover, silicene is more compatible with modern semiconductor technology based on silicon than graphene.¹³ Therefore, silicene is believed to have a bright future¹⁴ and shows potential applications in nanoscale devices.¹⁵

In actual experiment and application of monolayer materials, defects play an important role and cannot be ignored. Vacancy defects are unavoidable during the fabrication of monolayers and some kinds of defects are even purposively introduced for specific applications in some cases.¹⁶ Vacancies in graphene and silicene have been studied,^17–26 i.e. it is reported that the different kind of point defects in low dimension materials can modify their thermal stability, local structure, band structures and gold deposition. Furthermore, silicene with double vacancies has a small band gap, while the silicene with single vacancy may transfer from semimetallic into metallic.²² Vacancy defects could also modify the adsorption of molecules on silicene. It is reported that silicene can be functionalized via creating meshes of a single vacancy, where specific molecules can be either dissociated or pinned at the vacancy site.²⁷ It is well known that hydrogenation, i.e. the adsorption of hydrogen atoms is an efficient method to modulate the electronic structures. Therefore, we expect that the silicon vacancy defects can depress the energy barrier for the dissociative adsorption of hydrogen molecule due to the enhanced activity of silicene. Thus, creating vacancy defects would be an alternative method to accelerate the hydrogenation of silicene.

In this article, we perform a comprehensive study on the dissociative adsorption of hydrogen molecules on pristine and defective silicene based on density functional theory (DFT) calculations to understand the effects of vacancy defects on the hydrogenation of silicene and corresponding electronic properties of silicene. The energy barriers for this reaction process are also discussed and the corresponding mechanism is discussed by analysing the electronic structures.

2. Calculation methods

All the DFT results were calculated by using the DMol³ module.²⁸ The exchange–correlation functional has adopted generalized gradient approximation (GGA) with the Perdew–Burke–Ernzerhof (PBE)²⁹ functional. The core electrons were replaced by a single effective potential, i.e. core treatment is implemented with DFT semicore pseudopotentials (DSPPs) for relativistic effects. The basis set has used a double numerical plus polarization (DNP) method. The convergence tolerance of energy is taken as 10⁻⁵ Ha (1 Ha = 27.21 eV), with a maximum allowed force and displacement of 0.002 Ha and 0.005 Å respectively. It was reported that the exchange–correlation functional has much smaller effect on the energy barriers of calculated reaction than that of adsorption energies.³⁰ In order to investigate the minimum energy pathway for the dissociative adsorption of hydrogen molecule on silicene, linear synchronous transition/quadratic synchronous transit (LST/QST)³¹ and nudged elastic band (NEB)³² tools in DMol³ were employed. To take into consideration of the effect of van der Waals forces, Grimme scheme³³ for DFT-D correction is used. Spin-polarized calculations are added for band structure observation. The simulation has used three-dimensional periodic boundary conditions with 5 × 5 × 1 supercell for all studied structures, where a vacuum value of 20 Å is employed above the silicene surface to avoid the interaction between silicene and its periodic images. All atoms have been relaxed to reach their most stable positions, where the k-point is set to 6 × 6 × 1 based on the previous reports.²² After geometry optimization, a finer k-point grid of 15 × 15 × 1 is applied for calculating the density of states (DOS) for the silicene systems.

3. Results and discussions

The pristine silicene is shown in Fig. 1(a), while Fig. 1(b)–(d) show the totally relaxed silicene sheet with single and double silicon vacancy defects. For silicene with single vacancy (SV) defect, there are two possible reconstructed structures, SV(55|66) and SV2(5|3),^22,23 where SV(55|66) includes a sp³-hybridized central silicon atom and SV2(5|3) has three dangling atoms. In addition, the data in bracket indicate the number of silicon atoms for each member ring in the defective region indicated by different colours, as shown in Fig. 1. Previous studies indicate that SV2(5|3) is metastable in energy, which could transform to SV (55|66) by surmounting a small energy barrier of 0.91 eV.²⁶ Therefore, only SV(55|66) is studied for the dissociative adsorption of H₂ in the following as shown in Fig. 1(b). The double vacancy (DV) in silicene is formed by removing two neighbouring silicon atoms without dangling bond appearing after structure reconstruction. Fig. 1(c) shows the DV1(5|8|5) structure, where the two pentagons are connective by an octagon. When rotating one bond in the octagon, DV1(5|8|5) structure will transform into DV2(555|777) structure with three pentagons and three heptagons, as shown in Fig. 1(d), which has lower formation energy of 0.86 eV.²² The transformation energy barrier from DV1 to DV2 is up to 1.20 eV,²² therefore both configurations for double vacancy are possible existed and taken into account for the dissociative adsorption of H₂ molecules in this work.

Fig. 1 The most stable configurations of pristine silicene (a), defective silicene with single atom vacancy SV (b), double atom vacancy DV1 (c), and double atom vacancy DV2 (d). The red colour indicates the pristine region of silicene, while the green and purple colours indicate the defective region of silicene.

The band structures of pristine and defective silicene are also studied. The band structure of pristine silicene is shown in Fig. 1(a), where the π and π* bands show linearly character at the Fermi level, demonstrating the high charge mobility of about 10⁶ m s⁻¹.⁷ For silicene, it is known that the point defects can affect the electronic structure significantly.²² The defective silicene with SV shows metallic character with two energy bands crossing the Fermi level [see Fig. 1(b)], where the SV defect in silicene is nonmagnetic, which is different from that in graphene,^34,35 because the sp³-hybrid silicon atom in the center of 55|66 defective region forms fourfold bonds, while the defective graphene has dangling bonds at vacancies.^34,35 The band structure for silicene with DV1 defect is shown in Fig. 1(c), where an indirect band gap of 0.16 eV is opened. The silicene with DV2 defect shows semimetallic behaviour with zero band gap at Γ point instead of K point [see Fig. 1(d)], where the linear dispersion at the Fermi energy is destroyed. Note that at 0.45 eV [see Fig. 1(b)] and 0.41 eV [see Fig. 1(d)] above the Fermi level, two bands show linearly character at K point with a small band gap for silicene with SV and DV2, respectively, which are the Dirac cones, and they are shifted upwards when the presence of vacancies in silicene.

We first take pristine silicene as an example to study the dissociative adsorption process of hydrogen molecules on silicene as shown in Fig. 2. There are four possible positions for the adsorption of hydrogen molecule on pristine silicene [see Fig. 2(a)]: T and D represent the top position above the silicon atom in the up or down plane of the buckled silicene, B and H represent the bridge site of Si–Si bond and the hollow site of the six-member ring, respectively. After DFT + D calculations, the hydrogen molecule prefers to adsorb at the hollow site, which is also taken as the initial state (IS) of the hydrogenation process for silicene. For the atomic hydrogen adsorption, one hydrogen atom is first fixed at position 0. Three positions are possible for the adsorption of the second hydrogen atom, which are labeled as 1–3 as indicated in Fig. 2(b). After careful calculations, it is found that the hydrogenated silicene with the second hydrogen atom adsorbed at position 1 is the most stable configuration [FS, see Fig. 2(b)], where the bond length l_Si–H = 1.50 Å and l_Si–Si = 2.355 Å, similar to literature data l_Si–H = 1.50 Å, l_Si–Si = 2.36 Å (ref. 36) and l_Si–H = 1.519 Å, l_Si–Si = 2.359 Å (ref. 37) in silicane.

Fig. 2 (a) The most stable configuration of pristine silicene adsorbed with H₂ molecule, where the letters indicate the possible H₂ adsorption sites. (b) The most stable configuration of pristine silicene adsorbed with dissociative H atoms, where the numbers indicate the possible adsorption sites for the two H atoms. The yellow and cyan balls are respectively Si and H atoms in this and following figures.

For defective silicene, there are four, five and five possible adsorption positions for the hydrogen molecule (IS) as shown in Fig. 3(b)–(d). After careful examinations, the H₂ molecule prefers to adsorb at H2 site on silicene with SV defect, while that prefers to adsorb at H1 site on silicene with DV1 and DV2 defects. The adsorption energy of the H₂ molecule on both of pristine and defective silicene is summarized in Table 1. For atomic hydrogen adsorption (FS), when one hydrogen atom is fixed at position 0, there are one, three and three possible positions for the adsorption of the second hydrogen atom [see Fig. 3(b)–(d)]. After careful calculations, it is found that the second hydrogen atom prefers to adsorb at positions 1, 1 and 2 for the silicene with SV, DV1 and DV2 defections respectively.

Fig. 3 The reaction pathway of the dissociative adsorption of a H₂ molecule on pristine silicene (a), defective silicene with single atom vacancy SV (b), double atom vacancy DV1 (c), double atom vacancy DV2 (d).

Table 1 The adsorption energy E_ad, dissociative energy barrier E_bar, dissociative reaction energy E_r and Mulliken charge Q of H₂ molecule on pristine and defective silicene

	Pristine	SV	DV1	DV2
Ead (eV)	−0.113	−0.108	−0.125	−0.105
Ebar (eV)	2.23	0.71	0.99	1.14
Er (eV)	−0.16	−0.85	−0.84	−0.70
Q (e)	0.00	−0.05	−0.04	−0.07

---

The detailed reaction path for the dissociation of hydrogen molecule on silicene is then analyzed. Fig. 3(a) shows the minimum reaction path for the dissociation of a hydrogen molecule on pristine silicene. During this process, the energy barrier E_bar = E_TS − E_IS for dissociation is calculated to be 2.23 eV. Before the transition state (TS), the hydrogen molecule is split in two free hydrogen atoms and then one hydrogen atom binds with silicon atom on the upper layer while the other hydrogen atom keeps free at TS state. Furthermore, the silicon atom which will bind with the free hydrogen atom moves upwards but the H–H interaction is still strong, which avoids the formation of the Si–H bond. Thus, this process can be divided into two steps: the hydrogen molecule is split in two isolated atoms and one hydrogen atom binds with a silicon atom; subsequently, the other hydrogen atom binds with the corresponding silicon atom and the two bonded silicon atoms adjust their positions upwards. This process totally releases heat of 0.16 eV, the energy barrier is as large as 2.23 eV in the first step, and the following step releases heat energy of 2.39 eV. Therefore, the rate-limiting step is the first step due to the large energy needed in the dissociation.

For silicene with SV, the H₂ molecule lies at the 5-member ring and is vertical to the silicene surface as shown in Fig. 3(b). Before TS, the H₂ molecule moves downwards obviously and then it is dissociated into two free H atoms which both remain free without binding with Si atoms at TS. This could be explained by the larger distance between Si atom and H atom due to the existence of the defect. Meanwhile, distance between the two hydrogen atoms are still very short indicating strong H–H interactions. The two Si atoms that will bind with H atoms have both shifted upwards obviously even though one of the Si atoms is at the lower layer of silicene. After surmounting the energy barrier with 0.71 eV, the two H atoms bind with the corresponding silicon atoms (at the final state FS). One Si–H bond is parallel to the silicene plane, while the other Si–H bond is nearly perpendicular to the silicene plane. This is different from the two repellent H atoms on the product in Fig. 3(a) due to the adjacent position. The dissociation energy barrier is much lower than the critical energy barrier 0.91 eV,³⁸ below this barrier, a reaction can go through at ambient temperature. Therefore, the dissociation of hydrogen molecule on silicene with SV defect can go through smoothly at ambient temperature.

For silicene with DV1 defect, the H₂ molecule lies at the 8-member ring and is parallel to the silicene surface as shown in Fig. 3(c). After surmounting the dissociative energy barrier of 0.99 eV, two hydrogen atoms are bound with the corresponding silicon atoms and this reaction releases heat of 0.84 eV. The two Si–H bonds are almost perpendicular to the silicene surface, which is different from the two repellent H atoms on the pristine silicene in Fig. 3(a) due to the adjacent position. The dissociation energy barrier for H₂ molecule is slightly higher than the critical energy barrier 0.91 eV,³⁸ but much lower than that on pristine silicene (2.23 eV).

For silicene with DV2 defect, the H₂ molecule lies vertically at the 7-member ring as shown in Fig. 3(d). At TS state, the H–H bond is broken and there is no binding between H and Si atoms. After surmounting an energy barrier with 1.14 eV, the two hydrogen atoms bind with the corresponding silicon atoms and the Si–H bonds are almost perpendicular to the silicene surface with releasing heat of 0.70 eV. This energy barrier is slightly larger than that on silicene with DV1 defect, but much lower than that on pristine silicene (2.23 eV). Therefore, both DV1 and DV2 defect can significantly promote the dissociation of H₂ molecules on silicene. Therefore, the defects can significantly depress the dissociative energy barrier of H₂ molecule on silicene. The reversing energy barriers E′_bar (the energy difference between FS and TS) for the dissociation of H₂ molecule on defective silicene with SV, DV1 and DV2 are 1.56, 1.83, 1.84 eV, respectively, which are much larger than the critical energy barrier 0.91 eV,³⁸ indicating that the vacancy saturation can not be reverted and the defective silicene is not suitable for hydrogen reservoir.

To study the further performance of the defective silicene towards H₂ adsorption after vacancy saturation, the subsequent dissociative adsorption of the second H₂ molecule on defective silicene is shown in Fig. 4. We found that the second H₂ molecule prefers to adsorb on the alternative side of silicene after the vacancy saturation. For hydrogen molecule adsorption (IS), there are five, five and four possible adsorption positions as shown in Fig. 4(a)–(c). For atomic hydrogen adsorption (FS), when one hydrogen atom is fixed at position 0, there are three, three and five possible positions for the adsorption of the second hydrogen atom [see Fig. 4(a)–(c)]. After careful calculations, the most stable configurations before and after H₂ dissociations are shown in Fig. 4. The dissociative energy barriers for the second H₂ molecule are 1.21, 1.20 and 1.26 eV on silicene with SV, DV1 and DV2 defects, respectively, which are all larger than those for the dissociation of the first H₂ molecule and much higher than the critical energy barrier 0.91 eV,³⁸ indicating that the dissociation for the second H₂ molecule is difficult to occur at ambient temperature. Therefore, the following discussions are mainly focused on the dissociation for the first H₂ molecule on silicene.

Fig. 4 The reaction pathway of the dissociative adsorption of the second H₂ molecule on defective silicene with single atom vacancy SV (a), double atom vacancy DV1 (b), double atom vacancy DV2 (c). The white balls are the H atoms from the second adsorbed H₂ molecule.

The Mayer bond order³⁹ values of the first adsorbed H₂ molecule on silicene are also calculated and compared with the free H₂ molecule. The calculated bond order for free H₂ molecule is 1.0, while those for the adsorbed H₂ molecule on pristine silicene and defective silicene with SV, DV1 and DV2 are 0.9845, 0.9919, 0.9809 and 0.9924 eV, respectively. Therefore, the H–H bond slightly decreases and remains single bond character for the adsorbed H₂ molecule. After the dissociative adsorption of H₂ molecule on silicene, the Si–H bond orders for pristine silicene and defective silicene with SV, DV1 and DV2 are 0.8891, 0.8568, 0.8664 and 0.8691 eV, respectively, which indicates that single covalent bond (with bond order value larger than 0.5) is formed between Si and H atoms. The average Mulliken charges for the Si/H atoms are −0.030/0.019, −0.091/0.059, −0.037/0.016 and −0.057/0.015 e on pristine silicene and silicene with SV, DV1 and DV2 defects, respectively, which indicates the small charge transfer between Si and H atoms and further confirms the covalent character of Si–H bond.

The mechanism for depressing the dissociative energy barrier of hydrogen molecule on silicene with vacancy defects can be understood through analysing the partial density of states (PDOS) of the TS configurations with different defects as shown in Fig. 5. The PDOS of s orbital of the two split hydrogen atoms and the p orbital of the two silicon atoms binding with the two hydrogen atoms at TS are shown. It can be seen that the overlapping area of the two lines at Fermi level on defective silicene [see Fig. 5(b)–(d)] becomes much weaker compared with that at the pristine silicene [see Fig. 5(a)], indicating that the interaction of Si–H band is significant weakened on defective silicene. It is known that lower energy barrier is closely related to the weaker interaction nearby the Fermi level.⁴⁰ Meanwhile, the interaction in the area from about −10 eV to −4 eV is enhanced on defective silicene [see Fig. 5(b)–(d)]. Therefore, more interactions happen in the low energy area which could also lead to lower energy barrier on defective silicene.

Fig. 5 The PDOS of the two H atoms and the two corresponding Si atoms at TS on pristine silicene (a), defective silicene with single vacancy SV (b), double vacancy DV1 (c), double vacancy DV2 (d). The black curve is the PDOS of 2p orbital of the two silicon atoms and the blue one indicates 1s orbital of the two hydrogen atoms. The dotted line stands for the Fermi level.

After hydrogenation, the band structures of pristine and defective silicene are also studied. The band structure of hydrogenated pristine silicene is shown in Fig. 6(a), where linear dispersion is destroyed at Fermi energy and a direct band gap of 0.03 eV is obtained. The hydrogenated silicene with SV defect changes from metallic state to semiconducting state with direct band gap of 0.06 eV [see Fig. 6(b)]. For hydrogenated silicene with DV1 defect in Fig. 6(c), the band gap changes from indirect to direct and decreases from 0.16 to 0.05 eV. The hydrogenated silicene with DV2 defect shows semiconductor behaviour with band gap of 0.07 eV [see Fig. 6(d)]. Therefore, the dissociative adsorption of H₂ molecule can significantly tune the band structures of pristine and defective silicene.

Fig. 6 The band structure of silicene after the dissociation of the H₂ molecule: pristine silicene (a), defective silicene with single atom vacancy SV (b), double atom vacancy DV1 (c), and double atom vacancy DV2 (d).

4. Conclusion

In conclusion, the dissociative adsorption of hydrogen molecule on pristine and defective silicene is studied by using DFT calculations. The results show that the dissociative energy barrier for hydrogen molecule on two dimensional silicene can be significantly depressed with all kinds of the studied defects. Particularly, the energy barrier of 0.71 eV for the dissociation of H₂ molecule on single vacancy defect SV1(55|66) is smaller than the critical energy barrier of 0.91 eV, which indicates that the dissociative adsorption of H₂ molecule can proceed smoothly at ambient temperature. In addition, the band gaps of the pristine and defective silicene can be tuned through hydrogenation significantly. Our work provides meaningful insights for silicene applied in electronic devices.

Acknowledgements

This work was supported by the Fundamental Research Funds for the Central Universities (Grant No. 2014B13414 and 2015B01914), China Postdoctoral Science Foundation (2015M571652), National 973 Plan Project (2015CB057803) and Natural Key Foundation of Jiangsu Province (BK2011025). ZA acknowledges the financial supports from “100 talents” program of Guangdong University of Technology, “1000 plan” for young professionals program of Chinese Government.

References

1. C. C. Liu, W. Feng and Y. Yao, Phys. Rev. Lett., 2011, 107, 076802.
2. P. Vogt, P. D. Padova, C. Quaresima, J. Avila, E. Frantzeskakis, M. C. Asensio, A. Resta, B. Ealet and G. L. Lay, Phys. Rev. Lett., 2012, 108, 155501.
3. B. J. Feng, Z. J. Ding, S. Meng, Y. G. Yao, X. Y. Cheng, P. He, L. Chen and K. H. Wu, Nano Lett., 2012, 12, 3507–3511.
4. A. Fleurence, R. Friedlein, T. Ozaki, H. Kawai, Y. Wang and Y. Yamada-Takamura, Phys. Rev. Lett., 2012, 108, 245501.
5. L. Meng, Y. L. Wang, L. Z. Zhang, S. X. Du, R. T. Wu, L. F. Li, Y. Zhang, G. Li, H. T. Zhou, W. A. Hofer and H. J. Gao, Nano Lett., 2013, 13, 685–690.
6. M. R. Tchalala, H. Enriquez, A. J. Mayne, A. Kara, S. Roth, M. G. Silly, A. Bendounan, F. Sirotti, T. Greber, B. Aufray, G. Dujardin, M. A. Ali and H. Oughaddou, Appl. Phys. Lett., 2013, 102, 083107.
7. S. Cahangirov, M. Topsakal, E. Aktürk, H. Şahin and S. Ciraci, Phys. Rev. Lett., 2009, 102, 236804.
8. Z. Y. Ni, Q. H. Liu, K. C. Tang, J. X. Zheng, J. Zhou, R. Qin, Z. X. Gao, D. P. Yu and J. Lu, Nano Lett., 2012, 12, 113–118.
9. N. Gao, W. T. Zheng and Q. Jiang, Phys. Chem. Chem. Phys., 2012, 14, 257–261.
10. N. Gao, J. C. Li and Q. Jiang, Chem. Phys. Lett., 2014, 529, 222–226.
11. Q. G. Jiang, J. F. Zhang, Z. M. Ao and Y. P. Wu, Sci. Rep., 2015, 5, 15734.
12. Q. G. Jiang, J. F. Zhang, Z. M. Ao and Y. P. Wu, J. Mater. Chem. C, 2015, 3, 3954–3959.
13. A. Kara, H. Enriquez, A. P. Seitsonen, L. C. Lew Yan Voon, S. Vizzini, B. Aufray and H. Oughaddou, Surf. Sci. Rep., 2012, 67, 1–18.
14. G. Brumfiel, Nature, 2013, 495, 152–153.
15. L. Tao, E. Cinquanta, D. Chiappe, C. Grazianetti, M. Fanciulli, M. Dubey, A. Molle and D. Akinwande, Nat. Nanotechnol., 2015, 10, 227–231.
16. F. Banhart, J. Kotakoski and A. V. Krasheninnikov, ACS Nano, 2011, 5, 26–41.
17. A. Hashimoto, K. Suenaga, A. Gloter, K. Urita and S. Iijima, Nature, 2004, 430, 870–873.
18. A. V. Krasheninnikov and F. Banhart, Nat. Mater., 2007, 6, 723–733.
19. Y. Kim, J. Ihm, E. Yoon and G. D. Lee, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 84, 075445.
20. S. Huh, J. Park, K. S. Kim, B. H. Hong and S. B. Kim, ACS Nano, 2011, 5, 3639–3644.
21. S. U. Yu, B. Park, Y. Cho, S. Hyun, J. K. Kim and K. S. Kim, ACS Nano, 2014, 8, 8662–8668.
22. J. Gao, J. Zhang, H. Liu, Q. Zhang and J. Zhao, Nanoscale, 2013, 5, 9785–9792.
23. V. O. Özçelik, H. H. Gurel and S. Ciraci, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 88, 045440.
24. H. Sahin, J. Sivek, S. Li, B. Partoens and F. M. Peeters, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 88, 045434.
25. G. R. Berdiyorov and F. M. Peeters, RSC Adv., 2014, 4, 1133–1137.
26. S. Li, Y. F. Wu, Y. Tu, Y. H. Wang, T. Jiang, W. Liu and Y. H. Zhao, Sci. Rep., 2015, 5, 7881.
27. H. H. Gurel, V. O. Ozcelik and S. Ciraci, J. Phys. Chem. C, 2014, 118, 27574–27582.
28. B. Delley, J. Chem. Phys., 2000, 113, 7756–7764.
29. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865.
30. A. Roldán, J. M. Ricart and F. Illas, Theor. Chem. Acc., 2009, 123, 119–126.
31. T. A. Halgren and W. N. Lipscomb, Chem. Phys. Lett., 1977, 49, 225–232.
32. G. Henkelman and H. Jonsson, J. Chem. Phys., 2000, 113, 9978–9985.
33. S. Grimme, J. Comput. Chem., 2006, 27, 1787–1789.
34. P. O. Lehtinen, A. S. Foster, Y. C. Ma, A. V. Krasheninnikov and R. M. Nieminen, Phys. Rev. Lett., 2004, 93, 187202.
35. Q. G. Jiang, Z. M. Ao, W. T. Zheng, S. Li and Q. Jiang, Phys. Chem. Chem. Phys., 2013, 15, 21016.
36. M. Z. Hosssain, H. S. Kato and M. Kawai, J. Am. Chem. Soc., 2005, 127, 15030–15031.
37. T. H. Osborn, A. A. Farajian, O. V. Pupysheva, R. S. Aga and L. C. Lew Yan Voon, Chem. Phys. Lett., 2011, 511, 101–105.
38. D. C. Young, Computational Chemistry: A Practical Guide for Applying Techniques to Real World Problems, Wiley, New York, 2001.
39. I. Mayer, Int. J. Quantum Chem., 1984, 26, 151–154.
40. Z. M. Ao and F. M. Peeters, Appl. Phys. Lett., 2010, 96, 253106.

---