Plain Language Summary
During the propagation of earthquakes, the creation of new fractures produces damage in the surrounding rock at all scales. However, the mechanism of nanoscale damage production remain enigmatic. For brittle materials, crack propagation may evolve from straight to a pattern that includes oscillations above a threshold rupture speed. Using α-quartz as a representative material for crystalline rocks, we show that above 15 % of the Rayleigh wave speed, crack propagation involves oscillations and microbranching events during which two cracks propagate simultaneously until one branch dominates. The aborted branch and the propagating branch produce rock fragments with nanoscale dimensions.
1 Introduction
The propagation of dynamic fractures during earthquakes produces damage along faults and in the surrounding rock that may lead to rock pulverization under high strain rates (e.g., Dor et al., 2006; Doan et al., 2009). Pulverized rocks are characterized by multi-scale fractures (Rempe et al., 2013) and rock fragments with sizes ranging from sub-micrometers to centimeters, often having power-law size distributions (Sammis et al., 1987; Keulen et al., 2007; Muto et al., 2015). The observation of asymmetric damage zones in crustal faults may occur due to the propagation of ruptures along a bi-material interface, producing high tensile stresses on one side of the fault (e.g., Andrews, 2005; Ben-Zion and Shi, 2005; Reches and Dewers, 2005; Di Toro and Pennacchioni, 2005; Petley-Ragan et al., 2019). During dynamic rupture, the rock splits into fragments with a wide range of sizes. Previous work has suggested that some materials have a minimum fragment size achieved through grinding, i.e., the grinding limit (Kendall, 1978). For quartz, the grinding limit is 0.9 µm (Prasher, 1987); however, quartz grains with radius of 15 nm have been observed in experiments and naturally produced gouges (Keulen et al., 2007). The observation of rock fragments with dimensions below the grinding limit prompts a reexamination of the mechanical origin of these fragments. Here, we use molecular dynamics simulations to identify a mechanism that produces nanoscale damage caused by dynamic rupture and due to out-of-plane unstable crack propagation under tensile loading conditions.
A dynamic fracture may develop as a straight crack that propagates in a stable manner along a single fracture plane (e.g., Bouchbinder et al., 2014), or it may become unstable when the crack tip oscillates in thez -direction, which may lead to microbranching (e.g., Bleyer and Molinari, 2017). The type of propagation varies with rupture speed. At low speeds, fracture surface roughness may develop with three morphologies: mirror, mist, or hackle (Fineberg et al., 1991; Buehler, 2008). For mirrors, the crack propagates at low speeds creating a perfect cleavage, and as the speed increases, the crack surface changes from slightly rough (mist) to significantly rough (hackle). The mirror-mist-hackle regimes have been observed in experiments (e.g., Fineberg et al., 1991), as well as in models containing a perfect atomic lattice (e.g., Abraham et al., 1997; Buehler and Gao, 2006). Dynamic instabilities at the crack tip result in rupture propagation speeds below the Rayleigh wave speed, C R, the rupture speed predicted for dynamic straight cracks by linear elastic fracture mechanics (Fineberg et al., 1991). The transition between the stable and unstable propagation regimes occurs when the crack tip speed exceeds a critical value that has been measured experimentally in several materials: 0.44C R for silica glass (Sharon and Fineberg, 1998), 0.44C R for polymethyl methacrylate (Fineberg et al., 1991; Sharon and Fineberg, 1996), and 0.34C R for neo-Hookean brittle gels (Livne et al., 2005). Above these critical rupture speeds, the crack deviates from a linear propagation and oscillate or branch into two or more cracks.
At the nanoscale, dynamic fractures nucleate when the breaking of an interatomic bond starts an avalanche of subsequent bond breakage, resulting in the propagation of a crack across length scales. Crack instabilities have been observed all the way down to the nanoscale (e.g., Zhou et al., 1996). Molecular dynamics simulations of a two-dimensional homogenous brittle solid with linear elastic properties have shown that crack branching occurs when the rupture speed exceeds a critical value, similar to experimental observations. The first such study reported a critical rupture speed of 0.35C R(Zhou et al., 1996), and a later study indicated a value of 0.73C R (Buehler and Gao, 2006). The difference between these critical speeds can be attributed to the shape of the interatomic potential used in the simulations. Zhou et al. (1996) used a Morse pair potential while Buehler and Gao (2006) used a smoothed biharmonic potential, which allows a hyperelastic zone size and the cohesive force between interacting atoms to be tuned. Bueher and Gao (2006) found that hyperelasticity and the ratio between the distance of maximum interatomic attraction and the equilibrium interatomic distance controls the critical rupture speed.
Here, we simulate dynamic tensile fracture and damage generation at the molecular level in α-quartz, an abundant mineral in the Earth’s continental crust, to quantify how rupture speed controls crack branching and damage production. Because this mineral requires a more detailed molecular potential than the Morse pair potential or the biharmonic potential used in previous studies (Zhou et al., 1996; Buehler and Gao, 2006), we use the Vashishta potential (Vashishta et al., 1990), which has been designed to reproduce the properties of crystalline quartz (Broughton et al., 1997). Our results quantify a mechanism of strain energy dissipation and damage production at the nanoscale throughout the tensile region of a fault zone during coseismic loading.
2 Methods
2.1 Simulation setup and procedure
The setup consists of a three-dimensional block of α-quartz with a pre-cracked notch structure (Figure 1a). The block is cut along three crystallographic planes, and the orientation of the system (Figure 1b) is such that the xy plane (z -cut) corresponds to the (\(0\ 0\ 0\ 1\)) crystallographic plane, the yz plane (x -cut) is the (\(\overset{\overline{}}{2}\ 1\ 1\ 0\)) crystallographic plane and the xz plane (y -cut) is the (\(0\ \overset{\overline{}}{1}\ 1\ 0\)) crystallographic plane. The system is periodic in the y -direction, while the boundaries at top, bottom and both ends in the x -direction are fixed in order to apply a specific stress state.