Geophys. J. Int.
(2012)
189,
1797–1806
doi: 10.1111/j.1365-246X.2012.05464.x
GJI Seismology
A new paradigm for simulating pulse-like ruptures: the pulse
energy equation
Ahmed E. Elbanna
1
and Thomas H. Heaton
2
1
Department of Physics, University of California at Santa Barbara, Broida Hall, Santa Barbara, CA
93106
, USA. E-mail: ahmedettaf@gmail.com
2
Seismological laboratory, California Institute of Technology,
1200
E. California Blvd., MS
252-21
, Pasadena, CA
91125-2100
,USA
Accepted 2012 March 19. Received 2012 March 14; in original form 2011 September 14
SUMMARY
We investigate the chaotic behaviour of slip pulses that propagate in a spring block slider
model with velocity weakening friction by numerically solving a computationally intensive
set of
n
coupled non-linear equations, where
n
is the number of blocks. We observe that the
system evolves into a spatially heterogeneous pre-stress after the occurrence of a sufficient
number of events. We observe that, although the spatiotemporal evolution of the amplitude of
a slip pulse in a single event is surprisingly complex, the geometric description of the pulses
is simple and self-similar with respect to the size of the pulse. This observation allows us to
write an energy balance equation that describes the evolution of the pulse as it propagates
through the known pre-stress. The equation predicts the evolution of individual ruptures and
reduces the computational time dramatically. The long-time solution of the equation reveals
its multiscale nature and its potential to match many of the long-time statistics of the original
system, but with a much shorter computational time.
Key words:
Fractals and multi-fractals; Non-linear differential equations; Self-organization;
Friction; Earthquake dynamics.
1 INTRODUCTION
Seismic inversions have indicated that most earthquake ruptures propagate in a pulse-like mode (Heaton 1990); each point on the fault slips
for only a short time compared to the total event duration. This is to be contrasted with the crack-like mode of rupture, in which the slipping
region continues to expand until the rupture reaches the outermost boundaries of the fault. Once the rupture surface stops growing, arrest
waves propagate from the border of the rupture and cause the cessation of rupture at the interior points. In the crack-like mode, each point on
the fault continues to slip for a duration that is comparable to the overall event time. A possible mechanism for the fast healing of slip behind
the rupture front in the pulse-like mode of rupture is the strong dependence of the fault friction on the slip rate. This has been confirmed in
many experimental and numerical studies (Heaton 1990; Perrin
et al.
1995; Zheng & Rice 1998; Lapusta 2001; Lu
et al.
2007; Ampuero &
Ben-Zion 2008; Beeler
et al.
2008; Noda
et al.
2009).
A major problem in the field of computational earthquake mechanics is that the available computational resources are still insufficient
for simulating dynamic failure of a 3-D continuum with strong velocity weakening friction (Aagaard & Heaton 2008). With this limitation,
there is an increased interest in developing physically based approximate models that can reliably simulate the main rupture variables, such
as the final slip, within a reasonable computational time. For systems with pulse-like ruptures, we ask the following question: is it possible to
exploit the localized nature of the slip pulse in designing computational methods that can solve for some of the macroscopic rupture variables,
such as the final slip, without conducting the full dynamic simulations? An affirmative answer to this question is useful in many ways. For
instance, since many of the rupture variables are thought to be correlated with each other, for example, the final slip, the duration of slip at
a point (sometimes called the rise time), the pulse maximum slip rate,
...
, etc., knowledge of the final slip would facilitate estimating these
other variables. Another useful application would be in exploring the long-time evolution of the pre-stress and event size distributions in
different fault models.
Since we are unable to numerically simulate the long-time behaviour of a continuum that slips on an interface with strong velocity
weakening friction, we explore the much simpler problem a simple 1-D spring block slider model. The model has the benefit that it can
be tuned to produce slip pulse failures with dynamic complexity (Fig. 1) and its long-time behaviour can be explored through numerical
simulation. The spring block slider model is a discrete model where all direct elastic interactions are short range and there is no radiated
wavefield (Burridge & Knopoff 1967; Carlson & Langer 1989; Erickson
et al.
2008; Elbanna 2011). When the model fails in slip pulses, the
pulses carry energy and momentum from one part of the model to another. That is, slip pulses are an additional form of interaction.
C
2012 The Authors
1797
Geophysical Journal International
C
2012 RAS
Geophysical Jour
nal Inter
national
1798
A.E. Elbanna and T.H. Heaton
Figure 1.
The spring block slider model. A chain of
N
blocks of identical masses (
m
) interconnected by coil springs of stiffness (
k
c
) is driven by a loading plate
moving at a very low velocity (
v
). The blocks feel the effect of the loading plate through a series of leaf springs of stiffness (
k
l
). Any block is stuck to the ground
as long as the total elastic force acting on it is less that the static friction threshold. Once the static friction is exceeded, a sliding block experie
nces a dynamic
friction force
that varies inversely with the block sliding velocity. The blocks are subjected to an initial stress distribution that is constant everywhere but is
larger than the static friction threshold for the middle three blocks to initiate rupture. As a block moves, it transfers stresses into neighbouring b
locks causing
them to subsequently move if the static friction level is exceeded. When all the blocks that have already slipped stop, we designate this as the complet
ion of
one event.
Figure 2.
Contours of slip rate in one of the events generated by the spring block model showing a nearly constant rupture speed and a pulse-like character.
The changing colour of the slip pulse indicates the average particle velocity in the pulse, which correlates with the slip left behind by the pulse.
It is known that the spring block model with pure rate-dependent friction cannot be a complete frictional description in a model of a
smooth continuum fault (Rice 1983, 1993, 2006; Cochard & Madariaga 1996; Ben-Zion & Rice 1997; Rice
et al.
2001). However, the discrete
formulation of our problem has several advantages; we can explore the long-time evolution of systems with multiscale ruptures and we can
also investigate the effect of the strong non-linear feedback between the friction and slip rate. In many systems, from turbulent fluids to stock
markets, the existence of positive feedback leads to instability, localization and chaos. The spring block slider model presents a handy tool to
explore the implications of this feedback and to create new techniques that can be applied to the more complicated continuum problems.
The dynamics of the spring block slider models have been extensively studied (Carlson & Langer 1989; Carlson
et al.
1991; Shaw
et al.
1992; Elbanna 2011). The most important results relevant to us here are: (i) events propagate in a pulse-like mode with a nearly constant
rupture speed (Fig. 2); and (ii) pulses are approximately self-similar in shape, with their amplitude and width changing as they propagate but
their Gaussian—or nearly triangular-shape is generally preserved (Fig. 3).
The pulse width is found from the ensemble of simulations to be typically 1–5 per cent of the total event length. This extreme localization
of the pulse, as well as its approximate self-similar nature, suggests that the system’s dynamics are governed by the interplay between only a
few variables that characterize the slip pulse and the pre-event stress distribution (also known as the pre-stress).
2 ENERGY BALANCE FOR PROPAGATING SLIP PULSES IN THE SPRING
BLOCK MODEL
Energy conservation has to hold for the spring block model. The relevant forms of energy that contribute to the energy balance are: the
potential energy of the linear coil and leaf springs (these are derivable from the pre-stress), the frictional dissipation and the pulse kinetic
energy.
As the pulse propagates, we can identify two regions in the system: the region that has slipped and is now stuck, and the region that is
currently slipping (the slip pulse). This is shown schematically in Fig. 4. As the pulse propagates, it changes the potential energy of the system
by deforming the coil and leaf springs. Part of this change in potential energy is lost as frictional heat and part goes into the kinetic energy of
the slip pulse. In our model, the kinetic energy of the pulse is readily calculated as the sum of the kinetic energy of the blocks that are moving
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
The pulse energy equation
1799
Figure 3.
Pulse profiles at different instants of time for the event shown in Fig. 3 showing growth and decay in a nearly self-similar pattern.
Figure 4.
A schematic diagram showing the propagation of a slip pulse. As the pulse propagates, the portion of the fault behind the healing tip is no longer
slipping; the change in potential energy of the springs in this region can be fully expressed in terms of the final slip distribution and the pre-stress.
At any
given instant, the total energy of the slip pulse is calculated as the sum of the change in potential energy of the springs in slip pulse plus the sum of the
kinetic
energies of the blocks moving within the pulse.
within the pulse. (The kinetic energy of the block is half the block mass times its velocity squared.) In other words, for all instants in time
Total dissipated frictional work
=
Change in potential energy of the system
−
Kinetic energy of the pulse
.
(1)
By a slight rearrangement of terms, the above statement could be written as
Total energy of the pulse
=
Change in potential energy of the already ruptured portion
−
Total frictional work
.
(2)
The total pulse energy is the sum of the pulse kinetic energy previously defined and the change in potential energy of the springs, leaf and
coil, in the region occupied by the pulse.
Eq. (2) has the advantage that behind the pulse healing tip, the system has reached its final state and consequently the change in potential
energy in that portion can be exactly described in terms of the final slip distribution and the pre-stress. All the time-dependent variables are
now grouped in the pulse and frictional energy terms. Interestingly, those quantities correlate with the amplitude of the slip pulse. Except for
very small slip values, the total pulse energy can be approximated by
E
p
=
γ
D
ξ
where
D
is the slip. This is verified in Fig. 5 for one of
the spring block slider events. Using a log–log plot, the parameters
ξ
and
γ
are readily obtained:
ξ
is the slope of the straight segment and
log
γ
is the value of its intercept at zero abscissa. A similar approach can be used for the frictional work term. The total frictional dissipation
at a point is given by
W
f
=
̄
τ
f
D
where ̄
τ
f
is the slip-averaged frictional stress at that point (see Appendix). For the same event previously
considered, Fig. 6 shows that ̄
τ
f
can be approximated as ̄
τ
f
=
α
D
β
,
where
β
and log
α
are the slope and the intercept of the fitting straight
line, respectively. These simple approximations leave us then with just two variables to work with: the pre-stress and the final slip.
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
1800
A.E. Elbanna and T.H. Heaton
Figure 5.
Dependence of the total pulse energy on pulse slip. Except for very small slips, the pulse energy depends in a power-law fashion on the pulse slip
and the power-law form is depicted in the rectangle.
Figure 6.
Dependence of the slip averaged frictional force per block on the block slip.
In the Appendix, we derive the equation that relates the final slip distribution to the pre-event stress state. Here, we just mention the
equation in its final form
d
D
d
x
=
−
(
γξ
D
ξ
−
1
−
σ
oc
)
+
√
(
γξ
D
ξ
−
1
−
σ
oc
)
2
−
2
k
c
(
1
2
k
l
D
2
−
σ
ol
D
+
α
D
β
+
1
)
k
c
,
(3)
where
k
c
is the stiffness of the coil spring,
k
l
is the stiffness of the leaf spring,
σ
oc
is the pre-stress in the coil spring and
σ
ol
is the pre-stress
in the leaf spring. The scaling coefficients and exponents
α, β, γ , ξ
are empirically determined from the numerical simulations, as discussed
earlier, and are generally functions of the system parameters (see Appendix).
In solving eq. (3), we use the same values of springs stiffness as those used in the spring block model and we discretize it with
x
=
1
to match the model discretization. We then get a discrete map in the form
D
i
+
1
=
D
i
+
g
(
D
i
,σ
oc
i
,σ
ol
i
)
.
(4)
The discrete map (4) estimates the slip at the point (
i
+
1) given the stress and slip at point (
i
). Since the equation is first order, we need to
know the slip value at one point to get it started. For the test purposes, we choose this point to be at the location when the nucleation phase is
over and a single pulse is fully developed. The results of the equation simulation, along with the slip solution generated by the spring block
model, are shown in Fig. 7.
The equation performed well in the pulse-like region, matching the slip profile both qualitatively and quantitatively. The equation also
took just 1/10 000 the time required to simulate this event using the full spring block model. Two problems remain, however, if we are planning
to extend this approach to simulate sequence of events systematically. First, how are we going to model the nucleation phase where the rupture
is crack-like and our equation does not apply? And second, how are we going to handle bilateral events?
We decided that it would be best to simulate the nucleation phase accurately using the spring block model and then once the pulse is fully
developed, we can use the equation to predict the slip outside the nucleation region. In this way, we would be able to get the correct slip initial
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
The pulse energy equation
1801
Figure 7.
(Colour Online): Testing the pulse energy equation. (Top panel): Pre-stress distribution for the test case. (Bottom panel): The solution of the prop
osed
differential equation (in blue) mimics the slip solution from the spring block model in the pulse-like region (between blocks 200 and 420). The initia
l condition
for the equation was taken as the slip value at the point where the nucleation phase is over and the pulse has fully developed. The nucleation phase is han
dled
separately as discussed in the text.
Figure 8.
Flow chart for simulating a sequence of events using the hybrid model.
conditions for the application of the pulse energy equation as well as to simulate bilateral events. Our model in this case becomes hybrid. The
model is hybrid in the sense that a microscopic model, the full spring block model, is used to simulate nucleation whereas a macroscopic
model, the energy pulse equation, is used to simulate further propagation.
With the hybrid model, the evolution of the system can be tracked systematically as outlined in Fig. 8. Although the full spring block
model is used in simulating the nucleation phase, the computational time is two orders of magnitude smaller than the time spent in solving
the full set of the equations of motion for the whole rupture.
We test the long-time performance of pulse energy equation by generating 1000 events from the hybrid model and the full spring block
model. Although the two solution techniques produce solutions that evolve differently, Fig. 9 shows the general statistical properties of the
heterogeneous pre-stress for the 1001th event are similar.
Fig. 10 shows the corresponding results for the events length distribution. Although there are quantitative differences, the distribution
generated by the hybrid model is qualitatively similar to that generated by the spring block model; in particular, there is a hump in the
distribution corresponding to the increased probability of generating a large event for this choice of system parameters.
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
1802
A.E. Elbanna and T.H. Heaton
Figure 9.
(Colour Online): The statistical properties of the solution generated by the hybrid model compared to that generated by the full Spring block model.
(a) The amplitude of the Fourier transform of the pre-stress for the 1001th event generated by the hybrid model (blue) and the full spring block model (r
ed)
showing similar spectral slopes. (b) The probability distribution of the pre-stress for the 1001th event for both the hybrid model (red) and the sprin
g block
model (blue) showing similar trends. The probability distribution is generated by binning the stresses into bins of equal sizes (10 MPa) and count the
number
of blocks belonging to each bin.
Figure 10.
Event length distribution generated by both the hybrid model (dots) and the full spring block model (blue).
3 DISCUSSION
Designing computationally efficient physically based earthquake simulators is a central problem in seismology and earthquake physics.
Although a detailed dynamical modelling is essential for understanding the mechanics of earthquake ruptures on the level of individual
events, there are many applications in which we are interested in the statistics of the system more than in the details of single ruptures. Such
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
The pulse energy equation
1803
applications include, for example, seismic hazard analysis, statistics of seismicity and ground motion prediction. For such problems, we might
not be interested in the detailed dynamic evolution of the rupture but a physically based kinematic description may be sufficient.
In this manuscript, we have introduced a new paradigm for simulating earthquake slip for pulse-like ruptures that is based on the energy
budget of the propagating slip pulse. The proposed method is inspired by the localization of the slip pulse and its approximate self-similarity.
At any instant of time during a pulse-like rupture, the slip pulse only occupies a very limited portion of the fault. Moreover, the pulse width
and amplitude may grow or shrink during propagation depending on the conditions, but the pulse retains approximately its shape. These two
observations, localizations and the approximate self-similarity, suggested that only a small number of degrees of freedom may be used to
describe the evolution of the pulse.
We have found that a key to formulate a coarse-grained model for pulse-like ruptures is energy balance. In the spring block model studied
here, where all the dynamic changes are localized to the pulse zone, the energy balance in this case takes the simple form: change in potential
energy of the system at any given instant – cumulative frictional work dissipated by the pulse up to this instant
=
The current kinetic energy
of the pulse. By writing down the exact change in potential energy in terms of the final slip and pre-stress, and by approximating the pulse
energy and the frictional work as power laws in the pulse slip with power-law parameters empirically derived from the numerical simulation
of the spring block model, we were able to formulate a differential equation that relates the final slip with the pre-event stress. The solution
of the equation revealed many interesting features, but showed some limitations as well.
The equation performed well on the level of individual events. It predicted the basic features of the final slip distribution including slip
amplitudes and total rupture length. The saving in computational time was at least two orders of magnitude. However, with the equation
designed for pulse-like ruptures, it could not be used in simulating the nucleation phase of the rupture which is usually crack-like. To overcome
this problem, we used the spring block model to accurately simulate the nucleation phase and once the rupture assumed the pulse-like shape
we switched to the equation for the prediction of the final slip in the rest of the rupture. By using this hybrid approach (i) we accurately
computed the nucleation, (ii) we got the initial slip value necessary for starting the equation right and (iii) we could simulate bilateral ruptures
as two independent ruptures propagating opposite to each other. We were also able to use this hybrid model to simulate sequence of ruptures
in a systematic way.
The hybrid model performed satisfactorily when used to simulate the long-time statistics. By performing a limited search in the scaling
parameters space, we were able to find a combination of scaling parameters that led to long-time pre-stress statistics and large events size
distribution similar to the corresponding statistics from the spring block model. The saving in computational time was almost two orders of
magnitude. The success of our approach suggests that knowing the details of the dynamic propagation, manifested for example in the exact
values of the pulses energy or details of friction, might not be very important for approximating the long-time statistics. This feature has
analogies in many complex systems where the macroscopic long-term behaviour might not be very sensitive to the underlying microscopic
dynamics.
The proposed pulse energy equation is based on principles that are shared with more complicated continuum models, such as pulse
localization, energy balance and approximate self-similarity of the pulse. Extension of this approach to the continuum has to account for the
existence of long-range interactions carried by the radiated wavefield. In the spring block model, long-range wavefield is absent since all
interactions are short-ranged or done through the pulse propagation. We acknowledge that accounting for radiated energy will not be an easy
task. Nonetheless, we hypothesize that the more localized the pulse is, the more sensitive it gets to the local pre-stress conditions than to the
fluctuations carried by the radiated wavefield. Accordingly, we think it is still possible to extend the current approach to the continuum by
adding a term in the energy balance equation that accounts for the radiated energy flux.
Finally, it is possible to extend this approach to other types of friction law as long as localized slip pulses are generated. Possible
candidates are rate and state friction laws and slip weakening laws (if slip pulses are formed via the propagation of arrest waves). The energy
balance is still valid in these cases. The pulse energy and the frictional work scaling laws will probably take a different form, however, and
may not be simple power laws.
4 CONCLUSIONS
In this manuscript, we have presented a new paradigm for simulating slip in pulse-like ruptures using an energy balance approach. The key
ingredients for this paradigm are the spatial localization of the slip pulse and its approximate self-similarity during propagation. The energy
budget of the pulse as it propagates along the fault dictates that at any instant of time, the change in potential energy of the system minus
the cumulative frictional dissipation should be balanced by the kinetic energy. Using this energy balance approach, and by approximating
the kinetic energy and the frictional work by power laws in terms of the pulse slip, we were able to formulate, within a spring block model
framework, a differential equation that can be used to predict the final slip in an event given the pre-stress distribution. The equation replicated
the main features of slip in the pulse-like region of ruptures, including slip amplitudes and total rupture length, to a very good degree and
in a much shorter computational time. The non-pulse-like portions of the rupture, for example, nucleation regions, are handled through the
solution of the fully dynamic spring block model. By combining the spring block model for nucleation and the equation for further propagation
(the hybrid model), we were able to match some of the long-time statistics of the solution of the spring block model, such as the pre-stress
spectrum and large events size distribution, with two orders of magnitude savings in computational time. A similar approach might be applied
to continuum models, but in this case, long-range interactions due to the radiated wavefield have to be included.
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
1804
A.E. Elbanna and T.H. Heaton
ACKNOWLEDGMENTS
We thank Nadia Lapusta, Jean-Paul Ampuero and Deborah Smith for their comments and we acknowledge the collaboration with Jerrold
Marsden. This research was supported by the National Science Foundation, the Southern California Research Center, and the Gordon and
Betty Moore Foundation (the Caltech Tectonics Observatory).
REFERENCES
Aagaard, B.T. & Heaton, T.H., 2008. Constraining fault constitutive
behavior with slip heterogeneity,
J. geophys. Res.,
113,
B04301,
doi:10.1029/2006JB004793.
Ampuero, J.-P. & Ben-Zion, Y., 2008. Cracks, pulses and macro-
scopic asymmetry of dynamic rupture on a bimaterial interface
with velocity-weakening friction,
Geophys. J. Int.,
173
(2), 674–692,
doi:10.1111/j.1365-246X.2008.03736.x.
Beeler, N.M., Tullis, T.E. & Goldsby, D.L., 2008. Constitutive relationships
and physical basis of fault strength due to flash heating,
J. geophys. Res.,
113,
doi:10.1029/2007JB004988.
Ben-Zion, Y. & Rice, J.R., 1997. Dynamic simulations of slip on a smooth
fault in an elastic solid,
J. geophys. Res.,
102,
17 771–17 784.
Burridge, R. & Knopoff, L., 1967. Model and theoretical seismicity,
Bull.
seism. Soc. Am.,
57,
341–371.
Carlson, J.M. & Langer, J.S., 1989. Mechanical model of an earthquake
fault,
Phys. Rev. A,
40,
6470–6484.
Carlson, J.M., Langer, J.S., Shaw, B.E. & Tang, C., 1991. Intrinsic properties
of a Burridge-Knopoff model of an earthquake fault,
Phys. Rev. A,
44,
884–897.
Cochard, A. & Madariaga, R., 1996. Complexity of seismicity due to highly
rate dependent friction,
J. geophys. Res.,
101,
25 321–25 336.
Erickson, B., Birnir, B. & Lavallee, D., 2008. A model for aperiodicity in
earthquakes,
Nonlinear Process. Geophys,
15,
1–12.
Elbanna, A., 2011. Pulselike ruptures on strong velocity-weakening fric-
tional interfaces: dynamics and implications,
PhD thesis,
Caltech,
Pasadena, CA.
Heaton, T.H., 1990. Evidence for and implications of self-healing
pulses of slip in earthquake rupture,
Phys. Earth planet. Inter.,
64,
1–20.
Lapusta, N., 2001. Elastodynamic analysis of sliding with rate and state
friction,
PhD thesis,
Harvard University, Cambridge, MA.
Lu, X., Lapusta, N. & Rosakis, A.J., 2007. Pulse-like and crack-like rup-
tures in experiments mimicking crustal earthquakes,
Proc. Natl. Acad.
Sci. U.S.A.,
doi:p.10.1073/pnas.0704268104.
Noda, H., Dunham, E.M. & Rice, J.R., 2009. Earthquake ruptures with ther-
mal weakening and the operation of major faults at low overall stress
levels,
J. geophys. Res.,
114,
B07302, doi:10.1029/2008JB006143.
Perrin, G., Rice, J.R. & Zheng, G., 1995. Self-healing slip pulse on a fric-
tional surface,
J. Mech. Phys. Solids,
43,
1461–1495.
Rice, J.R., 1983. Constitutive relations for fault slip and earthquake insta-
bilities,
Pure appl. Geophys.,
121,
443–475.
Rice, J.R., 1993. Spatio-temporal complexity of slip on a fault,
J. geophys.
Res.,
98
(B6), 9885–9907.
Rice, J.R., 2006. Heating and weakening of faults during earthquake slip,
J. geophys. Res.,
111,
doi:10.1029/2005JB004006.
Rice, J.R., Lapusta, N. & Ranjith, K., 2001. Rate and state dependent fric-
tion and the stability of sliding between elastically deformable solids,
J. Mech. Phys. Solids,
49,
1865–1898.
Shaw, B.E., Carlson, J.M. & Langer, J.S., 1992. Patterns of seis-
mic activity preceding large earthquakes,
J. geophys. Res.,
97,
479–
488.
Zheng, G. & Rice, J.R., 1998. Conditions under which velocity-weakening
friction allows a self-healing versus a cracklike mode of rupture,
Bull.
seism. Soc. Am.,
88,
1466–1483.
APPENDIX
In this appendix, we construct the differential equation governing the energy transportation in the spring block slider model via the slip pulse.
Our starting point will be eq. (2) where we stated the energy balance principle for the pulse. As we discussed in the text, the formulation in
eq. (2) allows us to express the change in potential energy entirely in terms of the final slip and the pre-stress in the stuck region. This is done
as follows.
Change in potential energy of the system in the healed region
=
change in potential energy of the leaf springs (
PE
leaf
)
+
Change in
potential energy of the coil spring (
PE
coil
)
.
(A1)
The change in potential energy in the leaf springs per unit length is given by
PE
leaf
=−
k
l
u
o
D
−
1
2
k
l
D
2
=
σ
ol
D
−
1
2
k
l
D
2
,
(A2)
where
D
is the final slip and
σ
ol
is the pre-stress in the leaf springs.
The change in potential energy in the coil springs per unit length is similarly given by
PE
coil
=−
k
c
∂
u
o
∂
x
∂
D
∂
x
−
1
2
k
c
(
∂
D
∂
x
)
2
=
σ
oc
∂
D
∂
x
−
1
2
k
c
(
∂
D
∂
x
)
2
,
(A3)
where
σ
oc
is the pre-stress in the coil springs and
x
is the spatial coordinate. These expressions of the change in potential energy for the leaf
and coil springs are exact in the continuum limit of a 1-D elastic system.
On the other hand, the total energy of the pulse can be written in principle in the following functional form:
E
p
=
G
(
u
,
̇
u
,∂
u
/
∂
x
)
,
(A4)
where
u
and
̇
u
are the current displacement and slip rates of the blocks within the pulse region. We can write down the exact form of
G
;itis
the sum of kinetic energies of the blocks forming the slip pulse plus the sum of the change in potential energy of the leaf and coil springs
within the pulse region. However, this will not be of much use as it adds more unknowns to the energy balance in eq. (2) without bringing an
equivalent number of equations. Rather, we will pursue a different route. We will seek an approximation for
G
in terms of the final slip value
at the healing tip of the pulse (
D
). One expects that the bigger the pulse gets the more energy it has and the more slip it carries and hence it is
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
The pulse energy equation
1805
reasonable to assume some sort of correlation between the pulse energy and pulse slip. Accordingly, we may write
E
p
=
E
p
(
D
)
.
(A5)
Unfortunately, we cannot write down a closed-form expression for
E
p
in terms of
D
as this will depend on the detailed dynamics of
slip pulse propagation which in turn will depend on the system parameters (e.g. springs stiffness and rate of frictional weakening) and the
statistical properties of the pre-stress. Nonetheless, we can derive an empirical relation from the results of the numerical simulations of the
spring block model. An example of such relation was shown in Fig. 5. Except for very small slips (
D
<
0.03), the total pulse energy varies
with the pulse slip, as a power law given by
E
p
=
10
2
.
9
D
0
.
4
.
(A6)
Having approximated the total pulse energy, the next step is to find an analogous approximation for the frictional work term. The cumulative
frictional work from the initiation of rupture and up to the pulse current position is given by
W
F
=
∫
x
0
w
f
d
x
=
∫
x
0
d
x
∫
t
0
̇
uf
(
u
,
̇
u
)d
t
,
(A7)
where
w
f
is the frictional work per unit length and
f
(
u
,
̇
u
) is the friction functional form. The integral
∫
t
0
̇
uf
(
u
,
̇
u
)d
t
can be evaluated if the
time histories of
u
and
̇
u
are given. Since this requires a detailed dynamic modelling, evaluating the correct value of the integral will defy our
model reduction purposes. Instead, similar to what we have done in the total pulse energy case, we will seek an approximation. Note that the
frictional work density, or equivalently the time-dependent integral, can be expressed as
w
f
=
̄
τ
D
,
(A8)
where ̄
τ
is the slip averaged frictional work per unit slip ̄
τ
=
∫
t
0
̇
uf
(
u
,
̇
u
)d
t
/.
D
. As a first-order approximation for ̄
τ
, we may write
̄
τ
=
̄
τ
(
D
)
.
(A9)
The relation between ̄
τ
and
D
was shown in Fig. 6. It is well approximated by a power law with some scatter. The power law is given by
̄
τ
=
10
−
0
.
21
D
−
0
.
53
.
(A10)
Now we are in a position to formulate the sought energy balance equation in a mathematical form. From eq. (2) in the text, we may write
U
−
W
F
=
E
p
,
(A11)
where
U
the total change in potential energy in the healed region is,
W
F
is the cumulative frictional work in the healed region and
E
p
is the
total pulse energy.
However,
U
=
∫
x
0
u
d
x
=
∫
x
0
[
PE
coil
+
PE
leaf
]
d
x
.
(A12)
where
u
is the change in total potential energy per unit length in the healed region,
PE
coil
and
PE
leaf
have been previously defined.
W
F
=
∫
x
0
w
f
d
x
.
(A13)
Hence, one can write the energy balance equation in the following form:
∫
x
0
u
−
w
f
d
x
=
E
p
.
(A14)
Taking derivatives of both sides, the fundamental theorem of calculus then yields
u
−
w
f
=
d
E
p
d
x
.
(A15)
Substituting from (A2), (A3), (A6) and (A10) in (A15), we arrive at the following first-order second-degree differential equation:
−
1
2
k
c
(
∂
D
∂
x
)
2
+
σ
oc
∂
D
∂
x
−
1
2
k
l
D
2
+
σ
ol
D
−
α
D
β
D
=
d
d
x
(
γ
D
ξ
)
.
(A16)
In eq. (A16), we have used the power-law approximations for the total pulse energy and frictional work per unit length:
E
p
=
γ
D
ξ
and
̄
τ
=
α
D
β
. Simplifying and rearranging, we get
1
2
k
c
(
d
D
d
x
)
2
+
(
γξ
D
ξ
−
1
−
σ
oc
)
d
D
d
x
+
(
1
2
k
l
D
2
−
σ
ol
D
+
α
D
β
+
1
)
=
0
.
(A17)
One can readily recognize that the third term in eq. (A17) yields the condition for a pulse propagating with uniform slip. For the case of
uniform slip, d
D
/.
d
x
vanishes and we are left with
1
2
k
l
D
2
−
σ
ol
D
+
α
D
β
+
1
=
0, which states that the change in potential energy
1
2
k
l
D
2
−
σ
ol
D
balances the frictional dissipation
α
D
β
+
1
pointwise. This is the condition for steady propagation.
Factoring eq. (A17) gives
d
D
d
x
=
−
(
γξ
D
ξ
−
1
−
σ
oc
)
±
√
(
γξ
D
ξ
−
1
−
σ
oc
)
2
−
2
k
c
(
1
2
k
l
D
2
−
σ
ol
D
+
α
D
β
+
1
)
k
c
.
(A18)
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS
1806
A.E. Elbanna and T.H. Heaton
To determine which sign of the radical is to be taken, we once again consider the limiting case of steady propagation. In this limit d
D
/.
d
x
is zero and as we just showed,
1
2
k
l
D
2
−
σ
ol
D
+
α
D
β
+
1
vanishes. This is only possible if we take the positive sign. Hence,
d
D
d
x
=
−
(
γξ
D
ξ
−
1
−
σ
oc
)
+
√
(
γξ
D
ξ
−
1
−
σ
oc
)
2
−
2
k
c
(
1
2
k
l
D
2
−
σ
ol
D
+
α
D
β
+
1
)
k
c
.
(A19)
This is the final form of the pulse energy balance equation. The derivation of the equation is exact in the continuum limit of 1-D elastic
system. The source of approximation lies in the functional form of the dependence of the total pulse energy and frictional work on the pulse
slip. In the system under consideration, this dependence is found to be a power law as a direct consequence of the approximate self-similarity of
the slip pulse. The scaling coefficients and exponents
α, β, γ , ξ
are empirically determined from the numerical simulations, and are generally
functions of the system parameters as previously discussed.
C
2012 The Authors,
GJI
,
189,
1797–1806
Geophysical Journal International
C
2012 RAS