Graphdiyne family-tunable solution to shock resistance

The excellent mechanical properties of Graphdiyne (GDY) family has enabled it as an appealing candidate in the field of impact protection. In this in silico study, Monolayer GDY nanosheets of different morphology including GDY, GY-3, GY-4, GY-5 and GY-6 are assessed under hypervelocity impacts (from 1 to 6 km s−1). Tracking the deformation mechanisms under impacts as well as the Probability density function based on atomic Von Mises stress distribution, the length of acetylenic chain clearly alters ductile behavior as well as the energy dissipation/delocalization rate of GDY family during the impact. Results also suggest the penetration energy is not only determined by the energy delocalization rate but also sensitive to impact velocity for nanosheet with various acetylenic chain length. GY-5 with a much lower energy delocalization rate presents a close penetration energy comparing with GDY at a low impact at ∼2.0 km s−1, its superior ductility granted by long acetylenic chain not only dissipates kinetic energy of projectile via deformation, but also extends time for acceleration during the contact with projectile. Considering the impact resist performance of GDY family in terms of Specific penetration energy, GY-5 with the perfect balance between material density, ductility and Young’s modulus makes it the superior anti-ballistic material for impact velocity at <5 km s−1. For impact velocity >5 km s−1, it induces severer local deformation, and leaves no time for a well-developed distributed pattern as observed in a lower impact velocity scenario. As such, extensive elastic deformation of the nanosheet is not captured under impact, nanosheets with shorter acetylenic chains and hence greater material strength demonstrates superior impact resist. This study provides a fundamental understanding of the deformation and penetration mechanisms of monolayer GDY nanosheets under impact, which is crucial in order to facilitate their emerging applications for impact protection.


Introduction
Two-dimensional (2D) carbon nanosheet are regarded as a group of the most prominent materials for their excellent electronic, optical properties, ultrahigh mechanical properties and thermal conductivity [1]. Graphene is a representative of 2D carbon allotrope with honeycomb structure [2]. Substituting the carbon-carbon sp 2 bond in graphene by sp-sp 2 -hybridized bonds (i.e. -C≡C-), graphdiyne (GDY) with highly conjugated π system with non-zero band gap is obtained [3]. Consider the length of acetylenic chain length in GDY family, they can be classified as GDY, GY-3, GY-4, GY-5 GY-6 etc. Since its proposal in 1997 [4] extensive researches have been conducted for the realization of GDY family. The state-of-art technology provides two major approaches for the synthesis of GDY family, namely, dry chemical routes [5][6][7] and wet chemical routes [8][9][10][11] where the wet chemical routes are more cost effective and scalable for large-area GDY which enables a wide range of applications including electrocatalyst for ammonia production [12], electrode of high-performance batteries [13], photoelectrochemical water splitting cell [14-1614-16], ammonia molecules detector/sensor [17], high electro-mechanical transduction efficiency actuators [18]. Compared with graphene, GDYs are considered superior in electrical conductivity [6,19], which makes them a better building block of electric device such as highly efficient electron transport perovskite solar cell [20]. Unlike graphene, GDY family has non-zero band gap, which allows the applications as semiconductor, and it has been realized as high-performance diodes [21].
Along with all the applications, mechanical properties of GDY family are the fundamental requirement for their implementation, which have been extensively studied both numerically and experimentally under quasistatic loading condition. In terms of in-plane stiffness, results from Density Function Theory (DFT) and molecular dynamics (MD) approach agree well with each other, which is about half of graphene (∼166 N m −1 ) [22][23][24] and the recorded yield strength of GDY in zigzag and armchair direction are 45 GPa and 29 GPa respectively. In silico studies also suggest the length of acetylenic chain in GDY family plays a crucial role in determination of their mechanical properties, Young's modulus of GDY monolayers are inversely correlated with the length of acetylenic groups, whereas fracture strain extended with the acetylenic linkage [25,26]. Recently, Atomic force microscope indentation presents the Young's modulus of GYD synthesized via wet chemical route to be ∼218.5 GPa which is about half of the MD simulation (∼513.8 GPa) , this discrepancy is suggested by the vacancy defects in the nanosheet [27].
GDY family show excellent mechanical properties under quasi-static loading condition comparing with conventional engineering materials [28][29][30][31]. Therefore, with its low density and structure tunability, it is of great interest to investigate the application of GDY nanosheet networks in the field of impact protection such as combat armor and orbital debris protective shield for spacecraft [32,33]. Considering GDYs' counterpartgraphene, miniaturized ballistic experiment demonstrates the specific penetration energy of multilayer graphene is 10 times higher than that of conventional metal [34]. The followed up in silico studies suggest the impact resistance of 2D carbon members are strongly correlated with the density of mono-vacancies [35] and the axial-wave/cone-wave propagation pattern in graphene sheet under hyper velocity impact allows the graphene sheets to transfer more momentum per unit area and hence provide better ballistic protection [36]. In addition, our previous work suggests that the impact resist mechanism works differently for carbon nanosheet under various impact velocities, a higher impact velocity induces a severer local deformation, and there is no time for a well-developed stress distribution pattern as observed in a lower impact velocity scenario [37,38]. The versatile mechanical and anti-impact properties of GDYs for different morphologies are still necessary to be invested. In this paper, MD method is utilized to investigate the mechanical behaviour of GDYs with different acetylenic chain lengths under various impact loadings for the first time. The impact resistance performance of GDYs is compared with other carbonadoes nanosheets, via accessing (specific) penetration energy, stress propagation, stress distribution performance as well as the overall energy transfer during the impact process.

MD setups
In this work, the fracture behavior and mechanical performance of GDY family with different acetylenic chain lengths are mainly investigated, including GDY, GY-3, GY-4, GY-5 and GY-6. GDY family under hypervelocity impact are assessed through a series of large-scale MD simulations, performed based on open-source package LAMMPS [39]. Spherical diamond projectile (radius ∼2.5 Å) with an initial velocity up to 6 km s −1 is constructed to impact the center of Monolayer square nanosheets (side length ∼500 Å) with fixed boundaries. The armchair edge and the zigzag edge of nanosheet samples are in x and y direction respectively (see figure 1).
Consider the effect of thermal fluctuations during the impact process, environment temperature of 10 K is adopted. Adaptive intermolecular reactive empirical bond order (AIREBO) potential [40,41] with C-C cut-off distance of 2.0 Å [42][43][44] which well-represents the elastic properties of 2D graphene allotropes is employed to describe the C-C bonding energy as well as the bond breaking phenomenon of GDY family under impact loading. To differentiate the binding energy of C-C within diamond projectile and nanosheets, Tersoff potential is adopted for diamond [45] and the interactions between the projectile and nanosheets are defined by a Morse potential [46]. In the energy minimization stage of the simulation, nanosheets are relaxed using the conjugate gradient algorithm. Then, the Nose-Hoover thermostat [47] is employed to equilibrate the whole system at 10 K (NVT ensemble) for an extended duration of 2000 fs. In order to capture detailed deformation/bond break phenomena during the hyper-velocity impact, a small time step of 0.05 fs is selected and the Velocity Verlet algorithm is applied to solve the integrated the equations of motion [48]. In this work, in order to reproduce the energy conversion from kinetic energy to potential energy, thermostat is not applied during the whole process and non-periodic boundary condition is adopted.

Atomic stress
During the impact process, the atomic stress is estimated by virial stress P ab , which is presented as [49] p ab i is the atomic stress correlated with atom i. v i and W stands for the effective volume of the ith atom and the system volume respectively. m i and v i represent the mass and velocity of the ith atom. F ij and r ij serve as the force and length between atoms i and j, respectively, where the cartesian components of the corresponding variable are denoted by indices a and b.
Nanosheets in this work are assumed as continuum media with a thickness equals to the interlayer distance of graphite which is ∼3.35 Å. On the basis of the atomic virial stress (equation (1)), Von Mises stress s VM in the nanosheet (equation (2)) is tracked to visualize the complex stress distribution and transfer during the impact. Notably, atomic volume might alters the magnitude of the stress, however, it won't change the trends of the results which are stressed in this study 3. Results and discussion

Penetration energy under different impact velocity
Initially, impact performance of GDY family nanosheets under the impact velocities from 1 to 6 km s −1 is focused. The total energy variation of the GDY nanosheet (ΔE tot,GDY ) and the diamond projectile (ΔE tot,GYP ) are presented in figure 2. Here the total energy E tot of the system refers to the sum of kinetic energy and potential energy of all components in the system. Theoretically, the total energy change in projectile (ΔE tot,GDYP ) equals to the change in GDY nanosheet (ΔE tot,GDY ) for a closed system, which is confirmed by figure 2, as ΔE tot,GDYP overlaps with ΔE tot,GDY before perforation. However, after perforation, ΔE tot,GDYP maintains at a constant level as projectile travels in vacuum space, whereas ΔE tot,GDY increases slightly before it stabilizes at certain value. Such phenomenon is supposed as resulted from the fixed boundary conditions [37]. In this regard, the stabilized total energy change in projectile is considered as the penetration energy (E p ), which is ∼983.66 eV for GDY. Notably, graphene nanosheet shows a penetration energy (ΔE tot,GP ) ∼2118.68 eV [37], which is larger compared to GDY nanosheet. Considering the different atomic configurations, the penetration energy per atom is about 0.014 eV and 0.017 eV for tested GDY and graphene sample, respectively, which are close to each other. In this study, though there is significant amount of energy transfer between the projectile and GDY nanosheet within a short contact period, potential energy of the projectile shows marginal change, signifying ignorable deformation during impact. After perforation, both kinetic energy and potential energy of the project remain unchanged due to the vacuum condition ( figure 2(b)).
To investigate the impact resistance of GDY family, penetration energy of different samples under various impact velocity amplitudes are examined ( figure 3). Generally, greater impact velocity possesses higher penetration energy, besides, the penetration energy is inversely correlated with the length of acetylenic linkages i.e. GDY nanosheet with shorter acetylenic linkages trends to present a higher penetration energy and a higher threshold velocity for perforation. Individually, penetration energy of GDY, GY-3 and GY-4 present an overall parabolic shape, while GY-5 and GY-6 demonstrate more complicated fluctuation patterns which all present a local maximum at impact velocity ∼1.6 km s −1 and a local minimum at ∼3 km s −1 , indicating three identical penetration energy in the whole impact velocity domain e.g., GY-5 presents similar penetration energy of ∼860 eV at impact velocity of 1.4, 2.0 and 5.0 km s −1 . Interestingly, despite significant acetylenic linkage length difference between GDY and GY-5, they present close penetration energy at relatively low impact velocity ∼2.0 km s −1 .

Morphology influences
It is of great interest to explain how acetylenic chain length alters impact performance of GDY family under various impact velocity and explore the mechanism behind the penetration energy fluctuation phenomena, thus, GDY nanosheet deformation at different velocity amplitude are examined first. Figure 4 illustrates the atomic configuration of GDY at different stages of deformation process. The highest planer stress locates in the centre of conical impact area and decays from the centre to its boundary (figures 3(a), (c)). Different from the pure zigzag kicking fracture mechanisms in graphene [35,36], the cracks in GDYs propagate along both armchair and zigzag directions and their finally failure shape varies under different projectile velocity (figures 4(b), (d)). Specifically, the number of cracks increase with impact velocity, the projectile velocity of 2 km s −1 contributes to an irregular quadrilateral shape, while hexagon shape deformation is formed by a higher velocity of 6 km s −1 (figures 4(b), (d)). It is also noted that more break of bonds is observed as significantly more dangling parts are knocked out of the nanosheet under impact velocity of 6 km s −1 compare with impact velocity of 2 km s −1 (inset of figure 4(c)).  In addition, GDY experiences large deformation (with impact area radius around 53.22 Å and out-of-plane deformation about 27.53 Å, figure 4(a)) before the crack initiation and more kinetic energy of the projectile is dissipated in elastic deformation stage under an impact of 2 km s −1 compare with higher velocity impact. As impact velocity increases, deformation size prior to bond declines, which also brings down the penetration energy. Increase impact velocity above 3 km s −1 potential energy becomes less dominant as deformation prior to bond break is becoming insignificant. For impact velocity of 6 km s −1 , it leads to deformation area similar to the strike surface of the projectile (with the radius of ∼19.73 Å) and smaller out-of-planer deformation (∼11.89 Å) onset of crack propagation. In this case, penetration energy is mainly contributed by kinetic energy gain of the impact region. Reasons above explains the parabolic shape in GDY penetration energy, which agrees well with previous researches [37,38].
The elastic deformation varies size under various impact velocities, which is related with material stress/ elastic wave propagation velocity. For membrane with a higher elastic wave propagation velocity, it dissipates energy at a faster rate. Therefore, a better local stress distribution is expected during the impact, which delays the arrival of yield stress locally and helps the structure to stretch more prior to bond break. Elastic wave propagation velocity /r = v E s is a material constant. Hence, increase the impact velocity, less significant deformation is observed.
The elastic wave propagation velocity varies with different GDY configurations, so Probability density function (PDF) of the Von Mises atomic stress distribution is utilized to compare the stress wave propagation/ distribution in different GDYs under the same impact velocity of 2 km s −1 at the same simulation time of 2.3 ps (figure 5). It is clear that GDY with the lowest probability density peak possess the most evenly distributed stress pattern among all the nanosheets, followed by GY-3, GY-4, GY-5, GY-5, which agrees with the order of penetration energy. This phenomenon indicates the length of acetylenic chain is inversely correlated with the elastic wave propagation velocity.
GDY possesses superior elastic wave propagation velocity as well as significantly higher Young's modulus comparing with GY-5 [25]. Due to the two factors penetration energy gap between GDY and GY-5 enlarges as the impact velocity rises. However, the penetration energy of the two materials draws very close at ∼2.0 km s −1 impact. To unveil the mechanism behind the observation as well as the penetration energy fluctuation in GY-5, the atomic configurations of GY-5 under various impact velocities are studied. Owning 4 more repeating acetylene groups compare with GDY, the critical penetration velocity of GY-5 reduced to 1.4 km s −1 , poses a significant in-plane and out-of-plane deformation which are 101.51 in radius and 54.55 Å respectively before bond break ( figure 6(a)). Increase the impact velocity to 2.0 km s −1 , in-plane deformation radius and out-of-plane deformation of GY-5 reduced to 60.23 Å and 36.53 Å respectively (figure 6(c)). Compare figures 4(a) and 6(c), GY-5 deforms more in both directions prior to bond break, indicating a higher material ductility property under impact velocity 2.0 km s −1 , despite the difference in stress propagation velocity. It's superior ductility granted by long acetylenic chain [25,27] not only dissipates kinetic energy of projectile via deformation, but also offers extended time for the contact the with projectile which allows more directly transferred kinetic energy and bring explanation to the close penetration energy. Further increase the impact induces more severe local deformation, and there is no time for the stress to develop a well distributed pattern and extensive elastic deformation as observed in a lower impact velocity scenario. Thus, for an impact velocity of 5 km s −1 , the impact region melted immediately into more debris as the projectile approaches the nanosheet, comparing with the lower velocity impact cases (figures 6(b), (d), (f)).

Energy break-down
Notably, there exists three points which share identical penetration energy for GY-5 and GY-6 under different impact velocities, and the above analysis presents that GY-5 experiences close penetration energy under different projectile velocities of 1.4 km s −1 , 2.0 km s −1 and 5 km s −1 . To unveil the mechanism behind such observation, an Energy break-down analysis for GY-5 is carried out. Figure 7(a) illustrates a representative energy change process of GY-5 under an impact velocity of 2.0 km s −1 . The red line represents the rectified magnitude of GY-5 potential energy with respect to time, which is donated as Pe absr . Similarly, the blue line stands for the rectified kinetic energy Ke absr at different time instant. The rectification in both cases make sure all lines share same minimum value of 0. Generally, three stages can be identified as discussed below. Stage I records the energy transfer between projectile and nanosheet, it starts with the initiation of the simulation. Initial Pe absr and Ke absr energy of GY-5 maintains a constant level before the collision. Projectile reaches the nanosheet at ∼1.1 ps, Pe absr declines continuously due to deformation in nanosheet, meanwhile, projectile transfers its kinetic energy to GY-5 nanosheet via direct impact. Stage II describes the energy transfer within the nanosheet. Initiation of stage II is recognized by the peaked Ke absr at ∼3.4 ps, i.e. at this instant, projectile perforates the nanosheet and energy transfer between the two objects is minimized. After the local climax of Ke absr at ∼3.4, the variation in kinetic energy and potential energy are wellmatched. Stage III is signified by the stabilization of both lines in the diagram at ∼5.8 ps.
Considering the energy transfer in the system is a dynamic process, magnitude of potential energy and kinetic energy of nanosheet fluctuates with respect to time. In order to analyses the ballistic resist mechanism of GY5 under impact velocities 1.4, 2.0 and 5.0 km s −1 , the perforation instant is select for comparison. As illustrated in figure 7(b), though three impact velocities share a similar penetration energy, coherent decline in potential energy change (ΔPe) are presented despite a higher impact velocity generates more bond break and debris as observed in figures 4 and 6 (also see supplementary information is available online at stacks.iop.org/ MRX/7/115602/mmedia). As the impact velocity increases, ΔPe in GDY declines mainly due to the shrink in deformation prior to bond break. Meanwhile GY-5 nanosheet kinetic energy change (ΔKe) gains with a decreasing rate considering the ΔKe bars in figure 7(b). This phenomenon indicates the kinetic energy transfer efficiency (via direct impact) declines with respect to the impact velocity, which is related with the energy dissipation/delocalization capability of the material and the reduced contact time with nanosheet. These results indicate the penetration energy of GDY family is not only determined by the acetylenic chain length, but also sensitive to the impact velocity, which bring insights for the design of ballistic protection equipment targeting specific application scenario.

Impact resistance
Acetylene group not only modifies the overall strength and ductility of sample nanosheets but also varies their density. To ensure a more quantitatively comparison, specific penetration energy of GDY family under various velocity amplitudes are assessed to define the impact resistance. Theoretically, the penetration energy can be estimated from when the ratio between the projectile diameter D is much larger than the thickness h, i.e., D/h?1 [34]. Here, the first term refers to the minimum inelastic energy transferred to the sample and the second term represents the contribution from other energy dissipation mechanisms (e.g., bond break). A s represents the strike face area ( p = A r , s 2 r is the radius of the projectile); and v stands for the impact velocity. However, in some cases of our study (i.e. GY-5), the out of plane deformation is even more significant than diameter of the projectile ( figure 6(a), (b)). Consider bond breaks and debris created during the impact both elastic and inelastic deformation contribute significant amount to the total penetration energy. To ensure a fair comparison for all studied samples, A s is replaced by A GDY which equals the nanosheet sample area. Consider the different morphology and density, the specific or gravimetric penetration energy is adopted, which is calculated as and can be further expressed as *Apparently, E d * is a figure of merit to evaluate the impact energy delocalization ability of a material as more sample mass beyond the projectile contributes to the energy dissipation. As illustrated in figure 8, taking material density into consideration, the overall ranking of each material varies comparing with figure 3. For impact velocity >3 km s −1 , E P * of GDY family generally share a similar shift as that of graphene from the material-independent energy dissipation baseline (i.e., / v 2 2 ), especially. Such result is explained by the insignificant nanosheet deformation prior to bond break as impact velocity increases, in such instance, E P * is dominated by the kinetic energy directly transferred to atoms in strike face area. For an impact velocity smaller than 3 km s −1 , the penetration energy is found to deviate from the parabolic relationship with the impact velocity, which is resulted from the combination of nanosheet deformation and kinetic energy gain from direct impact. Particularly, GY-5 is observed to possess a much higher E P * for an impact velocity <3.0 km s −1 . Recall the penetration energy in figure 3, GY-5 shows comparable penetration energy with GDY at impact velocity ∼2.0 km s −1 , however, the density of GY-5 is ∼67.8% of GDY, thus GY-5 clearly demonstrates better performance in terms of E . P * Though GY-6 possesses a longer acetylenic chain, longer acetylenic chains mean a lower density as well as a lower Young's modulus which drags down energy delocalization efficiency, thus impact resistance in terms of E P * declines compare with GY-5. However, GY-6 still possesses an outstanding impact resistance comparing with other counter parts for impact velocity <3 km s −1 . Increasing impact velocity above 3 km s −1 , the performance of GY-6 declines quickly, mainly due to its poor energy delocalization capability among the tested samples. The perfect balance between material density, ductility and Young's modulus makes GY-5 the best anti-ballistic material for impact velocity <6 km s −1 . Notably, impact velocity ∼1.6 km s −1 seems critical for all tested carbonaceous materials, as all materials clearly present a local maximum around the velocity, despite it is not the minimum perforation velocity for all cases.

Conclusion
In summary, the fracture behaviour of nanosheet from GDY family namely GDY, GY-3, GY-4 GY-5 and GY-6 under various hypervelocity impacts are explored. The length of acetylenic chain alters elastic deformation size of the nanosheets from GDY family i.e. a longer acetylenic chain possesses a larger deformed area in the impact region prior to bond break, hence a more ductile behaviour is observed. Tracking the Probability density function based on atomic Von Mises stress distribution, nanosheet of GDY family with shorter acetylenic chain possesses higher Young's modulus to density ratio thus a better stress distribution pattern is captured due to superior energy dissipation/delocalization rate, which presents a positive correlation with the overall penetration energy ranking. Penetration energy is not only determined by the energy delocalization rate of the material but also sensitive to the impact velocity magnitude. Despite the difference in energy delocalization capability, GY-5 presents a close penetration energy comparing with GDY at low impact ∼2.0 km s −1 . G −1 Y-5's superior ductility granted by long acetylenic chain not only dissipates kinetic energy of projectile via deformation, but also offers extended time for the contact the with projectile which allows more time for the acceleration of nanosheet. Considering the impact resist performance of GDY family in terms of Specific penetration energy, GY-5 with the perfect balance between material density, ductility and Young's modulus makes it the best anti-ballistic material for impact velocity <6 km s −1 . For a higher impact velocity ∼6 km s −1 , it induces severer local deformation, and there is no time for a well-developed distributed pattern as observed in a lower impact velocity scenario. As such, extensive elastic deformation of the nanosheet does not occur under a high impact velocity and nanosheet with shorter acetylenic chain comes at top in terms of impact resist. This study provides a fundamental understanding of the deformation mechanisms for monolayer of GDY family under a hypervelocity impact, which should shed lights on the design of GDY-based textile targeting specific bullet-proof application, or shielding structure for aerospace systems considering the impact velocity of projectile.