An Innovative Look at Buckling of Cracked Ghraphene Nano Sheets Using Novel Extended Molecular Mechanics Method

In this paper, we take a creative look at the effect of cracking on the buckling of graphene sheets. We introduced a significant improvement in molecular mechanical calculations that more accurately examine the behavior of atomic bonds during the instability of nanostructures in the presence of defects. A molecular mechanics method based on the modified couple stress (MCS) theory to consider the size effect on the graphene plate is introduced. Primary cracks with different lengths, orientations, and numbers were examined. For each crack, buckling parameters such as critical load and buckling modes of graphene sheets were investigated. This model has better computational performance and less computational compactness than other models. Our studies show that cracks perpendicular to the load direction have a greater effect on buckling parameters, including modes and critical buckling load than cracks aligned to the critical load.

Increasing the crack length will further reduce the critical buckling load. The effect of the number of cracks in two directions in buckling parameters was investigated. Increasing the number of cracks always decreases a critical buckling load. Increasing the number of cracks affects the first buckling mode more than others. The study of the eccentricity of cracks shows that with increasing the eccentricity distance, in cracks perpendicular to the loading direction, the critical buckling load will be further reduced. However, for cracks aligned with the load direction, increasing the eccentric distance has less effect on reducing the critical buckling load.

Keywords: Size Effect, Novel Molecular Mechanics, Graphene Sheets, Buckling, Crack

Graphene as a nanostructure has been considered in many studies due to its many applications. In many applications, graphene may be exposed to loads and become unstable and fail. There are researches have analyzed the instability and properties of graphene and nanotubes [1-3]. The buckling of graphene has been investigated by different methods. For example, Genoese et al. [4] used the molecular mechanics method to investigate buckling. They investigated the nonlinear buckling of a graphene plate. They calculated multiple interactions between atoms from the molecular mechanics model by introducing the Morse potential functio0n. The effect of the orientation of the atoms on the graphene plate has been studied. By studying the equilibrium equations of the atomic structure through the arc length strategy, they studied the critical behaviors of graphene under compressive and shear loading. By introducing two-dimensional potentials, they modeled interatomic interactions.

Furthermore, Sakhaee-Pour [5] investigated buckling in elastic monolayer graphene sheets using structural molecular mechanics. He found that the elastic buckling per SLGS width varies nonlinearly concerning the side length. It also showed that in an equivalent geometry, the buckling load of the zigzag sheet is larger than the armature method. The statistical nonlinear regression models were used to develop the prediction equations. These equations provide a quick tool for calculating buckling force while their results remain 5% different from atomic modeling data.

The elastic instability of bilayer graphene sheets is studied with Chandra et al. [6]. Different boundary conditions and ratios are considered for bilayer and monolayer graphene sheets. Double-layer graphene sheets have higher buckling strength compared to single-layer sheets. The results showed that the higher stiffness of bilayer graphene sheets leads to a higher critical buckling load. The critical buckling load increases significantly with the increasing number of graphene layers.

The buckling of concentrated carbon nanotubes and the elastic properties were investigated with Eltaher et al. [7]. An equivalent energy model with a mechanical structural approach was presented to describe the atomic interactions and chemical energy. The bonds between two adjacent carbon atoms are considered as an isotropic beam with a uniform circular cross-section.

As another paper, Korobeynikov et al. [8] studied the bonding interactions between carbon atoms in single-layer graphene sheets by developing geometric parameters and other characteristics of the Bernoulli-Euler beam. In their research, they tried to present a new version of the molecular mechanics method. Using the modified molecular mechanics, the behavior of the graphene plate during vibration and buckling is investigated. They used the DREIDING force field.

In all cases, graphene was considered perfect. In many applications, graphene is subjected to conditions in which several bonds are broken. Some researchers have studied the effect of defects on graphene [9, 10]. For example, Dewapriya et al. [11] performed molecular dynamics simulations to investigate the interaction between boron nitride atoms and graphene with an edge crack in graphene. They examined the complex nature of the stress state at the crack tip. The simulation results investigated the effect of the underlying crystal structure of two-dimensional materials on crack stress fields and illustrated their fracture toughness. In their study, they examined cases such as the presence of gaps near the tip of the crack. Their results show that the crack can be stopped by placing the symmetrical holes at the nanoscale. For the first time, their research showed the growth of multistage graphene cracks and crack repair under high levels of strain applied in the presence of reciprocal atomic holes, which opens a new possibility for graphene strain engineering with crack control.

Zhigong et al. [12] studied the role of geometric distortion in predicting the mechanical response of materials in small dimensions. The Griffith strength of graphene was examined in the presence of distortion using molecular simulation. Their results showed that in the first mode of failure, distortion has no significant effect on reducing tensile strength. While for the second failure mode, the decrease was very significant. Tapia et al. [13] studied graphene with voids using the atomic finite element method. They showed that when the gap under the strain gradient is larger, it causes a greater reduction in shear modulus and Poisson's ratio. They determined that the relationship between the elastic modulus and the number of vacancies was relatively linear.

The possibility of predicting fracture and crack propagation in graphene are investigated by Theodosiou et al. [14] using molecular mechanics. Using the ABAQUS code library, they applied an extended method of molecular mechanics model. Graphene failure was investigated in all three failure modes. Their results showed that in the first mode, graphene is more flexible and more sensitive than in the third mode. Parashar et al. [15] also used a finite element atomic model to determine the refractive properties in the first mode of single-layer graphene. They estimated the energy release rate in the cracked graphene plate. The model has less computational compactness than others due to the lower number of degrees of freedom. Using the Griffith criterion in linear elastic fracture mechanics, Xin-Liang et al. [16] attempted to calculate the graphene fracture toughness with initial crack velocity. They also investigated the elastic properties of cracked graphene. The truss finite element models and a two-dimensional finite element plate are used by Tserpes et al. [17] to obtain the critical buckling stress for a graphene plate. Then, using the molecular mechanics model, they investigated the effect of vacancies on the critical buckling load and buckling states. They found that as the number of gaps in the graphene plate increased, the critical buckling stress would be greatly reduced.

As stated in the research background, different methods have been used to investigate the buckling on the graphene plate. Among these methods, the molecular mechanics approach can predict the buckling behavior of graphene sheets to an acceptable extent. But the classical molecular mechanics method have not considered the effect of nanoscale size in calculations. In the present study, considering the size effect at the nanoscale, we presented a new version of the molecular mechanics method. The following shows that this version has more computational accuracy than others.

On the other hand, none of the previous studies has investigated the effect of cracking on the buckling parameters of graphene sheets. We used the newly introduced method to study the instability of the graphene plate in the presence of different cracks. In the following, the effect of crack on buckling parameters such as critical buckling load and buckling modes are estimated. First, the graphene with a crack in the center of the sheet is investigated. Then, the effect of the distance between the cracks is investigated. In the next step, the number of cracks are increased and the effect of 3, 5, 7, and 9 cracks on the buckling parameters of the graphene plate are studied. Finally, the cracks are placed in the graphene horizontally and vertically at a distance from the center and investigate the effect of eccentricity on the buckling modes and the critical buckling load. Our results show that as the number of cracks increases, the critical load of the cracked graphene plate will decrease. Our studies show that cracks perpendicular to the load axis have a greater impact on the critical buckling load. Increasing the crack length also makes the graphene plate easier to buckle and reduces the critical buckling load. Also, for crack that is aligned to the load axis, the critical load does not change much at shorter crack lengths. But by increasing the crack length more than 65% length of the graphene plate, the critical load will suddenly decrease. But for a crack perpendicular to the load, increasing the tributary distance from the center causes the critical buckling load to be significantly reduced, even during the crack constant.

In the modified molecular mechanics method and the spatial structure of carbon atoms, atomic bonds in graphene sheets are considered as elastic beams. This beam element has loads such as flexural tension and torsional loads (Figures 1-a, 1-b, 1-c)

The bending moment will be created in the beams, due to the structure of the carbon rings. Using the modified couple stress theory, a new demonstration of the bending of atomic bonds is presented. At this way, the size effect into the equations of molecular mechanics is assumed. It is named the new equations as MCS-MM. By introducing the deformation created in the beams as a function of the length scale parameter, parameters such as tensile strength, bending, and torsion are calculated. First, the total energy of the beam despite all three changes is investigated. Beam strain energy is:

According to the theory of classical molecular structural mechanics, the strain energy of a uniform beam of length L subjected to pure axial force is:

By equating the energy in Equation 2 with stretching energy equation [18], , tension stiffness kr calculated as:

Two methods are used to calculate the strain energy due to bending in the structure of Figure 2-b. Then, with equal strain energy, the bending stiffness of the beam are calculated. The first way is to obtain the chemical strain energy in exchange for changing the angles shown in the figure 2 for a bond. The chemical strain energy according to the figure 2, is:

Due to the small angles, the relation within the bond angle change is, Therefore, the flexural energy of a structure according to Figure 2-b is obtained as follows:

In the second method, the strain energy due to bending in the beam assembly is obtained using the strain energy distortion based on the modified couple stress theory. In a linear elastic material, using the modified coupling stress theory, the strain energy density contains the classical stresses and strains of the coupling. According to reference [19], the strain energy for deformation in the Ω region is:

Where, the stress tensor, σ, strain tensor, ε, deviatoric part of the couple stress tensor, m, and symmetric curvature tensor, χ, are, respectively, defined by:

With, λ and μ being Lame’s constants, l a material length scale parameter, u the displacement vector and θ the rotation vector. To consider the size effect in calculating strain energy, the coupling part of the stresses and strains due to the curvature of the beam are calculated. Considering dv = Adx and using the equations 7 in equation 1, the strain energy of the beam will be equal to:

Where, For pure bending, the moment and deflection of the beam were related as:

Then, the deflection of “w”, calculated as:

The strain energy of a uniform beam considering deflection “w” from equation 10 leads to:

By applying the zero-angle change condition at point B, the relation between force and moment on the beam, calculated as:

According to the deformation diagram created in the beam in Figure 2, the following equation can be introduced between the small angle change and the created displacement at the end of the beam:

Finally, the strain energy for bending calculated as:

By calculating the strain energy of both methods, it is observed that the bond energy of the chemical method is equivalent to the energy of the MCS-MM method. Finally, the elastic strain energy for the set shown in Figure 2 will be:

The flexural stiffness of the beam, regarding the modified couple stresses, is computed as:

The size effect in the modified couple stress theory is included in the study of the torsional behavior of the beam. Torsional torque at the end of the beam is defined by considering the modified couple stress parameters as follows [20]:

Where, is the resultant torque in classical Cauchy elasticity, and Tμ is couple part of torque. The strain energy due to torsion is:

The torsional stiffness calculated as:

One of the reasons for considering the modified coupling stress theory is considering the size effect in determining the beam properties as an alternative to bonds between carbon atoms. Using stiffness constants related to each type of deformation in the beam is:

The length of the length scale parameter is present in the beam diameter determination equation and indicates the presence of the size effect in determining the beam properties. For the Bernoulli-Euler beam, the Poisson effect is small and negligible, and the equations are simplified [21]. The Young's modulus of the beam is also obtained using the stiffness parameters as follows:

The shear modulus, as a variable of the stiffness constants and the length scale parameter considering equations 3, 16, 19, is calculated as follows:

The numerical values of these stiffness constants, k_r,k_θ,k_∅, for the carbon-carbon bond are derived from reference [8]. To estimate the length scale parameter, l, the strain energy of graphene plates for different values of "l" must be calculated. As shown in the beam properties equations, the geometric and elastic properties of the beam depend on the l-parameter. For each of the proposed values of l, the amount of strain energy using the new method is obtained. On the other hand, one of the best and most realistic analyzes is the molecular dynamics method. In molecular dynamics models, intermolecular interactions are appropriately considered. To verify the results of the new model, the results compared with molecular dynamics simulations. The amount of strain energy from the molecular mechanics method are selected that had the least deviation from the molecular dynamics results and the corresponding length scale parameter as the appropriate length scale parameter for that particular problem.

The molecular dynamics simulation was applied to estimate the behavior of the graphene sheets under various conditions. The AIREBO potential are used to define the atomic interactions between carbon atoms in graphene. The AIREBO potential can represent by a sum over pairwise interactions, including covalent bonding interactions, Lennard-Jones potential terms, and torsion interactions. AIREBO potential describes as bellow [22]:

The molecular dynamics method is used on the platform of the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) package [23]. The molecular dynamics simulations have been done with LAMMPS. All geometries and boundary conditions of graphene sheets are in accordance with those prepared for molecular mechanics models.

The non-equilibrium MD simulations are performed with a velocity-Verlet integration algorithm in a micro-canonical ensemble system. The system temperature to a constant value is controlled using the Berendsen thermal scheme. The energy minimization for MD models is performed using the conjugate gradient algorithm to relax the system. This process is set in 2000 steps with a time step of 1 femtosecond. The set of graphene atom velocities is first produced at a given ambient temperature using the Gaussian distribution. The ambient temperature of the system is maintained at the temperature (300 K) via the Nose-Hoover thermostat. In the study of elastic properties and buckling behavior of graphene sheet, a one-directional load is applied to one end of the graphene sheet, while keeping the other end fixed. Various factors such as thermal conditions, time step, and boundary conditions significantly affect the results obtained from MD simulations.

Using molecular dynamics modeling in Lammps software [23] and using the AIREBO potential function [22], the strain energy in graphene is estimated. Graphene is under axial load. Using the method introduced in this research (MCS-MM), and using the following equation, the strain energy in the single-layer graphene structure is calculated:

To estimate the length scale parameter in the modified coupling stress theory, at first the energy deflection from the MCS-MM method and the MD method must be calculated using the following equation:

In this regard, using the results of molecular dynamics, strain energy is obtained. Also, using the method introduced in this paper, the strain energy is recalculated (Figure 3) As seen, the length scale parameter has a direct effect on the strain energy deviation between the two methods. Reducing this parameter reduces energy dissipation. At l = 0.15d, this amount of energy deviation reaches its minimum value and then increases.

The best length scale parameters for the beam elements are obtained using the presented method.

In this way, the results of the present method approach to the results of molecular dynamics by correctly selecting the length scale parameter and increase the accuracy of the calculations introduced in this article. The set value of the length scale parameter has been proposed as the most appropriate choice for predicting the elastic behavior of graphene.

Considering graphene plate as a discrete structure of covalently bonded carbon atoms, a molecular mechanical model including beams has been developed as an alternative to bonding between atoms through the finite element code ABAQUS. The bonds between the carbon atoms are defined as the elastic elements of the Bernoulli-Euler beam using the parameters obtained from the MCS-MM method.

The buckling for the graphene plate with different dimensions is shown in Figure 4. The results of the MCS-MM method are evaluated together with the results of the classical molecular mechanics method and molecular dynamics simulation method. Comparison between MCS-MM and Classical-MM methods shows that the method introduced in this paper predicts less critical load for graphene plate with width equal to 8 nm. The results of this method are very close to the results of molecular dynamics simulations. Molecular dynamics simulations provide accurate results due to the accurate consideration of interatomic interactions. The closeness of the results of the MCS-MM method to molecular dynamics shows the accuracy of the method introduced in this paper.

Buckling modes for pristine graphene are shown for the first five mode (Figure 5) In the following sections, we will use these buckling modes to compare the effect of cracking on the deformations created in each mode

In this section, the graphene sheets are considered with a crack in the center of it. According to Figure 6, on the left figure, a horizontal crack is located in the direction of the armchair. The graphene plate has 2L length, 2w wide, and 2a crack length. In the right figure, the crack is located in the vertical and zigzag directions. The buckling of the graphene sheet is checked by loading in the zigzag direction. The edges of the graphene plate are closed on armchair directions and free on the zigzag directions.

Figure 7 consists of 2 parts. The first part of Figure 7 shows the first five modes of buckling in the presence of the horizontal crack in the armchair direction. By comparing the buckling modes in the graphene plate without cracks and the graphene plate with horizontal cracks, it is observed that the presence of cracks has caused a noticeable change in the deformation created while buckling in the graphene plate. In the second part of figure 7, the buckling modes are shown on the graphene plate despite the vertical cracks in the zigzag direction.

Comparison of graphene buckling modes with horizontal and vertical cracks shows that in the first mode, in both cases, the crack edges are vertical, and in the first mode, they are separated from each other. Therefore, in the first buckling mode, the crack edges are spaced vertically from each other, and this behavior is equivalent to the first failure mode. This also happens in the third and fifth modes. Therefore, the buckling modes in a graphene plate with the boundary conditions determined in this research cause cracks to open in the first mode.

The figure 8 shows the critical open buckling changes for the graphene plate with horizontal and vertical cracks. As can be seen, when the load applied to the plane is perpendicular to the crack direction, the critical buckling load will decrease significantly with increasing crack length. However, when the applied load to the plate is in the same line as the crack in the plate, increasing the crack length will not cause much change in the critical buckling load. The buckling load for a crack length equal to 0 is equal to the value obtained for the crack-free plate, which obtained in the previous section and is validated with those obtained in molecular dynamics simulations.

The buckling in the graphene plate investigated with two cracks at a distance d from each other and crack length as 2a (Figure 9) In different cases, the distance between the cracks is increased and considered the distance "d" as 4.9 Ao and 63.9 Ao for horizontal cracks and 5.6 Ao to 56.7 Ao for vertical crack. The effects of two cracks in the graphene plate are investigated on the critical buckling load and the shape of the buckling modes. Also, the effect of distance between cracks on these two buckling parameters is studied.

Figure 10 shows the buckling modes in the case of two symmetrical cracks with a certain distance from the loading axis. The first part of Figure 10 is related to cracks perpendicular to the load axis, and the second part is related to cracks parallel to the load axis. As you can see, the presence of two symmetrical cracks in the graphene plate will cause deformation in different buckling modes compared to what is shown in Figure 5 for a crack. This difference is because the presence of two cracks will create a free boundary condition in the middle of the plate that is not affected by other parts of the graphene atomic structure. Due to the presence of cracks, these free parts cause the load of the free edges to be buckled separately, and this will somewhat disturb the classic shape of the buckling modes, as you can see in Figure 10. But because the two cracks escape at close distances from each other, this deformation of the buckling modes is closer to what considered in the previous cases.

As before, for the two-crack mode, the buckling load changes despite the distance d between the two cracks are shown in Figure 11. Compared to Figure 8, the presence of two symmetrical cracks in the graphene plate than in the case of a single crack can further reduce the critical force required for buckling. This change is greater if the cracks are perpendicular to the load axis. The presence of two cracks has been able to create more freedom in the atomic structure of graphene, and for this reason, it needs less force to reach the buckling state. In this case, the larger the crack length, the greater the freedom of the buckling structure and the greater the critical buckling load. But, where the cracks are parallel to the applied load, there is no more change in the critical buckling load.

The buckling behavior of the graphene plate by increasing the distance between two cracks are investigated. Where the cracks are located at distances 34 Ao horizontally and 31 Ao vertically, the deformation of graphene in different modes can be seen in Figure 12. As the distance between cracks perpendicular to the loading direction increases, it has been shown that the buckling load decreases further with increasing crack length. As the length of the crack increases, at greater distances from the two cracks, the graphene plate becomes a few pieces with free boundary conditions. The deformations created in different modes are significantly different from those shown for graphene without the crack in Figure 5.

In the last step, increase the distance between the cracks and bring them closer to the edges of the graphene plate (54 Ao horizontally and 48 Ao vertically) You can see the buckling modes in Figure 13. This figure depicts changes in the critical buckling load relative to the crack length. By increasing the distance between the cracks, the critical buckling load will decrease. This change is more noticeable than before. In this step, in addition to examining parameters such as buckling modes and critical buckling load, the eigenvalues in different buckling modes are compared for the three distances which is defined between cracks to see the effect of the change in distance between cracks on Eigen-modes. The buckling modes are placed on the graphene plate.

The buckling load changes for different crack distances versus the crack length are plotted in the figure 14. The results show that the smaller the distance between the cracks, the lower the buckling load. As the distance between the cracks increases and the cracks approach the outer edges, the effect of the crack distance on the buckling will decrease. However, increasing the crack length compared to the distance of the cracks from each other has a greater effect on the buckling parameters of the graphene plate.

Figure 15 shows the changes in the critical buckling load relative to the crack length. The crack aligned with the applied load has less effect on reducing the critical buckling load. In general, the buckling load decreases with increasing crack length.

1st Mode of Buckling

By increasing the number of cracks, the behavior of the graphene plate under buckling conditions is investigated. For this purpose, the number of cracks are increased. The graphene plate with 3, 5, 7, and 9 cracks are investigated, respectively (Figure 16)

Figure 17 shows the deformation of graphene with three cracks in the first five buckling modes. Buckling modes are shown for two types of cracks, horizontal and vertical.

2nd Mode of Buckling

Figure 18 shows the second mode of buckling in the presence of different numbers of cracks. In Figure 18-a, the number of cracks increases from 3 to 9 from left to right. The higher the number of cracks perpendicular to the load, the more flexible the graphene plate. Section b shows the graphene plate with cracks aligned with the load. Increasing the number of cracks also makes the graphene plate more flexible.

Investigation of figures 17 and 18 show that increasing the number of cracks in the graphene sheet causes a significant part of the deformation created in each buckling mode to be absorbed by the middle of the sheet. In this way, increasing the number of cracks causes the deformations to be far from the standard state of that buckling mode. As a reason for this, it can be considered that increasing the number of cracks increases the free boundaries in the graphene plate. In this way, graphene will have more flexibility under one load.

The critical loads for each of the above cases are obtained. Figure 19 shows the critical buckling load variation in terms of crack length. For cracks perpendicular to the load, the critical buckling decreases steadily with increasing crack length. Also, increasing the number of cracks in the graphene plate causes this reduction of the critical buckling load to occur more rapidly.

Figure 20 shows the critical buckling load changes to crack in the direction of the applied load. It is observed that at shorter lengths of crack, fewer changes occur in critical buckling load. When the crack length reaches about 60% of the plate length, increasing the crack length will have a much greater impact on the critical buckling load. But always more cracks will reduce the critical buckling load with more slope.

In the next step, the effect of the number of cracks in the graphene plate on the buckling in the 2nd to 5th modes is estimated. Figure 21 shows the changes in the Eigen values of the graphene plate for the second buckling mode despite cracks of different lengths. The same effect that cracking has on the critical buckling load on the graphene plate can be seen in specific values for other buckling modes. The presence of cracks always causes the free edges in the graphene to increase, and in each mode, the Eigenvalue decreases accordingly.

Figures 22 to 24 show changes of eigenvalues of the 3rd to 5th modes at different crack lengths. According to the presented results, same as the crack length, the number of cracks also affects the specific values. The higher the number of cracks, the greater the decrease in specific values in different modes. The graphs of the change in eigenvalues for different modes per change in crack length and the number of different cracks are compared. Studies show that the specific values in the first mode are always more affected by the number and crack length than other buckling modes.

In the following, the effect of the eccentricity of crack on the buckling behavior of graphene plates is investigated. As the first case, a vertical crack is placed at a certain distance from the loading axis (Figure 25) For each distance, the crack length is increased and calculate the critical buckling load is calculated.

Figure 26 represent the deformation of graphene for cracks of length 1 nm and three different eccentric distances. It is observed that with increasing distance from the center, the symmetric distribution of strain energy is disturbed.

Investigating the results of the critical buckling load shows that over a given crack length, increasing the eccentric distance reduces the effect of crack presence (Figure 27) The closer the crack is to the loading axis, the greater the critical buckling load. From the crack length a = 0.65w, the increase in crack distance with the loading axis is reversed. It means that the critical buckling load increases with increasing crack distance to the loading axis.

Next, the eccentricity of a horizontal crack in the graphene sheet is investigated. The crack geometry is shown on the graphene plate on one side of the loading axis in Figure 28. Figure 29 investigated two issues.

First, it shows the deformation due to buckling in different modes for cracks with different distances from the loading axis. On the other hand, it shows the changes in the strain energy distribution at different positions and different distances in the graphene plate. As shown, the eccentricity of the horizontal crack affects the crack deformation in different modes.

For example, in the fifth buckling mode, as the crack grows to the right edge, the deformation at this edge is more affected. Also in the first mode, which is known as the critical mode, this can be seen clearly. It happens because by increasing the eccentricity of the crack, the length of the crack on one side is increased. It creates a larger free edge. Therefore, the boundaries of the graphene plate will be larger and, as before, it will increase its flexibility. Then, less force is needed to create buckling.

The critical buckling load changes for this type of crack presence in the graphene plate and different distances from the crack center to the loading axis are shown in Figure 30. As illustrated, by increasing the distance of the center of the horizontal crack relative to the load axis, the critical buckling load is significantly reduced.

According to the results presented in this study, to increase the critical load of buckling and consequently increase the resistance of graphene plate to buckling, it must prevent the multiple cracks in the graphene plate. The crack growth should be prevented in graphene and prevent the growth of cracks, such as creating holes in the crack path or other methods. The results show that the asymmetry in the crack distribution in the graphene plate causes the critical load and buckling states to be more affected and the critical load value to be less.

This paper tried to investigate the effect of crack on buckling parameters using a modified molecular mechanics method based on modified coupling stress theory. The initial cracks created in the graphene plate were examined based on size, orientation, and number. The buckling of the monolayer graphene plate was studied in different modes and with different cracks. This method has better computational performance than other methods of nanostructure’s analysis. The MCS-MM method can be used in comparison with other methods with desirable accuracy.

At first, the crack in the center of the graphene sheet is placed and applied pressure from two directions. The results showed that increasing the crack length can reduce the critical buckling load. For a graphene sheet with two cracks perpendicular to the load axis, with large distances between the cracks, the critical buckling load remains almost constant. The distance between two cracks has less effect on the critical buckling load than the crack length. But the smaller the distance between the cracks, the less the critical buckling load.

This study shows that with a large number of cracks, the graphene sheet becomes more flexible and the critical buckling load is reduced. By increasing the number of cracks to 3, 5, 7, and 9, the critical buckling load will be less in the higher number of cracks. That is, increasing the number of cracks reduces the critical load. Its reason can be explained by the fact that as the number of cracks in the graphene plate increases, the number of free boundaries increases. The presence of open borders in the middle of the graphene plate causes uncontrolled deformations to occur at these borders and the graphene plate deform more and more easily. Therefore, it will require less force to achieve an unstable state in graphene. Then, the critical buckling load is reduced.

The parallel cracks to the applied load will have less impact on the buckling of the graphene plate. It can be seen in the changes in the critical buckling load and buckling modes. Investigation of buckling modes versus the number of different cracks illustrated that increasing the length or number of cracks has a direct effect on the Eigenvalues of different buckling modes and decreases them. But always the number and size of cracks had a greater effect on the first mode of buckling than others.

This study shows that to prevent buckling in graphene sheets, the number of cracks should be reduced, and crack propagation in graphene should be prevented. Some ways must be used that can both minimize the number of cracks and stop the growth of cracks. Also, the study of the presence of cracks outside the center shows that for coaxial cracks with loading direction, any amount of asymmetry in the crack distribution in graphene increases the critical buckling load and leads to a longer life of graphene. However, increasing the eccentricity for cracks perpendicular to the loading direction increases the critical buckling load.