Strength of graphene with curvilinear grain boundaries

Journal Pre-proof Strength of graphene with curvilinear grain boundaries Sankha Mukherjee, Robert Alicandri, Chandra Veer Singh PII:

S0008-6223(19)31190-X

DOI:

https://doi.org/10.1016/j.carbon.2019.11.058

Reference:

CARBON 14815

To appear in:

Carbon

Received Date: 5 September 2019 Revised Date:

4 November 2019

Accepted Date: 19 November 2019

Please cite this article as: S. Mukherjee, R. Alicandri, C.V. Singh, Strength of graphene with curvilinear grain boundaries, Carbon (2019), doi: https://doi.org/10.1016/j.carbon.2019.11.058. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.

Strength of Graphene with Curvilinear Grain Boundaries Sankha Mukherjeea,#, Robert Alicandria, # and Chandra Veer Singha,b* a

Department of Materials Science and Engineering, University of Toronto, 184

College Street, Suite 140, Toronto, ON M5S 3E4, Canada. b

Department of Mechanical and Industrial Engineering, University of Toronto, 5

King's College Road, Toronto M5S 3G8, Canada. *

corresponding author email: [email protected] (Chandra Veer Singh)

#

These authors contributed equally to this work.

Abstract Understanding the mechanical properties of grain boundaries (GB) in graphene is critical to commercial realization of large-area graphene in a vast array of applications. Previous pursuits on this topic mainly deal with tilt GBs, overlooking rotational GBs. Rotational GBs are composed of pentagon-heptagon pairs and are observed in mass-produced graphene. Herein, molecular dynamics simulations were conducted to determine the effect of GB curvature on the mechanical properties of graphene. While there was minimal effect on the Young’s Modulus (less than 5%), a Weibull distribution-based analysis revealed weakening of strength and strain to failure with decreasing grain boundary curvature (up to 23 GPa). Mechanical failure in rotational GBs was observed to occur by catastrophic failure of C-C bonds parallel to the tensile axis either away from or on the GB, differing by the initiating failure event depending on GB curvature. At higher GB curvatures, initial failure occurred in the bulk away from the GB. At lower GB curvatures, however, failure was initiated by a pseudo-plastic bond rotation, which propagated to catastrophic failure under increasing applied stress. Mitigation of crack propagation at lower curvatures was observed by the formation of monoatomic carbon chains, prolonging material survival under an applied stress.

2

1.

Introduction

Large-scale application of graphene-based structures will depend on the success of costeffective and fast synthesis routes[1]. Recent efforts for producing atomically thin graphene heavily rely on techniques such as chemical vapour deposition [2, 3] (CVD) on either Cu or Ni substrates, thermal graphitization of SiC [4-7], and liquid-phase exfoliation of graphite [8]. While these methods enable one to produce large quantities of graphene, in general, the quality of the final product is inferior to mechanically exfoliated graphene owing to the presence of a host of lattice defects such as point defects, dislocations and GBs [9, 10]. For example, there exists a general consensus on the fact that the presence of point defects adversely affect the mechanical properties [11, 12] and carrier mobilities of graphene [12, 13]. On the other hand, the role of GBs in determining the mechanical strength of graphene is not clear due to the existence of conflicting reports. Atomistic simulations of GBs with different misorientations in polycrystalline graphene have shown varying mechanical properties and failure mechanisms[14]. While there are reports which claim a detrimental effect of GBs on the breaking strength of graphene [2, 15, 16], several other reports also show insignificant deterioration or even improvement in strength caused by GBs [15, 17]. Interestingly, in a recent article, using nanoindentation experiments on CVD grown graphene, Xu et. al. demonstrated that graphene samples containing suitable combination of isolated grains, (i.e. in the absence of triple junctions) posses breaking strengths larger than those of pristine samples [18]. Furthermore, several theoretical calculations reported that symmetric GBs can lead to an increase in overall strength of graphene and additionally induce magnetism[15, 19, 20]. These conflicting results make GBs particularly interesting and therefore investigating defect arrangements along GBs in graphene and understanding their effects on GB cracking is essential for engineering applications of graphene. In this context, high-resolution tunneling electron microscopy of CVD or thermally produced graphene has shown not only the polycrystalline nature of the samples, but also the presence of rotational GBs, otherwise closed loops of adjacent pentagon-heptagon defect pairs. The ubiquity of these rotational GBs as noted in Ref. [6, 9, 21] necessitates further research to systematically understand their impact on the mechanical performance of graphene. A rotational grain

3

boundary, more commonly known as a ‘flower defect’ in graphene is a GB that ends by itself, host various vacancies or interstitial reconfigurations, without creating any unsaturated bonds due to a net Burger’s vector of zero. Furthermore, the overall effect of the varying orientation of the component dislocations throughout the circumference of the rotational grain boundaries in contrast to the constant, but repeating orientation of dislocations in tilt grain boundaries is not clear. GB loops posses interesting physical [22], thermal [23, 24] and electronic properties [25] due to their unique shape and structure. For example, Krasavin and Osipov demonstrated that in the long-wavelength limit the phonon mean free paths for flower defects vary as ω-3, whereas for open boundaries of any shape it varies as ω-1 [23]. There exists an agreement in the literature that the density and distribution of the pentagon-heptagon defects along the grain boundary have a significant effect on the mechanical strength of the graphene sheet through interaction (and possible cancellation) of adjacent stress fields—leading to the counter-intuitive result that higher defect densities can strengthen the material relative to less defective systems [15, 20, 26]. The relatively short range of defect stress fields has also been invoked for explaining a Hall-Petch behavior for polycrystalline graphene such that grain size is inversely proportional to the failure stress [27, 28]. Decreases in the tensile strain rate and increases in system temperature have also been shown to weaken the graphene sheet due to increased sampling and activation of failure events respectively [11, 29]. Previous studies studied the energetics associated with the growth of flower defects [30]; their creation via dislocation nucleation [9]; their kinetics of growth via dislocation motion [22] and subsequent annihilation via radiation-induced bond rotation [31]. However, although technologically important and scientifically interesting, the effects of this class of GBs on graphene’s mechanical properties or fracture behavior has not yet been fully explored. For example, the only available results report the basic mechanical properties such as Young’s Modulus and strength, and correlation between increasing atomic pre-stress along the boundary for low index GBs [32]. However, effects of flower defect size on the mechanical properties and failure mechanisms are still unexplored. Furthermore, polycrystalline graphene is known to exhibit significant statistical fluctuations in strength, which follow ‘weakest-link’ statistics especially for samples with high defect densities [16]. Whether graphene samples with rotational GBs follow any statistics is a looming question. If they do, whether statistical analysis of tensile

4

survival probabilities can be used to obtain estimates of mean strength is a further question still to be answered. It is the goal of this work to summarize how the evolution of these defect-systems under tensile loading is related to defect size as an effective representation of grain boundary curvature. Relative to other studies on the mechanical properties of defective graphene, this work offers a unique perspective in addressing the role of interfacial configuration in curvilinear and finite GBs more akin to those found experimentally than the straight and infinite GBs of previous studies.14-17 We attempt to answer the questions raised above by systematically conducting molecular dynamics simulations to isolate these effects by providing a parameterized study of the mechanical behavior and failure mechanisms of rotational GBs in graphene. A combination of rate dependent continuum modelling with Weibull statistics was employed to conclusively assess the statistical variations in strength for flower defects with varying curvature. Furthermore, atomistic origins of failure were investigated which showed curvature dependence at the nanoscale. These results offer important insights into the statistical strength and failure mechanisms of rotational GBs in graphene. 2. Simulation Method Classical

molecular

dynamics

simulations

were

performed

using

the

Large-scale

Atomic/Molecular Massively Parallel Simulator (LAMMPS) [33]. The Adaptive Intermolecular Reactive Empirical Bond Order (AIREBO), as implemented in LAMMPS, was used for all molecular dynamic’s simulations with the modification of the minimum bond cut-off radius of 1.92 Å while the maximum bond cut-off radius of 2.0 Å was kept unchanged [20, 34]. The thickness of the graphene sheet was taken as 3.35 Å as in previous studies [11, 26, 28, 29]. The time-step for all the simulations were chosen to be 1 fs, consistent with previous studies [11, 35]. Uniaxial tensile tests were performed in three steps: first, the energies of the graphene samples containing flower defects were minimized to allow for proper coordination of the grain boundary atoms. Secondly, the resulting structures were relaxed at 300 K for 20 ps in an isobaric-isothermal (NPT) ensemble. Finally, the grain boundary systems were deformed in the armchair direction until complete failure. Stresses on the defect systems while under deformation were computed as the time-averaged virial stress on each atom in the system. Stresses acting on

5

every atom was evaluated after every 500 timesteps and a spatial average was computed to obtain the estimate of total stress acting on the structure. Tensile simulations were conducted for a large range of strain rates, 5×107/s to 109/s, to carry out strain rate sensitivity analysis. The Young’s Modulus was calculated as the slope of the stress-strain curve in the linear regime, i.e., in our case within the first 5% strain. The ultimate tensile strength was taken as the maximum stress value obtained on a given curve. Tensile simulations were parameterized with respect to graphene sample sizes (12 nm to 80 nm) to allow for extrapolation of material behavior at realistic experimental conditions. Using the guidelines and notation presented in reference[30], a series of rotational GBs of the form C6(m,n)|m=n were constructed by rotating a hexagonal patch of graphene relative to the bulk lattice, originally oriented with the armchair direction in the x-axis. As in reference [30], the angle of rotation, 30°, was set by the geometrical parameters m and n using the formula

×

60° . The size of the rotated patch was radially dependent on the geometric parameters m and n, such that the number of radial hexagons in each hexagonal patch was determined through the parameter x=m+n-1. This parameter determines the notation used in the remainder of this work for describing the constructed rotational GBs in the form of C6(x). The curvatures of the GBs were varied by tuning the x-parameter. Furthermore, for any value of x more than 1, it was necessary to truncate the rotated crystal to minimize pre-strain along the grain boundary by allowing for the saturation of the boundary by pairs of 5|7 dislocations, for details refer to Figure S1. Otherwise, without truncating the rotated patch, the grain boundary would be composed of alternating sequences of grossly deformed, in other words energetically unfavourable, sequences of hexagons, heptagons, and pentagons. Note that for values of x of the form 4a+1 (where a is a positive integer) truncation of the rotated patch would leave the truncated atoms in their original position in the bulk, while for values of x of the form 4a+3, truncated atoms are rotated first then removed from the lattice. Molecular dynamic snapshots were generated through the freely-available OVITO software [36]. Examples of constructed rotational GBs with C atoms colored according to their potential energies is shown in Figure 1.

6

Figure 1- (a-d) Atomic configurations of various constructed rotational GBs, C6(x), where x=1,3,5,7, with x being the number of radial hexagons (rounding down) contained in a rotated hexagonal patch. Atoms are colour coded according to their potential energies with yellow mapped to lowest value. Depicted snapshot taken at the end of relaxation at 300 K.

To ensure sp2-coordination of atoms along the rotational grain boundary, each topology was first equilibrated with atomic motion out-of-plane disabled, followed by allowing the graphene sheet to buckle out-of-plane once each atom was properly coordinated. This two-step approach was required in that under-coordinated atoms would effectively repel themselves out-of-plane during minimization and never become coordinated otherwise. The latter buckling of the sheet was initiated by offsetting boundary atoms in a sinusoidal pattern with an amplitude of 0.1 Å. This offset is validated by the lowering of system energy upon buckling where the relaxed flat planar

7

geometry can be seen to represent a local energy minimum with the relaxed buckled geometry the global energy minimum. Within the limits of maintaining a ratio of simulation cell width (L) to defect diameter (D) of 10,15 and 20 for a given topology, the grain boundary system was placed in the approximate centre of the simulation cell. The simulation cell itself was made to be of as equal length as possible in both the armchair and zigzag directions as to minimize overlap of stress fields between periodic images of the system. The choice of diameter-to-simulation cell width was qualitatively a ratio that prevents periodic image interaction while also minimizing the simulation size for the sake of computational efficiency. 3. Results and Discussion

3.1 Quasi-intrinsic strength and survival analysis Mechanical failure of graphene at finite temperatures is a thermally activated process and depends on the initial velocity distribution ascribed to the constituent atoms [37]. To understand the stochasticity associated with the quasi-static (although at a very high strain rate) tensile MD simulations of flower defects of more than 1700 samples were studied. These samples were deformed uniaxially using different initial velocity distributions, strain rates, L/D ratios and grain boundary curvatures. At finite temperatures, disclination dipoles generated by the rotational GBs created local wrinkles in the graphene sheet, causing variations in the thickness direction (Figure S2), also seen in polycrystalline graphene[38]. Furthermore, the extent of out-of-plane displacement increased from C6(1) to C6(7). During tensile loading, initially these corrugations get stretched, reducing such fluctuations. In Figure 2, representative stress versus strain data in the armchair direction for samples C6(1), C6(3), C6(5), and C6(7) with a L/D ratio of 10. For comparison, the tensile response of pristine graphene in the same direction is also presented. Here, the maximum value of stress for each case is identified as the tensile strength of the sample. Different breaking points represented by triangles are indicated to show the stochasticity in strength. In all the structures, stress-strain responses showed signatures of nonlinear elasticity and failure occurred instantaneously with a sudden drop in the stress, showing signatures of brittle failure. All the samples showed the similar stress/strain curves as pristine graphene pulled in the armchair direction. This is evident from the construction of the grain boundary structure with armchair-oriented graphene being the bulk lattice, but this also demonstrates the minimal

8

impact of the rotational grain boundary on the Young’s Modulus of the graphene sheet at the tested grain boundary concentrations. Note that the modulus of any given grain boundary, calculated as the slope of the stress-strain curve at 5% strain, varies by no more than 5% of the Young’s Modulus for armchair-oriented graphene of 895 GPa as tabulated in Table 1 below. These values agree with the literature values presented in reference [12] for polycrystalline graphene.

Figure 2: Typical stress-strain plots for rotational GBs of varying curvature for a L/D ratio of 10. Tensile data from 5 computational runs performed at 300K and a strain rate of 109 s-1 are shown for each topology.

Rotational GBs, under all conditions tested, exhibited a weakening of the graphene sheet in terms of strain at failure as shown in Figure 3. This weakening was relatively unaffected by an applied strain rate; however, the strain required for failure is slightly reduced for smaller strain rates. This finding will define the scope of future work for energetics calculations, as all seemingly key failure modes have been identified. Furthermore, in terms of grain boundary curvature, the weakening of the rotational grain boundary trends increasingly with decreasing grain boundary curvature (i.e., larger radii boundaries). In terms of grain boundary curvature, the weakening of the rotational grain boundary trends increasingly with decreasing grain boundary curvature (i.e., larger radii boundaries).

9

Figure 3: Strain to failure GBs of varying curvature. Data plotted from tensile runs at 300K and a strain rate of (a)108 s-1 and (b) 109 s-1. Error bars show 95% confidence intervals. Error bars may be within marker. Table 1: Calculated Young’s Moduli at 300K.

Topology

Armchair-Oriented

C6(1)

C6(3)

C6(5)

C6(7)

E [GPa]

895 ± 9

912 ± 4

819 ± 16

884 ± 4

864 ± 8

The possibility of failure and fracture strengths of brittle materials containing defects of random distributions of size, and orientation is traditionally predicted using the Weibull distribution [39, 40]. Assuming GBs to be non-interacting defects, Shekhawat and Ritchie combined Weibull distribution with Kramer’s rate theory to obtain expression for survival probability for polycrystalline monolayer graphene at a fixed temperature as a function of strain rate ( ε& ), grain size (µ ) and sample size (L) [16], given by

S (σ | L, µ , ε& ) = e

−

L2ε&0  σ −σ 0    µ 2ε&  ν 

m

,

(1)

where σ is the stress, σ0 and v are rescaled location and scale parameters, ε&0 is a reference strain rate, and m is the Weibull modulus. Herein, larger values of m indicate greater variability in material behavior, and that the final failure of the material is a weak function of the defects present or that the density of defects present in the material is insignificant. The location parameter (σ0) is an indicator for the lower bound of strength.

10

S(σ) for all the rotational GBs were found to be a weak function of sample size (L), as an example S(σ) profiles for C6(3) and C6(5) samples are shown in Figure 4 as a function of sample size. Furthermore, as shown in Figure 5, the slope of S(σ) was also found to be sensitive to applied strain rates, similar to plasticity associated with dislocation nucleation [41]. For all the flower defects the likelihood of survival decreased with strain rate due to more time available for nucleation events which eventually lead to failure. In order to account for these effects, we rewrite the survival probability given by Shekhawat and Ritchie in the following form,

S (σ | µ , ε&0 , ε& ) = e

−

ε&0  σ −σ 0    ε&  ϕ 

m

,

(2)

2

 µ m Where ϕ = ν   . The mean breaking strength ( σ Mean ) and its variance ( σ V ) for Eq. (2) are L

given by, 1

σ Mean

 ε&  m  1  = σ 0 + ϕ   Γ 1 +  ,  ε&0   m  2

 ε&  m   2   1  σ V = ϕ   Γ  1 +  −  Γ  1 +    ε&0    m    m  

(3)

2

 , 

(4)

where Γ(⋅) is the gamma function. For a detailed derivation behind Eq. 3&4, refer to Section 3 in the Supplementary Material. Following the approach from reference [16], Figure 6 shows the survival probability of C6(3),C6(5), and C6(7) rotational GBs under strain rates of 5×107/s, 108/s, 5×108/s and 109/s, and the fits to Eq. 2 are indicated by black lines. In the same figure, the symbols represent stress at failure σ, where the damage tolerance of rotational GBs, being the difference between the onset and terminal values of stress at failure. Survival probability C6(1) is shown in Figure S3, for clarity. It can be clearly observed that the slope of S(σ) increases with decreasing grain boundary curvatures. The magnitudes of σ0 and m for C6(1), C6(3), C6(5), and C6(7) boundaries for 108/s and 109/s strain rates studied are presented in Table 2. The fitted values for the Weibull modulus (m) range increased from 2.5 to 18 with an increase in the radius of the rotational grain boundaries, similar in magnitude (i.e., to the value of 10) determined in

11

reference [23]. Furthermore, the magnitude of σ0 also decreased monotonously with decreasing curvature of the GBs. It is well known that the strength of 2D materials strongly depend on defect concentration [42, 43]. Interestingly, current results suggest that for identical levels of defect concentrations, different flower defects could have varied levels of detrimental effects on the strength of graphene. From the viewpoint of Weibull statistics which describes the survival probability of a disordered brittle material, the strength of material depends on the statistics of the weakest region in the sample. That is if a sample is considered to compose of a set of noninteracting volumes the distribution related to the strength of the samples will be determined by the statistical response of the weakest regions. For such systems, the probability that the defects with a stress dependent barrier height of ∆E(σ), endure a stress increase of dσ, is given by the product of the individual probabilities of all the defects. Therefore, the overall probability of survival decreases with increasing number of defects regardless of their volumetric concentration. From a materials standpoint, for flower defects with a larger disordered volume provide increased number of nucleation sites for critical flaw generation, eventually reducing the mean strength of the material. Therefore, although for a given L/D ratio, different samples containing different flower defects had comparable concentrations of 5-7 defects, the volume of defective atoms increased from C6(1) to C6(7). This behavior is observed in Figure 6 wherein C6(1) and C6(3) shows a degenerate distribution for S(σ) unlike C6(5) and C6(7) which show a positive right skewed distribution. Furthermore, while the σ Mean values for C6(1), and C6(3) structures were almost similar to those of pristine graphene (95.6 GPa at a strain rate of 109/s), those for C6(5), and C6(7) were decreased by up to 23 GPa. This is similar in magnitude to the weakening observed by Shekhawat and Ritchie16 for polycrystalline graphene, however in the present study, such weakening can be attributed to the presence of only a single flower defect rather than a grain boundary network. Determining if this trend holds for tilt angle GBs with an infinite radius of curvature is a key element for future work.

12

Figure 4: Red and blue symbols represent the Survival probabilities of C6(3) and C6(5) containing graphene samples under tensile loads at 300K for strain rates of 108/s. Squares represent a L/D ratio of 20, whereas for triangular and circular symbols the magnitude of L/D are 15 and 10, respectively.

Figure 5: Survival probabilities of C6(5) containing samples as a function of strain rate for a L/D ratio of 10.

13

Figure 6: Survival probability of C6(3), C6(5), and C6(7) containing graphene samples under tensile loads at 300 K for strain rates of, (a) 109/s, and (b) 5×108/s, (c) 108/s and (d) 5×107/s. Symbols represent data obtained from MD simulations and black lines are fits using Eq. 2. Table 2: Different parameters obtaining by fitting MD data to Eq. 2.

C6(1)

C6(3)

C6(5)

C6(7)

108/s

109/s

108/s

109/s

108/s

109/s

108/s

109/s

m

18.4

9.3

7.3

5.4

2.8

3.5

2.5

3.2

σ0 (GPa)

92.3

94.2

93.5

94.4

78.3

82.7

68

73.1

RMSE

0.95

0.99

0.99

0.94

0.98

0.99

0.99

0.99

Figure 7a (a complimentary log-log plot between − log( S | σ ) and (σ − σ 0 ) / ϕ )) shows that for the parameters obtained using Eq. 2, all the data obtained from MD simulations collapse in a similar fashion irrespective of strain rates, grain boundary curvatures, and sample sizes. Such a result

14

indicates that Weibull distribution, specially Eq. 2 can capture the survival probability of graphene samples with flower defects. The intercept or slope of these lines may vary depending 0

σ

on the values of

or ϕ . Furthermore, a lack of linearity could be observed towards the

extremities, which implies a departure of the stress distribution from Weibull statistics. However, this is a typical phenomenon since the stresses toward the tail have greater tendency to follow a family of extreme value distributions [44]. Furthermore, the slope of the data points for C6(1), C6(3) are smaller than those for C6(5), C6(7), highlighted in Figure 7a, which is indicative of the presence of different failure mechanisms for these classes of grain boundary structures. These mechanisms are discussed in detail in next section. Figure 3b presents the variation of σ Mean for C6(1), C6(3), C6(5), and C6(7) containing samples as a function of strain-rate obtained using Eq. 3, also in the same figure are superimposed average strength values obtained from our MD simulations which show reasonable agreement. Note that Eq. 3 does not contain free parameters. It can be observed, that while the magnitudes of σ Mean do not deteriorate significantly for C6(1) or C6(3) samples over several order of magnitudes of strain rate, however, weaken steeply for C6(5) and C6(7) samples. For example, regardless of sample sizes, compared to pristine graphene, the strength C6(5) and C6(7) samples can decrease by more than 13 GPa and 23 GPa for a strain rate of 109/s. The variance is mean strength, i.e., σ V (calculated using Eq. 4) as a function of applied strain rate is shown in Figure S4 which indicates that variance of σ Mean increases nonlinearly with strain rate. Furthermore, the rate of increase in σ V increased from C6(1) structures to C6(7), with σ V reaching up to 1.6 GPa for C6(7) GBs which is small compared to the strength of C6(7). These results reflect the effect of Weibull modulus (Eq. 4) which deteriorates in small curvature GBs, a decrease in the magnitude of m results in the increase of σ V at higher strain rates.

15

(

Figure 7: (a) Complimentary log-log plot between log ( − log( S | σ ) ) and log (σ − σ 0 ) / ϕ

) for C (1), C (3), C (5), and C (7) 6

6

6

6

containing graphene samples under tensile loads at 300 K. Open circles represent a strain rate of 5×107/s, solid squares 108/s, open triangles 5×108/s, and open squares 109/s. Color code: red, C6(1); black, C6(3), blue, C6(5); and green, C6(7). (b) Mean failure stress for different rotational GBs as a function of strain rate obtained using Eq. 3. Symbols represent average failure stresses obtained from MD simulations. Color code: red, C6(1); black, C6(3), blue, C6(5); and green, C6(7).

Similar to Han et al., we considered atomic virial stress along the loading direction (i.e., σXX) in different rotational defects (i.e., near 5-7 defects forming the boundaries) at the end of thermal relaxation at 300 K, prior to the application of any load. The per-atom stress values for different flower defects are shown in Figure 8. Examination of the stress fields for a given grain boundary prior to tensile deformation indicated that defect-induced stress decays rapidly to a far-field value. Furthermore, These pre-stress levels are consistent with previously observed stresses in graphene with tilt boundaries obtained from both MD simulations and continuum mechanics based disclination dipole models [20, 45]. It can be seen that without prior loading, stress concentrations are formed in the grain boundaries composed of 5-7 disclination dipoles. Atoms with such stress concentrations (i.e., high levels σXX) either break or rotate and thereby act as nucleation sites for the formation of critical flaws under tensile loading. In graphene samples with symmetric tilt boundaries, per-atom prestresses in 5-7 disclinations are known to alleviate each other owing to their systematic alignment [45]. In contrast, for asymmetric tilt boundaries such effects are less pronounced owing to the random arrangement of disclination dipoles. Similarly, as rotational GBs are also similar to asymmetric tilt GBs with continuously varying tilt angle, aforementioned stress cancellation effect is less pronounced. As a result, similar to asymmetric tilt boundaries, for rotational defects, the relationship between initial prestress and overall strength is less obvious. For example, in Figure 9 the maximum per-atomic pre-stress ( max ) is plotted for different grain boundaries in all the rotational grain boundaries. For all the σ XX strain-rates, the magnitude of σMean in general decreased with maximum pre-stress and grain boundary length, except for C6(3).

16

Figure 8: Stress field (

along the X-direction of (a) C6(1), (b) C6(3), (c) C6(5), and (d) C6(7) graphene samples with a L/D

ratio of 10. The atomic snapshots were created at the end of relaxation at 300K. Atoms are colour coded according to the stress acting on them with blue mapped to the lowest value.

max Figure 9: Dependence of mean tensile strength (σMean) on initial maximum tensile pre-stress ( σ XX ) as a function of strain rate in

C6(1), C6(3), C6(5), and C6(7).

17

3.2 Progressive Failure Mechanisms Two failure mechanisms for rotational GBs were observed from the atomic snapshots, (a) failure of C-C bond parallel to the tensile axis away from a rotational grain boundary and (b) failure of C-C bond parallel to the tensile axis along a rotational grain boundary. The first mechanism was observed in the high-curvature (small-radius) C6(1) and C6(3) structures while the second mechanism was seen in all larger-radius (lower-curvature) boundaries. This distinction in failure mechanism is reflected in the survival probability (Figure 7a) and σ Mean (Table 2) values for these defects. For either mechanism, catastrophic failure perpendicular to the tensile axis was observed upon bond cleavage, forming an elliptical crack with its major axis perpendicular to the applied loading. For high-curvature boundaries, failure occurred at tensile loadings comparable to that of pristine graphene. Furthermore, void/crack formation and simultaneous propagation did not relate to the flower defects. This behavior is analogous to flaw insensitive fracture of graphene below critical defect sizes [46]. In contrast, low-curvature boundaries only exhibited catastrophic failure after rotation and propagation of a radially oriented crack initiated from a failed bond rotation until a parallel C-C bond at the edge of the grain boundary is broken. Recently, Fox et al. performed extensive MD simulations to study the effect of orientation of the grain boundary relative to the tensile axis to highlight critical bonds that lead to the failure [47] in graphene at 300 K. Even though varying orientations were considered in their study, bond rotation still was not observed, an observation quite different from the current study. On the other hand, Terdalkar et al. had shown that in cracked single crystal graphene fracture paths evolve through alternating sequence of bond rotation of rupture at the cracktip, using nudged elastic band simulations [48]. The details of this bond rotation-assisted failure of low curvature flower defects are discussed in the following sections. It is important to note here that similar bond-rotation was previously observed to occur in tilt angle GBs [49] when tested at much higher temperatures, i.e., 1500 K. Furthermore, catastrophic crack propagation in low curvature boundaries was found to be prevented through the formation of monoatomic carbon chains, which limit separation of the fracture surfaces, was also observed in previous works [26, 28]. In this regard, material failure becomes a balance between formation of monoatomic carbon chains and crack propagation. Figure 10 below depicts this balance via a series of failure snapshots for a C6(7) rotational grain

18

boundary. This varying failure mechanisms observed can be understood as dependant on two factors.

Figure 10: Atomic snapshots of crack propagation in a C6(7) rotational grain boundary at 300 K for a strain rate of 108/s. Atoms are colour coded according to potential energies with yellow mapped to highest value. The tensile axis is parallel to the armchair direction of the bulk graphene lattice in all cases. Sub-figures (a) through (c) show propagation of a crack initiated from an incomplete bond rotation over approximately 52 ps, green arrows represent loading direction. The rotating bond is composed of the two-atom chain with coordination 2 in sub-figure (b) and partially completes its rotation in sub-figure (c). Note the formation of monoatomic carbon chains (with coordination 2) in sub-figures (b) through (c)—where the tensile effect of the chain keeps the adjacent lattice areas in a state of compression or minimal tension.

Availability of fracture initiation sites due to bond rotation. Lower-curvature defects exhibited a

pseudo-plastic behavior of bond rotation as a means of tolerating higher tensile stresses by reducing the applied moment on the bond by moving from an orientation perpendicular to the tensile axis to a parallel one. In order to obtain further insights, grain boundary formation E GB − nµc GB energies ( E GB , form ) for each rotational grain boundary was calculated, given by: E form = π Dt where E GB is the energy of the defective system, n is number of atoms in the defective system, µ c is the per-atom cohesive energy of a C atom in pristine graphene, and t is the thickness of graphene sheet [50]. The resulting grain boundary formation energies are plotted in Figure 11 with error bars showing standard deviation in a 95% confidence interval. A general trend can be noted in the increase in E GB form with decreasing grain boundary curvature and can be qualitatively correlated with the initial onset of bond rotation during tensile formation. In this manner, C6(1) rotational grain boundaries do not exhibit bond rotation prior to catastrophic fracture due to the unavailability of high-energy sites along the boundary to undergo rotation. In contrast, all lowercurvature boundaries exhibit bond rotation at progressively lower tensile strains due to the higher-energy (metastable) nature of said boundaries. This would imply that the residual stress fields (Figure 8) introduced by lower-curvature rotational boundaries lead to the activation of bond rotation as a failure mechanism compared to failure of parallel-oriented C-C bonds. We

19

also note that although the magnitude of E GB form for C6(3) is comparable to those of C6(5) and C6(7), bond rotation was not observed. This could be due to flaw indifferent fracture of graphene below critical size scales [46].

Figure 11: Grain boundary formation energy for different flower defects after a thermal relaxation at 300 K. Error bars represent standard deviation in energies obtained from 20 different simulations.

Available curvature for crack propagation. Noting that attempted bond rotation was observed

at one of four sites oriented 45°-off the tensile axis; stable crack propagation was limited until the cleavage of a C-C bond parallel to the tensile axis. For lower-curvature structures, the distance between observed bond rotation sites and a parallel C-C bond along the boundary increased as a function of grain boundary radius. Correspondingly, the available sites for formation of monoatomic carbon chains increased, enabling material survival at higher loadings. Further variability in material survival was observed due to the activation of all four bond rotation sites with the level of stress required for catastrophic failure from an incomplete bond rotation increasing with the number of attempted bond rotations. This relationship is displayed schematically in Figure S5. In the case of the high-curvature boundaries, the tight curvature limits the number of energetically favourable sites for monoatomic chain formation, leaving cleavage of a parallel C-C bond as the more probable failure mechanism, thereby exhibiting the higher failure strengths of an armchair-oriented graphene lattice. In combination, these factors explain the observed trends regarding the weakening of rotational GBs as a function of curvature as well as the ability of rotational GBs at lower curvatures to fail in a non-catastrophic manner as

20

a balance between available failure sites and the ability to propagate that failure. Further snapshots of either failure mechanism as well as a successful bond rotation are shown in Figure 12 and Figure 13 below respectively.

Figure 12: Observed failure mechanisms in rotational GBs. The rotational grain boundary is highlighted using a black dashed curve and instances of incomplete bond rotations a green circle. The tensile axis (green arrow) is parallel to the armchair direction of the bulk graphene lattice in all cases. (a-b) shows failure of a parallel-oriented C-C bond away from a C6(1) rotational grain boundary at 300K, failure initiation away from the flower defect is indicated using a yellow arrow. (c-d) shows failure induced from a micro-crack from an incomplete bond rotation in a C6(7) rotational grain boundary at 300K. Note that all four observed sites for bond rotation are active in this snapshot.

21

Figure 13: (a-c) Atomic snapshots of a successful bond rotation in a C6(5) rotational grain boundary, loading direction is indicated by a green arrow. The rotating bond is denoted by a thick black bar. This sequence was taken from a run performed at 300 K. Note the out-of-plane deformation of the rotating bond in the second sub-figure. Stress in the x-direction is reduced by a factor of 2 from approximately 30 GPa to 15 GPa due to the bond rotation.

4. Conclusions

In general, the effect of GBs on the strength of metals is understood in terms of the Hall-Petch effect. Contrastingly, in graphene, high-angle GBs possess strengths comparable to that of pristine samples, whereas low-angle boundaries are much weaker. Earlier reports (both experimental and theoretical) on the mechanical properties of GBs in graphene focused on tilt GBs, ignoring flower defects, which are observed in CVD graphene. As a result, our knowledge about how these defects influence strength and fracture mechanisms is lacking. It follows that the net effect of the varying orientation of the component dislocations throughout the circumference of the rotational grain boundaries in contrast to the constant, but repeating orientation of adjacent dislocations in tilt grain boundaries is effectively unknown.

Herein, molecular dynamics

simulations were conducted to understand the mechanical properties of rotational GBs, specifically the effect of their curvatures on the modulus, strength and deformation mechanisms. Four different flower defect structures were studied with varying sizes of the graphene sheets. Uniaxial tensile simulations were conducted to assess their mechanical properties which indicated that these defects were ineffective in the elastic regime, as the Young’s Modulus of the samples only differed by no more than 5% regardless of their curvature. However, weakening in strength and strain to failure was seen to increase with decreasing grain boundary curvature. A continuum elastic theory-based approach in combination with Weibull statistics was used to analyze MD data. Eq. 3&4 provides scaling relations for the prediction the mean strength of flower defects in graphene and beyond. Further analysis revealed that the presence of a single C6(5) and C6(7) could deteriorate the strength of graphene sheets by more than 15 GPa and 20 GPa (i.e., more than 20% deterioration) compared to pristine graphene. Two separate failure mechanisms have also been determined: (a) failure of C-C bond parallel to tensile axis away from rotational grain boundary, and (b) failure of C-C bonds parallel to tensile axis along rotational grain boundary. Stability of failure by either mechanism was observed as a balance between formation of monoatomic carbon chains and crack propagation. Formation of said monoatomic chains prolonged survival under an applied stress, and in some cases, allowed for toleration of stresses on the level survivable by a non-defective sheet. Additionally, material

22

survival was also prolonged by an observed pseudo-plastic response of bond rotation at ambient temperatures. It remains for future work to characterize the activation energy for bond rotation along a rotational grain boundary as a function of grain boundary curvature. These results will be invaluable to future grain-boundary engineering in graphene by not only characterizing the effect of a commonly seen class of defects on the material’s mechanical properties, but also by describing the effect of a key design variable (grain boundary curvature) in terms of available failure mechanisms. Author Contributions: All the authors contributed towards preparing the manuscript and have

given approval to the current version. Acknowledgement: The authors acknowledge funding by the Natural Sciences and Engineering

Research Council of Canada (NSERC) through the Discovery Grant and undergraduate scholarship programs, the Hart Professorship, and the University of Toronto. We also acknowledge Compute Canada for providing computing resources at the SciNet, CalculQuebec, and Westgrid consortia. The authors also thank Prof. Sandip Barui from the University of South Alabama for important discussions on Weibull distributions and extreme value statistics.

References [1] [2] [3] [4]

C. Zhu et al., "Supercapacitors based on three-dimensional hierarchical graphene aerogels with periodic macropores," Nano letters, vol. 16, no. 6, pp. 3448-3456, 2016. P. Y. Huang et al., "Grains and grain boundaries in single-layer graphene atomic patchwork quilts," Nature, vol. 469, no. 7330, p. 389, 2011. H. J. Park et al., "Growth and properties of chemically modified graphene," physica status solidi (b), vol. 247, no. 11-12, pp. 2915-2919, 2010. B. An, S. Fukuyama, and K. Yokogawa, "Graphitization of 6H–SiC (0001̄) Surface by Scanning Tunneling Microscopy," Japanese journal of applied physics, vol. 41, no. 7S, p. 4890, 2002.

23

[5] [6]

[7]

[8] [9] [10] [11] [12] [13] [14]

[15] [16] [17]

[18] [19] [20]

[21]

[22]

[23]

G. G. Jernigan et al., "Comparison of epitaxial graphene on Si-face and C-face 4H SiC formed by ultrahigh vacuum and RF furnace production," Nano letters, vol. 9, no. 7, pp. 2605-2609, 2009. N. Guisinger, G. M. Rutter, J. Crain, C. Heiliger, P. First, and J. A. Stroscio, "Atomic-scale investigation of graphene formation on 6 H-Si C (0001)," Journal of Vacuum Science & Technology A: Vacuum, Surfaces, and Films, vol. 26, no. 4, pp. 932-937, 2008. Y. Cui, H. Zhang, W. Chen, Z. Yang, and Q. Cai, "Structural evolution of flower defects and effects on the electronic structures of epitaxial graphene," The Journal of Physical Chemistry C, vol. 121, no. 28, pp. 15282-15287, 2017. Y. Hernandez et al., "High-yield production of graphene by liquid-phase exfoliation of graphite," Nature nanotechnology, vol. 3, no. 9, p. 563, 2008. S. Kurasch et al., "Atom-by-atom observation of grain boundary migration in graphene," Nano letters, vol. 12, no. 6, pp. 3168-3173, 2012. A. P. Kauling et al., "The Worldwide graphene flake production," Advanced Materials, vol. 30, no. 44, p. 1803784, 2018. M. Daly and C. Veer Singh, "A kinematic study of energy barriers for crack formation in graphene tilt boundaries," Journal of Applied Physics, vol. 115, no. 22, p. 223513, 2014. A. Zandiatashbar et al., "Effect of defects on the intrinsic strength and stiffness of graphene," Nature communications, vol. 5, p. 3186, 2014. J.-H. Chen, L. Li, W. G. Cullen, E. D. Williams, and M. S. Fuhrer, "Tunable Kondo effect in graphene with defects," Nature Physics, vol. 7, no. 7, p. 535, 2011. H. Zhang, Z. Duan, X. Zhang, C. Liu, J. Zhang, and J. Zhao, "Strength and fracture behavior of graphene grain boundaries: effects of temperature, inflection, and symmetry from molecular dynamics," Physical Chemistry Chemical Physics, 10.1039/C3CP44716B vol. 15, no. 28, pp. 11794-11799, 2013. R. Grantab, V. B. Shenoy, and R. S. Ruoff, "Anomalous strength characteristics of tilt grain boundaries in graphene," Science, vol. 330, no. 6006, pp. 946-948, 2010. A. Shekhawat and R. O. Ritchie, "Toughness and strength of nanocrystalline graphene," Nature Communications, Article vol. 7, p. 10546, 01/28/online 2016. H. I. Rasool, C. Ophus, W. S. Klug, A. Zettl, and J. K. Gimzewski, "Measurement of the intrinsic strength of crystalline and polycrystalline graphene," Nature communications, vol. 4, p. 2811, 2013. J. Xu, G. Yuan, Q. Zhu, J. Wang, S. Tang, and L. Gao, "Enhancing the strength of graphene by a denser grain boundary," ACS nano, vol. 12, no. 5, pp. 4529-4535, 2018. M. P. López-Sancho, F. de Juan, and M. A. Vozmediano, "Magnetic moments in the presence of topological defects in graphene," Physical Review B, vol. 79, no. 7, p. 075413, 2009. Y. Wei, J. Wu, H. Yin, X. Shi, R. Yang, and M. Dresselhaus, "The nature of strength enhancement and weakening by pentagon–heptagon defects in graphene," Nature materials, vol. 11, no. 9, p. 759, 2012. C. Gong, K. He, Q. Chen, A. W. Robertson, and J. H. Warner, "In situ high temperature atomic level studies of large closed grain boundary loops in graphene," ACS nano, vol. 10, no. 10, pp. 9165-9173, 2016. F. A. Lavergne, A. Curran, D. G. Aarts, and R. P. Dullens, "Dislocation-controlled formation and kinetics of grain boundary loops in two-dimensional crystals," Proceedings of the National Academy of Sciences, vol. 115, no. 27, pp. 6922-6927, 2018. S. Krasavin and V. Osipov, "Impact of grain boundary characteristics on thermal transport in polycrystalline graphene: Analytical results," Journal of Applied Physics, vol. 125, no. 8, p. 084301, 2019.

24

[24]

[25] [26]

[27] [28] [29] [30] [31]

[32] [33] [34]

[35]

[36]

[37] [38] [39]

[40] [41] [42] [43] [44]

N. Khosravian, M. Samani, G. Loh, G. Chen, D. Baillargeat, and B. Tay, "Effects of a grain boundary loop on the thermal conductivity of graphene: A molecular dynamics study," Computational Materials Science, vol. 79, pp. 132-135, 2013. L. Linhart, J. Burgdörfer, and F. Libisch, "Accurate modeling of defects in graphene transport calculations," Physical Review B, vol. 97, no. 3, p. 035430, 2018. T.-H. Liu, C.-W. Pao, and C.-C. Chang, "Effects of dislocation densities and distributions on graphene grain boundary failure strengths from atomistic simulations," Carbon, vol. 50, no. 10, pp. 3465-3472, 2012. J. Kotakoski and J. C. Meyer, "Mechanical properties of polycrystalline graphene based on a realistic atomistic model," Physical Review B, vol. 85, no. 19, p. 195447, 2012. Z. Song, V. I. Artyukhov, B. I. Yakobson, and Z. Xu, "Pseudo Hall–Petch strength reduction in polycrystalline graphene," Nano letters, vol. 13, no. 4, pp. 1829-1833, 2013. Y. Sun, F. Ma, D. Ma, K. Xu, and P. K. Chu, "Stress-induced annihilation of Stone–Wales defects in graphene nanoribbons," Journal of Physics D: Applied Physics, vol. 45, no. 30, p. 305303, 2012. E. Cockayne, G. M. Rutter, N. P. Guisinger, J. N. Crain, P. N. First, and J. A. Stroscio, "Grain boundary loops in graphene," Physical Review B, vol. 83, no. 19, p. 195425, 05/12/ 2011. O. Lehtinen, S. Kurasch, A. V. Krasheninnikov, and U. Kaiser, "Atomic scale study of the life cycle of a dislocation in graphene from birth to annihilation," Nature Communications, Article vol. 4, p. 2098, 06/28/online 2013. T. Zhang, Z. Yang, T. L. Ho, N. Zhu, and D. Zhu, "17 Investigation on Mechanical Behavior of Single-Layer Graphene with Grain Boundary Loops," 2016. S. Plimpton, "Fast parallel algorithms for short-range molecular dynamics," Journal of computational physics, vol. 117, no. 1, pp. 1-19, 1995. S. J. Stuart, A. B. Tutein, and J. A. Harrison, "A reactive potential for hydrocarbons with intermolecular interactions," The Journal of chemical physics, vol. 112, no. 14, pp. 6472-6486, 2000. W. Wang, S. Li, J. Min, C. Yi, Y. Zhan, and M. Li, "Nanoindentation experiments for single-layer rectangular graphene films: a molecular dynamics study," Nanoscale research letters, vol. 9, no. 1, p. 41, 2014. A. Stukowski, "Visualization and analysis of atomistic simulation data with OVITO–the Open Visualization Tool," Modelling and Simulation in Materials Science and Engineering, vol. 18, no. 1, p. 015012, 2009/12/15 2009. D. Akinwande et al., "A review on mechanics and mechanical properties of 2D materials— Graphene and beyond," Extreme Mechanics Letters, vol. 13, pp. 42-77, 2017. O. V. Yazyev and S. G. Louie, "Topological defects in graphene: Dislocations and grain boundaries," Physical Review B, vol. 81, no. 19, p. 195420, 2010. C.-Y. Hui, S. Phoenix, and D. Shia, "The single-filament-composite test: a new statistical theory for estimating the interfacial shear strength and Weibull parameters for fiber strength," Composites science and technology, vol. 57, no. 12, pp. 1707-1725, 1998. W. Weibull, "A statistical distribution function of wide applicability," Journal of applied mechanics, vol. 18, no. 3, pp. 293-297, 1951. T. Zhu, J. Li, A. Samanta, A. Leach, and K. Gall, "Temperature and strain-rate dependence of surface dislocation nucleation," Physical Review Letters, vol. 100, no. 2, p. 025502, 2008. R. E. Roman and S. W. Cranford, "Defect sensitivity and Weibull strength analysis of monolayer silicene," Mechanics of Materials, vol. 133, pp. 13-25, 2019. G. Lopez-Polin, J. Gomez-Herrero, and C. Gómez-Navarro, "Confining crack propagation in defective graphene," Nano letters, vol. 15, no. 3, pp. 2050-2054, 2015. D. R. Cox, Analysis of survival data. Chapman and Hall/CRC, 2018.

25

[45] [46] [47] [48] [49] [50]

J. Han, S. Ryu, D. Sohn, and S. Im, "Mechanical strength characteristics of asymmetric tilt grain boundaries in graphene," Carbon, vol. 68, pp. 250-257, 2014. T. Zhang, X. Li, S. Kadkhodaei, and H. Gao, "Flaw insensitive fracture in nanocrystalline graphene," Nano letters, vol. 12, no. 9, pp. 4605-4610, 2012. A. Fox, U. Ray, and T. Li, "Strength of graphene grain boundaries under arbitrary in-plane tension," Carbon, vol. 142, pp. 388-400, 2019. S. S. Terdalkar, S. Huang, H. Yuan, J. J. Rencis, T. Zhu, and S. Zhang, "Nanoscale fracture in graphene," Chemical Physics Letters, vol. 494, no. 4-6, pp. 218-222, 2010. J. Zhang, J. Zhao, and J. Lu, "Intrinsic strength and failure behaviors of graphene grain boundaries," ACS nano, vol. 6, no. 3, pp. 2704-2711, 2012. H. Liu, Y. Lin, and S. Luo, "Grain boundary energy and grain size dependences of thermal conductivity of polycrystalline graphene," The Journal of Physical Chemistry C, vol. 118, no. 42, pp. 24797-24802, 2014.

26

Declaration of interests ☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ☒The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.