Dynamic antiplane self-similar crack with distance-weakening friction: an analytical solution for source parameter estimation

Shiro Hirano$^1$ & Hiromichi Itou$^2$
$^1$Dept. Physical Science, Ritsumeikan Univ., Japan; $^2$Dept. Mathematics, Tokyo Univ. of Science, Japan.


Author's website (discussions via e-mail are welcome)
This work has been published as the following OpenAccess paper:
Hirano, S., & Itou, H. (2020). Parameter interdependence of dynamic self-similar crack with distance-weakening friction, Geophys. J. Int., 223(3):1584$-$1596, doi:10.1093/gji/ggaa392


Movie 1: Slip-rate (black) and slip (blue) functions on a dynamic anti-plane self-similar crack analytically derived by Kostrov (top) and this study (bottom).

Slip rate function and energy release rate of a dynamic self-similar crack model have contributed to a theoretical understanding and estimation of on-fault energetics. Although Kostrov's self-similar singular crack model provided an unbounded(=unphysical) slip-rate function characterized by square-root singularity, his analytical work revealed a relationship between rupture velocity and energy release rate.

As an extension of the self-similar crack model, we introduce some weakening zone of friction behind crack tips to bound peak slip rate and clarify some relationships among physical parameters, including stress state, process zone size, rupture velocity, peak slip rate, and energy release rate.

Our analytical solution within the framework of linear elastic fracture mechanics tells us that:

  1. the rupture velocity is determined when a fracture energy gradient and initial stress are given.
  2. the rupture velocity is restricted below the shear wave speed if peak slip rate < 5.7 m/s holds, as previously implied by a numerical off-fault damage model, meaning that our model can mimic off-fault inelastic energy dissipation.
  3. a convenient approximation holds to estimate fracture energy based on seismically observable parameters: rupture velocity and peak slip rate.

1. Introduction

Classical self-similar crack solution
Multi-scale inversion analysis of M2$-$M6-class events [Uchide+ 2010], the onset of ground displacement [Meier+ 2016], and AE event analysis based on rock experiments [Yoshimitsu+ 2007] have revealed that seismic moment evolution follows a unified power law: $\displaystyle \mu \int_\text{fault} D(x,t) \, dx \propto t^3, $ where $\mu$ is the rigidity, and $D(x,t)$ is slip amount on the fault at position $x$ and time $t \, (>0)$, which implies the self-similarity: \begin{align} V(x,t) = V(X), \label{SelfSimilarity} \end{align} where $V := \partial_t D$ is slip velocity, $X := \dfrac{x}{ct}$ is the scaled position, and $c$ is rupture velocity.
For 2-D and 3-D elastic body, Kostrov [1964] derived an analytical solution to the self-similar crack model with constant rupture velocity ratio $k := \dfrac{c}{\beta}$ ($\beta$ is S-wave velocity) and stress drop $\Delta \sigma$ as \begin{align} V(X) = \frac{2 c \Delta \sigma}{\mu I(k) \sqrt{1-|X|^2}} \label{Kostrov} \end{align} (see also the top of Movie 1) and energy release rate as a function of $k$. Especially, $I(k) = E(\sqrt{1-k^2})$ in case of mode-III rupture, where $E(\cdot)$ is the complete elliptic integral of the second kind.
Problem after Kostrov's UNPHYSICAL solution \eqref{Kostrov}
Slip rate is unbounded due to square-root singularity at the tip, i.e., $\displaystyle \lim_{|X| \to 1}V(X) = \infty$. But this should be bounded and $\displaystyle \lim_{|X| \to 1}V(X) = 0$ must hold.
How can we consider some physical parameters: peak stress, peak slip rate, process zone size, and critical slip-weakening distance?
Aim of this study: derivation of an analytical solution under physically reasonable condition, $V(\pm 1)=0$
Self-similarly extending process zones should be introduced behind crack tips to restrict slip rate.
We may obtain a bounded slip rate function, and some relationships among stress state, peak slip rate, process zone size, rupture velocity, and energy release rate.

2. Mathematical Modeling

  1. For an antiplane case, on-fault ($z=0$) and fault-parallel shear traction perturbation based on Aki & Richards' [2002] representation theorem is: \begin{align} T(x,t) = - \frac{\mu \beta}{2 \pi} \lim_{z \to 0} \frac{\partial^2}{\partial z^2} \int_0^t \left( \int_{-vt}^{+vt} \frac{D(x',t')}{\sqrt{\beta^2 (t-t')^2 - (x-x')^2-z^2}} dx' \right) dt'. \label{TofD} \end{align}
  2. We assume symmetry (i.e., $V(-x,t) = V(x,t)$). Eqs.\eqref{SelfSimilarity}\eqref{TofD} with $ X := \dfrac{x}{t}, \, X' := \dfrac{x'}{t'}, $ results in \begin{align} \frac{d T(X)}{d X} = \frac{\mu}{2 \beta} \frac{\sqrt{\beta^2-X^2}}{X} \frac{d}{dX} \int_{-1}^{+1} \frac{V(X')}{X - X'} \frac{dX'}{\pi}, \label{IntEq} \end{align} and its general solution, including an arbitrary constant $C$, is [Tricomi 1985] \begin{align} \hspace{-2em} V(X) = \frac{1}{\sqrt{1-X^2}} \left( C - \frac{2\beta}{\mu} \int_{-1}^{+1} \frac{\sqrt{1-X'^2}}{X-X'} \left( \int_0^{X'}\frac{\varphi T'(\varphi) \, d\varphi}{\sqrt{k^2-\varphi^2}} \right) \frac{dX'}{\pi} \right). \label{GeneralSolution} \end{align}
  3. We assume linear distance-weakening friction shown in Fig. 1.
    Fig. 1. A distance-weakening friction model.
    Finally, given stress ratio $S:= \dfrac{\sigma_b}{\Delta \sigma} - 1 = \dfrac{\sigma_E}{\Delta \sigma}$ [Andrews 1976], we get the following solution: \begin{align} \frac{\mu V(X)}{\Delta \sigma \beta} = \frac{1+S}{1 - X_0} \frac{2}{\pi} \left( \frac{s(k X)}{k} L\left(\frac{s(X_0) s(k X)}{s(X) s(k X_0)}\right) - \frac{s(k X_0)}{k} L\left(\frac{s(X_0)}{s(X)}\right) - s(X) L\left(\frac{k s(X_0)}{s(k X_0)}\right) \right), \label{Solution} \end{align} where $L(x) = \log\left|\dfrac{1-x}{1+x}\right|$ and $s(x) = \sqrt{1-x^2}$.
  4. RHS of eq.\eqref{GeneralSolution} must be zero as $X \to \pm 1$ to satisfy $V(\pm 1)=0$. Therefore, \begin{align} C = \lim_{X \to \pm 1} \frac{2\beta}{\mu} \int_{-1}^{+1} \frac{\sqrt{1-X'^2}}{X-X'} \left( \int_0^{X'}\frac{\varphi T'(\varphi) d\varphi}{\sqrt{k^2-\varphi^2}} \right) \frac{d X'}{\pi} \label{C} \end{align} must hold. Moreover, eq.\eqref{C} and $T(X)|_{X=0} = -\Delta \sigma$ yields \begin{align} \frac{1}{1+S} = \frac{E(s(k))}{\pi k^2 (1 - X_0)} \left( \frac{k^{-1}-k}{2} L\left(\frac{k s(X_0)}{s(k X_0)}\right) + s(X_0) s(k X_0) \right), \label{RkS} \end{align} which means that $k$, $S$, and relative process zone length $R := 1 - X_0$ depend on each other.
Fig. 2. Schematic illustration of crack-like and pulse-like slip rate.
Characteristics of self-similar crack-like rupture mode
Process zone size and energy release rate are proportional to time [Hirano 2020].
We found that $R$ depends on $S$ and $k$.
Our analytical solution consists of eqs.\eqref{Solution}\eqref{RkS} enables us to investigate the interdependence of parameters: $S$, $k$, peak slip rate, and energy release rate.
In case of steady-state pulse-like rupture mode
A counterpart of the self-similar crack
Process zone size and energy release rate are constant, meaning that any value of $0 < k < 1$ is possible.
More simplified comparing to eq.\eqref{IntEq}:
$$ T(X) = \displaystyle \frac{\mu \sqrt{1-k^2}}{2 c} \int_{-1}^{+1} \frac{V(X')}{X' - X} \frac{dX'}{\pi}, $$ where $X := x-ct$, that can be solved in a straightforward manner based on the theory of singular integral equation [Muskhelishvili 1953].
A bounded slip-rate function has been obtained for the model, in which $S$ depends only on $R$ (i.e., independent of $k$) [Rice+ 2005].

3. Results

Fig. 3: $S = S(k, R)$ based on eq.\eqref{RkS} to fulfill $V(\pm 1)=0$.
Relationship among $k$, $R$, and $S$
Given each fixed value of $S \le 2.625$, the minimum value of relative rupture velocity $k$ exists and increases as $S$ decreases (white crosses in Fig. 3).
Given each fixed value of $S$ and $k$, process zone size $R$ can be a two-valued function (e.g., $R \sim 0.2$ and $0.85$ for $k=0.6$ and $S=2$, as white circles in Fig. 3). Physical meaning, stability, and reasonability of the larger $R$ values (e.g., $>0.5$) are still veiled. But we assume that $R < 0.3$ in the following analyses.
Given each fixed value of $S$, process zone size $R (< 0.3)$ shrinks as rupture accelerates. However, its rate differs from the Lorentz contraction (i.e., $R \propto \sqrt{1-k^2}$) of pulse-like rupture.
Time- and slip-weakening friction
Given self-similar distance-weakening friction (Fig.4a) is equivalent to self-similar time- and slip-weakening (Fig.4b,c) friction.
Self-similar Energy release rate
The energy release rate is numerically calculated as the area below the slip-weakening curve (Fig.5a), that is, $$ \dfrac{\mu}{\Delta \sigma^2} \dfrac{\mathcal{G}}{x} = \dfrac{1}{k} \displaystyle \int_{X_0}^1 \dfrac{T(X) + \Delta \sigma}{\Delta \sigma} \dfrac{\mu V(X)}{\Delta \sigma \beta} \dfrac{dX}{X^2}.$$
Fig. 4: Distribution of slip rate (solid), slip (broken; right axis), and traction change(short-broken) in (a) Space, and (b) Time, and (c) slip-stress curve for $S = 2$ and $k = 0.9$. The grey rectangles in (a) and (b) represent the process zone and weakening time, respectively. Curves for $c t = 0.25, 0.5, 0.75$ (grey), and $1.0$ (black) are plotted in (a), while curves for $x = 0.25, 0.5, 0.75$ (grey), and $1.0$ (black) are plotted in (c).
Fig. 5: (a) Normalized energy release rate equivalent to area below the slip-weakening curve in Fig. 4c and (b) Normalized peak slip rate picked from slip-rate function in Fig. 4a as a function of $k$ and $S$. The left-bottom blank region signifies the absence of a solution satisfying $V(\pm 1) = 0$.

4. Discussion

Determination of rupture velocity
Fig. 5a shows that normalized energy release rate, or energy release rate gradient, $\dfrac{\mu}{\Delta \sigma^2} \dfrac{\mathcal{G}}{x}$ is a decreasing function of $k$. However, less sensitive to $S$. Thus, rupture velocity would be determined automatically if $S$ and $\dfrac{\mathcal{G}}{x}$ are given, which is different from the result of steady-state pulse-like rupture mode that rupture velocity can be arbitrary even if $S$ and $\mathcal{G}$ are given [Rice+ 2005].
Fig. 5b shows a relationship among $k$, $S$, and peak slip rate $V_\text{peak}$, meaning that rupture velocity is automatically determined if $S$ and $V_\text{peak}$ are given.
Andrews [2005] numerically demonstrated that artificial restriction of $V_\text{peak}$ yields macroscopically the same energy dissipation and limitation of rupture velocity as of off-fault damage model even under linear elasticity. For example, parameters assumed by Andrews (Table 1) result in $\dfrac{\mu V_\textrm{peak}}{\Delta \sigma \beta} = 4.6$ and suggest $k \sim 0.63$ for $S=3$ and $k \sim 0.75$ for $S=2$, etc (Fig. 5b). Hence, our analytical result may approximate Andrews' numerical model.
$V_\textrm{peak}$ might be a threshold that energy dissipation due to off-fault damage irrupts and if so, would be estimated by experimental studies of rocks in the future.
Conclusion: under the fixed value of $S$, rupture velocity is determined if the energy release rate gradient OR peak slip rate is given.
Table 1: Parameters assumed by Andrews [2005].
$V_\textrm{peak}$$\Delta \sigma$$\beta$$\rho$$\mu \, (= \rho \beta^2)$
5.7 m/s10 MPa3,000 m/s2,700 kg/m$^3$24.3 GPa
Estimation of $\mathcal{G}/x$ based on seismological observables
Fig. 6. Distribution of $D_c'/D_c$, showing that the proxy of characteristic slip distance, $D_c'$, works well.
We confirmed that following approximations proposed by some prior studies hold even in our model:
  • $\mathcal{G} \sim \frac12 \sigma_b D_c$ (Fig. 4c).
  • $D_c' \sim 0.6 D_c$, where $D_c'$ is slip distance in which $V_\textrm{peak}$ appears [Mikumo+ 2003] (Fig. 6).
  • $\dfrac{\mu V_\textrm{peak}}{\sigma_b \beta} \sim \dfrac{k}{\sqrt{1-k^2}}$ (Fig. 7a), based on transrational semi-infinite crack model [Ida 1972] and numerical solution for spontaneous rupture propagation with velocity-weakning friction [Gabriel+ 2013].
Combination of the above approximations yields the following: \begin{align} \mathcal{G} \sim \frac{\mu \sqrt{1-k^2} V_\textrm{peak} D_c'}{c}, \label{ApproxG} \end{align} where the right-hand side consists only of seismological observables (rupture velocity, peak slip rate, and $D_c'$) and may contribute to the seismological estimation of energy release rate. The ratio of both side in the approximation \eqref{ApproxG} is ranging from $0.5$ to $1.2$ for various values of $S$ (Fig. 7b), which is sufficient considering the accuracy of seismological observations.
In the case of the steady-state pulse-like rupture model, Rice+ [2005] provided an approximation of fracture energy based on some seismological observables. Now we have a different one for the self-similar crack model, a counterpart of the pulse-like model.
Fig. 7. Accuracy of approximations for (a) $V_\textrm{peak}$ and (b) $\mathcal{G}/x$. The grey lines are for every 0.25 interval of $S$-value.

5. Conclusion


This work was supported by JSPS and RFBR under the Japan–Russia Research Cooperative Program (project no. J19–721), JSPS KAKENHI (grant no. 17H02857), and the Research Institute for Mathematical Sciences, an International Joint Usage/Research Center located in Kyoto University. We are grateful for insightful and constructive comments by Eiichi Fukuyama, Pierre Romanet, and Steven Day throghout the review process.