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.