Investigation of the Elastic Properties of Graphenylene Using Molecular Dynamics Simulations

Molecular dynamics simulations are used here to study the mechanical behavior of graphenylene under uni-directional and bi-directional loadings. The effects of nanosheet chirality and size on Young’s modulus of graphenylene are investigated. Compared to graphene, graphenylene possess a smaller elastic modulus. It is shown that for large armchair and zigzag graphenylenes, the effect of nanosheet size on the mechanical properties can be neglected. It is observed that increasing temperature results in decreasing Young’s modulus of graphenylene. Besides, fracture of graphenylene occurs at large strains. Moreover, it is represented that for small graphenylenes, bulk modulus is significantly sensitive to the size variation. However, this sensitivity disappears for large nanosheets.


Introduction
Since the discovery of graphene by Novoselov et al. in 2004 1 , this carbon allotrope has attracted the attention of many researchers.This can be attributed to its great physical, chemical, and mechanical properties [1][2][3] .Due to these great properties, graphene has found potential applications in several areas such as nanoelectronics, sensors, transistors, batteries, as well as polymer composites [4][5][6] .
Recently, Song et al. 22 studied the properties of graphenylene as the first example of a non-delocalized sp 2 -carbon structure.They showed that graphenylene has periodic pores with the diameters of 3.2 Å. Besides, it was claimed that it is a semiconductor with a narrow direct band gap which leads to its applicability in electronic devices and efficient hydrogen separation.Using ab initio density functional theory (DFT) calculations, Yu 23 illustrated that graphenylene can store lithium with great density of energy.He showed that graphenylene can adsorb a lithium atom stronger than pristine graphene.Hankel and Searles 28 showed that the lithium storage capacity of graphenylene is larger than graphite.Andrew T. Koch et al. 29 used DFT and tight-binding methods to investigate the graphenylene based nanotubes.
Mechanical properties of graphenylene are studied here by employing molecular dynamics (MD) simulations.To this end, the effect of some parameters such as chirality and size on Young's and bulk moduli of graphenylene are investigated.Besides, the fracture process of armchair and zigzag nanosheets under uni-directional and bi-directional tensions are explored.

Models
The structure of a sample graphenylene is shown in Figure 1.To obtain Young's modulus of graphenylene, the loads are applied in two directions, namely armchair and zigzag, which are shown schematically in Figure 1.As it is seen in Fig 1(b) and (c), for zigzag graphenylene, more C-C bonds are directed along the loading direction (the ellipses have been used to show the directed bonds along the loading direction).

Simulation method
To obtain the mechanical properties of graphenylene, MD simulations are performed by LAMMPS MD code 30,31 .
To define the interactions between the carbon atoms in the structure of graphenylene, Adaptive Intermolecular Reactive Empirical Bond Order (AIREBO) potential function 32 is employed.This potential function is the improved version of Brenner's second-generation reactive empirical bond order potential 33 in which the Lennard-Jones term has been added to include van der Waals interactions.Moreover, torsional term is included to consider torsions of σ bonds.So it is anticipated to give more accurate results than Brenner's

Results and discussion
To study the effect of atomic structure on the mechanical properties of graphenylene, armchair and zigzag nanosheets are modeled here.Moreover, considering different side lengths (a) and widths (b), the effect of graphenylene geometry on the mechanical properties are studied.Figure 2 shows armchair and zigzag nanosheets with the corresponding geometrical parameters.The loads are applied on the upper and lower edges of nanosheets.In each loading step, the boundary atoms are displaced and then the systems are let to relax for 2ps (2000 steps with the time step of 1fs).So, the stress would be distributed all over the nanosheets and the quasistatic condition of the tensile tests are simulated.All of the simulations are performed under NVT (constant number of molecules, constant volume, and constant temperature) ensemble.Besides, to ensure the stability of the simulations, Nose-Hoover thermostat is employed here.

Young's modulus
To obtain the mechanical properties of graphenylene, the lower sides of the nanosheets are fixed and upper sides are displaced incrementally.In each step, the strain and the system energy are recorded.As a unique value has not been reported for the thickness of nanostructures [34][35][36][37] , the in-plane Young's modulus is used here which is obtained as 38 : where S E , 0 A and ε are strain energy of system, initial unit cell area and strain, respectively.The recorded strain energies of armchair and zigzag graphenylenes with the dimensions of 30×50 Å 2 are shown in Figs. 3 and 4 against strain.The drops of the curves show the bond breakages.As it is seen, the fracture occurs in different steps.In each step, some of the bonds break.This fact can be verified with the graphs of temperature versus strain which are shown in Figs. 5 and 6.As it is seen, each bond breakage results in a pick in the temperature of the systems.
The gradual breakage of the armchair and zigzag nanosheets are also represented in Figs.7 and 8.It is seen that in each fracture step, some of the bonds are broken down.Finally, by connecting the cut bonds, a line of atoms is formed and stretched.For both of the nanosheets, the boundary effect is clearly seen.For armchair graphenylene, this effect results in breakage near the end of nanosheet.However, for zigzag graphenylene, although the some of the bonds near the left boundary are broken, the final fracture occurs at the middle of nanosheet.Moreover, comparing the figures at different steps, one can conclude that at each step, by breaking some bonds, some other new bonds are formed.The large fracture strain of graphenylene can be attributed to the formation of these new bonds.
Figure 9 depicts in-plane Young's moduli of armchair and zigzag graphenylenes with the aspect ratio of 1.5 versus nanosheet side length.The sizes of the considered nanosheets are given in Table 1.Since the dimensions can not be completely equal, the closest dimensions are considered.As it is seen at a same size, zigzag graphenylenes have larger elastic modulus than armchair ones.However, the difference can be neglected for side lengths larger than 50Å.Besides, increasing side length of graphenylene at constant aspect ratio result  Represented in Figure 10 are in-plane Young's moduli of armchair and zigzag graphenylene with the side length of 50 Å versus graphenylene aspect ratio.The exact sizes of the considered nanosheets are given in Table 2.As the previous figure, it is seen that increasing aspect ratio of the nanosheets results in decreasing Young's modulus.Both of the effects of aspect ratio and chirality (being armchair and zigzag) diminishes for the aspect ratios larger than 2.
Young's moduli of 30×50 Å 2 armchair and zigzag graphenylenes are obtained as 68.31 N/m and 79.42 N/m, respectively, which are approximately 25% of Young's modulus of graphene (272 N/m).The ultimate strength of armchair and zigzag 30×50Å 2 graphenylenes are 429.43nNand 1140.59nN,respectively.Comparing these values with the ultimate strength of graphene (1070 nN), it is seen that while the ultimate strength of armchair graphenylene is 40% of graphene, zigzag graphenyle has a larger ultimate strength than graphene.Besides, the large difference between the ultimate strength of armchair and zigzag nanosheet is interesting.The yield and ultimate strains of 30×50 Å 2 armchair graphenylene are obtained as 14.43% and 94.46%, respectively.These values are 15.43% and 101.26% for zigzag graphenylene.Note to the large ultimate strains of armchair and zigzag graphenylene which were shown schematically in Figs.7 and 8. Poisson's ratio of armchair and zigzag graphenylenes were obtained as 55% and 73%, respectively which are smaller than Poisson's ratio of graphene (approximately 88%).
The effect of temperature on the elastic modulus of the grapheneylene is given in Figure 11 for a 40×60 Å 2 armchair nanosheet.Unlike Figs. 10 and 11 in which the dependence of the grapheneylene elastic modulus to the nanosheet side length and aspect ratio were investigated, no uniform behavior is observed by increasing the temperature.Therefore, a linear interpolation is used to view the general trace of the graphenylene elastic modulus by time variation.It is observed that nanosheets possess smaller Young's modulus at larger temperature.It can    ( ) be said that increasing the temperature results in higher kinetic energies of the atoms of nanosheets which results in smaller strength of nanosheets.

Bulk modulus
Bulk moduli of graphenylene are computed by bi-directional loading on nanosheets.Here, one horizontal edge and one vertical edge are constrained and the loads are applied in the form of displacement to the opposite edges.The following relation can be used to obtain the nanosheet bulk modulus: in which A is the nanosheet area.The configurations of a 50×50 Å 2 graphenylene under bi-directional loading at different timesteps are given in Figure 12.It is seen that the fracture initiates at the corners of the nanosheet and then propagates linearly parallel to the graphenylene boundaries.Linear propagation of fracture leads to breaking some bonds and formation some new bonds.Interestingly, bond breakage and formation of new bonds happen in the center of the nanosheet which is not in the fracture region.Besides, the fracture happens perpendicular to the zigzag direction and near the graphenylene boundary.
The computed bulk moduli of square graphenylenes with different side lengths are given in Figure 13.Table 3 represents the exact sizes of the selected nanosheets.It is seen that the effect of side length on the graphenylene bulk moduli is significant at the small side lengths.However, after 100Å, the graph take a linear form and the effect of graphenylene size on the bulk modulus can be neglected.

Conclusions
Employing MD simulations, mechanical properties of graphenylene were computed.Applying unidirectional displacement on the armchair and zigzag graphenylene, Young's moduli of nanosheets with different sizes were obtained.It was shown that at a same size, zigzag graphenylenes possess larger Young's moduli than armchair.The effect of chirality was disappeared for large graphenylenes.Tracking the fracture process of armchair and zigzag graphenylenes, it was seen that the fracture happens at different steps in each of them some bonds are destroyed and some new ones are created.The boundary effect was observed for all unidirectional and bi-directional loading of graphenylenes.Finally, the bulk moduli of square graphenylenes were computed.

Figure 1 :
Figure 1: (a) Schematic of graphenylene with the armchair and zigzag direction, (b) Loading direction and C-C bonds detected along the loading direction for armchair and (c) zigzag nanosheets

Figure 2 :
Figure 2: (a) Armchair and (b) zigzag graphenylene with the corresponding geometrical parameters

Figure 9 :
Figure 9: In-plane Young's modulus of armchair and zigzag graphenylenes with the aspect ratio of 1.5 versus graphenylene side length.

Figure 10 :
Figure 10: In-plane Young's modulus of armchair and zigzag graphenylenes with the side length of 50 Å versus graphenylene aspect ratio.

Figure 11 :
Figure 11: Effect of temperature on the elastic modulus of the grapheneylene for a 40×60 Å 2 armchair nanosheet.

Table 1 :
Dimensions of armchair and zigzag nanosheets with the aspect ratio of 1.5

Table 2 :
Dimensions of armchair and zigzag nanosheets with the aspect side length of 50 Å