Mechanical Performance of Graphene-Based Artificial Nacres under Impact Loads: A Coarse-Grained Molecular Dynamic Study

Inspired by the hierarchical structure and outstanding mechanical performance of biological nacre, we propose a similar multi-layered graphene–polyethylene nanocomposite as a possible lightweight material with energy-absorbing characteristics. Through coarse-grained molecular dynamics simulations, we study the mechanical performance of the nanocomposite under spall loading. Results indicate that the polymer phase can serve as a cushion upon impact, which substantially decreases maximum contact forces and thus inhibits the breakage of covalent bonds in the graphene flakes. In addition, as the overlap distance in graphene layers increases, the energy absorption capacity of the model increases. Furthermore, the polymer phase can serve as a shield upon impact to protect the graphene phase from aggregation. The dependence of mechanical response on the size of impactors is also explored. Results indicate that the maximum contact force during the impact depends on the external surface area of impactors rather than the density of impactors and that the energy absorption for all model impactors is very similar. Overall, our findings can provide a systematic understanding of the mechanical responses on graphene–polyethylene nanocomposites under spall loads.


Introduction
The design of lightweight materials with energy-absorbing characteristics for blast and impact applications remains a big challenge for civil, aerospace, and military applications. However, nature provides us with plenty of inspirations. For example nacre [1], bone [2], and wood [3] are good cases of materials which are light-weight and exhibit high mechanical performance in terms of stiffness, toughness, and strength, stimulating the development of bio-inspired materials [4][5][6][7][8][9][10].
Nacre is also known as mother of pearl, and is a biological material with exceptional mechanical performance. Composed of 95 vol. % of layered aragonite platelets and 5 vol. % organic materials which serve as glue, nacre exhibits much higher fracture toughness (10 MPa·m 1/2 ) than pure aragonite (0.25 MP·m 1/2 ) [11]. Moreover, the strength of hydrated nacre (80 MPa) is still comparable to that of aragonite (160 MPa) [11]. Both stiff and soft components in a layered structure are necessary to achieve the high strength and toughness found in nacre. Although the highly mineralized layers contribute to its high strength, nacre would be quite brittle without the ability to dissipate strain. To achieve this, the interlayer shearing of the organic glue layers generates inelastic deformation which can remove high local stress [12], resulting in such high toughness.
Inspired by the extraordinary strength and toughness of nacre, plenty of efforts have been devoted to synthesize similar biomimetic materials. One successful example is graphene-based artificial nacre nanocomposites. There are plenty of publications adopting this strategy to fabricate graphene-based composites which show enhanced performance in terms of Young's modulus [13][14][15][16], glass transition temperature [16], and especially fracture strength [13,14,[17][18][19] and toughness [17][18][19]. Despite the successful synthesis of artificial nacres by adopting a cross-linking strategy to improve performance, the corresponding mechanism of how these added cross-links influence the nanoscale behavior of the polymer-graphene nanocomposites remains largely unexplored. In our recent work [6], the effect of interfacial cross-links on the mechanical performance of graphene-based artificial nacre has been systematically studied through coarse-grained molecular dynamics simulations, and the underlying molecular mechanism has been discussed. However, the performance of graphene-based artificial nacre under blast loads remains largely unexplored. Compared with full atomistic molecular dynamics, coarse-grained molecular dynamics (CGMD) is an ideal method to study the mechanical properties of the polymer-graphene nanocomposites due to its high computational efficiency and the wide range of spatial and temporal scales involved in the nanocomposite system. Therefore, coarse-grained molecular dynamics simulations are adopted in this paper to investigate the mechanical performance of graphene-based artificial nacre under impact loads.

Materials and Methods
As shown in Figure 1a, two different sets of graphene assembly are studied in the simulations: with and without the presence of polymers. For the first set-the graphene layers, typically composed of seven pieces of graphene flakes-are placed layer by layer in a staggered manner, connected by polyethylene (PE) glues with a thickness 1 nm. In order to create the computational model for the PE in between graphene layers in the nanocomposite, a self-avoiding walk (SAW) [20] algorithm is adopted. For the second set, the polymer phase is removed and only graphene flake layers are considered. Impactors are constructed using a face-centered cubic lattice with a lattice constant of 0.54 nm. To explore the dependence of impact response on the size and shape of impactors, three different impactors, named Im 1 , Im 2 , and Im 3 , are studied as shown in Figure 1b. Among them, two are solid cylinders with radiuses 1.9 nm and 3.8 nm, respectively, and the other one is a hollow cylinder with an outside radius 3.8 nm and an inside radius 3.29 nm. All the impactors have the same mass, and for all simulations, the initial velocity of the impactor is fixed to be a constant 50 Å/ps. The justification of our choice regarding the initial velocity of the impactor and the effect of that on the mechanical responses of targets can be found in the supporting information (see Figures S1 and S2 and the relevant discussions in the supporting information). Note that inside impactors Im 1 and Im 3 , the atoms have the same atomic mass, which is four times of that inside the impactor Im 2 . To explore the dependence of impact response on the size of the graphene layers, the overlap distance L OL varies from 4 to 24 nm, while the size of the simulation box (48 nm × 8.4 nm × 29 nm) remains the same as shown in Figure 1c,d.
To improve the simulation efficiency and thus reach a relatively large temporal and spatial scale, coarse-grained models are adopted to describe the interactions for both graphene and polyethylene. The justification of the potential force field can be found in the supporting information (see Figure S3 and relevant discussions in the supporting information). For polymer, a simple coarse-grained chain model is adopted to reproduce the atomistic structure of polyethylene [21,22]. Each bead represents a CH 2 group in the middle or a CH 3 group at the end of the polymer chain. Adjacent beads are connected by harmonic bonds. Three-body bending interactions are represented by a harmonic function of the cosine of the bond angle θ. Four-body torsion deformations are described by a third-order polynomial function of the cosine of the dihedral angle ϕ. Non-bonded interactions are captured by the LJ 12-6 potential with a cutoff of 2.5σ, where σ is the interatomic distance at which the LJ potential equals to zero. The corresponding forms and parameters of the potential field are listed in Table 1. With respect to graphene sheets, we follow the coarse-grained model from Ruiz et al. [23] to describe the relevant structure and interactions. In this model, the hexagonal lattice of graphene remains, in which each bead represents four atoms of the atomistic lattice. Interactions between adjacent beads are captured by the Morse bond potential. Three-body bending deformations are described by a harmonic function of the bond angle θ. Four-body torsion interactions are represented by a harmonic-type dihedral function of the dihedral angle ϕ. With respect to the non-bonded interactions, a LJ 12-6 potential is adopted with a cutoff of 2.5σ. The justification of the choice regarding cutoff distance can be found in the supporting information (see Figure S4 and relevant discussions in the supporting information). The corresponding forms and parameters of the model are listed in Table 2. Note that the non-bonded interactions between polymer beads and graphene beads are determined by the Lorentz-Berthelot mixing rule [24]. The impactors were fixed as rigid bodies, and the interactions between the impactors and the graphene-polyethylene composites are also described using LJ potential. Corresponding parameters like ε and σ are the same as those used for interactions between graphene and polyethylene. Note that the above potential function is truncated and shifted to zero at 6 √ 2σ so that only repulsive forces exist between impactors and graphene-polyethylene composites. represented by a harmonic-type dihedral function of the dihedral angle φ. With respect to the nonbonded interactions, a LJ 12-6 potential is adopted with a cutoff of 2.5σ. The justification of the choice regarding cutoff distance can be found in the supporting information (see Figure S4 and relevant discussions in the supporting information). The corresponding forms and parameters of the model are listed in Table 2. Note that the non-bonded interactions between polymer beads and graphene beads are determined by the Lorentz-Berthelot mixing rule [24]. The impactors were fixed as rigid bodies, and the interactions between the impactors and the graphene-polyethylene composites are also described using LJ potential. Corresponding parameters like ε and σ are the same as those used for interactions between graphene and polyethylene. Note that the above potential function is truncated and shifted to zero at √2 σ so that only repulsive forces exist between impactors and graphene-polyethylene composites. (a) Schematic view of the initial setup of the simulation (red represents graphene fibers, blue represents polymer glues, and green represents impactors); (b) Impactors used in the simulations are named as Im1, Im2, and Im3 (all with the same mass and the length of each impactor is 8.4 nm, the same size as the out-of-plane dimension of the sample); (c) The samples to be tested in the impact simulations are named as S1, S2, S3, S4, and S5 (with polymer glues). The length of the sample is fixed at 48 nm while the overlap distance varies from 4 to 24 nm; (d) The samples to be tested in the impact simulations are named as S1', S2', S3', S4' and S5' (without polymer glues). Table 1. Functional forms and parameters of force field for the CG polyethylene (PE). Non-bonded = 4ε σ − σ = 57 K, σ = 4.28 Å, is Boltzmann's constant, is the distance between beads i and j Figure 1. Geometrical configurations of computational models. (a) Schematic view of the initial set-up of the simulation (red represents graphene fibers, blue represents polymer glues, and green represents impactors); (b) Impactors used in the simulations are named as Im 1 , Im 2 , and Im 3 (all with the same mass and the length of each impactor is 8.4 nm, the same size as the out-of-plane dimension of the sample); (c) The samples to be tested in the impact simulations are named as S 1 , S 2 , S 3 , S 4 , and S 5 (with polymer glues). The length of the sample is fixed at 48 nm while the overlap distance L OL varies from 4 to 24 nm; (d) The samples to be tested in the impact simulations are named as S 1 ', S 2 ', S 3 ', S 4 ' and S 5 ' (without polymer glues). Table 1. Functional forms and parameters of force field for the CG polyethylene (PE).

Type of interaction Potential form and parameters
Bond ε kB = 57 K, σ = 4.28 Å, k B is Boltzmann's constant, r ij is the distance between beads i and j Table 2. Functional forms and parameters of force field for the CG Graphene.

Type of interaction Potential form and parameters
Bond ε LJ = 0.82 kcal·mol −1 , σ = 3.46 Å, r ij is the distance between beads i and j.
To perform coarse-grained molecular dynamics simulations, the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) [25] was used. Before impact simulations, samples were relaxed under the canonical (NVT) ensemble for 400 ps at temperature 200 K in order to make the model physically reasonable. Periodic boundary conditions were adopted along x and y direction, while the z direction had a free boundary condition to prevent self-interference of the material along the direction of impact. Unless otherwise stated, time step is 1 femtosecond for all simulation runs to maintain a good balance between computational efficiency and accuracy. Once the initial models were generated, annealing was performed where the samples were heated up to 500 K over 1400 picoseconds using the NVT ensemble and then cooled down to 200 K in another 1400 picoseconds using the isothermal-isobaric (NPT) ensemble, where the pressure was set equal to 100 kPa. Finally, impactors were placed 10 nm above the top surface of the target, and the initial velocities of the impactors were given as 5 nm/ps. Subsequently, impactors were released and the microcanonical (NVE) ensemble was used during impact simulations.

Effects of Polymer Glues
In this section, the effect of polyethylene polymers on mechanical responses is investigated during the impact simulations; for these cases, the impactor Im 1 was used. In all scenarios, the initial kinetic energy for the sample was fixed at 18.253 keV (2.92 fJ). Figure 2 shows the results for samples with (S 1 ) and without (S 1 ') polymers when the overlap distance was 24 nm. The discussion regarding temperature change can be found in the supporting information (see Figure S5 and the relevant discussions in the supporting information). According to Figure 2a, the interacting forces between the composites and the impactor remain zero at the beginning due to there being no contact between composites and impactors. At around 2 ps, the interacting forces between the impactor and the sample S 1 undergo a sharp spike upon impact and eventually decay to zero. Similarly, at around 2 ps, the interacting forces between the impactor and the sample S 1 ' increase dramatically and then decrease to zero. Compared with the latter case, the maximum force for S 1 is much smaller due to the buffering effect of polymers during impact. As shown in Figure 2b, after a short period of silence, the potential energy undergoes an intensive increase for both S 1 and S 1 '. Subsequently, the potential energy decreases a small amount due to the breaking of covalent bonds in graphene sheets. After that, it increases moderately, and finally enters a plateau. Note that the final change of potential energy for sample S 1 is bigger than that for S 1 '. To further clarify the origin of the above differences, several snapshots during the impact are singled out and presented in Figure 3. The dynamic process of bond-breaking inside the target can be found in the supporting information (see Figure S6 in the supporting information for details). According to Figure 3a, the contact zone of sample S 1 was highly densified as the impactor hit. Subsequently, sample S 1 bent and the polymers on the top layer detached from their adjacent graphene layers. Later on, some of the flakes on the top layer were broken. In addition, due to the intensive movement of polymers, a void occurred inside the sample. Finally, the bottom layer of graphene bent, slithered, and detached from the whole sample. For S 1 ' (the sample without polymers), the dynamics were similar as seen in Figure 3b. However, unlike S 1 , the impactor penetrated through the top layer, breaking it into two parts. Thus, from Figures 2 and 3, it can be seen that polyethylene polymers can serve as a buffer, decreasing the maximum contact forces and thus protecting the embedded graphene layers from fracture. supporting information for details). According to Figure 3a, the contact zone of sample S1 was highly densified as the impactor hit. Subsequently, sample S1 bent and the polymers on the top layer detached from their adjacent graphene layers. Later on, some of the flakes on the top layer were broken. In addition, due to the intensive movement of polymers, a void occurred inside the sample. Finally, the bottom layer of graphene bent, slithered, and detached from the whole sample. For S1' (the sample without polymers), the dynamics were similar as seen in Figure 3b. However, unlike S1, the impactor penetrated through the top layer, breaking it into two parts. Thus, from Figures 2 and 3, it can be seen that polyethylene polymers can serve as a buffer, decreasing the maximum contact forces and thus protecting the embedded graphene layers from fracture.   supporting information for details). According to Figure 3a, the contact zone of sample S1 was highly densified as the impactor hit. Subsequently, sample S1 bent and the polymers on the top layer detached from their adjacent graphene layers. Later on, some of the flakes on the top layer were broken. In addition, due to the intensive movement of polymers, a void occurred inside the sample. Finally, the bottom layer of graphene bent, slithered, and detached from the whole sample. For S1' (the sample without polymers), the dynamics were similar as seen in Figure 3b. However, unlike S1, the impactor penetrated through the top layer, breaking it into two parts. Thus, from Figures 2 and  3, it can be seen that polyethylene polymers can serve as a buffer, decreasing the maximum contact forces and thus protecting the embedded graphene layers from fracture.   To verify the above conclusion, multiple simulation scenarios were considered by varying the overlap distance L OL . Figure 4a shows the maximum contact force during the impact, indicating that the maximum force on samples with polymers is much smaller than those counterparts without polymers for all overlap distances, L OL , covered. In addition, the maximum force only fluctuates to a small extent, as the overlap distance varies for samples both with and without polymers, indicating the negligible effect of overlap distance on maximum contact forces. When the impactor hits the target with polymer, it first contacts with polymers. The polymers are compressed intensively and the velocity of the impactor is decreased. Therefore, when the impactor starts to contact with the graphene layers, the interacting forces increase moderately compared with the pure graphene sample. In contrast, for samples without polymers, the impactor directly hit the graphene layers, resulting in dramatic increase of contact forces. In Figure 4b, the maximum potential change as marked in Figure 2b is shown as a function of the overlap distance. For the samples with polymers, the maximum potential energy for the time of impact, Pe 1 , is smaller than that of the final state of the sample, Pe 2 , for all overlap distances covered. In addition, as the overlap distance increases, the maximum potential energy, both Pe 1 and Pe 2 , increases due to the strengthening connections inside the sample. For samples without polymers, Pe 1 is slightly bigger than Pe 2 regardless of changing the overlap distance; for the samples with polymers, this trend is reversed. Due to the absence of polymers, the bonds of graphene flakes are more susceptible to the impact, which can be confirmed by the snapshots displayed in Figure 3. The large number of broken bonds destroys the integrity of the material system, degrading its energy absorption capacity through irreversible deformation. In summary, the polyethylene polymer phase serves as cushion when encountering ballistic loads. To verify the above conclusion, multiple simulation scenarios were considered by varying the overlap distance LOL. Figure 4a shows the maximum contact force during the impact, indicating that the maximum force on samples with polymers is much smaller than those counterparts without polymers for all overlap distances, LOL, covered. In addition, the maximum force only fluctuates to a small extent, as the overlap distance varies for samples both with and without polymers, indicating the negligible effect of overlap distance on maximum contact forces. When the impactor hits the target with polymer, it first contacts with polymers. The polymers are compressed intensively and the velocity of the impactor is decreased. Therefore, when the impactor starts to contact with the graphene layers, the interacting forces increase moderately compared with the pure graphene sample. In contrast, for samples without polymers, the impactor directly hit the graphene layers, resulting in dramatic increase of contact forces. In Figure 4b, the maximum potential change as marked in Figure 2b is shown as a function of the overlap distance. For the samples with polymers, the maximum potential energy for the time of impact, Pe1, is smaller than that of the final state of the sample, Pe2, for all overlap distances covered. In addition, as the overlap distance increases, the maximum potential energy, both Pe1 and Pe2, increases due to the strengthening connections inside the sample. For samples without polymers, Pe1 is slightly bigger than Pe2 regardless of changing the overlap distance; for the samples with polymers, this trend is reversed. Due to the absence of polymers, the bonds of graphene flakes are more susceptible to the impact, which can be confirmed by the snapshots displayed in Figure 3. The large number of broken bonds destroys the integrity of the material system, degrading its energy absorption capacity through irreversible deformation. In summary, the polyethylene polymer phase serves as cushion when encountering ballistic loads.

Effects of Impactor Size
In this section, the effect of the impactor size on the mechanical response during impact is investigated. For these impact simulations, three different impactors with the same mass are used, as shown in Figure 1b. In addition, the initial kinetic energy of the impactors are fixed at 18.253 keV (2.92 fJ) in all scenarios. Contact forces during the impact are shown in Figure 5 for both sample S5 and S5', indicating that the forces only exist for a very short period-around 1 ps in our simulations. At the beginning, the force remains zero until the impactor hits the sample. Subsequently, the force increases intensively to the peak and then decays rapidly to zero. In addition, the maximum forces generated by impactor Im1 are smaller than those generated by impactors Im2 and Im3. Since the outer radii for impactors Im2 and Im3 are both 3.8 nm (twice that of impactor Im1), the contact areas between the impactors (Im2 and Im3) and the samples are much larger than that between the impactor (Im1) and the samples, resulting in larger contact forces. In Figure 6, the change in potential energy is shown

Effects of Impactor Size
In this section, the effect of the impactor size on the mechanical response during impact is investigated. For these impact simulations, three different impactors with the same mass are used, as shown in Figure 1b. In addition, the initial kinetic energy of the impactors are fixed at 18.253 keV (2.92 fJ) in all scenarios. Contact forces during the impact are shown in Figure 5 for both sample S 5 and S 5 ', indicating that the forces only exist for a very short period-around 1 ps in our simulations. At the beginning, the force remains zero until the impactor hits the sample. Subsequently, the force increases intensively to the peak and then decays rapidly to zero. In addition, the maximum forces generated by impactor Im 1 are smaller than those generated by impactors Im 2 and Im 3 . Since the outer radii for impactors Im 2 and Im 3 are both 3.8 nm (twice that of impactor Im 1 ), the contact areas between the impactors (Im 2 and Im 3 ) and the samples are much larger than that between the impactor (Im 1 ) and the samples, resulting in larger contact forces. In Figure 6, the change in potential energy is shown as a function of time. As seen in the previous section, after a short period of fluctuation, the potential energy increases fiercely to its first peak, decreases slightly, and then increases moderately to a plateau value. Interestingly, although different impactors are used to perform the impact simulations, the differences among the potential energy change using different impactors are negligible. Figure 7 shows the dynamic process of the impact simulation (the dynamic process of bond-breaking can be found in Figures S7-S9 in the supporting information), indicating that all the impactors penetrate through the graphene-polyethylene composites regardless of the size of the impactors. In addition, from a comparison of Figure 7a-c with Figure 7d-f, the degree of damage is very similar even though impactors with different sizes are used. Interestingly, it is worth noting that in the samples with polymers, some of the graphene grains remain intact, even though the grains in the middle are seriously damaged. However, according to Figure 7d-f, not only are the grains in the middle damaged, but the other grains also aggregate together. This implies that the polymers can serve as a shield for the embedded graphene grains in the composite structure.
Polymers 2017, 9,134 7 of 11 as a function of time. As seen in the previous section, after a short period of fluctuation, the potential energy increases fiercely to its first peak, decreases slightly, and then increases moderately to a plateau value. Interestingly, although different impactors are used to perform the impact simulations, the differences among the potential energy change using different impactors are negligible. Figure 7 shows the dynamic process of the impact simulation (the dynamic process of bond-breaking can be found in Figures S7-S9 in the supporting information), indicating that all the impactors penetrate through the graphene-polyethylene composites regardless of the size of the impactors. In addition, from a comparison of Figure 7a-c with Figure 7d-f, the degree of damage is very similar even though impactors with different sizes are used. Interestingly, it is worth noting that in the samples with polymers, some of the graphene grains remain intact, even though the grains in the middle are seriously damaged. However, according to Figure 7d-f, not only are the grains in the middle damaged, but the other grains also aggregate together. This implies that the polymers can serve as a shield for the embedded graphene grains in the composite structure.   as a function of time. As seen in the previous section, after a short period of fluctuation, the potential energy increases fiercely to its first peak, decreases slightly, and then increases moderately to a plateau value. Interestingly, although different impactors are used to perform the impact simulations, the differences among the potential energy change using different impactors are negligible. Figure 7 shows the dynamic process of the impact simulation (the dynamic process of bond-breaking can be found in Figures S7-S9 in the supporting information), indicating that all the impactors penetrate through the graphene-polyethylene composites regardless of the size of the impactors. In addition, from a comparison of Figure 7a-c with Figure 7d-f, the degree of damage is very similar even though impactors with different sizes are used. Interestingly, it is worth noting that in the samples with polymers, some of the graphene grains remain intact, even though the grains in the middle are seriously damaged. However, according to Figure 7d-f, not only are the grains in the middle damaged, but the other grains also aggregate together. This implies that the polymers can serve as a shield for the embedded graphene grains in the composite structure.   Simulations are performed to identify the role of impactor size in the mechanical responses of samples with different overlap distance, and results are shown in Figure 8. As expected, the maximum forces during the impact fluctuate as the overlap distance increases. Moreover, the maximum forces for impactors Im2 and Im3 always remain very similar, and both are much greater than the forces caused by impactor Im1. Although the density of Im3 is much greater than that of Im2, the external (contact) surface area of these two impactors are the same, and both are greater than for impactor Im1. This leads to the differences in maximum forces shown in Figure 8a. Furthermore, due to the similar force trajectories between Im2 and Im3, the maximum potential energy of the initial impact, Pe1, for impactor Im2 remains close to that for impactor Im3 as the overlap distance changes in Figure 8b. In addition, the maximum potential energy Pe1 from impactor Im1 is much smaller than those for impactors Im2 and Im3. However, in the final state, the differences among scenarios using different impactors in terms of maximum potential energy Pe2 are very small, as shown in Figure 8c. Although impactors have different size, they do have the same mass. In addition, they have the same initial velocity. Therefore, their impact energies are exactly the same. Another important factor influencing the potential energy change is the effective contact region between the impactor and the target. For example, the inner region of Im2, has no interaction with the target, and thus the effective contact region only includes the outer region. Therefore, the effective contact region of Im2 is pretty close to that of Im3. As a result, the responses under impact are pretty close to each other for Im2 and Im3, as shown in Figure 8. In summary, it might suggest that the velocity and mass of the impactor play a more important role in evaluating the energy capacity of the targets than the size and/or density of the impactor. Simulations are performed to identify the role of impactor size in the mechanical responses of samples with different overlap distance, and results are shown in Figure 8. As expected, the maximum forces during the impact fluctuate as the overlap distance increases. Moreover, the maximum forces for impactors Im 2 and Im 3 always remain very similar, and both are much greater than the forces caused by impactor Im 1 . Although the density of Im 3 is much greater than that of Im 2 , the external (contact) surface area of these two impactors are the same, and both are greater than for impactor Im 1 . This leads to the differences in maximum forces shown in Figure 8a. Furthermore, due to the similar force trajectories between Im 2 and Im 3 , the maximum potential energy of the initial impact, Pe 1 , for impactor Im 2 remains close to that for impactor Im 3 as the overlap distance changes in Figure 8b. In addition, the maximum potential energy Pe 1 from impactor Im 1 is much smaller than those for impactors Im 2 and Im 3 . However, in the final state, the differences among scenarios using different impactors in terms of maximum potential energy Pe 2 are very small, as shown in Figure 8c. Although impactors have different size, they do have the same mass. In addition, they have the same initial velocity. Therefore, their impact energies are exactly the same. Another important factor influencing the potential energy change is the effective contact region between the impactor and the target. For example, the inner region of Im 2 , has no interaction with the target, and thus the effective contact region only includes the outer region. Therefore, the effective contact region of Im 2 is pretty close to that of Im 3 . As a result, the responses under impact are pretty close to each other for Im 2 and Im 3 , as shown in Figure 8. In summary, it might suggest that the velocity and mass of the impactor play a more important role in evaluating the energy capacity of the targets than the size and/or density of the impactor.

Conclusions
In this paper, coarse-grained molecular dynamics simulations are performed to investigate the mechanical responses of graphene-polyethylene nanocomposites upon spall-like impact. Results indicate that the polymer phase can serve as a cushion upon impact, which substantially decreases maximum contact forces and thus inhibits the breakage of covalent bonds in the graphene flakes. In addition, as the overlap distance in graphene layers increases, the energy absorption capacity of the model increases. Furthermore, the polymer phase can serve as a shield upon impact to protect the graphene phase from aggregation. The dependence of mechanical response on the size of impactors is also explored. Results indicate that the maximum contact force during the impact depends on the external surface area of impactors rather than the density of impactors, and that the energy absorption for all model impactors is very similar. Overall, our findings can provide a systematic understanding of the mechanical responses on graphene-polyethylene nanocomposites under spall loads.
Supplementary Materials: The following are available online at www.mdpi.com/2073-4360/9/4/134/s1. Figure  S1: Mechanical responses under different impact loads (achieved through varying the initial velocity of the impactor) (a) Reaction forces (b) Potential energy change (The impactor is Im1 and the target is S1). Figure S2: Snapshots of the whole system at 30ps under different impact loads (The impactor is Im1 and the target is S1). Figure S3: Schematic view of the coarse-graining strategy for graphene (a) Coarse-grain lattices overlaid over the atomistic structure (white) (b) Illustration of the contributions of the coarse-grained force field. Coarse-grain lattices are colored by blue while the bonded interactions are highlighted ball-stick representations in red. Note that non-bonded interactions are represented by virtue lines in red. Figure S4: The effect of cutoff distance of force field on responses of target (a) Reaction forces (b) Potential energy change (The impactor is Im1 and the target is S1). Figure S5: Temperature change during the impact simulations. Figure S6: Dynamic process of bond breaking (a) Im1 S1 (b) Im1 S1' (Note that in Figure S2a polymer chains are removed for clarity). Figure S7: Dynamic process of bond breaking (a) Im1 S5 (b) Im1 S5' (Note that in Figure S3a polymer chains are removed for clarity). Figure S8: Dynamic process of bond breaking (a) Im2 S5 (b) Im2 S5'. Figure S9: Dynamic process of bond breaking (a) Im3 S5 (b)Im=3 S5'.

Conclusions
In this paper, coarse-grained molecular dynamics simulations are performed to investigate the mechanical responses of graphene-polyethylene nanocomposites upon spall-like impact. Results indicate that the polymer phase can serve as a cushion upon impact, which substantially decreases maximum contact forces and thus inhibits the breakage of covalent bonds in the graphene flakes. In addition, as the overlap distance in graphene layers increases, the energy absorption capacity of the model increases. Furthermore, the polymer phase can serve as a shield upon impact to protect the graphene phase from aggregation. The dependence of mechanical response on the size of impactors is also explored. Results indicate that the maximum contact force during the impact depends on the external surface area of impactors rather than the density of impactors, and that the energy absorption for all model impactors is very similar. Overall, our findings can provide a systematic understanding of the mechanical responses on graphene-polyethylene nanocomposites under spall loads.
Supplementary Materials: The following are available online at www.mdpi.com/2073-4360/9/4/134/s1. Figure S1: Mechanical responses under different impact loads (achieved through varying the initial velocity of the impactor) (a) Reaction forces (b) Potential energy change (The impactor is Im 1 and the target is S 1 ). Figure S2: Snapshots of the whole system at 30ps under different impact loads (The impactor is Im 1 and the target is S 1 ). Figure S3: Schematic view of the coarse-graining strategy for graphene (a) Coarse-grain lattices overlaid over the atomistic structure (white) (b) Illustration of the contributions of the coarse-grained force field. Coarse-grain lattices are colored by blue while the bonded interactions are highlighted ball-stick representations in red. Note that non-bonded interactions are represented by virtue lines in red. Figure S4: The effect of cutoff distance of force field on responses of target (a) Reaction forces (b) Potential energy change (The impactor is Im 1 and the target is S 1 ). Figure S5: Temperature change during the impact simulations. Figure S6: Dynamic process of bond breaking (a) Im 1 S 1 (b) Im 1 S 1 ' (Note that in Figure S2a polymer chains are removed for clarity). Figure S7: Dynamic process of bond breaking (a) Im 1 S 5 (b) Im 1 S 5 ' (Note that in Figure S3a polymer chains are removed for clarity). Figure S8: Dynamic process of bond breaking (a) Im 2 S 5 (b) Im 2 S 5 '. Figure S9: Dynamic process of bond breaking (a) Im 3 S 5 (b) Im 3 S 5 '.

Conflicts of Interest:
The authors declare no conflict of interest.