Classical molecular dynamics simulation on the dynamical properties of H 2 on silicene layer

This study investigates the diffusion of hydrogen molecule physisorbed on the surface of silicene nanoribbon (SiNR)using the classical molecular dynamic (MD) simulation in LAMMPS (Large-scale Atomic/Molecular Massively Parallel Simulator). The interactions between silicon atoms are modeled using the modified Tersoff potential, the Adaptive Intermolecular Reactive Empirical Bond Order (AIREBO) potential for hydrogen – hydrogen interaction and the Lennard – Jones potential for the physisorbed H2 on SiNR. By varying the temperatures (60 K 130 K), we observed thattheΔxdisplacement of H2 on the surface SiNR shows a Brownian motionon a Lennard-Jones potential and a Gaussian probability distribution can be plotted describing the diffusion of H2. Thecalculated mean square displacement (MSD) wasapproximately increasing in timeand theactivation energybarrier for diffusion has been found to be 43.23meV.


Introduction
The search for novel materials is widespread nowadays in response to the advancing world of science and technology.As the field of research goes nanotechnology, many materials, those that are experimentally impossible in the past, have been synthesized at present.These materials are able to adapt and are technologically fit for what is required of them.
One of these materials is graphene [1], a two-dimensional (2-D) honeycomb structure material composed of carbon atoms, which exhibits properties that are very unique and exceptional.Since the production of graphene [1][2], many single-layered materials have been synthesized and theorized over the years.Some of the recent researches investigate the interaction of hydrogen on the surface of these 2-D materials.Due to the increasing demand for energy, most of the studies focuses on finding an efficient hydrogen storage material [3][4][5][6] which requires a binding energy of 0.2 ~ 0.6 eV.Applications to microelectronics using these materials are also considered since these materials are semimetals and can have a tunable band gap when saturated with hydrogen.
Studies on the interaction of the hydrogen atom on the surface of graphene were conducted through MD simulation and density functional theory (DFT) [5].Furthermore, an investigationusing classical MD simulation on the static and dynamical properties of molecular hydrogen on the surface of graphenepresented the calculation of the activation energies for the diffusion and desorption of H 2 [6].
One of the recent discoveries is silicene [7], an allotrope of silicon with a pure sp 2 hybridization, one-atom thick, honeycomb structure material similar to that of graphene.Since it has a graphenelike structure, it is thought to possess properties that resemble that of graphene.Since silicene is a silicon-based material, it can easily be integrated with the present silicon industry.Unlike graphene with a planar structure, silicene structure is slightly buckled and due to this buckling, silicene exhibits properties that are unique and different from graphene [8].The interaction of hydrogen on silicene and other 2-D materials has attracted many since hydrogenation is thought to affect the properties of the material like band gap tuning and gas sensing.While most of the studies relating to the H 2 interaction focuses mainly on graphene [7][8][9], only few have considered the hydrogen interaction on silicene [10] and most of it uses the first principle approach, although has higher accuracy as compared to classical approaches but a computationally intensive simulation.
In the present work,we investigated the interaction of H 2 molecule on the surface of the silicene nanoribbon (SiNR diffusion barrier energy has been found to be 43.23 meV.The interaction of molecular hydrogenon silicene is an on-going study at present and this study has its intent to pave more ways, shedding new light for moretechnological applications of silicene.

Computational Method
Dynamical properties of H 2 on silicene were investigated using classical molecular dynamics (MD) simulation [11] through LAMMPS (Large-scale Atomic/Molecular Massively Parallel Simulator).The interactions between silicon atoms were modeled using the modified Tersoff potential [12] and the interaction between the two hydrogen atoms were modeled using the AIREBO potential [6].The hydrogen molecule is considered to be physisorbed on the surface of silicene and thus modeled with Lennard-Jones potential.The Lennard -Jones potential for the diffusing H 2 on silicene is given by the equation [13], wherer is the distance between H 2 and silicene, r c is the cut off distance, ߳ is the measurement on how strongly the H 2 and silicene attract, and the distance where there is no intermolecular potential is given byߪ = 2.5 Å.
The simulation box for the silicene has a periodic boundary condition.The dimension of SiNR is 50 Å × 50 Å, a square nanoribbon which contains 418 atoms which lie flat on the simulation box and placed at z = 5Å due to its buckled structure.The hydrogen atoms were placed 2.454 Å [13] at the top of the hollow site as shownin Fig. 1.The distance of H 2 from the surface of silicene in the hollow site is shorter when compared to the other adsorption sites.This is due to the stronger binding energy on the other sites as compared to the hollow site [14].The simulation uses an NVT ensemble and a Brendsen thermostat to regulate the temperature of the system during the simulations.A timestep of 0.2 fs was used to investigate the dynamical properties of H 2 on silicene during the simulation.In this study, low temperatures were used ranging from 60K to 130K, with a total number of 1 × 10 5 steps to equilibrate the system and another 1 × 10 8 steps to take the averages.For low temperatures, simulation times were increased to minimize statistical errors [6].Since the boundaries on every side of the simulation box are periodic, the adsorption, diffusion and desorption of H 2 from the surface of silicene were observed and were tracked during the simulation.The displacement of the H 2 atoms in three dimensions was tracked.
The calculation of the diffusion coefficient, D at any given temperature is calculated through the mean square displacement (MSD),〈‫ݎ‬ 2 〉 [15][16], where the mean square displacement in two dimensions calculated in every temperature is given as [15], The coordinates x and y corresponds to the center of mass of the H 2 molecule.The activation energy, E a was calculated using the relationship of the temperature and the diffusion coefficient, D through the empirical formula: [15][16] The plot of the logarithm of the diffusion coefficient versus the inverse of the temperature gives a linear fit from which the activation energy [17] is calculated.
The probability of finding the H 2 moleculefrom its original position on the SiNR surface is calculated using the Gaussian probability density distribution function [18], employing the Δx displacement as the random variable representative to the H 2 coordinate, (5) whereߪis the standard deviation and ߤ is the mean.
02022-p.2In Figure 2, the Δxdisplacement, a representative of the H 2 coordinate, shows the Brownian motion behavior on a Lennard-Jones potential as H 2 moves on the SiNR surface with varying temperatures.The probability of finding the H 2 molecule from its original position is shown in Figure 3. Through the Δx displacement, a Gaussian probability distribution from equation ( 5) can be plotted describing a diffusing molecule.At lower temperatures, there is a high probability that H 2 can found at its original position,and with increasing temperatures, there is a low probability that H 2 can be found and diffuses randomly away from its original position.This only demonstrates that the energy of the H 2 molecule has increased with temperature causing it to be more diffusive on the surface.In Table 1, the diffusion coefficients were calculated from the MSD.The plot of the logarithm of the calculated diffusion coefficient versus the reciprocal of the temperature is shown in Figure 5.This shows the temperature dependence of the diffusion of the H 2 on the surface of SiNR.Based on equation (4), the activation energy barrier for diffusion which was calculated from the slope of the graph is43.23meV with the calculated preexponential D 0 = 1.942 × 10 5 Å 2 /ps.The activation energy barrier for the diffusion of H 2 on SiNR is greater compared to the diffusion barrier calculated for H 2 on graphene [6] which has the value of 9.8 meV.Factors affecting the diffusivity of H 2 on these materials are also considered.The structure of silicene is buckled, unlike graphene with a planar structure, affects the behavior of H 2 on its surface.This buckling breaks the symmetry of the silicene causing ripples on the surface [19].The difference in the binding energy of both materials may affect the diffusivity of H 2 .Silicene has shown a greater binding energy of 94 meV [20] than that of graphene, 25 meV per H 2 molecule.Based on our calculations, H 2 tends to be more stable on silicene than in graphene.H 2 is less likely to be diffused on the SiNR surface perhaps due to its buckled structure which may affect the interaction by increasing the activation energy barrier for diffusion through the change in the parameters of the Lennard-Jones potential.This calculation may be of help to some possible applications of silicene like hydrogen sensor since it has higher sensitivity on H 2 as compared to other gases or a hydrogen separator [21].This could also be seen in the anomalous increase of the thermal conductivity of silicene as the number of hydrogen increases on the surface of silicene.[22] Since the production of a free-standing silicene is still experimentally impossible, it is now a challenge as to how will we unlock more properties of silicene that will support its production and application when the time comes when all the impossibilities become a possibility.

Summary
The diffusion of H 2 molecule on the surface of silicene was investigated using the classical molecular dynamics simulation in LAMMPS.Various temperatures ranging from 60K-130K were introduced during the simulation to record of the Δx displacement of H 2 on the SiNR surface.The Δx displacement shows the Brownian motion on a Lennard-Jones potential.Through the Δxdisplacement, a Gaussian probability distribution can be plotted describing the diffusion of the H 2 molecule.At lower temperatures, there is a high probability that H 2 can found at its original position, and with increasing temperatures, there is a low probability that H 2 can be found and diffuses randomly away from its original position.A linear fit is plotted from the average MSD showing that it increases in time.The diffusion coefficients were calculated andthe activation energy barrier for diffusion has been found to be 43.23meV.With this investigation, a need arises as to till more ground in the study of H 2 interaction on silicene for its future applications considering that the production of a free-standing silicene is still experimentally impossible at present.

Figure 1 .
Figure 1.Top view of the silicene structure (yellow) topped with the hydrogen molecule (green) on the hollow site (left) and the distance between silicene and H 2 molecule (right).

Figure 2 .
Figure 2.The Δxdisplacement of H 2 on SiNR with a time step of 0.2 fs per simulation run.

Figure 3 .
Figure 3.Gaussian distribution calculations of the displacement of H 2 at various temperatures.The mean square displacement (MSD) for H 2 on the surface of SiNR was recorded during the simulation for the calculation of the diffusion

Figure 4 .
Figure 4.The mean square displacement of H 2 on SiNR for the calculation of the average diffusion coefficient of various temperatures.

Figure 5 .
Figure 5.The temperature dependence of the diffusion coefficient of H 2 on SiNR.

Table 1 .
The data of the diffusion coefficient for the hydrogen diffusion on SiNR at various temperatures.